Fast Localization of Small Inhomogeneities from Far-Field Pattern Data in the Limited-Aperture Inverse Scattering Problem

: In this study, we consider a sampling-type algorithm for the fast localization of small electromagnetic inhomogeneities from measured far-ﬁeld pattern data in the limited-aperture inverse scattering problem. For this purpose, we designed an indicator function based on the structure of left- and right-singular vectors of a multistatic response matrix, the elements of which were measured far-ﬁeld pattern data. We then rigorously investigated the mathematical structure of the indicator function in terms of purely dielectric permittivity and magnetic permeability contrast cases by establishing a relationship with an inﬁnite series of Bessel functions of an integer order of the ﬁrst kind and a range of incident and observation directions before exploring various intrinsic properties of the algorithm, including its feasibility and limitations. Simulation results with synthetic data corrupted by random noise are presented to support the theoretical results.


Introduction
In this study, we consider the fast localization of a set of small electromagnetic inhomogeneities embedded in a homogeneous space from far-field pattern data measured over a limited-aperture configuration. The inhomogeneities are characterized by either a contrast of dielectric permittivity in relation to the exterior space, a contrast of magnetic permeability, or both contrasts at a fixed angular frequency ω = 2π f . Throughout the paper, we address the mathematical treatment of the scattering of time-harmonic electromagnetic waves from thin infinitely long cylindrical obstacles. The problem in question is a two-dimensional issue and is correspondingly associated with a scalar scattering problem for Transverse Magnetic (TM) and Transverse Electric (TE) waves in a cylindrical waveguide [1,2]. In the process, we reference various studies [3][4][5][6][7][8][9][10] related to the application of the limited-aperture inverse scattering problem.
The first part of this paper is focused on designing specific indicator functions for permittivity and/or permeability contrast cases to detect the location of small electromagnetic inhomogeneities from the constructed Multistatic Response (MSR) matrix, the elements of which are far-field pattern data for various incident fields [11]. The initial assumption was that the elements of the MSR matrix can be represented by an asymptotic field formulation in the presence of inhomogeneities [12]. As such, based on the singular-value decomposition of the MSR matrix, we generated an appropriate test vector consisting of the incident field at each search point before we used the orthonormal property of the left-and right-singular vectors of the MSR matrix and present a method for designing the indicator functions for localizing the inhomogeneities. In terms of the related works on sampling-type imaging techniques, we refer to [13][14][15][16] for subspace migration, [17][18][19][20] for the linear sampling method, [21][22][23][24] for the direct sampling method, and [25][26][27][28] for the orthogonality sampling method. An extension to the real-world problem can be found in [29][30][31] for subspace migration, [32][33][34][35] for the direct sampling method, and [36][37][38][39][40] for the orthogonality sampling method.
In the second part of the paper, we establish a new mathematical theory for the structure of the indicator functions. Using the asymptotic expansion formula and the uniform convergence property of the so-called Jacobi-Anger expansion formula [2], we demonstrate that the indicator function can be expressed by an infinite series of Bessel functions of an integer order of the first kind and a range of incident and observation directions. Meanwhile, in terms of the dielectric permittivity contrast case, we examine the notion that the main factor of the indicator function is the Bessel function of a zero order of the first kind, while in terms of the magnetic permeability contrast case, by adopting specific test vectors consisting of a standard basis for R 2 , we focus on how the main factor of the indicator function is the Bessel functions of a zero order and two of the first kind. The attendant established theory indicates that the imaging resolution depends on the selection of the applied frequency and that the performance is highly dependent on the selection of the range of incident and observation directions. Moreover, it also provides the condition for the range of both the incident and the observation directions to guarantee rigorous results.
In the last part of the paper, we perform various numerical simulations using noisy data generated via the Foldy-Lax formulation [41] to demonstrate the feasibility and limitations of the designed indicator function. Specifically, we address the fact that while the method is fast and stable, it involves certain difficulties, such as the discrimination of appropriate nonzero singular values and the low-quality imaging results obtained when the range of incident or observation directions is narrow. Hence, while only the location or outline shape of the inhomogeneities can be retrieved via the designed algorithm, the results can be regarded as a good initial representation of the iterative-based techniques [42][43][44][45][46][47][48][49][50].
This paper is organized as follows. In Section 2, we outline the basic concept of the direct scattering problem, as well as the far-field pattern formula, including in terms of the asymptotic expansion formula. Then, in Section 3, we design the sampling-type indicator functions for the permittivity and permeability contrast cases, investigate the mathematical structure of the indicator functions by establishing relationships with the Bessel functions of the integer order of the first kind and the range of incident and observation directions, and explore various properties of the indicator functions. Following this, in Section 4, the numerical simulation results with synthetic data for the purely dielectric permittivity and magnetic contrast cases are presented to support the investigated theoretical results and the examined properties of the indicator functions.

Direct Scattering Problem and Far-Field Pattern
In this section, we briefly introduce the two-dimensional direct scattering problem and the asymptotic expansion formula for the far-field pattern in the presence of a set of small electromagnetic inhomogeneities Σ = {Σ s : s = 1, 2, · · · , S}. Throughout the paper, we assume that all Σ s are small balls and are well separated from each other, which can be expressed as follows: where r s denotes the location, α s characterizes the size, and S 1 is the unit circle centered at the origin. Here, we let Ω ⊂ R 2 be the region of interest and set Σ s ⊂ Ω for all s = 1, 2, · · · , S. Meanwhile, we assumed that all materials are characterized by their value of dielectric permittivity and magnetic permeability at a given angular frequency ω. Let ε 0 and µ 0 denote the value of the permittivity and permeability of R 2 , respectively. Analogously, we can use ε s and µ s to represent those of Σ s and can correspondingly introduce the following piecewise constants of permittivity and permeability: ε(r) = ε s for r ∈ Σ s ε 0 for r ∈ R 2 \Σ and µ(r) = µ s for r ∈ Σ s µ 0 for r ∈ R 2 \Σ, Following this, we can denote k = 2π/λ as the background wavenumber that satisfies k 2 = ω 2 ε 0 µ 0 and the following relationship among the well-separated inhomogeneities: for s, s = 1, 2, · · · , S and s = s , where λ denotes the positive wavelength.
Here, we can consider the plane-wave illumination, that is the incident field with a propagation direction θ ∈ S 1 is given by u inc (r, θ) = e ikθ·r . As such, let u tot (r, θ) be the time-harmonic total field that satisfies the following: with transmission conditions at the boundaries of Σ s . Let u scat (r, θ) be the corresponding scattered field that satisfies u scat (r, θ) = u tot (r, θ) − u inc (r, θ) and the Sommerfeld radiation condition, uniformly in all directions ϑ = r/|r|. Following [11], the unknown scattered field u scat (r, θ) can be written in terms of the following single-layer potential with unknown density function ϕ: 0 is the Hankel function of a zero order of the first kind. The measurement data are the so-called far-field pattern u ∞ (ϑ, θ) of u scat (r, θ). Based on [2], u ∞ (ϑ, θ) can be defined in terms of S 1 , and this satisfies: uniformly on ϑ = r/|r|, as |r| −→ ∞. It should be noted that u ∞ (ϑ, θ) can be written as follows: Remark 1. According to the research [2,11], methods for solving the inverse scattering problem can be classified as either quantitative or qualitative ones. Quantitative methods are used to retrieve the parameter distributions ε(r) or µ(r) in a region of interest. However, it is very hard to solve the problem directly due to the ill-posedness and nonlinearity. Thus, most algorithms are based on iterative schemes. Instead, qualitative methods are used to identify the location r s . In this paper, as already mentioned in the Introduction, we consider a sampling-type imaging technique classified as a qualitative method to retrieve r s from measurement data u ∞ (ϑ, θ).

Introduction and Analysis of Indicator Functions
We can introduce the indicator function in terms of the limited-view inverse scattering problem. To this end, we can let S 1 obs = ϑ p : p = 1, 2, · · · , P and S 1 inc = θ q : q = 1, 2, · · · , Q be the sets of observation and incident directions, respectively, which are connected, appropriate subsets of S 1 (see Figure 1). Here, ϑ p and θ q could be set as: respectively. Following this, we can consider the following MSR matrix, the elements of which are far-field pattern data: We could then proceed with the design and analysis in terms of three specific cases: a dielectric contrast-only case (ε(r) = ε 0 and µ(r) = µ 0 ), a magnetic contrast-only case (ε(r) = ε 0 and µ(r) = µ 0 ), and a case using both contrasts (ε(r) = ε 0 and µ(r) = µ 0 ).

Permittivity Contrast Case
First, we must focus on the dielectric permittivity contrast-only case. In this case, we assumed that the total number of observation and incident directions was larger than the total number of anomalies, i.e., P, Q > S, and ε 1 ≥ ε 2 ≥ · · · ≥ ε S . Since the complete form of ϕ in (4) is unknown, we cannot directly use far-field pattern data for identifying Σ s . Instead, according to [11], the far-field pattern could be represented as follows: and correspondingly, K could be represented as follows: where: e ikθ 1 ·r 1 e ikθ 2 ·r 1 · · · e −ikθ Q ·r 1 e ikθ 1 ·r 2 e ikθ 2 ·r 2 · · · e −ikθ Q ·r 2 . . . . . . . . . . . . and: Based on this representation, we designed an indicator function as follows. First, we conducted the singular-value decomposition (SVD) of K such that: where the superscript * represents the Hermitian operator and τ s ∈ R, U s ∈ C P×1 , and V s ∈ C Q×1 are the nonzero singular value and the left-and right-singular vectors of K, respectively. By comparing Equations (7) and (8), we could introduce the following unit vectors to test the orthonormality relationship with the singular vectors U s and V s : for r ∈ Ω, Then , the following orthonormality relationships hold: where U, V = U · V (refer to, e.g., [13]). Based on the above relationships, we could then introduce the following sampling-type indicator function F ε : Ω −→ [0, 1]; for r ∈ Ω, As such, the value of F ε (r) will become the largest one 1 at r = r s ∈ Σ s , s = 1, 2, · · · , S, and the small ones (between 0 and 1) at r ∈ R 2 \Σ, which made it possible to identify the locations of all Σ s . To better understand the further properties of the indicator function, including the feasibility, the effect on the range of incident and observation directions, and any fundamental limitations, we investigated the mathematical structure of F ε (r), with the results outlined below. Theorem 1. Let r − r s = |r − r s |[cos φ s , sin φ s ] T . If ϑ and θ are small, F ε (r) can be represented as follows: where J n denotes the Bessel function of the integer order n and: Proof. See Appendix A.

Permeability Contrast Case
In terms of the magnetic permeability contrast-only case, we assumed that the total number of observation and incident directions were larger than twice the total number of anomalies, i.e., P, Q > 2S, and µ 1 ≤ µ 2 ≤ · · · ≤ µ S . In this case, since the complete form of the far-field pattern formula is still unknown, we use the following asymptotic expansion formula for the far-field pattern: where M(Σ s ) is a 2 × 2 matrix, the elements of which are 2µ 0 /(µ s + µ 0 ); refer to [11]. Correspondingly, K could be represented as follows: where: To design the indicator function, as noted above, we conducted the SVD of K such that: Then , by comparing Equations (15) and (16), we introduced the following unit vectors to test the orthonormality relationship with the singular vectors U s and V s : for r ∈ Ω, Then, the following orthonormality conditions hold for s = 1, 2, · · · , S (see, e.g., [13]): Based on the above relationships, we could introduce the following sampling-type indicator function F µ : As such, the value of F µ (r) will become the largest one 1 at r = r s ∈ Σ s , s = 1, 2, · · · , S, and the small ones (between 0 and 1) at r ∈ R 2 \Σ, which made it possible to identify the locations of all Σ s . To better understand the further properties of the indicator function, including the feasibility, the effect on the range of incident and observation directions, and any fundamental limitations, we investigated the mathematical structure of F µ (r). Here, we derived the following useful identity: Lemma 1. For sufficiently large P and Q, ϑ = [cos ϑ, sin ϑ] T ∈ S 1 obs , θ = [cos θ, sin θ] T ∈ S 1 inc , and r = |r|(cos φ, sin φ) ∈ R 2 , the following relationship holds: where: Now, we can explore the structure of F µ (r).
If ϑ and θ are small, F µ (r) can be represented as follows: where: and: Proof. See Appendix C.

The Case of Both Permittivity and Permeability Contrasts
Here, we consider the case involving both the dielectric permittivity and the magnetic permeability contrasts. In this case, we assumed that the total number of observation and incident directions was larger than the total number of anomalies, i.e., P, Q > 3S. Then, according to (6) and (14), the elements of K could be represented as follows: and correspondingly, based on Equations (7) and (15), K could be represented as follows: Throughout the processes outlined in Sections 3.1 and 3.2, we introduce the following indicator function F ε,µ : Ω −→ [0, 1]; for r ∈ Ω, inc (r), V 2s . (24) Then, by combining the results of Theorems 1 and 2, we could obtain the following result: If ϑ and θ are small, F ε,µ (r) can be represented as follows: where inc (r) are as defined in Theorems 1 and 2.

Properties of the Indicator Functions
Based on Theorems 1 and 2, we could examine certain properties of the indicator functions, which can be summarized as follows: Property 1 (Feasibility). Given that J 0 (0) = 0, we can observe that J 0 (k|r − r s |) = 1 when r = r s ∈ Σ s (see Figure 2). Since the maximum value of indicator functions is equal to one, the term J 0 (k|r − r s |) contributed to the imaging for localization of the inhomogeneities in both the permittivity and the permeability contrasts. Here, it should be noted that since J n (0) = 0 for n ∈ N (see Figure 2) and the terms Λ inc (r) contain J n (k|r − r s |), n ∈ N, their values are close to zero when r = r s ∈ Σ s , meaning they did not contribute to the imaging for the localization of inhomogeneities. Moreover, due to the oscillating property of the Bessel function, these terms had a deleterious effect on the imaging (i.e., disturbed the localization) as they caused the generation of many artifacts. Meanwhile, the oscillation pattern depended on the value of k, i.e., the applied frequency. Here, if a high frequency were applied, a high-resolution imaging result could be obtained, albeit that it would contain several artifacts. In contrast, when a low frequency is applied, the appearance of artifacts would be reduced, but a low-resolution image would be obtained (see Figure 3).   inc (r) are inversely proportional to the values ϑ P − ϑ 1 and θ Q − θ 1 , and disturbing terms are significantly dependent on the range of the incident and observation directions. More specifically, if the range of both the incident and observation directions is sufficiently wide, the indicator functions F ε (r), F µ (r), and F ε,µ (r) will be dominated by the factor J 0 (k|r − r s |). Therefore, it can be expected that a good result will be obtained. However, if the range of the incident or observation directions is narrow, the effect of the factor J 0 (k|r − r s |) becomes negligible and the quality of the imaging result will, therefore, likely be poor).
Property 3 (Limitation: dependency on the total number of directions). Based on the design of the indicator functions, if there exists S inhomogeneities, at least S and 2S incident and observation directions will be required for the permittivity and permeability contrast cases, respectively. If the range of these directions is narrow, this condition cannot be satisfied, meaning some inhomogeneities could not be identified.
Property 4 (Least range of directions). By using a tedious calculation, it is possible to examine the idea that the terms Λ inc (r) disappear when using a full-view setting ϑ P − ϑ 1 = θ Q − θ 1 = 2π. However, the full-view problem is not considered in this paper. When 0 < ϑ P − ϑ 1 < 2π and 0 < θ Q − θ 1 < 2π, the following condition can be found: ϑ 1 = θ 1 = φ s and ϑ P = ϑ Q = π + φ s for all s = 1, 2, · · · , S, or equivalently: This condition tells us that if the range of both the incident and the observation directions is equal to or wider than π, good results can be expected. This is supported by the results presented in Figures 4c and 5c.
Throughout this study, the measurement data u ∞ (ϑ p , θ q ) of the MSR matrix were chosen as the z-component of the electric (permittivity contrast case) and magnetic (permeability contrast cast) fields (see, e.g., [1]) and were generated by solving the Foldy-Lax formulation [41] to prevent any inverse crime and 20 dB Gaussian random noise being added to the generated data. To discriminate the nonzero singular values of K, a 0.1threshold technique (counting the number of singular values τ s that satisfies τ s /τ 1 ≥ 0.1) was adopted (see, e.g., [51]). −π/6 −π/6 −π/6 −π/3 −π/3 −π/3 Q 6 6 6 11 11 11 Example 1 (Permittivity contrast case). Figure 4 shows the distribution of the normalized singular values of K and the imaging results of F ε (r) when the permeabilities equated to µ s ≡ µ 0 and the permittivities of Σ s were ε 1 = 5, ε 2 = 3, and ε 3 = 2. Based on the distribution of the normalized singular values, it was feasible to select three nonzero singular values to define the indicator function, with this number equal to the total number of inhomogeneities. However, when the range of both the incident and the observation directions was narrow, it proved to be extremely difficult to identify the true location of Σ s due to the appearance of a blurring effect in the neighborhood of r s (see Figures 4(b) and 4a). Here, it was interesting to ascertain whether, based on the results shown in Figure 4c, it was possible to obtain a good result when the range of observation directions was π, despite the narrow range of the incident directions. However, while the range of the incident direction was wider than that shown in Figure 4c, blurring effects still emerged in the map of F ε (r) when the range of the observation directions was narrower than π (see Figure 4e).  Example 2 (Permeability contrast case). Figure 5 shows the distribution of the normalized singular values of K and the imaging results of F µ (r) when the permittivity equated to ε s ≡ ε 0 and the permittivities of Σ s were µ 1 = 5, µ 2 = 3, and µ 3 = 2. Based on the distribution of the normalized singular values, unlike in the permittivity contrast case, it proved difficult to select the appropriate number of nonzero singular values to define the indicator function, with this number not equal to twice the total number of inhomogeneities. Correspondingly, the location of certain inhomogeneities (specifically, r 1 ∈ Σ 1 , the permeability of which was larger than that of the others) could not be identified via the map of F µ (r). However, similar to the permittivity contrast case, it was possible to obtain a good result when the range of observation directions was π (see Figure 5c), while blurring effects still emerged in the neighborhood of r s when the range of observation directions was narrower than π (see Figure 5d,e).

Example 3 (Effects on nonuniform directions).
Throughout this paper, we assume that the sets of observation S 1 obs and incident S 1 inc directions are uniformly distributed. Here, we examine the effect of the distribution of observation and incident directions. Figure 6 shows the distribution of the normalized singular values of K and the imaging results of F ε (r) with the same simulation configuration of Example 1, except the application of a nonuniform observation and incident directions. Notice that opposite Figure 4b, it is natural to select two singular values, and the result is poor; refer to Figure 6b. Moreover, although the range of observation (Figure 6c) and incident directions (Figure 6e) was considered to be wide enough, the obtained imaging results were poorer than the ones by using uniform incident directions. A similar phenomenon for the permeability contrast case can be examined in Figure 7.   Example 4 (Comparison with the orthogonality sampling method). For the final example, we applied the orthogonality sampling method to compare the imaging performances. Figure 8 shows the imaging results through the indicator function F OSM (r) of the orthogonality sampling method with the same simulation configuration of Example 1. The traditional indicator function of the orthogonality sampling method with multiple incident directions was given by [28] (Section 3.3 ). By comparing the result in Figure 4, we can examine that it is hard to identify all inhomogeneities through the map of F OSM (r), and only the location r 1 , whose permittivity is larger than others, can be recognized. Thus, based on this result, we can say that the designed technique guarantees good results.

Conclusions
Based on the structure of the left-and right-singular vectors of the MSR matrix and the asymptotic expansion formula in the presence of small electromagnetic inhomogeneities, we designed a sampling-type indicator function for identifying the location of small inhomogeneities in the limited-aperture inverse scattering problem. To demonstrate the feasibility and the limitations and to explore various the properties of the designed indicator function, we proved that the indicator function can be expressed by an infinite series of Bessel functions of an integer order and the range of the incident and observation directions. The numerical simulation results with noisy data supported the theoretical results and demonstrated the explored properties. In this paper, the identification of small inhomogeneities was addressed, and extending the investigation to the shape identification of arbitrarily shaped targets, such as extended objects or arc-like cracks, will undoubtedly be an interesting research topic. Moreover, the analysis of the indicator function when the sets of observation and incident directions are uniformly distributed will be the forthcoming work. Finally, applications to real-world microwave imaging [30,35,52] and extending matters to include the three-dimensional problem [53][54][55] will also prove to be interesting approaches.

Acknowledgments:
The author would like to acknowledge the anonymous referees for their valuable comments, which helped to enhance the quality of the paper.

Conflicts of Interest:
The author declares no conflict of interest.

Appendix A. Proof of Theorem 1
From (10), the following relationships hold: Given that ϑ and θ are small, the following approximation holds uniformly (see, e.g., [15]) for ϑ p , ϑ ∈ S 1 obs and r = |r|[cos φ, sin φ] T : 1 P P ∑ p=1 e ikϑ p ·r ≈ 1 Thus, we could derive the following: obs (r) is given by (13). This was also the case for θ q , θ ∈ S 1 inc and r = |r|[cos φ, sin φ] T , where, since the following approximation holds uniformly, we could derive the following: inc (r) is given by (13). Given that τ s ∈ R, e iα s e iβ s = e i(α s +β s ) = e i arg(τ s ) = 1, we could obtain: which led us to the structure (12).