Assessing the Capabilities of a New Linear Inversion Method for Quantitative Microwave Imaging

1Department of Electrical, Electronics and Computer Engineering (DIEEI), University of Catania, Viale A. Doria 6, 95126 Catania, Italy 2Department of Information Engineering, Infrastructures and Sustainable Energy (DIIES), University Mediterranea of Reggio Calabria, Via Graziella, Localita Feo di Vito, 89124 Reggio di Calabria, Italy 3Institute for Electromagnetic Sensing of the Environment (IREA), National Research Council of Italy (CNR), Via Diocleziano 328, 80124 Napoli, Italy


Introduction
The development of reliable solution procedures in electromagnetic inverse scattering is an important topic as it would open the way to effective microwave imaging techniques in many different areas ranging from biomedical to subsurface imaging, nondestructive test, and civil engineering.However, it is well known that the possibility to achieve quantitative results, that is, characterizing both the geometrical and the electromagnetic features of a scattering system, is a tough challenge from both a theoretical and a computational point of view.In fact, microwave imaging entails the solution of a nonlinear and ill-posed problem [1,2].
In particular, ill-posedness requires the adoption of proper regularization strategies [2] that lead to an approximation of the actual unknown, unless a priori information is exploited [3][4][5][6][7][8][9][10]. On the other hand, nonlinearity is usually tackled by formulating the problem in terms of an optimization task, which is computationally heavy or even unfeasible within a global optimization framework.As a result the solution is pursued by means of local optimization strategies, which are prone to the "false solutions" problem [11].
As an alternative, linearized approaches seem to be very attractive from an applicative point of view, since they show inherent advantages in handling both ill-posedness and nonlinearity of the problem.Unfortunately, they can be safely applied only in a limited range.This is the case of the first order Born approximation (BA) that can be exploited to achieve quantitative reconstructions only for weak scatterers, that is, objects which are small with respect to the probing wavelength and whose electromagnetic features are similar to those of the embedding medium [12].As a result, the method applicability is severely limited in any actual instances.This notwithstanding, BA can be applied in all those cases when a qualitative imaging, that is, position, size, and possibly shape of an unknown scattering system, is the goal of the diagnostic survey.

International Journal of Antennas and Propagation
Recently, a linearized approach has been proposed [13], which is based on a novel field approximation that implicitly takes into account the unknown targets.As such, this approach makes the imaging of nonweak targets viable within a linear inversion framework.This interesting capability is achieved by means of a "smart" preprocessing of the measured data, which allows rearranging the original incident fields in such a way to "condition" the scattering phenomenon in a predictable way, which is giving rise to total fields which can be approximately foreseen (and which differ from the "simple" incident field used in typical linearizations).
The recombined incident fields can be seen as a set of new experiments, which we refer to as virtual experiments, as they are achieved without performing additional measurements.While this framework has been applied in different ways to develop linear and nonlinear, as well as approximated and nonapproximated, approaches [13][14][15], in this paper, we investigate, for the first time, the range of validity of the approximation underlying the method in [13].To this end, we consider the canonical case of a homogeneous circular cylinder and perform a parametric numerical study meant to investigate the method's imaging capabilities with respect to the electrical size and the electromagnetic properties of the scatterers.The analysis allows determining the range of applicability of the inversion procedure and shows that its range of validity is significantly larger with respect to the BA.Moreover, both the method and the validity criterion are further assessed against noncanonical targets, confirming the capability of outperforming the BA in imaging nonweak scatterers.
The paper is structured as follows.In Section 2, the mathematical formulation of the imaging problem and the main issues are briefly recalled.In Section 3, the imaging strategy is briefly recalled and outlined.Sections 4 and 5 address the numerical analysis.Discussion of the results and conclusions follow.Throughout the paper, the time factor exp{} is assumed and dropped.

Statement of the Problem
We assume a 2D geometry with the electric field polarized along the invariance axis, that is, the -axis.Let  denote the investigation domain that hosts the cross sections of one or more dielectric scatterers with support Σ and complex dielectric permittivity ε (, ) =   () −   ()/( 0 ),  being the working angular frequency and  = (, ) denoting the generic point in .We consider a homogeneous background medium with complex permittivity ε () =   −   /( 0 ).Both the scatterers and the background are assumed to be nonmagnetic.The region  is probed by means of incident fields impinging from several directions radiated by filamentary currents positioned all around . Accordingly, the field is assumed to be measured on a curve Γ.
Under the above assumptions, the integral "data" and "state" equations governing the scattering phenomenon read, respectively, as wherein   and   denote the incident and the total field, respectively,   denote the th impinging direction , and   = (  , ) denotes the th position on the measurement curve Γ.Moreover, () = ε ()/ε  − 1 is the contrast function, and (  , ) = −(/4) 2 0 (  |  − |) is the Green's function for the considered homogeneous background,   being the (complex) wavenumber of the host medium and  2 0 the second kind zeroth order Hankel function.Accordingly, the inverse scattering problem is cast as the retrieval of the unknown function  for given incident fields in  and measured total (and possibly incident) fields on Γ.
Equations (1) represent the scalar 2D inverse scattering problem which is nonlinear due to the functional relationship between data and unknowns.As a result, linearized models like the Born approximation, may be introduced to avoid this drawback of the problem.However, it is well known how the latter completely neglects the effect of the scattering system on the field, thus leading to unavoidable deterioration of the imaging when it is applied for nonweak scatterers.As a result, no quantitative information about the scatterer's nature can be inferred.Iterative approaches, like the Born [16] and the distorted Born iterative method [17,18], have been introduced to possibly restore the effect of the scatterers on the field solving a cascade of linearized problems.However, these approaches may suffer from lack of convergence, as well as the difficulty in choosing a stopping criterion as tradeoff between reconstruction accuracy and computational burden.

A (Scatterer-Oriented) Approximation and Method's Implementation
In this section we recall the method proposed in [13].In particular, we defer the reader to the original paper for theoretical issues, while addressing mainly guidelines for the method's implementation.We start from the observation that the linearity of Maxwell's equations, and then of the scattering phenomena, suggests that a superposition of the incident fields gives rise to a scattered field which is nothing but the same superposition, that is, with the same "weighting coefficients," of the original scattered fields.Saying it in formulas, we consider a superposition of the incident fields given by which, interacting with the scatterers, will give rise to the scattered field: wherein   =   −   is a noisy version of the scattered field.
Indeed, according to the above, the basic idea underlying the virtual experiments' framework is finding a linear combination of the original data problem (the measured field), so that the scattering phenomenon behaves in such a way that the internal field can be suitably foreseen.
Accordingly, the first step of the method is concerned with the solution of a "design equation."Among many different possibilities, we consider the far field equation, which is the basic equation of the linear sampling method (LSM), an imaging approach developed to recover the shape of an unknown scattering system [19].The discretized version of this equation is cast for each point   of an arbitrary grid that samples the imaged domain, that is, where E  is the  ×  scattered field data matrix as arranged in a multiview-multistatic configuration,  is the dimensional vector unknown, and g is the -dimensional vector containing the values of the Green's function at the receiver positions   .Note that this equation is nothing but the linear combination expressed by means of (3).Equation ( 4) is an ill-posed one that can be effectively solved via singular value decomposition (SVD) of the matrixoperator E  and Tikhonov regularization [2] to find the following solution: wherein k  and u  are the right and left singular vectors of E  , respectively,  = min(, ),   denotes the th singular value of E  , and ⟨ , ⟩ denotes the standard scalar product in  2 .Finally  is the Tikhonov regularization parameter that, according to [13], can be chosen once for all the sampling points.
Once the design equation is solved, a first outcome of the imaging approach is the qualitative imaging pursued by the LSM.As a matter of fact, the energy of the solution, that is, its  2 -norm, assumes low values for sampling points belonging to the scatter's support and high values elsewhere (with respect to its overall dynamics over the sampling grid) [19].Such a behavior can be exploited by considering the following support indicator: which continuously varies assuming the lowest values in those sampling points belonging to the scatterer support Σ. Equation ( 6) is the common way to exploit solution of (4) [19].This result is not the only one actually brought by the LSM solution.Indeed, in the following, we are going to describe how the LSM can be exploited to conveniently afford the inverse scattering problem.In this respect, the second step of the imaging procedure amounts to select a number   of sampling points among those belonging to the retrieved scatterer's support, namely, the "pivot points," and to use the pertaining solution   to achieve the linear combinations (2) and (3).By doing so, the multiview-multistatic standard scattering experiments are converted into a set of multiple virtual (or synthetic) experiments which do not require an additional measurement process since they are built only via software processing.
The advantage to deal with such a (recombined) version of the inverse scattering problem is related to the possibility of introducing an effective approximation for the internal field, which reads Equation ( 7) assumes the total field in  to be the sum of the virtual incident field and the field radiated by a filamentary current centered in the pivot point   .In particular, the second term at the right hand side of (7) represents a lowpass filtered version of the elementary source field, which avoids singularity of the (approximated) secondary field for  =   [13].This approximation for the total field is based on the assumption that expression (7) can be considered not only on Γ, where it is actually enforced via the solution of the "design equation" ( 4), but also in the whole domain  since the target, regardless of its geometrical and electromagnetic nature, is forced to react like a point scatterer.Without going into further details on the nature of the "virtual experiments," we only recall that this circumstance arises when the contrast source  =    is mostly focused in a neighborhood of the pivot point [20].
Once the (virtual) data equation has been linearized through field approximation (7), one can take profit from the singular value decomposition (SVD) to solve the following linear problem: where V is the matrix-operator which rules the scattering phenomena for a given set of incident fields and contrast function .A regularized solution of problem ( 8) can be International Journal of Antennas and Propagation  found by means of truncated SVD (TSVD) scheme, which yields to the following explicit solution: wherein y  and w  are the right and left singular vectors of V, respectively, and   are its singular values.The truncation index  can be easily set using the Picard's plot as suggested in [13].

Numerical Analysis: Range of Validity
In this section we apply the method outlined in the previous section to achieve quantitative imaging of circular cylinders having different radius and dielectric permittivity; see Table 1.The total field over the measurement circumference Γ has been simulated by means of a full-wave forward solver based on the finite element method (COMSOL Multiphysics) and the incident field has been subtracted from it to extract the scattered field data.In Figure 1, the layout of the simulation setup for a free space homogeneous background is sketched.
In particular, the target has been placed at the center of the imaging domain.Moreover, in order to collect in a nonredundant way as much information as possible [21], we take  = 2Re[  ]  / √ 2, wherein   is the side of the (square) investigated domain which is assumed to embed the target.The data have been generated considering an equal number ( = ) of transmitters and receivers evenly spaced on a measurement circumference Γ of radius  ≃ 2  .In all the cases, according to the Richmond rule [22], we have considered   = 42×42 number of cells for the discretization of the imaging domain.To evaluate the imaging results we have considered the reconstruction error defined as where  is the actual contrast profile and χ is the estimated one.In particular, the outcome of the imaging procedure has been quantitatively appraised by means of the absolute difference Δ err between the reconstruction error achieved with field approximation (7) and the reconstruction error obtained with the solution of the "ideal" data equation, that is, when ( 8) is solved assuming that the total field is exactly known.It is worth noting that although the processing of the ideal data equation is meaningless in any practical instance, it provides in our analysis an immediate benchmark for the method, being the "best result" achievable from regularized solution of the data equation.As a result, considering Δ err allows comparing the reconstruction error introduced by approximation ( 7) with that achieved in the best possible  case [13].Figure 2 shows Δ err for different values of electrical dimension and refractive index of the circular cylinder.
It is possible to see that, for increasing values of the refractive index √  and electrical dimensions of the scatterer, the accuracy of the reconstruction becomes worse.In particular, for the considered values of the diameter and permittivity, we have experienced that the reconstruction error for the ideal data equation keeps almost constant (≤12.5%)while the error achieved in the case of field approximation (7) increases with the electrical dimension and refractive index of the scatterer.It is also worth noting that its slope rate is larger for increasing dimension of the scatterer rather than for increasing refractive index.
In order to quantitatively appraise the validity of the method, we have considered the following index: with   and  being the wavenumber and diameter of the homogeneous lossless scatterer, respectively.Note that for lossy scatterers the index can be generalized by considering |  |.In Table 1 the values of  for the considered circular scatterers are reported.In particular, we have experienced that the largest index still corresponding to a satisfactory reconstruction (err ≤ 20%) is almost 7.70.We consider this value as limit for the method's applicability.Interestingly, if we convert and compare this value with the limit found in [12] to analyze the validity of the Born approximation (√  ), our index results to be 1.22  , which is more than three times larger than the one in the case of the BA (0.35  ).As a result, the range of the approximation provided by the method at hand is much larger than the standard scattering linearization based on the BA.As already stated in the previous sections, this is due to the peculiarity of the virtual experiments to take implicitly into account the nature of the scattering system and then to provide an effective "scatterer-oriented" field approximation.
In order to corroborate the above analysis, we have addressed another example dealt with a scatterer of different shape.In particular, we have considered a lossless kite scatterer described by the following parametric equation: ( = 0.025,  = 0.01625,  = 0.0375) with a leading dimension of about 1  ; see Figure 3  result is in agreement with the validation analysis as method's performances become critic when  ≃ 7.70.Accordingly, the reconstruction errors are 21% and 53%, respectively.

Numerical Analysis: Comparison with the Born Approximation
In this section we report some numerical examples concerned with nonhomogeneous scatterers in order to compare the performances of the imaging approach at hand with those pursued by means of the BA.The actual permittivity profiles of the considered examples are shown in Figure 4. Also in this case the imaging results have been appraised by means of (10).For these examples the validity index  cannot be adopted as defined in (11), being the scatterers nonhomogeneous.As a result, for these scatterers the index has been evaluated considering an effective permittivity to take into account the nonhomogeneous nature of the scatterer,  Σ being the number of the pixels belonging to the target's support.
In the first example the unknown profile is a stepwise constant nonhomogeneous rectangular scatterer with permittivity  = 2 ÷ 1.4 and conductivity  = 50 ÷ 10 mS/m, Figure 4(a).A number of transmitters and receivers equal to 21 and   = 42 × 42 cells for the discretization of the imaging domain are considered.
The second example deals with a target made of two nonconcentric lossless circles with different permittivity values ( = 1.8 ÷ 1.4); see Figure 4(b).In this case, we have considered  =  = 21 and   = 42 × 42 cells.
In the last example two identical squares and one circular shaped scatterers with  = 2,  = 1.4,and  = 1.6, respectively, and the same conductivity  = 0.1 S/m are considered; see Figure 4(c).The number of transmitters and receivers is equal to 26 and the imaged domain is discretized into 50 × 50 cells.
In The comparison between the inversion results clearly shows that the approximation based on the virtual experiments confirms its good capability in appraising the dielectric nature of nonweak scattering systems as far as both the real and the imaginary part of the recovered contrast profile are concerned.As a matter of fact the reconstruction errors are 20%, 10%, and 38% and 67%, 64%, and 65% for the virtual experiments based approximation and BA, respectively.It is worth to underline that, although in the third example the reconstruction error is rather large, in any case the reconstruction error for the "ideal" data equation achieves 36%.This is due to the peculiar features of the scattering system made of three very close disjoint scatterers which are very small in terms of the probing wavelength.

Conclusion
In this paper we have investigated the range of validity of a new linear approach for microwave imaging which relies on the emerging framework of the "virtual scattering experiments."The imaging results clearly show that the proposed method can be applied in a range much wider than the standard first order Born approximation while not increasing either the complexity of the overall solution procedure or its computational burden.To this aim, we have addressed a numerical analysis in order to show the capability of the method in tackling quantitative imaging of nonweak targets.Interestingly, the method can be applied also in aspect limited configurations wherein the design of the virtual experiments is expected to be more challenging [13,15].Such a numerical study can be useful to improve the reliability of the method in developing effective microwave imaging diagnostic tools.

Figure 1 :
Figure 1: Sketch of the measurement setup used in COMSOL to simulate the total field data on the measurement curve Γ.

Figure 2 :
Figure 2: Validation analysis: (a) plot of Δ err for the parameters' values reported inTable 1 and (b) cut of Δ err for different values of   .

Figure 3 :
Figure 3: Kite example: (a) layout of the reference profile; (b) support estimation via LSM indicator and the selected pivot points (  = 17) and (c) real part of the retrieved contrast profile for   = 1.5;(d) support estimation via LSM indicator (  = 17) and (e) real part of the retrieved contrast profile for   = 2.The truncation index in the TSVD inversion scheme is  = 73 and  = 52, respectively.
(a).Moreover, we have considered two different values for the target's permittivity; that is,   = 1.5 and   = 2.The first step of the imaging strategy allows recovering in both the cases the support of the scatterer and choosing equispaced pivot points to devise the virtual experiments; see Figures3(b)-3(d).On the other hand, the quantitative reconstruction is satisfactory only for the kite with the smaller refractive index.Indeed for the latter  = 7.40, while  = 7.85 for the kite with   = 2

Figure 5 :
Figure 5: From top to bottom: imaging results for the three examples.From left to right: real and imaginary part of the retrieved contrast profile by means of the proposed approach ( = 65,  = 62, and  = 132) and the BA ( = 76,  = 78, and  = 64).
Figures 4(d)-4(f) the indicator Υ and the pivot points chosen to devise the virtual experiments are depicted for the three different examples.For each example, in Figure 5 are shown the outcome of the TSVD based inversion by means of the proposed approach and the reconstruction performed via BA.
Table 1 and (b) cut of Δ err for different values of   .

Table 1 :
Values of  for different dimensions and wavenumber of the cylindrical circular scatterers.The values corresponding to an unsatisfactory reconstruction error are boldfaced.