Abstract
An approach to reconstruct buried objects is proposed. It is based on the integral equations of the electromagnetic inverse scattering problem, written in terms of the Green’s function for half-space geometries. The full nonlinearity of the problem is exploited in order to inspect strong scatterers. After discretization of the continuous model, the resulting equations are solved in a regularization sense by means of a two-step inexact Newton algorithm. The capabilities and limitations of the method are evaluated by means of some numerical simulations.
1. Introduction
Subsurface imaging is an important subject in several applicative areas, including seismic and geophysical prospecting, nondestructive, testing and medical diagnostics [1–39]. Electromagnetic techniques are widely applied to face this problem. Ground penetrating radar [40] is the basic instrumentation able to retrieve discontinuities in the lower subspace (the ground). Processing of GPR data can allow an improvement in the detection and localization capabilities of the system. In some cases, neural network-based approaches have been adopted [41]. Moreover, inverse-scattering-based techniques have been proposed, too. By using these methods, the measured values of the field scattered by the discontinuities present in the lower half space are “inverted” in order to obtain the spatial distributions of dielectric parameters (e.g., the dielectric permittivity and the electric conductivity) inside a fixed inspection domain. Approaches of this kind have been proposed, for example, for the detection of water leaking from pipes [42]. If an accurate model of the buried scatterers is available, it would be possible to retrieve information about the buried objects in terms of their dielectric parameters.
The present authors developed in [43] an inverse-scattering-based method for buried object imaging under the second-order Born approximation (SOBA) [21, 44–48]. By using the SOBA, the reconstruction of buried objects has been obtained with a better accuracy with respect to linearized approximations. In particular, due to the nonlinear nature of the scattering equations, in which the electric field inside the inspected region is an unknown quantity as well as the spatial distributions of the dielectric properties of the medium, the SOBA allows to approximately retrieve the field distribution, too.
However, in order to detect strong discontinuities, the exact nonlinear inverse scattering formulation must be considered. This approach has been faced in [49, 50], in which strong scatterers have been reconstructed in free space under tomographic imaging conditions. The approach proposed in [49], which is based on the use an inexact Newton method, has been experimentally tested in [51]. The same method has been successfully applied in [52] within the contrast source formulation (under free-space conditions).
In the present paper, the approach developed in [43] for buried object detection under the SOBA is extended to treat strong scatterers. In this case, the exact nonlinear inverse scattering formulation, in terms of the half-space Green’s function, is taken into account. In particular, starting from the integral equations for the scattered field for cylindrical objects, slides of the spatial distributions of the dielectric parameters are reconstructed by numerically solving these equations in a regularization sense.
The paper is organized as follows. In Section 2 the mathematical formulation of the approach is reported. In Section 3, the application of the inexact Newton method is described. Numerical results are reported in Section 4, and conclusions are drawn in Section 5.
2. Mathematical Formulation
Let us consider a cylindrical object buried in a homogeneous-half-space medium (Figure 1). A 2D transverse-magnetic (TM) configuration is assumed. The dielectric properties of the buried object are assumed to be independent of the coordinate. The transmitting and receiving antennas are located in the transverse plane and generate a -polarized electric field in both the half spaces (this field is independent of the axial coordinate, too).
The dielectric characteristics of the investigation domain S, in the lower half space, are described by the so-called object function, which is defined as where is the complex dielectric permittivity at point r, and denotes the dielectric permittivity of the soil.
The investigation area is illuminated by line-current sources, with angular frequency ω (in the following, the time dependence is assumed and omitted). A mixed illumination configuration (partially borehole and partially from the top), which allows direct comparisons with the results reported in [43], is exploited. Under such assumptions, also the scattered electric field is -polarized and independent of . Consequently, only the transverse plane can be considered for the imaging procedure.
The scattered electric field, , , is collected for any incident wave in measurement points , which constitute the so-called observation domain . The dielectric properties of can be related to the scattered field by means of the following set of nonlinear integral equations: where , and , are the data and state operators, where is the half-space Green’s function [43] (superscripts in the operators recall where they are calculated, whereas subscripts recall where they are integrated).
Equations (2) can be rewritten in a compact operator form as where is an array containing the values of the object function and of the total electric field inside the investigation area for any illumination. In (3), is an array containing the values of the scattered electric field acquired by the measurement probes in the upper half space and the values of the incident electric field inside the investigation area. The half-space operator can be expressed, according to (2), as As it is well known, the inverse problem (4) is non-linear and ill-posed. Consequently, nonlinear regularizing methods must be used.
3. Regularized Solution of (3)
Due to the mentioned ill posedness, (3) must be solved in a regularized sense. The inexact Newton regularizing scheme developed in [49, 50] is applied here. This method includes two nested loops, in which a linearization (with respect to the current solution guess ) of the operator is constructed (outer loop) and solved by means of the truncated Landweber method (inner loop). Exactly as in [49], the algorithm performs the following steps.(1) For (k denotes the outer iteration index), a starting guess is chosen (initialization phase).(2) The following linearized version of equation is constructed where , , and is the Fréchet derivative of at point . A regularized solution to (5) is obtained by applying an iterative truncated Landweber algorithm, for which we have where is the adjoint operator of .(3) The current solution is updated by using the following relation: (4) A predefined stopping rule is finally applied. If this criterion is not satisfied, a fixed maximum number of outer iterations, , is performed. In both cases, at any iteration, the procedure can be stopped, otherwise steps – are repeated.
It should be noted that step is actually the inner loop of the algorithm, whereas steps – constitute the outer loop.
4. Numerical Results
The developed method has been tested by means of numerical simulations. The working frequency has been set equal to 300 MHz (i.e., similar to that used in ground penetrating radars). A homogeneous lossy soil characterized by a relative dielectric permittivity and an electric conductivity has been considered (in order to simulate a dry sandy soil [53]). Similar to [43], a mixed measurement configuration is assumed, where reflection and cross-borehole operating modes are mixed together in order to acquire as much information as possible. Correspondingly, the antennas are located at points (Cartesian coordinates) where is the free-space wavelength in the upper medium. A subset of these antennas act both as transmitter (TX) and receiver (RX). When an antenna operates in TX mode, all the remaining ones collect the scattered electric field (i.e., they work in RX mode). Consequently, for every illumination, field samples are collected. In particular, in the cases reported in the following, antennas are used. They are located at positions denoted by the indexes given by .
The investigation area is a square domain of dimensions , whose upper side coincides with the air-ground interface. In order to numerically solve the equations involved in the inverse-scattering problem, such domain has been discretized into square subdomains.
The scattered field values used as input data for the considered simulations have been numerically obtained by using a computational code based on the method of moments [54]. In order to avoid the so-called inverse crime, a finer discretization (such as the ratio of the cell dimensions in the direct and inverse problems is not an integer number) has been used for the computation of the synthetic data. Moreover, in order to simulate a more realistic scenario, a Gaussian noise with zero mean value and variance corresponding to a signal-to-noise ratio has also been added to the computed data.
In all the simulations reported in this section, the initialization phase has been always performed by assuming an empty investigation domain in the lower half space. Consequently, in this phase, the object function is equal to zero, and the total electric field is assumed equal to the one transmitted in the lower half space.
The first reported case concerns the reconstruction of a single circular cylinder of radius centered at . The dielectric properties of the cylinder are and . The parameters of the inexact Newton algorithm are the following: maximum number of inner iterations, ; maximum number of outer iterations, . Figure 2 reports the behavior of two error parameters for different values of the number of outer iterations of the iterative procedure. In particular, the two parameters are defined as with and being the actual and reconstructed permittivity distributions (at the current iteration) and and the measured and computed values of the electric field data (at the current iteration, too). As can be seen from Figure 2, the plot of the “residual” on the data (the error parameter ) decreases monotonically as expected. Moreover, the plot of the solution “residual” (the error parameter ) shows a typical semiconvergence behavior, that is, after a certain value of the iteration number the reconstruction quality degrades [55, 56].
The final image of the distribution of the relative dielectric permittivity inside the investigation region (corresponding to the optimum iteration number) is shown in Figure 3. As can be seen, the method is able to correctly retrieve the dielectric profile of the investigation area. In particular, both the shape of the buried inclusion and the value of its dielectric permittivity are correctly identified.
For the present case, an analysis of the effects of the noise has also been performed. In particular, the mean relative error on the reconstruction, defined as the -norm has been calculated. In (10), and denote the values of the original and reconstructed complex permittivity inside the th discretization subdomain of the investigation region. Moreover, the optimal value of the outer iteration have been evaluated, too. The obtained results are reported in Table 1. As expected, as the noise level increases, the reconstruction “quality” decreases. However, even in presence of strong noise levels, the method is able to yield rather good reconstructions (e.g., the error parameter is equal to about 8% error for ).
In the second set of simulations, the reconstruction of a more complex configuration has been considered. In particular, the presence of a double-layer circular cylinder has been assumed. The outer layer has radius and is characterized by a relative dielectric permittivity and an electric conductivity S/m, whereas the inner one has a radius and is characterized by a relative dielectric permittivity and an electric conductivity S/m. The object is centered at . The parameters of the inexact Newton method are the same as the ones adopted in the first case. The reconstructed distribution of the relative dielectric permittivity (final image) is shown in Figure 4. As can be seen, the object is correctly located, and the values of the dielectric properties are identified with a quite good accuracy.
Finally, the reconstruction of two separate circular cylinders have been considered. The radii of the two cylinders are and . They are buried with centers located at points and . The dielectric properties of the two scatterers are the following: , , and S/m. The signal-to-noise ratio is . The parameters of the iterative inversion algorithm are the same as in the previous cases. The reconstructed distribution of the relative dielectric permittivity is shown in Figure 5. This image corresponds to the optimal solution. As can be seen, also in this quite complex case, the method is able to correctly reconstruct the investigated scene. In particular, the two objects are clearly identified and separated. Moreover, both shapes and dielectric properties of the two targets are estimated with quite good accuracies (although for the second cylinder the dielectric permittivity is slightly underestimated). In order to better assess the reconstruction accuracy, Figure 6 gives the profiles of the retrieved distribution of Figure 5 along some cuts parallel to the coordinate axes.
(a)
(b)
(c)
(d)
5. Conclusion
In this paper, a previously proposed approach for the imaging of buried objects has been extended to inspect strong scatterers, for which the exact nonlinear equations of the inverse scattering problem must be considered. The developed method is based on a two-step inexact Newton method, which allows us to obtain a regularized solution of the (discretized) integral equations relating the samples of scattered electric field collected in the upper half space to the distributions of the dielectric parameters inside the test domain located in the ground. The reported preliminary results of numerical simulations seem to indicate the validity of the proposed approach, which needs further evaluation assessments, including an experimental validation, which will be the subject of future research activity in this area.