A SAP-DoA Method for the Localization of Two Buried Objects

A localization technique for buriedmetallic and dielectric objects is proposed and tested. An array of isotropic antennas investigates a scenario with cylindrical targets buried in a dielectric soil.The targets are in the near field of the array, and a Sub-Array Processing (SAP) approach is adopted: the array is partitioned into subarrays, and Direction of Arrival (DoA) algorithms are used to process the electromagnetic field received by each subarray and estimate the dominant arrival direction of the signal. By triangulating all the estimated DoAs, a crossing pattern is obtained. It is filtered by a Poisson-based procedure and subsequently elaborated by the k-means clustering method in order to distinguish between targets and background, estimate the number of targets, and find their position. Several simulations have been performed to compare different DoA algorithms and to test the localization method in the presence of two buried cylinders. Different values of the permittivity of the involved dielectric materials have been considered; the target positions and size have also been varied. The proposed procedure can be useful for ground-penetrating radar applications, near-surface probing, and for the detection and localization of defects in a hosting medium.


Introduction
Localization of buried objects is a very important problem in remote sensing, geophysical prospecting, landmines reclamation, civil engineering applications, and near-surface probing [1][2][3].This has inspired researchers worldwide in the development of techniques for the solution of the inverse electromagnetic scattering problem by buried obstacles.
In this paper, we consider the localization of multiple cylindrical targets embedded in a dielectric half-space.The targets are dielectric or perfectly conducting objects with circular section, and the problem is two-dimensional.Various approaches can be applied to solve this problem by using approximated techniques [4][5][6][7], or adopting microwaveholography algorithms [8], or else exploiting diffraction tomography techniques [9].Moreover, stochastic methods have been proposed, such as the Support Vector Machine (SVM) [10,11] and others [12][13][14].The use of neural networks has been also introduced, showing very good performances [15].The interest in innovative methods is still strong, and new techniques have been recently developed [16].In this context, we here extend the technique proposed in [17,18] for the detection of one cylinder, by introducing into the localization procedure a clustering method for the investigation of scenarios with multiple objects.
It is worth noting that our localization procedure collects the DoAs estimated by the subarrays considering them as semi-infinite lines, starting from the centre of each subarray and extending to the half-space containing the cylinders, so that, in the case of metallic objects the lines are virtually prolonged also through the targets, and the procedure focuses on the geometrical problem of calculating the interceptions by triangulating all the lines.Direction of Arrival (DoA) techniques [19,20] enable an antenna array to estimate the number of incident signals and their directions of arrival.Moreover, a Subarray Processing (SAP) [21] approach can be adopted for the detection of targets lying in the near field of an antenna array.In particular, the receiver array can be partitioned in several subarrays, such that the field scattered by the targets can be assumed to be locally planar at each subarray [17,18].Then, by applying DoA estimation algorithms, it is possible to predict International Journal of Antennas and Propagation the dominant direction of the incoming signal at each subarray.In this paper, we apply and compare several DoA algorithms: nonparametric (Bartlett [22], Capon [23], Linear Prediction [24], and Maximum Entropy [25]), subspace based (Minimum Norm [26], Pisarenko Harmonic Decomposition (PHD) [27], MUltiple SIgnal Classification (MU.SI.C.) [28], Root MU.SI.C. [29], and Estimation of Signal Parameters via Rotational Invariance Technique (ESPRIT) [30]), and maximum likelihood (ML) based methods (Deterministic and Stochastic ML [31,32] and Weighted Subspace Fitting (WSF) [33]).We did not consider Root MU.SI.C., ESPRIT, and WSF algorithms before [17,18], since these algorithms are particularly suitable for the localization of multiple objects.
By triangulating all the DoAs estimated by the subarrays, we obtain a pattern of crossings.It can be filtered by using a Poisson-based procedure [17,18] in order to remove a noisy background of unwanted crossings.Subsequently, we estimate the number of targets with the Data Spectroscopic [34] algorithm.Finally, we process the pattern by using the -means clustering algorithm [35] and estimate the object position.The localization procedure has been implemented in MATLAB.
The proposed approach might be employed for an arbitrary number of targets; however, in this paper we focus our attention on the simple case of two buried objects.
To obtain input data for our simulations, we solve the forward scattering problem by using an in-house software implementing the rigorous Cylindrical-Wave Approach (CWA) [36,37].
The paper is organized as follows.The proposed localization technique is resumed in Section 2. In Section 3, numerical results are presented.We consider different permittivity values for the buried cylinders as well as for the hosting medium.The positions of the targets and their size are also varied.Conclusions and ideas for future work are discussed in Section 4.

The Localization Procedure
In this section, we briefly resume the procedure that we adopted for the localization of multiple buried cylinders.
The geometry of the problem is sketched in Figure 1.A known electromagnetic field illuminates a region with  cylindrical targets embedded in a uniform, homogeneous, unlossy half-space, with an array of sensors adhering to the medium in which the scatterers are buried in.The obstacles are in the near field of the whole array; the th object has axis in (  ,   ) and radius   ( = 1, . . ., ).The array of receivers has spacing  and is composed by  subarrays of  sensors.
For each subarray, if () is the signal impinging on the th ( = 1, . . ., ) sensor at time , the sensor output can be written as follows: where   () is the angular response of the sensor and   denotes its position; moreover  is the wavenumber in the medium hosting the objects.If   = ( − 1) and the sensors are isotropic, the steering vector assumes the following expression: where [⋅]  indicates transposition and () is a Vandermonde-like vector.For a uniform linear array or sensors, () is uniquely defined if and only if the condition |cos | ≤  is fulfilled, that is always satisfied with a spacing  ≤ /2.
A matrix of steering vectors can be built, in the presence of  <  signal sources (the buried scatterers), with  ×  independent columns.In (2),  is the vacuum wavenumber associated with the plane wave propagating through the medium in which the object is buried.If this medium is different from the one in which the array is placed, the inversion procedure neglects the refraction by the planar interface separating the media.However, when the array is very close to the interface (as happens in several ground penetrating radar applications), the numerical error that affects the DoA estimation due to the assumption of a homogeneous medium is very small.The simplest way to deal with the effects of the interface is to correct the DoAs by using the Snell formula.
We assume that an additive, zero-mean, uncorrelated Gaussian noise is present, with variance  2  .Therefore, the  ×  Hermitian correlation matrix   of the array can be defined as where () = ()+() is the vector of the received signal at each antenna plus noise,   is the source correlation matrix,   =  2   is the  ×  noise correlation matrix, being  the identity matrix, (⋅)  symbolizes the Hermitian transposition, and [⋅] represents the expected value operator.
The exact statistics for noise and signals are not a priori known; however, assuming the process is ergodic, the spatial correlation can be approximated by a time-averaged correlation, so that where the superscript symbol ̂denotes an estimation of the true value.
The pseudospectrum [20] is a function that gives an indication of the angles of arrival based upon maxima versus angle.The mean output power of the array, at the time instant , is given by the square magnitude of the array output.
Assuming that the components of  can be modeled as zeromean stationary processes, for a given weight vector  the mean output power of the array is The pseudospectrum function is strongly depending (through ) on , and the particular choice of the weight vector differentiates the various DoA estimation algorithms [19,20].By applying a DoA algorithm to each subarray, a set of arrival directions is obtained.By triangulating all of the estimated DoAs, a cloud of intersections is generated (see Figure 1), and it can be statistically analysed to give an estimation of the object centres [17,18].
Two kinds of regions can be distinguished in the crossing pattern: regions including a high number of crosses, that most probably locate the targets, and a background in which the amount of crosses is quite small and spread.These different classes of crossings are modelled using suitable spatial Poisson distributions [38].In particular, the whole region scanned by the array is divided into a number of nonoverlapping windows of a fixed size, with side lengths that can be suitably fixed depending on the size of the sought object.The number of crossings in each window is counted.A testing procedure is employed, with a decision threshold depending on the desired false-alarm probability  fa , to extract the window set belonging to the scatterers and to remove the background markers [17,18].
In clustering analysis, the goal is to understand the macroscopic structure and relationships among the objects, by considering the ways in which they are similar and dissimilar.Different clustering techniques can be distinguished, by whether they carve up the objects into disjoint clusters at a single level (partitional clustering) or give a complete hierarchical description of how objects are similar to each other (hierarchical clustering), using a dendrogram [39].A simple partitional clustering technique is -means [39]; it is successfully used in many application relevant to bioengineering [40], circuit theory [41], and signal processing [42].Given a set of points in a multidimensional space, a set of  cluster centres is randomly chosen, and each point in the dataset is allocated to the nearest cluster centre.The centroid of each allocated set of points is computed, and these centroids become the new cluster centres.The process is repeated until the positions of cluster centres converge.Since  is an input parameter of the algorithm, the number of clusters has to be known or guessed beforehand.It is also possible to use an additional clustering method just to estimate  and then apply the -means algorithm.To this aim, we implemented the Data Spectroscopic [34] algorithm, obtaining an exact prediction of  for all the simulated cases.The complexity of the -means algorithm is represented by (), where  is the desired number of clusters,  is the number of iterations, and  is the number of data objects.Finally, for each target we define the localization error as where ( x , ẑ ) are the estimated axes coordinates of the cylinder.In this way, the localization error is the ratio between the Euclidean distance of the estimated centre from the real one and the actual cylinder radius.When err  ≤ 1, the estimated centre of the th object is collocated in a point inside its transverse section; otherwise, if err  > 1, the estimation of the centre falls outside the object.As a reference threshold, we assume that a target can be considered correctly localized if the error is less than one.

Numerical Results
A uniform linear array of 51 equally spaced elements, with spacing  =  0 /4, is considered.A partition in subarrays with  = 3 elements is applied (the total number of subarrays is therefore  = 17).Two circular cross-section cylinders, with radii  1 and  2 , are located at distances ℎ 1 and ℎ 2 from the horizontal axis of the array; the cylinders can be perfectly conducting or dielectric, with relative permittivity  1 and  2 , respectively.They are buried in a dielectric half-space International Journal of Antennas and Propagation with relative permittivity   .The scenario is illuminated by a plane wave with TM polarization (electric field parallel to the cylinder axis).The noise variance is  2  = 0.25; the false-alarm probability is  fa = 10 −6 for conductive cylinders and  fa = 10 −3 for dielectric cylinders.Different values are considered for the radii and the burial depths of the targets, and for the permittivity of the involved materials.In the SAP procedure, square windows with side length  0 are used,  0 being the vacuum wavelength.
As mentioned in the Introduction, the input data of our simulations are obtained by solving the forward-scattering problem with an in-house software implementing the rigorous Cylindrical Wave Approach (CWA) [36,37].This is a rigorous and fast spectral-domain method that solves the scattering problem of perfectly conducting or dielectric cylindrical objects with circular cross-section, buried in a dielectric half-space or in a multilevel medium, as a finitethickness slab [43,44].The CWA can be applied for arbitrary values of permittivity, size and position of the targets; it can also be employed to study the scattering of an incident pulsed wave, with a rather general time-domain shape [45].Recently, the CWA has been extended to take into account the presence of roughness in the interface between air and soil and of losses in the ground; these novelties represent a significant advancement of the method, since they allow to model the environment hosting the cylinders in a more realistic way [46].As shown in [36] for perfectly conducting cylinders and in [37] for dielectric cylinders, obstacles of general shape can be simulated through the CWA with good results, by using a suitable set of small circular-section cylinders.
It is worth noting that CWA results take into account all the multiple-reflection and diffraction phenomena, thus also the refraction of the scattered field by the planar interface between air and ground; the localization procedure, instead, neglects such refraction phenomenon.However, if the array is very close to the interface (as in several GPR applications), the numerical error that affects the DoA estimation turns out to be very small, when the assumption that the array and the targets are in a homogeneous medium is made.
In the case presented in Figure 2, the cylinders are perfectly conducting, and their radii are  1 = 0.5 0 and  2 =  0 ; the cylinder axes are in ( 1 ,  1 ) and ( 2 ,  2 ), with  1 = −5 0 ,  2 varying in the interval [−2.5 0 , 5 0 ] with a  0 /8 step,  1 =  0 √  +  1 , and  2 =  0 √  +  2 .The relative permittivity of the ground is   = 4.The objects are both in the near field of the array and in the far field of the subarrays.In the (a), (c), and (e) plots, the localization error for the first cylinder is plotted as function of 2 2 / 0 ; in the (b), (d), and (f) plots, the same is reported for the second cylinder.The effect of the correlation among the signals backscattered by the cylinders is apparent in the left-column plots, where the estimation of the first cylinder is quite corrupted.The second-cylinder error, instead, is less affected by the signal correlation; in fact, even when the second cylinder is placed in the left half of the whole array ( 2 < 0), its sensing is mainly due to subarrays in the right half of the array, that do not have the other cylinder standing in the line of sight.When the cylinders are far enough, and the effect of the correlation is smaller, the localization results are good for both the cylinders.Not all the DoA estimation methods perform the same, in particular PHD and the Maximum Entropy seem to be less robust while the remaining methods are quite reliable.
In Figure 3, the same as in Figure 2 is reported for air cavities ( 1 =  2 = 1).Again, the signal correlation mainly affects the first-cylinder sensing.
In Figure 4, the same as in Figures 2 and 3 is reported, for dielectric cylinders with  1 =  2 = 9.The localization of dielectric targets reveals to be less precise than the localization of perfectly conducting and air cylinders.Also in this case, the effect of the correlated field is stronger for  2 < 0. The estimated position of the larger and shifting cylinder is quite precise.
The contemporary presence of a cavity ( 1 = 1) and a dielectric cylinder ( 2 = 9) is considered in Figure 5.The radii are  1 =  2 =  0 ; the other parameters are the same as in the previous figures.The general behaviour is similar to that of the previous cases, even if in this case the objects have the same size.This confirms that in the previous cases the second cylinder was better localized due to its position rather than its larger size.

Conclusions
In this paper, a procedure for the localization of a finite number of buried objects is presented.In particular, a Subarray Processing method, employing various Direction of Arrival (DoA) algorithms, is used for the localization of both dielectric and perfectly conductive cylindrical scatterers, and of air cavities, in the proximity of a linear array of receiving antennas.In particular, the following DoA algorithms are considered: Bartlett, Capon, Linear Prediction, Maximum Entropy, Minimum Norm, Pisarenko Harmonic Decomposition, Multiple Signal Classification (MU.SI.C.), Root MU.SI.C., Estimation of Signal Parameters via Rotational Invariance Technique (ESPRIT), Deterministic and Stochastic Maximum Likelihood, and Weighted Subspace Fitting (WSF).By triangulating all the estimated DoA, for each scenario a pattern of intersections is obtained, that subsequently is filtered by using a Poisson-based procedure, in order to remove a noisy background of unwanted crossings.Then, the number of targets is estimated by using the Data Spectroscopic algorithm.Finally, the crossing pattern is processed by using the -means clustering algorithms, and the object positions are estimated.
To obtain the input data, the forward-scattering problem is solved by using an in-house software implementing the Cylindrical-Wave Approach (CWA).Moreover, the presence of an additive Gaussian noise is taken into account.
Several configurations are considered, varying the constituting material of the cylinders, their position and size.The localization outcomes show an effective detection capability and a quite precise localization estimation in many cases.The main advantage of this procedure is the computational simplicity (making this approach suitable for an early-stage localization), while the most evident drawback is the worsening  of the estimation as the distance between the objects becomes smaller, or when the relative permittivity of the ground is high.We look forward to check the localization procedure in some laboratory tests and also to improve the localization procedure by taking into account the losses in the ground.We also plan to consider stochastic approaches as the Support Vector Regression [11], in order to evaluate their effectiveness when applied in conjunction with our formulation, as was done in [47] for the localization of a perfectly conducting cylinder.

Figure 1 :
Figure 1: Geometry of the problem with a sketch of (a) the DoA estimated by the subarrays; (b) the pattern of DoA crossings.

Figure 2 :Figure 3 :Figure 4 :
Figure 2: Localization of two perfectly conducting cylinders buried in a medium with   = 4; the position of the second cylinder is varied and the localization error is plotted.(a) Bartlett, Capon, Linear Prediction, Maximum Entropy, Minimum Norm and PHD estimation of the axis position for the first cylinder; (b) same as in (a), for the second cylinder; (c) MU.SI.C., root MU.SI.C., and ESPRIT estimation of the axis position for the first cylinder; (d) same as in (c), for the second cylinder; (e) WSF, DML, and SML estimation of the axis position for the first cylinder; (f) same as in (e), for the second cylinder.

Figure 5 :
Figure 5: Localization of a cavity ( 1 = 1) and a dielectric cylinder ( 2 = 9), both having a radius equal to  0 .All the other parameters are the same as in the previous figures.