Reconstruction of three-dimensional targets using frequency-diversity data

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers. L’archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d’enseignement et de recherche français ou étrangers, des laboratoires publics ou privés. Reconstruction of three-dimensional targets using frequency-diversity data Ting Zhang, Patrick C. Chaumet, Anne Sentenac, Kamal Belkebir


I. INTRODUCTION
In many applications, such as geophysical probing 1 or medical microwave imaging, 2 one aims at estimating the properties of unknown objects (basically their permittivity distribution) from their response (scattered fields) to a known electromagnetic illumination.][5] Unfortunately, in many cases, the minimization solver is trapped in a local minimum of the cost functional, leading to poor or even misleading estimations of the targets.0][11] First, the multifrequency data, by supplementing the spatial ones, reduce the appearance of artifacts due to the limited view configuration. 11econd, the good convergence of the minimization procedure at low frequency, when multiple scattering is negligible, can be used to stabilize outside a local minimum the solution obtained at higher frequency (when multiple scattering is present). 9,12The challenge of multi-frequency inversion is thus to combine the high resolution theoretically accessible at high frequency with the stability of the solution obtained at low-frequency and the reduction of limited view artifacts.
Most inversion techniques dealing with multiple frequency data use the frequency hopping procedure. 9,12It consists in inverting one single frequency dataset at a time, starting from the lowest one, by using, as initial guess, the final solution obtained with the previous single-frequency dataset.This approach has been thoroughly checked on two-dimensional synthetic and experimental data and was recently applied to three-dimensional experimental data.Despite its overall good performances, especially when used in conjunction with regularization, [13][14][15] the frequency hopping procedure may be ineffective or even harmful in configurations where the scattered field noise is not uniformly distributed over the frequencies.For example, when the target is buried in a cluttered medium, 8 the clutter can significantly deteriorate the inversion solution at some specific frequencies related to its spatial characteristics.Using this solution as initial guess for the subsequent hopping procedure, can corrupt dramatically the final result.To avoid this undesirable effect, reconstruction approaches using the full-frequency data simultaneously seem preferable. 16n this paper, we study the performances of a 3D inversion scheme involving the full-frequency fully-polarized data, hereafter named Multiple-Frequency inversion method (MF), and compare it to the Frequency Hopping procedure (FH).Both reconstruction algorithms are checked on the experimental three-dimensional microwave data of the Fresnel database. 17Contrary to most inversion schemes that have addressed this canonical problem, [13][14][15] our approach is conducted without any regularization procedure or a priori informations as the aim is to see the peformances of the mutli frequency procedure wihout be blurred with added informations.

II. STATEMENT OF THE PROBLEM
In this section, we describe the imaging configuration and the method used for simulating the field scattered by the estimation of the targets.The general geometry of the scattering problem under study is illustrated in Fig. 1.An unknown three-dimensional object is assumed to be entirely confined in a bounded box Ω ∈ R 3 which defines the investigating domain.The surrounding of the target is assumed to be homogeneous and non-magnetic, with relative permittivity ε b = ε 0 and permeability µ = µ 0 (ε 0 and µ 0 being the permittivity and permeability of vacuum, respectively).The investigating domain is illuminated successively by l = 1, . . ., L electromagnetic excitation E inc l, p , for p = 1, . . ., P operating frequencies f p .For each excitation l and operating frequency p the scattered fields f mes l, p are measured on a surface Γ located outside the investigating domain.The inverse problem is stated as finding the relative permittivity distribution ε(r) of the unknown object from the knowledge of the scattered data f mes l, p .The total field inside Ω and the scattered field observed by the receivers can be represented as, where χ(r ′ ) denotes the permittivity contrast which is assumed to be frequency independent, i.e. the constitutive material of both targets under test and the background medium are non-dispersive Equation ( 1), noted as the state or coupling equation also sometimes referred as the near-field equation, permits to calculate the total field E l, p at the position r of the investigating domain Ω.Kernels G Ω p and G Γ p , related to the Green tensor of the homogeneous background medium, describe the field relationship from L 3 (Ω) to L 3 (Ω) and from L 3 (Ω) to L 3 (Γ), respectively.Once the total field inside the scattering domain is obtained, the scattered field E sca l, p is calculated from the observation equation, Eq. (2).The forward scattering problem is solved numerically using the coupled dipole method (CDM) which was introduced by Purcell and Pennypacker in 1973. 18

III. INVERSE SCATTERING PROBLEM
In this section we sketch the hybrid inversion technique that has been used to retrieve the target relative permittivity from a single frequency dataset 19,20 and we show how it can be adapted to deal with multi-frequency data. 21IG. 1. Geometry of the scattering problem.The scatterers are confined in the scattering domain Ω.The scattered field is measured on the measurement surface Γ.

A. Single-frequency and multi-sources inversion procedure
We consider the imaging of targets supporting moderate multiple scattering.In other words, the field inside the investigation domain can not be approximated by the incident field.The inversion algorithms able to tackle this issue can be classified in two kinds: the linearized methods in which the total field at each iteration step is fixed to the value obtained with the previous estimation of the permittivity distribution [22][23][24] and the non-linearized methods in which the total field is minimized together with the sought permittivity. 25,26In this work, we consider the Hybrid method (HM) which combines the rapidity of the linearized approach and the robustness to noise of the non-linearized ones. 20t each iteration step n the permittivity distribution and the total field inside the investigating domain are updated by minimizing the cost functional where normalizing coefficients are defined as The subscripts Ω and Γ are included in the norm ∥.∥ and later in the inner product ⟨., .⟩ to indicate the domain of integration.h (1) l, n and h (2) l, n are residual errors computed from the near-field equation Eq. ( 1) and the observation equation Eq. ( 2), respectively.
More details on this inversion technique can be found in Ref. 7.

B. Multiple-frequency and multi-sources inversion procedure
When multiple-frequency data are accessible, a popular strategy consists in using the frequency hopping method. 9,27The inverse problem is solved sequentially on single frequency dataset from the low to the high frequencies.At each sequence, the initial estimate of the contrast distribution is given by the final result obtained at the previous sequence.For the lowest frequency inversion the initial estimate can be obtained thanks to the back-propagation procedure.Thus, P iterative single frequency inversion are achieved corresponding to the P frequency harmonic data.
Another possibility is to derive and minimize a cost functional involving the entire data.This approach, hereafter named the multiple-frequency inversion method, can be used as an alternative strategy to the frequency hopping procedure.For the sake of simplicity, the near-field and the observation equations Eqs. ( 1)-( 2) are rewritten using symbolic notation, E sca l, p = G Γ p χE l, p .where p = 1, . . ., P corresponds to the frequencies.In the multiple-frequency inversion the minimized cost functional at each iteration step is defined as where normalizing coefficients are given by Hence all the frequencies are involved in the cost function.The functions h (1) l, p, n and h (2) l, p, n are two residual errors defined as previously in Eqs. ( 6) and ( 7), i.e., the first one is the residual error with respect to the incident field in the investigating domain computed from Eq. ( 1) and the second residual error is the error on the scattered field computed from Eq. (2).
Two sequences related to the contrast and total field inside the investigating domain, χ n and E l, p, n , respectively, are built up according to the following recursive relations where ν l, p, n , ω l, p, n and d n are updating directions with respect to the total field E l, p, n and the contrast χ n , respectively, and κ l, p, n , β n are scalar coefficients.The updating directions ν l, p, n and d n are chosen as the standard Polak-Ribière conjugate-gradient directions, 28 while ω l, p, n is given by where Ẽl, p, n−1 represents the total field inside the investigating domain Ω, calculated from Eq. ( 1) with contrast χ n−1 .The scalar weight κ l, p, n and β n are chosen at each iteration step n so as to minimize the normalized cost functional mentioned in Eq. ( 9) We propose to use the a priori information that the real and imaginary parts of the relative complex permittivity are greater than unity and non-negative, respectively.Instead of retrieving a complex function χ n , two real auxiliary functions ξ n and η n are reconstructed such that the recursive relation with respect to contrast χ n Eq. ( 14) becomes As updating directions d n;ξ and d n;η , the authors take where g ξ and g η are the gradients of the cost functional F (ξ, η, E •, •, n ) with respect to ξ and η respectively, evaluated at the (n − 1)-th step, assuming that the total field inside the investigating domain does not change.These gradients are given by where the over bar denotes the complex conjugate, and (G Ω p ) † and (G Γ p ) † are the adjoint operators of G Ω p and G Γ p , respectively.The search directions ν l, p, n for the total field inside the test domain is similar to those chosen for the contrast functions ξ and η: where g l, p, n;E l, p is the gradient of the cost functional F (ξ, η, E l, p, n ) with respect to the field E l, p , evaluated at the (n − 1)-th step, assuming that ξ and η does not change: For the multiple-frequency measured data, we have L × P updating directions ω l, p, n and ν l, p, n .The cost function F n is a nonlinear expression with respect to 2 × L × P complex unknown (κ l, p, n;ν , κ l, p, n;ω ) and two real unknown ( β n;ξ , β n;η ).We now compare the performance of the frequency hopping procedure and the multiple-frequency inversion scheme on the experimental microwave data of the Fresnel database.

IV. EXPERIMENTAL RESULTS
We considered several targets supporting more or less multiple scattering and two imaging configurations, one being almost complete (i.e. the illumination and observation are performed all around the target) and the other one being significantly incomplete (the illumination and observation are performed on one side of the target only).

A. Geometry of the targets under test
Four different targets picked up from the Fresnel database are considered, see Figs.(c) same as (a) but the centers are located at (a/2, −a/2, 3a/2) and (a/2, 3a/2, 3a/2).(d) Two dielectric spheres in contact of relative permittivity ε = 2.6, radius r = 2.5 cm with centers located at (−r, 0, 0) and (r, 0, 0). are referred to target or object (A), (B), (C) and (D), respectively.In all the reconstructions, the targets are assumed to be confined in a bounded box Ω of volume size 12.5 × 12.5 × 12.5 cm 3 .The mesh size of the discretization of Ω is taken equal to d = 0.5 cm whatever the operating frequency.At the central frequency f 0 = 5.5 GHz, the lattice size is about λ 0 /10 which is small enough to provide an accurate solution of the forward scattering problem.Adapting the mesh size to the frequency would be an interesting option for diminishing the calculation time of the FH but is expected to have a marginal influence on the final reconstructions.In all the examples, sole the estimated real part of the permittivity is discussed, the imaginary part being always very small.To quantify the accuracy of the reconstruction, we define a contrast reconstruction error as, where χ actual is the permittivity contrast of the actual target while χ rec is the reconstructed one.
For both the FH and MF procedures, the iteration stops automatically when the variation of the cost functional residue is less than 0.01%, which means that the cost function has reached a plateau.

B. Imaging configuration
The experimental setup, depicted in Fig. 3, is described in detail in Ref. 29.The incident wave, radiated by a parabolic antenna, can be modeled by a linearly polarized plane wave propagating in the (x, y) plane with θ inc ranging from 0 • to 350 • by step 10 • (green line).The illumination frequency varies from 3 GHz up to 8 GHz, by step of 0.25 GHz.The detection angle θ sca varies from 0 • up to 320 • by step 40 • while φ sca varies from 30 • up to 150 • by step 15 • (red line).
In the complete configuration, the vectorial scattered fields, defined as the difference between the total field and the incident field, 30 are obtained at 81 different positions on a sphere surrounding the target for 36 angles of incidence regularly distributed along 2π in the x, y plane.In the incomplete FIG. 3. Sketch of the experimental setup. 17The illumination is performed in the (x, y) plane with θ inc varying from 0 • to 350 • step 10 • (green line).The observation angles θ diff varies from 0 • to 320  configuration, we consider a subset of these data corresponding to a single illumination with θ inc = 0 • and a detection in transmission at angles 90 • < θ sca < 270 • .The total noise of the experimental dataset in the complete and incomplete configuration is esti- ∥h (2) l, p ∥ 2 Γ where the simulated scattered field is obtained with the actual shape of the targets, see Table I.It is observed that the data are corrupted in the same way whatever the targets and configurations.A more precise analysis, not shown, indicates that the noise increases with the frequency and affects more the cross-polarized terms than the co-polarized ones. 17

Complete configuration
In this section are reported the reconstructions of the four targets using the frequency hopping (FH) approach and the multiple-frequency (MF) inversion method using the complete data.The reconstructions with both FH and MF techniques of the targets (A,B,C) are displayed in Fig. 4, Fig. 5, Fig. 6 respectively.Both techniques yield satisfactory results, the Frequency Hopping approach yielding slightly more artifacts.
The target (D), see Fig.

Incomplete configuration
We now consider the performances of FH and MF in the limited view configuration.In this case, there is only one illumination (instead of 81 in the complete configuration).
The reconstructions of targets (A,B,C) using both FH and MF are shown in Fig. 8, Fig. 9, Fig. 10 respectively.In all these examples, the MF reconstructions are significantly better than the FH ones.The latter are deteriorated by artifacts and over estimations of the permittivity contrast while the former remain close to that obtained in the complete configuration.
On the other hand, the reconstruction of the large target D failed with both the FH and MF methods (not shown).In this case, whatever the methods, the small number of data is not sufficient for estimating properly the large number of 'active' (i.e.non zero) target unknowns.Such a problem cannot be solved without any regularization.

D. Weighted multiple-frequency inversion procedure
Basically, in multiple frequency imaging, the low frequency data ensures the stability of the reconstruction but provides low resolution images while the high frequency data brings high resolution but may prevent the convergence of the algorithm.In this paragraph we try to ameliorate the Multiple Frequency scheme by balancing the contributions of low and high frequencies.In presence of strong multiple scattering (for objects large with respect to the central frequency), it may be advantageous to increase the weight of the low frequency data for ameliorating the reconstruction.The layout of the Weighted Multiple Frequency (WMF) inversion method is the same as that of the MF AIP Advances 4, 127135 (2014) FIG. 9. Same as Fig. 5 in the incomplete data configuration.FIG. 10.Same as Fig. 6, while under the incomplete data configuration.method, sole the cost functional is modified and reads, 16 Two normalizing coefficients are defined as Gradients with respect to ξ, η and E l, p read as gξ = −2Re We have tried different α (from one to three) with the effect of continuously reducing the role of the high frequency data in the cost functional.For targets (A), (B) and (C) the reconstructions were already satisfactory with α = 1 (due to the small size of the objects) and were made a smoother with increasing α.On the other hand, for target (D) which is large compared to the wavelength, we found that taking α = 3 ameliorated significantly the quality of the reconstruction with a reconstruction error Err χ =21%.Compared to Fig. 7, with Err χ =44%), the reconstructed permittivity distribution, Fig. 11 is almost homogeneous and the target shape is now retrieved.This example shows the possibility given by the Multiple Frequency method and the interest of weighting differently the high and low frequency data.

V. DISCUSSION AND CONCLUSION
To summarize our results, we report in Table II the error on the reconstructions, Err χ , for all the targets, imaging configurations and inversion techniques.In addition, we indicate the reconstruction errors obtained with a single set of data obtained at the highest frequency 8 GHz.Hereafter, the single frequency reconstruction is called SF.We empharize that all reconstructed results shown here are achieved without any regularization procedure.We observe that, whatever the targets and the imaging configurations, the frequency hopping reconstructions are not better (and even worse) than the reconstructions obtained with the highest frequency dataset except for the target D of larger size.This results points out the difficulty in using the frequency hopping procedure.A single corrupted dataset at an intermediate frequency being sufficient for trapping the technique in a local minimum, the FH is quite vulnerable to noise.On the other hand, the MF procedure is always better than the SF or FH techniques.It yields remarkably similar reconstructions when changing from the complete to the incomplete configuration.We recall that the incomplete configuration is performed with only one illumination, 36 observations (in transmission), and 21 frequencies, i. e. 81 times fewer data than in the complete configuration and the time required for the inversion is divided by thirty.This result stresses the robustness of MF to noise and its interest when dealing with objects that are larger than the central illumination wavelength.Finally, note that the computation times of the MH and FH are about the same (the MF time is two times larger than the FH one in the complete configuration but is smaller than the FH one in the incomplete configuration).However, the MF method requires more memory than the FH as the measured data at all frequencies are minimized together.In view of this study, we believe that accounting for all the frequencies simultaneously in the inversion procedure should be generally preferred to the frequency hopping technique when dealing with multiple frequency data.
Another advantage of the multifrequency approach is the possibility to weight differently the data at high and low frequencies.For targets large compared to the central wavelength, decreasing the weight of the high frequencies ameliorates significantly the quality of the reconstruction.This important point can provide acceptable reconstructions even in an incomplete configuration and without any regularization.We believe that for the next step, including the regularization procedure [13][14][15] into this multifrequency approach should be considerable for some more complicated configurations i.e., metallic targets of special shapes, in order to enhance the obtained reconstruction.This work has been financed by the ANR SURMITO.

FIG. 7 .
FIG. 7. Reconstructed permittivity of target (D), see Fig. 2 (d), under the complete configuration, in the (x, y) plane (a) and (c) for z = 0 cm, (b) and (d) for z = −2 cm.(a) and (b) frequency-hopping method.(c) and (d) multiple-frequency method.The multiple-frequency method reconstructs better the contact point between the spheres.

TABLE I .
Err r quantifies the noise corrupting the experimental data for each target and imaging configuration.

TABLE II .
Err χ for for different targets and for different inversion methods, under the complete configuration and the incomplete configuration.SFC: single frequency at 8 GHz under the complete configuration.SFIC: single frequency at 8 GHz under the incomplete configuration.
FHC: frequency-hopping approach under the complete configuration.MFC: multiple-frequency approach under the complete configuration.FHIC: frequency-hopping approach under the incomplete configuration.MFIC: multiple-frequency approach under the incomplete configuration.