Abstract
We proposed an improved solution to two problems. The first problem is caused by the sidelobe of the dominant scatterer masking a weak scatterer. The proposed solution is to suppress the dominant scatterer by modeling its electromagnetic effects as a secondary source or “extra dependent transmitter” in the measurement domain. The suppression of the domain scatterer reveals the presence of the weak scatterer based on exploitation of multipath effects. The second problem is linearizing the mathematical forward model in the measurement domain. Improving the quantity of the prediction, including multipath scattering effects (neglected under the Born approximation), allows us to solve the inverse problem. The multiple bounce (multipath) scattering effect is the interaction of more than one target in the scene. Modeling reflections from one target towards another as a transmitting dipole will add the multiple scattering effects to the scattering field and permit us to solve a linear inverse problem without sophisticated solutions of a nonlinear matrix in the forward model. Simulation results are presented to validate the concept.
1. Introduction
Radio frequency (RF) tomography is a mathematical process of 3D image reconstruction from a measurement using a multistatic distribution of transmitters and receivers. The geometric diversity of these elements increases the information in the measurements. The process of determining the permittivity and conductivity profile in the measurement domain, and, therefore, the shape of the target, from the scattered field measurements, is an inverse problem. To solve this problem, under conventional methods such as the Born approximation, we use the principles of linear scattering to determine a linear relationship between measured returns and target shape. The Born approximation is valid if the scatterer is small and does not interact strongly with other objects. However, strong scatterers within the domain may generate sidelobes masking weaker returns. This masking, in conjunction with multipath effects, may result in loss of features and subsequent failure to identify a target.
In this research, a novel method is proposed to increase overall image quality and extend the capabilities of RF tomography by modeling the strong scatterers in the measurement domain as dipoles that behave as secondary sources (transmitters). Unlike conventional methods, the dipole model reduces the effects of the sidelobes from the strong scatterers and exploits the multipath of multiple targets or complex shapes. The multipath phenomena contain more information about the targets permitting illumination in the shadowed region and an increase to the radar aperture length as shown in Figure 1. The electromagnetic characteristics for each modeled dipole are estimated by representing the cells in the measurement domain’s image. The eigenvalue and eigenvector from each cell represent the phase and magnitude for the modeled dipole and also the spatial orientation of the target. The process of modeling large scatterers as dipoles can be iterated, addressing one strong scatterer at a time. This method effectively suppresses the sidelobes and exploits the multipath within the measurement domain. Using the Born approximation, the linear relationship between the scattered fields and the target is updated for simplicity. With iterations, the “extra” dipole will account for the multipath effects, thus removing some limitations caused by the Born approximation. This concept has been successfully demonstrated in software (FEKO© by Altair). In addition, this work also presents an innovative conversion using a back-projection algorithm for multipath effects and modeling of an “additional” source or transmitter in the measurement domain. The result of implementing this method of modeling strong scatterers as dipoles successfully demonstrated an increase in the resolution and enhanced radar imagery. This research is focused on radar imaging applications. The general concept of radar imaging, such as multistatic synthetic aperture radar/inverse synthetic aperture radar (SAR/ISAR) can be viewed as radio frequency (RF) tomography. The theoretical analysis of multistatic radar imaging techniques shares the same framework. The common feature of multistatic radar imaging techniques is a wide distribution of transmitters and receivers. This provides more information about the shape and edges of the target. The formation of RF tomography images can be made within far-field or near-field regions. This imaging process utilizes detailed electromagnetic models of the measurement domain to estimate target echoes. The entire measurement domain is represented by a 3D “contrast” function, which has a nontrivial value at the boundaries between the target and the surrounding free space. Radar tomography has been proven to work in many different applications, from ground penetrating radar (GPR) to building penetration [1–3], to ISAR and even in passive modalities [4–7].
RF tomography has an inherent weakness with strong scatterers that mask a nearby weak scatterer, thereby decreasing overall image quality. This is due primarily to the creation of sidelobes by the strong scatterer. There are many current techniques to ameliorate the strong scatterer echoes resulting in an enhanced image. One method is simply to remove the dominant scatterer by subtracting its point spread function signal from the received signal [8, 9]. The most common algorithm, using this method to suppress the dominant scatterer is the CLEAN algorithm. Moreover, the adaptive removal of modeled strong echoes has been improved by using multilook techniques to improve the effectiveness of strong echoes removal [10]. Another technique models the dominant scatterer to describe uncertainty and estimate its causes and impact on the results [11]. Modeling the dominant interfering target (or all targets) as dipoles provides more information about the shape and edges of the target as presented in [12, 13]. The estimation of the dipole polarization using the sensitivities of eigenvalues and eigenvectors has been discussed in the literature [14]. However, these algorithms are linear models that only work in far-field regions.
In this research, an algorithm is proposed based on a unique electromagnetic model of the measurement domain. This approach overcomes the limitations of existing linear models. The contrast function used to represent the measurement domain consists of dipole scatterers of various strengths and polarization. The dominant scatterers in the measurement domain are viewed as secondary transmitters that exploit the effect of the multipath caused by more than one target or a complex shape target. As a result, solving the inverse problem of the secondary transmitter is achieved by suppressing the sidelobe of the dominant scatterers. The dominant dipoles are identified and represented as eigenvalue/eigenvector pairs to provide the characteristic electromagnetic behavior. Then their scattered field is estimated, and the effect of this field on the weak scatterers is suppressed.
Multiple or multipath scattering occurs when the targets reflect energy off other objects before being scattered back to the receiver. Multipath radar imaging (SAR, ISAR, and RF tomography) suffers from the effects of this multiple path scattering. In conventional radar imaging, the reflections from multiple scattering centers interfere with each other, leading to a loss in resolution. Also, the energy of the multipath reflection may appear at locations where no target exists, resulting in unwanted objects, a.k.a. “ghosts” on the image [15].
Another phenomenon is “masking,” caused by a strong sidelobe of the dominant scatterer within the scene. The strong sidelobe dominates and essentially masks the scattering from the lesser objects. Thus, scenes that have a single target are inherently better suited to produce accurate images and is in concert with the direct-path reflection assumption under the Born approximation. The CLEAN algorithm is used to suppress strong sidelobes within a scene of interest.
Our analysis applies to near-field and far-field regions for the multipath problem. Characterizing the dominant scatterer as a dipole gives the same response of near-field scattering as well as far-field. The problem with the near-field is that the point spread function of the target is varying due viewing angle. The operating regions (near-field and far-field) around the electromagnetic field can change the resolution of the tomographic image. The theoretical calculation of the resolution on far-field region is limited to a wavelength between 0.3536–0.5λ. In the near-field region, the resolution of the object function is smaller in half space where it is reported to be less than 0.1768λ [16]. Furthermore, any image resolution in far-field less than 0.25λ and 0.1λ for near-field will be considered as super-resolution phenomenon together with extracting weak scatterer features in the reconstructing object function. The next chapter explains resolution terms and the criteria for resolution definition.
2. RF Tomography Forward Model
Radar tomography views a target from many positions, due to the spatial distribution of transmitters and receivers. This geometric diversity increases the information obtained from the measurement domain. Multiple () dipole transmitters and () dipole receivers surround a region containing the targets to be imaged. Transmitter has location and polarization , and receiver has location and polarization . Figure 2 illustrates this configuration, showing the measurement domain and a single transmitter/receiver pair. At any given time, a transmitter radiates a known waveform while the remaining transmitters are left inactive. The distributed receivers collect the echoes from the targets. The position of each receiver is known up to a fraction of the wavelength, and all the receivers share the same reference signal [1]. The information obtained from each receiver/transmitter pair is stored after mitigating noise (and clutter).
Using this geometry, a simple forward and inverse model for the time-harmonic electric field can be produced. Note that time-domain waveforms can be expressed in time-harmonic fields using stretch processing or simply the Fourier transform. The 3D geometric view of the model is shown in Figure 2. The investigation domain of the target is assumed to be a homogeneous medium with estimated dielectric permittivity , conductivity , and permeability . The transmitters and receivers are located in free space. The measurement domain has an unknown dielectric permittivity and the conductivity profiles, where is a separate position vector.
The electric field response of each cell (the divided area of the measurement domain defined as pixel) in the measurement domain can be calculated by Green’s function [17, 18]. An accurate calculation of Green’s function will provide an improved image after the process. Green’s function determines the electric field response in free space of any point source at any chosen location [19]. In this research, a bistatic transmitter and receiver scenario was considered. In this scenario an isotropic source located at radiates an electromagnetic vector field in all directions into free space. This electromagnetic field is defined as the incident field, denoted by , with respect to the target. The region of interest for which an image is to be constructed is discretized into pixels. For a target located at , the electric scattering field is calculated theoretically based on the transmitted scattering field from the source. The scattered field is recorded by a receiver at , where is Green’s function that calculates the electric field response at distance . Assuming that the targets do not have any depolarization effects on the incident field, then and have the same polarization. Also, individual targets are assumed to be isotropic. These two common assumptions make it possible to simplify the vector wave equation into a scalar equation without changing its meaning, by using Born approximation to linearize the inverse problem of the relationship between the target contrast function and the scattering field as follows:
Both the Born approximation and the wavelength of the operation frequency limit the RF tomographic image. In the literature, many approximations were introduced to overcome these limitations. The iterative Born method has been proposed to provide a numerical solution of multiple inverses and forward problems in the time domain [20]. However, the convergence and the computational load remain issues for this numerical solution. The Distorted Born iterative method accelerated the convergence of the iterative Born method by updating Green’s function during each iteration [21]. Born approximation was extended where Habashy et al. proposed a more accurate algorithm in low-frequency cases [22]. Later, for more accurate results higher-order extended Born approximation was presented by Cui et al. [23]. In the proposed method of this study, the limitations of the Born approximation are removed by estimating the multipath effects and suppression of the sidelobe of the dominant scatterer. The Born approximation is linearly related to the contrast function, and (1) can be represented as a matrix multiplication. In this case, the relationship can be expressed aswhere multiplication by matrix implements the forward model. The operator has to be inverted to compute the unknown contrast function. Because is not square and is usually ill-conditioned, this is a difficult task. However, in the literature, many approaches are provided for inverting the linear operator [24–26]. One commonly used method is a simple inversion using the Hermitian matrix of the linear operator , which is effective for dealing with ill-conditioned operators.
The contrast function is obtained from the multiplication of the Hermitian linear operator as an inverse problem solution. The linear operator acts as a filter on the spatial harmonics of the reconstructed contrast function. This image is a vector image where each cell could be represented as a dyadic function of three dimensions, expressed as . The eigenvectors and eigenvalues for each cell are obtained, such thatIn this equation, is the cell’s vector representation. After the contrast function is obtained, its corresponding eigenvalues and eigenvectors are generated. These provide the required information to estimate the polarization and magnitude of the dominant scatterers. There are three eigenvalues for each cell. The largest eigenvalue contains the most important information about the polarization and the magnitude of the dominant scatterer. After detecting the center of first maximum scatterer field’s locations at where represented the order of dominant scatterers detecting in the measurement domain, the scattering field response is estimated for each dominant scatterer. The scattering field of the new transmitter illuminating the measurement domain is also updated asNow, the multipath caused by the dominant scatterer is considered in Green’s function calculation. The polarization of the dominant scatterer is not usually as the transmitters and receivers, so it can take any orientation in the measurement domain. The dominant scatterer electric field can be either calculated theoretically or simulated as forward scattering problem. The relationship between the contrast function and the scattering field in terms of the multipath interaction effects remains linear in (4) without considering the effects of filtering processing, but the accuracy of the mathematical model under the Born approximation has been improved since the multiple scattering effect is considered in the linear operator updated in the final stage.
The matrix form of the contrast function of the measurement domain was described as the inverse problem. A linear relationship under Born approximation is applied to obtain the contrast function from the electric field response of the transmitters and receivers at each location in the measurement domain. The contrast function of the measurement domain is represented as cells where each cell represents a contrast value at a specific location . The measurement domain has transmitters and receivers. The contrast function, representing a tomographic image, divides the measurement domain into individual pixel cells (th pixel). The scattered field is represented as a column vector of []. By applying the matrix form of (1),where is the location of the th pixel (cell) of number of pixels (cells) in the measurement domain. Equation (5) is the forward model of the linear relationship between the contrast function of the measurement domain and the measured scattering field under Born approximation. Figure 3 shows the orientation of the transmitters and the receivers can be expressed as [12]
The matrix form of (6) is represented as
The scattering field vector is , the linear operation is matrix, and the contrast function vector is . is the scattering electric field collected at position. is denoting the unknown dielectric permittivity element of the contrast vector at location as function of the target permittivity profile in permittivity background medium to create an image at final stage. is representing the contribution of a single transmitter and receiver to the pixel value; is calculated electric field response from the transmitter to a single point location at the measurement domain. The calculated electric field of a single receiver to a point location is denoted by . The sum of the multiplications of these quantities will be equivalent to the scattering electric field received by the transmitters.
Since (the multiplication of Green’s function values) is not a square matrix and ill-posed, a Hermitian transpose matrix can be employed to obtain the contrast vector as inverse problem solution in vector notation
Other inversion techniques are presented to deal with ill-posed condition matrix.
The multiplication of by in this case performs the convolution of the transmitted waveform within the contrast function. Figure 4 shows the first domain scatterer location identified and modeled as extra transmitter. The new transmitter in the measurement domain is estimated as delayed version of the original transmitted signal.where denoted the electromagnetics modeling of first dominant scatterer. The iterative process to locate all number of dominant scatterers takes place by updating the operator. Each step will add one dominant scatterer illuminating in the domain as new transmitter on the measurement domain. After locating all dominant scatterers in the measurement domain the operator will be asThe matrix of the scattering field vector is ; the relationship between the scattering filed and the contrast function of the measurement domain is kept linear by updating the operator which is a matrix of shown asThe final in (11) is representing the multiple scattering effects of dominant scatterers considered in the inverse problem under Born approximation.
3. Algorithm
The algorithm begins by forming a radio frequency tomographic image by inverting the linear operator under Born approximation to the scattered field generating the set of eigenvalues for each cell of the resulting reconstruction by (7). The maximum eigenvalue that is larger than a specified threshold is used to locate the dominant scatterer in the reconstruction array. This scatterer is then modeled as a dipole. Corresponding eigenanalysis of the dominant scatterer will provide orientation to model the dipole polarization [27]. This dipole can be treated as an extra transmitter in the measurement domain. The simulated field generated by this dipole acting alone is denoted as . This quantity is then processed to update the operator , while the linearity is kept under Born approximation. The extra modeled dipole scattered filed effect is added to vector. The forward model is then updated by the addition of the dipole, treating it as a new transmitter, leading to a modified operator by (11). The signal emitted by the dipole is a scaled and delayed version of the signal from the active transmitter. This leads to interference at the receivers, determined by the wavelength and path length differences involved. This updated operator is then used in the matrix approach inversion, applying it to the modified scattered field creating a reconstruction that has the dominant scatterer removed. As the block diagram shows in Figure 5, this process then repeats to identify and remove the next dominant scatterer from the previous set of transmitters, including the previously modeled dipoles. The process stops when there is no value more than the preset threshold or the noise floor value is reached.
This algorithm generates a series of reconstructions, and each subsequent reconstruction shows more of the weaker scatterers. This set of reconstructions is then combined into a single reconstruction of the measurement domain contrast function, showing both the strong and weak scatterers. Even after linearizing, the relationship to find the contrast object function is still a complex problem to invert. Moreover, the matrix operator is the mathematical model of the real world behavior where a good approximation of this model will lead to a good imaging system. Unfortunately, rectangular matrix operator is not possible to invert directly because it is ill-conditioned.
Conjugate Gradient is an algorithm used to solve a large system of linear equations based on the steepest descent method. The method is to simply represent , where is unknown vector (the object contrast function in this research), is a known vector (scattering field ), and is known (the matrix operator). Now, the steepest descent method starts with an initial guess. Then an iterative method seeks to optimize the direction until the minimum gradient is reached, while the gradient of the function is computed.
4. Simulations
Simulations using the computational electromagnetic software FEKO© package (Altair) were performed to obtain the scattered field for irregular complex target geometry scenario. The RF tomographic imaging principles were applied to reconstruct the contrast function of the measurement domain of this scenario to validate the proposed method. Preliminary simulation results, obtained by computing the diffracted fields with FEKO using Method of Moment (MoM), indicated that it is possible to reconstruct an image. The operating frequency is 8 GHz, and the measurement domain used is divided into an array of cells of dimension . The plane of the two-dimensional cut was simulated at . The simulation used 19 transmitters placed at the center of a 2 m radius circle within the measurement domain. 51 receivers were placed in a concentric circle with a radius of 0.4 m. A domain containing 12 cylinders, forming a complex target, was considered, as illustrated in -plan and 3D image in Figures 6 and 7, respectively. The cylinders were aligned with the -axis and are randomly spaced into an L-shape along the and directions. These cylinders have a radius of and a length of 0.05 m. A single weak scatterer is modeled as an additional cylinder, with a radius of , and located at coordinate (0.1 m, 0.112 m) in the -plane, as shown in Figure 7. The 12-cylinder complex target is larger as compared to the small single cylinder. Therefore, it acts as a dominant scatterer in the simulation and will mask the echo from the single weak scatterer. The forward model is applied to create simulated scattered field to begin the algorithm. The matrix operator consists of elements representing the predicted electric field response under Born approximation. The first application of the operator inversion leads to a reconstruction as shown in Figure 8. This corresponds to [], in which the weaker scatterer does not appear. After determining the cell location of the dominant scatterer, the eigenvalue and eigenvector analysis was carried out to determine its phase and magnitude. The algorithm then proceeded, as previously described, by detecting and replacing iteratively the 12 dominant scatterers, one at a time, with modeled dipoles that have the same electromagnetic effect. The power density distribution of the measurement domain in Figure 9 shows the reflected power of the PSF of the 12 cylinders where its sidelobe masks the weak scatterer. The weak scatterer is very thin wire with radius making its PSF very weak compared to the sidelobe of the dominant scatterers as shown for -axis cross range in Figure 10; another tomographic view from -axis shows no further information can be obtained whereas the weak scatterer is still buried under the sidelobe of the domain scatterer in Figure 11. The criterion is to measure the improvement of the resolution in extracting weak target features from very noisy measurement domain that is caused by multipath delay error. Exceeding the limits of super-resolution which conditioned the PSF of the dominant scatterer must be removed whereas the power of weak scatterer PSF must be preserved too to reach the super-resolution in radar imaging. The improvement in the image will be by maximizing SNR of weak scatterer. After suppressing all of the dominant scatterers in the complex target and adding their multiple scattering effects to the operator, the final image produced by the algorithm is shown in Figure 12. The weak scatterer that did not appear in Figure 8 is now clearly visible in Figure 12 after applying the proposed algorithm. Figure 13 shows the performance of the algorithm after the suppression of the dominant scatterers using the SNR criterion to determine the improvement in the resolution. The cross range figures showed the sidelobe of the dominant scatterer is gone and the scattering field of the weak scatterer can image the target over the measurement domain in Figures 14 and 15, respectively. The CLEAN algorithm is applied to compare the performance in near-field region. The general idea of CLEAN algorithm is to subtract the estimate of the dominant scatterer electromagnetics response from the image. The simulation was structured to form as image an irregular complex targets shape used before. The simulations were performed with the same parameters for the proposed algorithm results. After subtracting the PSF of 12 cylinders formed as complex target from the received scattered field, we reconstructed contrast function of the measurement domain as shown in Figure 16. The power distribution of the contrast function showed the weak scatterer PSF energy is subtracted too from the received scattering filed. Figure 17 shows the performance of the algorithm after applying the CLEAN algorithm. The cross ranges of the contrast function power density at - and -axis in Figures 18 and 19 showed the weak scatterer PSF is buried and attenuated.
4.1. Irregular Geometry SNR Analysis
The visibility of a weak scatterer is increased when we used MRL algorithm versus other approaches. To present the improvement in resolution, as mentioned in previous chapters, the resolution is improved via increased SNR. Here we selected an irregular geometry to analyze the SNR improvement obtained under the Born approximation, applying the CLEAN algorithm, and the MRL algorithm. The irregular geometry exhibits the highest noise (interference) as compared to other geometries. Table 1 presents the weak scatter signal, in both and , since the PSF of the dominant scatterer changes in azimuth in near-field region. The weak scatterer is located at 0.1 m in and 0.112 m in .
The SNR obtained via the MRL algorithm indicates a stable value of the intensity over the mean interference in both and . Moreover, the sidelobe level is reduced to maximize the SNR. The weak scatterer intensity under the Born approximation and CLEAN algorithm (in ) had a smaller SNR value than in because the sidelobe in is increased. Furthermore, the radiation pattern of the domain scatterer is masking the weak scatterer in . Comparing all three algorithms, we can say the MRL algorithm sufficiently increases the SNR in near-field region.
5. Extended Noisy Numerical Model
In this section we will extended the numerical model as the scatterers are extended. The simulation has been done using the same transmitters and receivers parameters of previous simulation at 8 GHz operating frequency. The measurement domain used is 50 × 50 cm divided into cells of size 0.1λ × 0.1λ. The two-dimensional tomographic image was simulated at . The domain has extended 12 cylinders placed in L-shape on the measurement domain creating complex target geometry. The extended 12 cylinders are equally spaced with distance . These cylinders have radius of and a length of 0.05 m. The weak scatterer has radius of (83% larger than the weak scatterer in previous numerical results) and located at coordinate (0.2 m, 0.223 m) in the -plane with length 0.1 m. The noise in the numerical model occurs from two sources. First source is the multipath effects (interaction between the scatterers); the second noise source was added to the scattering field vector with different mean and variance than the received scattering field. The initial tomographic image of contrast function is shown Figure 20. The measurement domain now is containing more power than the previous analysis. The extended complex targets are strong in some location and the MRL misses the cylinders buried under the noise, but after applying MRL algorithm to process the domain scatterers as dipoles, the weak scatterer is weakly revealed as shown in Figure 21. If the mean and the variance of the added noise are similar then the inversion algorithm diverges as sequence of noisy environments.
6. Conclusion
An algorithm to suppress the sidelobes due to dominant scatterers in tomographic radar images was proposed, examined, and studied. Accordingly, each dominant scatterer is replaced by an equivalent dipole source, and a new image is formed. The equivalents dipole is processed as extra transmitter in the measurement domain predicting the multipath effects, the multiple scattering between targets. The output no longer has dominant scatterers or their associated sidelobes. The results were positive showing an improved image with greater resolution of small objects and weak targets within the scene.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.