Analysis of Weighted Multifrequency MUSIC-Type Algorithm for Imaging of Arc-Like, Perfectly Conducting Cracks

Abstract

The main purpose of this paper is to investigate the structure of the weighted multifrequency multiple signal classification (MUSIC) type imaging function in order to improve the traditional MUSIC-type imaging. For this purpose, we devise a weighted multifrequency MUSIC-type imaging function and examine a relationship between weighted multifrequency MUSIC-type function and Bessel functions of integer order of the first kind. Some numerical results are demonstrated to support the survey.

1. Introduction

Inverse problem, which deals with the reconstruction of cracks or thin inclusions in homogeneous material (or space) with physical features different from space, is of interest in a wide range of fields such as physics, engineering, and image medical science which are closely related to human life; refer to [1–9]. That is why inverse problem has been established as one interesting research field. Compared to the early studies on inverse problem in which much research had been done theoretically, in recent studies, more practical and applicable approaches have been undertaken and the reconstructive way appropriate to each specific study field started to be investigated thanks to the development of computational science using not only computers but also mathematical theory. As we can see through a series of papers [10, 11], the reconstruction algorithm, based on the iterative scheme such as Newton’s method, has been mainly studied. Generally, in regard to algorithms using Newton’s method, in the case of the initial shape quite different from the unknown target, the reconstruction of material leads to failure with the nonconvergence or yields faulty shapes even after the iterative methods are conducted. Hence, in such an iterative method, several noniterative algorithms have been proposed as a way to find the shape of initial value close to that of the unknown target as quickly as possible.

The noniterative algorithms such as multiple signal classification (MUSIC), subspace migrations, topological derivative, and linear sampling method can contribute to yielding the appropriate image as an initial guess. Previous attempts to investigate MUSIC-type algorithm presented various experiments with the use of MUSIC-type algorithm. For instance, the use of MUSIC-type algorithm for eddy-current nondestructive evaluation of three-dimensional defects [12], and MUSIC-type algorithm designed for extended target, the boundary curves which have a five-leaf shape or big circle was presented [13]. In addition, MUSIC-type algorithm was introduced for locating small inclusions buried in a half space [2] and for detecting internal corrosion located in pipes [3]. Although the past phenomena about experimental results could not be theoretically explained because the mathematical structure about these algorithms was not verified, recent studies [14–19] managed to partially analyze the structure of some algorithms. On the basis of these studies, the present study examined the structure of algorithms to make improvements in imaging the defects. Therefore, this paper aims to improve traditional MUSIC-type imaging algorithm by weighting applied to each frequencies.

This paper is organized as follows. In Section 2, we discuss two-dimensional direct scattering problem in the presence of perfectly conducting crack and MUSIC-type algorithm. In Section 3, we introduce a weighted multifrequency MUSIC-type imaging algorithm and analyze its structure to confirm that it is an improved version of traditional MUSIC algorithm. In Section 4, we present several numerical experiments with noisy data. In Section 5, our conclusions are briefly presented.

2. Direct Scattering Problem and Single- and Multifrequency MUSIC-Type Algorithm

In this section, we simplify surveying the two-dimensional direct scattering problem for the existence of perfectly conducting cracks and the single- and multifrequency MUSIC algorithm. For more information, see [10, 20].

2.1. Direct Scattering Problem and MUSIC-Type Imaging Function

First, we consider the two-dimensional electromagnetic scattering by a perfectly conducting crack located in the homogeneous space ℝ². Throughout this paper, we assume that the crack Ω is a smooth, nonintersecting curve, and we represent Ω such that Ω={𝜉(𝑡)∶𝑡∈[−1,1]}, (1) where 𝜉∶[−1,1]→ℝ² is an injective piecewise smooth function. We consider only the transverse magnetic (TM) polarization case. Let us denote 𝑢_total to be the time-harmonic total field, which can be decomposed as 𝑢_total(𝐱)=𝑢_incident(𝐱)+𝑢_scatter(𝐱), (2) where 𝑢_incident(𝐱) = exp(𝑗𝜔𝜽⋅𝐱) is the given incident field with incident direction 𝜽∈𝕊¹ (unit circle) and 𝑢_scatter(𝐱) is the unknown scattered field that satisfies the Sommerfeld radiation condition lim_|𝐱|→∞|𝐱|𝜕𝑢_scatter(𝐱)𝜕|𝐱|−𝑗𝜔𝑢_scatter(𝐱)=0 (3) uniformly in all directions 𝐱=𝐱/|𝐱|. Now, the total field 𝑢_total satisfies the two-dimensional Helmholtz equation Δ𝑢_total(𝐱)+𝜔²𝑢_total(𝐱)=0inℝ²𝑢⧵Ω,_total(𝐱)=0onΩ (4) with a given positive frequency 𝜔. In the case that Ω is absent, incident field 𝑢_incident can also be a solution of (4).

The far-field pattern is defined as function 𝑢_far(𝐱,𝜽) that satisfies 𝑢_scatter(𝐱,𝜔)=exp𝑗𝑘₀|𝐱||𝐱|𝑢_far()1𝐱,𝜽+𝑜|𝐱|=exp𝑗𝜔|𝐱||𝐱|𝑢_far()1𝐱,𝜽+𝑜|𝐱| (5) as |𝐱|→∞ uniformly on 𝐱. Then, based on [21], the far-field pattern 𝑢_far(𝐱;𝜽) of the scattered field 𝑢_scatter(𝐱) can be expressed by the following equation: 𝑢_far(𝐱;𝜽)=−exp𝑗𝜋/4√8𝜋𝜔^Ωexp−𝑗𝜔𝐱⋅𝐲𝜑(𝐲;𝜽)𝑑𝐲, (6) where 𝜑(𝐲;𝜽) is an unknown density function (see [10]).

Second, we present the traditional MUSIC-type algorithm for imaging of perfectly conducting cracks. For the sake of simplicity, we exclude the constant −exp(𝜋/4)/√8𝜋𝑘 from formula (6). Based on [20, 22], we assume that the crack is divided into 𝑀 different segments of size of the order of half the wavelength 𝜆/2. Having in mind the Rayleigh resolution limit for far-field data, only one point at each segment is expected to contribute to the image space of the response matrix 𝕂 (i.e., see [20, 22, 23]). Each of these points, say 𝐲_𝑚, 𝑚=1,2,…,𝑀, will be imaged via the MUSIC-type algorithm. With this assumption, we perform the following singular value decomposition (SVD) of the multistatic response (MSR) matrix 𝕂=[𝑢_far(𝐱_𝑖;𝜽_𝑙)]^𝑁_𝑖,𝑙=1∈ℂ^𝑁×𝑁: 𝕂=𝕌𝕊𝕍^∗_𝑚=^𝑀_𝑚=1𝜏_𝑚𝐔_𝑚𝐕^∗_𝑚, (7) where superscript ∗ is the mark of Hermitian, 𝐔_𝑚 and 𝐕_𝑚∈ℂ^𝑁×1 are, respectively, the left- and right-singular vectors of 𝕂, and 𝜏_𝑚 denotes singular values that satisfy 𝜏₁≥𝜏₂≥⋯≥𝜏_𝑚>0,𝜏_𝑚=0,for𝑚≥𝑀+1. (8) If so, {𝐔₁,𝐔₂,…,𝐔_𝑀} are the basis for the signal and {𝐔_𝑀+1,𝐔_𝑀+2,…,𝐔_𝑁} span the null space of 𝕂, respectively. Therefore, one can define the projection operator onto the null subspace, 𝐏_noise∶ℂ^𝑁×1→ℂ^𝑁×1. This projection is given explicitly by 𝐏_noise∶=𝕀_𝑁−^𝑀_𝑚=1𝐔_𝑚𝐔^∗_𝑚, (9) where 𝕀_𝑁 denotes the 𝑁×𝑁 identity matrix. For any point 𝐳∈ℝ², we define a test vector 𝐟(𝐳,𝜔)∈ℂ^𝑁×1 as 𝐟()=𝐳,𝜔exp𝑗𝜔𝜽₁⋅𝐳,exp𝑗𝜔𝜽₂⋅𝐳,…,exp𝑗𝜔𝜽_𝑁⋅𝐳^𝑇. (10) Based on this, we can design a MUSIC-type imaging function 𝑊∶ℂ^𝑁×1→ℝ such that 𝔼(𝐳)=𝐏_noise(𝐟(𝐳,𝜔))⁻¹=1𝐏_noise(𝐟(.𝐳,𝜔)) (11) Then, the map of 𝔼(𝐳) will have peaks of large and small magnitudes at 𝐳∈Ω and 𝐳∈ℝ²⧵Ω, respectively.

2.2. Multifrequency MUSIC-Type Imaging Function

We design multifrequency MUSIC-type imaging function and try to describe its structure. First, we introduce a multifrequency MUSIC-type algorithm 𝔼_MF∶ℂ^𝑁×1→ℝ defined by 𝔼_MF()=1𝐳;𝑆𝑆^𝑆_𝑠=1𝐏_noise𝐟𝐳,𝜔_𝑠^2−1/2. (12) Then, we can introduce the following lemma. A more detailed derivation can be found in [16].

Lemma 1 (see [16]). Assume that 𝑘_𝑆 and 𝑆 are sufficiently large; then, 𝔼_𝑀𝐹(𝐳;𝑆)≈1𝑁1−^𝑀_𝑚=1Φ(𝐳−𝐲_𝑚;𝜔₁,𝜔_𝑆)^−1/2, (13) where function Φ(𝑥;𝜔₁,𝜔_𝑆) is defined as Φ𝑥;𝜔₁,𝜔_𝑆1∶=𝜔_𝑆−𝜔₁𝜔_𝑆𝐽₀𝜔_𝑆𝑥²+𝐽₁𝜔_𝑆𝑥²−𝜔₁𝐽₀𝜔₁𝑥²+𝐽₁𝜔₁𝑥²+^𝜔𝑆_𝜔1𝐽₁(𝜔𝑥)²𝑑𝜔. (14)

So we can recognize the mathematical structure of multifrequency MUSIC-type algorithm. However, the finite representation of ∫𝐽²₁(𝑥)𝑑𝑥 does not exist. Because of this term, although this can be negligible (see [15]), the map of 𝔼_MF(𝐳;𝑆) should generate unexpected points of small magnitudes. In order to solve this problem, the last term of [14] should be eliminated. For this, we suggest a weighted multifrequency MUSIC-type imaging algorithm in the upcoming section.

3. Weighted Multifrequency MUSIC-Type Algorithm and Its Structure

In order to propose the weighted multifrequency MUSIC-type imaging algorithm, we introduce the following lemma derived from [16].

Lemma 2 ([16, page 218]). For sufficiently large 𝑁 and 𝜔, the following relationship holds: 𝐏_{𝑛𝑜𝑖𝑠𝑒}(𝐟(𝐳,𝜔))≈√𝑁1−^𝑀_𝑚=1𝐽₀(𝜔|𝐳−𝐲_𝑚|)^21/2. (15)

With this, let us define an alternative projection operator weighted by applied frequency 𝐏_WN∶ℂ^𝑁×1→ℝ as 𝐏_WN(𝐟(𝐳,𝜔))=𝐏_noise√𝜔𝐟(𝐳,𝜔). (16) Then, the following result can be obtained.

Theorem 3. Assume that 𝑁 and 𝜔 are sufficiently large; then, 𝐏_WN(𝐟(𝐳,𝜔))≈√𝑁𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚^21/2. (17)

Proof. The following equations are satisfied by the definition of 𝐏_WN(𝐟(𝐳,𝜔)) and by Lemma 2: 𝐏_WN(𝐟(𝐳,𝜔))=𝐏_noise√𝜔(𝐟(=√𝐳,𝜔))𝜔𝐏_noise𝐟(≈√((𝐳,𝜔)))𝜔√𝑁1−^𝑀_𝑚=1𝐽₀𝜔𝐳−𝐲_𝑚^21/2=√𝑁𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚^21/2. (18)

Next, we introduce a weighted multifrequency MUSIC-type imaging function based on MUSIC-type imaging function 𝔼_WMF∶ℂ^𝑁×1→ℝ defined by 𝔼_𝑊𝑀𝐹()=1𝐳;𝑆𝑆^𝑆_𝑠=1𝐏_WN𝐟𝐳,𝜔_𝑠^2−1/2. (19) Then, we can obtain the structure of 𝔼_WMF(𝐳;𝑆).

Theorem 4. Assume that 𝑆 and 𝜔_𝑆 are sufficiently large; then, 𝔼_𝑊𝑀𝐹()≈𝐳;𝑆1𝑁𝜔_𝑆+𝜔₁2−^𝑀_𝑚=1Ψ𝐳−𝐲_𝑚;𝜔₁,𝜔_𝑆^−1/2, (20) where function Ψ(𝑥;𝜔₁,𝜔_𝑆) is defined as Ψ𝑥;𝜔₁,𝜔_𝑆1∶=𝜔_𝑆−𝜔₁𝜔²_𝑆2𝐽₀𝜔_𝑆𝑥²+𝐽₁𝜔_𝑆𝑥²−𝜔²₁2𝐽₀𝜔₁𝑥²+𝐽₁𝜔₁𝑥². (21)

Proof. By Theorem 3, we can calculate the following: 𝔼_WMF()=1𝐳;𝑆𝑆^𝑆_𝑠=1𝐏_WN𝐟𝐳,𝜔_𝑠^2−1/2≈1𝑆^𝑆_𝑠=1√𝑁𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚^21/22−1/2=1𝑆^𝑆_𝑠=1𝑁𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚^2−1/2=1𝑁^𝑆_𝑠=1𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚²1𝑆^−1/2. (22) Then, since 𝑆 is sufficiently large, we can observe that ^𝑆_𝑠=1𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚²1𝑆≈1𝜔_𝑆−𝜔₁^𝜔𝑆_𝜔1𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚²𝑑𝜔, (23) and applying an indefinite integral of the Bessel function (see [24, page 106]) 𝑥𝐽₀(𝑥)²𝑥𝑑𝑥=²2𝐽₀(𝑥)²+𝐽₁(𝑥)² (24) yields 1𝜔_𝑆−𝜔₁^𝜔𝑆_𝜔1𝜔𝐽₀𝜔𝐳−𝐲_𝑚²=1𝑑𝜔𝜔_𝑆−𝜔₁𝜔²_𝑆2𝐽₀𝜔_𝑆𝐳−𝐲_𝑚²+𝐽₁𝜔_𝑆𝐳−𝐲_𝑚²−𝜔²₁2𝐽₀𝜔₁𝐳−𝐲_𝑚²+𝐽₁𝜔₁𝐳−𝐲_𝑚²=Ψ𝜔𝐳−𝐲_𝑚;𝜔₁,𝜔_𝑆. (25)
Hence, we can obtain 1𝜔_𝑆−𝜔₁^𝜔𝑆_𝜔1𝜔−^𝑀_𝑚=1𝜔𝐽₀𝜔𝐳−𝐲_𝑚²=𝜔𝑑𝜔_𝑆+𝜔₁2−^𝑀_𝑚=1Ψ𝜔𝐳−𝐲_𝑚;𝜔₁,𝜔_𝑆. (26)
Therefore, 𝔼_WMF(𝐳;𝑆)≈1𝑁𝜔_𝑆+𝜔₁2−^𝑀_𝑚=1Ψ(𝐳−𝐲_𝑚;𝜔₁,𝜔_𝑆)^−1/2. (27) This completes the proof.

Looking at the results of Theorem 4, in contrast to the 𝔼_MF(𝐳;𝑆), the 𝔼_WMF(𝐳;𝑆) does not have the ∫𝐽²₁(𝑥)𝑑𝑥 term. Therefore, we expect that the imaging results of the 𝔼_WMF(𝐳;𝑆) will be better than 𝔼_MF(𝐳;𝑆). In the next section, numerical experiments will be presented to support this.

4. Numerical Experiments

In this section, some numerical examples are displayed in order to support our analysis in the previous section. Applied frequencies are of the form 𝜔_𝑠=2𝜋/𝜆_𝑠, where 𝜆_𝑠, 𝑠=1,2,…,𝑆(=10) is the given wavelength. The observation directions 𝜽_𝑛∈𝕊¹ are taken as 𝜽_𝑛=cos2𝜋𝑛𝑁,sin2𝜋𝑛𝑁^𝑇. (28)

For illustrating arc-like cracks, three Ω_𝑙 are chosen: Ω₁=1𝑠,2cos𝑠𝜋2+15sin𝑠𝜋2−110cos3𝑠𝜋2^𝑇∶[,Ω𝑠∈−1,1]₂=𝑠2sin2,sin𝑠^𝑇𝜋∶𝑠∈4,7𝜋4,Ω₃=Ω₃⁽¹⁾∪Ω₃⁽²⁾, (29) where Ω₃⁽¹⁾=1𝑠−5𝑠,−²2+35^𝑇−1∶𝑠∈2,12,Ω₃⁽²⁾=1𝑠+5,𝑠³+𝑠²−35^𝑇−1∶𝑠∈2,12. (30)

It is worth emphasizing that all the far-field data 𝑢_far of (6) are generated by the method introduced in [25, Chapter 3, 4]. After generating the data, a 20 dB white Gaussian random noise is added to the unperturbed data. In order to obtain the number of nonzero singular values 𝑀 for each frequency, a 0.1-threshold scheme (choosing first 𝑚 singular values 𝜏_𝑚 such that 𝜏_𝑚/𝜏₁≥0.1) is adopted. A more detailed discussion of thresholding can be found in [20, 22].

Figures 1 and 2 show the imaging results via multifrequency MUSIC and weighted multifrequency MUSIC algorithms for single crack Ω₁ and Ω₂, respectively. As we already mentioned, since the term ∫𝐽²₁(𝑥)𝑑𝑥 can be disregarded, it is very hard to compare the improvements via visual inspection of the reconstructions. However, based on Figure 3, we can examine that the proposed weighted multifrequency MUSIC algorithm successfully reduces these artifacts, so we can conclude that this is an improved version.

(a)  Map of 𝔼MF(𝐳;10)

(b) Map of 𝔼WMF(𝐳;10)

Figure 1 

Shape reconstruction of Ω1 via MUSIC algorithm.

(a) Map of 𝔼MF(𝐳;10)

(b) Map of 𝔼WMF(𝐳;10)

Figure 2 

Same as Figure 1 except the crack is Ω2.

(a) Map of 𝔼MF(𝐳;10)

(b) Graph of oscillating pattern

(c) Map of 𝔼MF(𝐳;10)

(d) Graph of oscillating pattern

Figure 3 

Blue- and red-colored lines are 𝔼WMF(𝐳;10) and 𝔼MF(𝐳;10), respectively, at 𝐳=[−0.8,𝑦]𝑇 ((a) and (b)) and 𝐳=[𝑥,0]𝑇 ((c) and (d)).

Figure 4 shows the imaging results via multifrequency MUSIC and weighted multifrequency MUSIC algorithms for multiple cracks Ω₃. Similar to the imaging of single crack, we can observe that weighted multifrequency MUSIC algorithm improves the traditional one, although it is hard to compare the improvements via visual inspection.

(a)  Map of 𝔼MF(𝐳;10)

(b) Map of 𝔼WMF(𝐳;10)

Figure 4 

Same as Figure 1 except the crack is Ω3.

Figure 5 shows the noise contribution in terms of SNR. In order to observe the effect of noise, 30 dB and 10 dB white Gaussian random noises are added to the unperturbed data. Based on these results, we can easily observe that both traditional and proposed MUSIC algorithms offer very good result when 30 dB noise is added. However, when 10 dB noise is added, the traditional MUSIC algorithm yields a poor result while the proposed algorithm yields an acceptable result.

(a) Map of 𝔼(𝐳)

(b) Map of 𝔼WMF(𝐳;10)

(c) Map of 𝔼(𝐳)

(d) Map of 𝔼WMF(𝐳;10)

Figure 5 

Influence of noise. 30 dB ((a) and (b)) and 10 dB ((c) and (d)) white Gaussian random is added.

Now, we consider the imaging of oscillating crack. For this, we consider the following cracks (see Figure 6): Ω₄=1𝑠,2𝑠²+110sin(4𝜋(𝑠+1))^𝑇[Ω∶𝑠∈−1,1]₅=1𝑠,2𝑠²+1120sin(20𝜋(𝑠+1))−100cos(15𝜋𝑠)^𝑇∶𝑠∈[−1,1]. (31)

(a) Crack Ω4

(b) Crack Ω5

Figure 6 

Shape of oscillating cracks Ω4 and Ω5.

Figure 7 shows the maps of 𝔼(𝐳) and 𝔼_WMF(𝐳;10) for the crack Ω₄. In this result, we can observe that both traditional and proposed MUSIC algorithms produce acceptable result, but the proposed algorithm successfully eliminates replicas.

(a) Map of 𝔼(𝐳)

(b) Map of 𝔼WMF(𝐳;10)

Figure 7 

Shape reconstruction of Ω4 via 𝔼(𝐳) (a) and 𝔼WMF(𝐳;10) (b).

Figure 8 shows the maps of 𝔼(𝐳) and 𝔼_WMF(𝐳;10) for highly oscillating crack Ω₅. Opposite to Figure 7, both traditional and proposed MUSIC algorithms yield poor result. This example shows the limitation of proposed algorithm.

(a) Map of 𝔼(𝐳)

(b) Map of 𝔼WMF(𝐳;10)

Figure 8 

Same as Figure 7 except the crack is Ω5.

5. Conclusion

Based on the structure of multifrequency MUSIC-type imaging function, we introduced a weighted multifrequency MUSIC-type imaging function. Through a careful analysis, a relationship between imaging function and Bessel function of the first kind of integer order is established, and we have confirmed that the proposed imaging algorithm is an improved version of the traditional one.

Although, the proposed algorithm produces very good results and improves the traditional MUSIC algorithm, it still needs some upgrade, for example, imaging of highly oscillating cracks. Development of this should be an interesting and remarkable research project. In this paper, we considered the MUSIC algorithm in full-view inverse scattering problem. Based on the result in [26], the MUSIC algorithm cannot be applied to limited-view problems but the reason is still unknown. Identifying the structure of the MUSIC algorithm in the limited-view inverse scattering problems will be the forthcoming work.

Acknowledgments

The author would like to express his thanks to Won-Kwang Park (Kookmin University) for his valuable discussions and MATLAB simulations to generate the forward data and MUSIC-type imaging function. The author would like to acknowledge two anonymous referees for their precious comments. This work was supported by the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education, Science and Technology (no. 2011-0007705) and the research program of Kookmin University in Korea.