Abstract
Many inverse scattering problems concern the detection and localisation of point-like scatterers which are sparsely enclosed within a prescribed investigation domain. Therefore, it looks like a good option to tackle the problem by applying reconstruction methods that are properly tailored for such a type of scatterers or that naturally enforce sparsity in the reconstructions. Accordingly, in this paper we compare the time reversal-MUSIC and the compressed sensing. The study develops through numerical examples and focuses on the role of noise in data and mutual coupling between the scatterers.
1. Introduction
The problem of detecting and localising point-like scatterers (i.e., small in terms of wavelength) is a classical inverse scattering problem which is relevant in a number of radar imaging applications. For instance, in GPR imaging it can occur in concrete pillar diagnostics in order to detect and localise the embedded rebar. This application is of high relevance in diagnostic operations not only to assess the structural stability of concrete structures [1, 2] but also to inspect them after disaster events (e.g., earthquakes) [3]. Another very important field of application concerns early breast cancer detection where the tumour actually can be considered to be a small scatterer [4].
In order to address the imaging problem for such a kind of scatterers, different methods can be taken from literature. Among them we mention the migration algorithms [5], the distributional approach [6], the time reversal [7], and the MUSIC [8]. In particular, the latter three methods have been properly conceived for point-like scatterers. More recently, further methods based on the compressed sensing theory have been exploited as well [2]. These methods basically enforce sparsity in the reconstructions and therefore appear suitable to deal with point-like scatterers, which represent a sparse scattering scene.
In this paper, we focus on the time reversal-MUSIC (TRM) and the compressed sensing (CS) approach [9] in order to get a comparison between a method that takes into account the scatterers’ features (TRM) and one which is founded on a sparse optimisation paradigm (CS) [2].
In the framework of radar imaging applications, CS approach is addressed in [10]. The comparison between CS and MUSIC algorithm is shown, for example, in DOA estimation [11–13], using multiple-input multiple-output (MIMO) radar in [14]. The hybrid method, called compressive-MUSIC, is introduced in [15] where the Multiple Measurement Vector (MMV) problem is addressed.
The electromagnetic scattering is modelled in terms of a nonlinear operator because of the mutual scattering occurring between the scatterers. Accordingly, tackling the inversion through the CS approach entails dealing with a nonlinear inverse scattering problem [9] which is fraught with reliability drawbacks (due to the false solutions) and burdened with computational complexity. Therefore, in order to avoid the drawbacks affecting optimisation algorithms implementing nonlinear inversions, here we choose to model the scattering phenomenon by a linear model whilst the CS approach is used. In particular, the linear model is established as follows. First, mutual scattering among different scatterers is neglected. Then, assuming that the scatterers’ sizes are small as compared to the probing wavelength, they are represented by delta distributions just supported over the scatterers’ centres. This is the so-called distributional model introduced in [6]. Of course, the corresponding reconstruction algorithm (being linear) offers relevant advantages in terms of the computational effectiveness and reliability (no local minima occur) but the role of the neglected mutual scattering must be studied.
As to the TRM algorithm, it works by projecting the so-called steering vector on the “noise subspace.” Now, the point is that such a noise subspace is not affected by mutual scattering. Therefore, TRM does not require a linear approximation for the scattering.
The analysis is developed for a 2D scalar scattering configuration. In particular, for the sake of simplicity, a homogeneous background medium is considered. Layered background can be easily addressed provided to exploit the pertinent Green function. The considered scenario, however, captures well the most important features against which we intend to compare the two methods. In particular, we are mainly interested in comparing TRM and CS against mutual scattering and noise. To this end, an ensemble of numerical experiments will be presented.
The plan of the paper is the following. In Sections 2 and 3, the TRM and CS inversion algorithms are briefly recalled along with the pertinent scattering models. Section 4 is devoted to the presentation of the numerical examples. Conclusions end the paper.
2. TR-MUSIC
We refer to a scalar 2D geometry (invariance occurs along the -axis) depicted in Figure 1. An ensemble of cylindrical scatterers, with cross sections much lower than the relevant wavelength, is probed by uniform filamentary current sources (directed along as well) whereas the scattered field is collected all around the scene. The point targets are located at positions , , within a background medium characterised by an outgoing scalar Green function , denoting the temporal frequency. source points, taken uniformly at , , over a circle enclosing the scattering scene and observation points at with are considered. We assume that transmitting and receiving positions coincide. However, for each transmitting antenna the scattered field is collected over all the receiving antennas so that a multistatic/multiview configuration is addressed. The collected scattered fields are arranged to form the so-called multistatic data matrix (MDM) whose th column represents the scattered field vector when the th source is firing. Accordingly,where represents the scattered field collected over all and is the column vector of the exciting currents.
The multistatic data matrix can be factorised as follows: represents the transpose operator, is the propagator matrix from the scatterers to the observation points whose th column is , and is the matrix that takes into account interactions between the scatterers.
The matrix can be expressed aswhere is diagonal matrix accounting for the objects scattering coefficients and is an matrix, with zero diagonal entries, which takes into account the propagation paths between the different objects.
The TRM algorithm is based on the singular value decomposition (SVD) of which consists of those vectors that verify : and are the th left and right singular vectors, respectively, and is the corresponding th singular value. The rank of is and, accordingly, possess singular values, , , different from zero. Moreover, the corresponding right singular functions, , , span the range of ; that is, . Hence, the dimension of the range of is equal to number of scatterers . The “noise” subspace is spanned by the right singular functions corresponding to singular values equal to zero; that is, . Clearly, it is orthogonal to so that . At the same time, the range of is spanned by the columns of (), so that any is orthogonal to . Accordingly, the positions of the point-like scattering objects can be detected by defining the following functional (pseudospectrum) as ranges within the investigated area [16] at prescribed discrete points:where is the Hermitian operator. The objects positions are then identified where the pseudospectrum reaches its maxima.
As is well known, in presence of noise, it is difficult to determine the exact number of targets and hence the noise subspace. This is because signal and noise singular values are no more very well separated (under ideal cases, noise singular values are theoretically equal to zero). A wrong estimation of the noise subspace leads to a detriment of the achievable performance. In order to cope with this problem, the number of targets (or equivalently, the number of signal singular values) is identified by exploiting the “Akaike Information Criterion” (AIC) [17]. AIC basically measures the relative cross entropy between the signal and noise singular values so that the number of signal singular values is identified by the minimum of such a measure.
3. Compressive Sensing
In order to present the CS approach, as we mentioned in the Introduction, a linear model for the scattering is first established. To this end, first multiple scattering between the scatterers is neglected. Accordingly, the scattered field reads aswhere is the scattering coefficient of the th object. However, the obtained model is not yet linear as the unknown positions appear embodied within the Green functions. On introducing the distributional unknown representation, that is, (which is the reflectivity function of the overall scenario), it yieldswhere now the unknown cylinders’ positions are represented as the support of distributions. Inverting (6) for is actually a linear inverse problem.
The discretised counterpart of (6) is given bywhere is the matrix that represents the discretised linear operator in (6) whose elements are and . Assuming that the investigation domain is pixelated in points, the matrix actually is a linear operator acting between and .
At this juncture, it is interesting to note that different ways data are exploited by the two algorithms. Indeed, while the CS exploits all the scattered field data (that actually ), the TRM works within as the different views of illumination are basically used to restore the rank of .
In order to invert (7), we consider a constrained optimisation method that exploits the a priori information about the sparse nature of the solution (i.e., a solution that contains the least number of nonzero elements) as a regularising tool. More in detail, the general regularisation problem can be cast aswith being the regularisation parameter.
According to the CS approach, optimisation in (8) for , is inherently convex and can be easily implemented by fast algorithms [9, 18] by addressing the equivalent problem
This procedure is known as the basis pursuit denoising (BPDN) problem. Analogously, it is possible to follow the LASSO (least absolute shrinkage and selection operator) approach [19]
Note that the sparseness of the reconstructed image is strictly controlled by the parameters and . For example, the lower the parameter the sparser the reconstruction, and vice versa. Therefore, the choice of is a crucial question that strongly affects reconstruction quality.
It is worth remarking that while in TRM the crucial point to address before imaging was the selection of the signal subspace dimension, the counterpart for CS is the choice of the regularisation parameters or (depending on the adopted minimisation strategy). To cope with this problem, an -curve method is used [20].
4. Numerical Results
The goal of this section is to compare the achievable performance by the two approaches at hand. To this end, a single frequency multistatic/multiview configuration is considered with the scattered field data collected over sources/observation positions deployed uniformly over a circle of radius ( being the wavelength) enclosing the scattering scene, as shown in Figure 2. Metallic point-like objects are considered as scatterers, which mimics, for example, metallic rebar embedded in a masonry structure. Moreover, as the antennas are very close to the pillar, for the sake of simplicity, the air-pillar interface is neglected and a homogeneous background medium with (i.e., resembling the one of concrete) is considered.
For each test case, we report the pseudospectrum (4) and the reflectivity (8) functions, respectively.
4.1. Setting TRM and CS
As we have anticipated in the previous sections, both the TRM and CS require defining some parameters (or a strategy to select them) before the imaging stage can be actually run. In particular, for TRM the noise subspace needs to be identified whereas for CS the regularisation parameter ( or ) should be fixed. As to the TRM, we have chosen to adopt the Akaike Criterion to select the most relevant singular values; instead, an L-curve approach is chosen to fix the regularisation parameter in the CS inversions.
In this first set of numerical experiments, we just check the appropriateness of such choices. To this end, a scattering scenario consisting of targets located at , , , , and is considered. Data are corrupted by an additive white Gaussian noise with a signal-to-noise ratio (SNR) equal to 20 dB.
Figures 3–6 refer to TRM. As expected (see Figure 3(a)), in absence of noise the singular values experience an abrupt decay which allows clearly discerning between the signal and the noise singular values. Actually, five scatterers can be clearly understood. When noise corrupts data (Figure 3(b)) the mentioned abrupt transition is very much smoothed so as to make it difficult to identify the noise subspace. In this case, for example, the visual inspection of the singular value behaviour can suggest thresholding the singular values at . This would lead to a detriment of the achievable resolution as is shown in Figure 4 where the target located at is missed at all.
(a)
(b)
We now switch to considering the AIC. The behaviour of AIC is shown in Figure 5. As can be seen, the minimum is got at , which is in correspondence with the actual number of scatterers. In this case, having succeeded in identifying the noise subspace, the corresponding reconstruction allows detecting and locating all the scatterers in the scene (see Figure 6).
Let us now turn to considering the CS. As said above, the L-curve criterion is used to determine the optimum value of and . In Figures 7(a) and 8 the curves in log-log scale are plotted as and vary within and with a step of 0.05, respectively. In both cases, those curves are L shaped and hence can be used to select the regularisation parameter. However, the BPDN criterion allows exploring a wider range of the L-curve with a narrower range of (cf. Figures 7(a) and 8). For this reason, in the following section, in order to decrease the computational cost, BPDN criterion is preferred. The rest of Figure 7 shows how the choice of the regularisation parameters is crucial. Actually, three reconstructions corresponding to three different values of (reported in Figure 7(a)) are shown. These figures allow appreciating the effectiveness of the strategy for choosing the regularisation parameter.
(a)
(b)
(c)
(d)
Reconstruction results reported in Figures 6 and 7 also show that both methods are effective in locating point-like scatterers when the targets are well resolved by the configuration, that is, when the distance from the target is larger than .
4.2. Reconstruction Results
Equipped with the criteria for the selection of the noise subspace for TRM and the regularisation parameter for CS, we now turn to addressing some reconstruction examples in order to more deeply compare the methods. We are particularly interested in scattering configurations where the targets are closely located to each other as for well-resolved scatterers we have already shown that both methods work well.
First of all, a glimpse of two-point resolution is given as this is often used to characterise the performance of an imaging method. So two targets are considered to be located along the -axis and at a distance of and , the latter being below the classical diffraction limits. Furthermore, data are corrupted by additive zero mean white Gaussian noise so as to have a dB and dB. Finally, both linear (LD) and nonlinear (i.e., accounting multiple scattering between the targets) data (NLD) are exploited. The corresponding reconstructions are shown in Figures 9 and 10.
As can be seen from Figure 9, TRM succeeds in locating both scatterers in all the considered cases. Also the effect of noise is clear since it entails a smoothing of the pseudospectrum when the SNR lowers. Note also that TRM reconstructions have been obtained by employing NLD. Reconstructions with LD are very much the same and have been omitted. Accordingly, it can be concluded that TRM is not affected by multiple scattering but only by the noise.
Figure 10 reports the reconstructions obtained by CS. In this case, both LD and NLD are considered. From such a figure, first, how sparse optimisation works can be appreciated as most of the reconstructions are set to zero. Moreover, it can be deduced that CS reconstructions are not too much affected by the level of noise. This is clear by comparing lines of different colours that refer to the same case (i.e., the same distance between the scatterers and the same kind of data, LD or NLD). Instead, reconstructions are very much affected by the mutual scattering (cf. dotted and solid lines). Of course, this is particularly true as the scatterers are more closely located (because multiple scattering becomes more relevant): in this case the two scatterers are not discerned. Therefore, it can be concluded that CS is more sensitive to model error than noise (note that this is just the opposite as for TRM) and the achievable two-point resolution is worse than the one that TRM allows obtaining.
Let us now pass to consider the case of more scatterers which on average are located more far apart than the previous case (for which we already know that CS fails). Scattering configurations consisting of three and four scatterers are considered. The three targets are located at , , and . The four-scatterer layout is obtained by adding a further target at .
The reconstructions obtained by TRM are reported in Figure 11. As can be seen, in the case of three scatterers TRM works in locating the scatterers even though (as expected) when the SNR lowers (for the case at hand the SNR reduces from 20 dB to 10 dB) the pseudospectrum smoothes so that its peaks have a lower dynamic as compared to the reconstruction noise pedestal. Note that the degradation of the achievable performance due to the noise is the same as shown in the previous figures where it can be more clearly appreciated. When the number of scatterers increases to four TRM fails to resolve the scatterers’ locations even for dB. This outcome is perfectly consistent with the results shown in [21]. Indeed, in such a paper it was shown that, for a 2D problem, TRM works as long as the number of scatterers enclosed within a circle of radius does not exceed , the latter being the so-called number of degrees of freedom (NDF) of the radiation operator defined on such a circle. Basically, the NDF corresponds to the maximum signal space dimension in the case of noisy data. Therefore, if the number of scatterers is higher than the NDF even a small amount of noise leads TRM to fail. This is what happened while passing from three to four scatterers. Actually, this example has been selected because it allows clearly highlighting TRM faults. Furthermore, it suggests inspection if CS can work better. The reconstructions obtained by CS are reported in Figure 12. As can be appreciated, in the case of linear data, CS works fairly well. Instead, when NLD are exploited CS does not work. Indeed, while for the three scatterers a delocalisation occurs for the case of four scatterers one of the scatterers is completely missed and a spurious target appears. Once again it can be concluded that while TRM performance is mainly limited by noise, CS instead is strongly affected by multiple scattering. Accordingly, in order to improve TRM performance a denoising scheme can be run before imaging. Instead, mutual scattering cannot be in general filtered out. As a limit case we consider noiseless NLD. In this circumstance, TRM works perfectly whereas CS outcome does not significantly change (see Figure 13).
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
(a)
(b)
5. Conclusions
In this paper the problem of detecting point-like targets has been addressed. The TRM and the CS approaches have been compared for a scalar 2D geometry where the scatterers were embedded within a homogeneous host medium. The study has been developed through some numerical results.
The outcome of the study is the following:(i)for well-resolved scatterers both the methods are effective;(ii)for closely located scatterers the TRM performance is limited by the available SNR whereas CS is strongly affected by multiple scattering between the targets;(iii)if the SNR is relatively high TRM works also for closely located scatterers and outperforms CS;(iv)TRM is faster than CS. This is because CS relies on two nested iterative procedures, one to build the L-curve and one, for each regularisation parameter, to solve the constrained optimisation problem. More in detail, by using a standard laptop PC, the examples presented in the paper required a fraction of second for TRM while CS runs at best in 5 minutes and at worst in 14 minutes.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This work was supported by the Italian Ministry of University and Research through the FIRB initiative under the project MICENEA (RBFR12A7CD) and the COST Action TD1301 MiMed.