Near-Field 3D Sparse SAR Direct Imaging with Irregular Samples

: Sparse imaging is widely used in synthetic aperture radar (SAR) imaging. Compared with the traditional matched ﬁltering (MF) methods, sparse SAR imaging can directly image the scattered points of a target and effectively reduce the sidelobes and clutter in irregular samples. However, in view of the large-scale computational complexity of sparse reconstruction with raw echo data, traditional sparse reconstruction algorithms often require huge computational expense. To solve the above problems, in this paper, we propose a 3D near-ﬁeld sparse SAR direct imaging algorithm for irregular trajectories, adopting a piece of preliminary information in the SAR image to update the dictionary matrix dimension, using the Gaussian iterative method, and optimizing the signal-processing techniques, which can achieve 3D sparse reconstruction in a more direct and rapid manner. The proposed algorithm was validated through simulations and empirical study of irregular scanning scenarios and compared with traditional MF and sparse reconstruction methods, and was shown to signiﬁcantly reduce the computation time and effectively preserve the complex information of the scenes to achieve high-resolution image reconstruction.


Introduction
Traditional near-field synthetic aperture radar (SAR) imaging systems and security scanners use high-precision motion control devices to build ideal regular synthetic arrays to achieve imaging [1,2], which is limited by the Nyquist sampling theorem.The amount of echo data acquired by applying this regime is large and the acquisition time is high.In contrast, a near-field imaging system with an irregular sampling method has the advantages of high flexibility, low cost, and various application scenarios.In addition, with the emergence of the new generation of millimeter-wave (mmWave) devices, near-field synthetic aperture radar (SAR) imaging has moved from indoor scenarios such as medical diagnostics [3,4], gesture recognition [5,6], nondestructive testing (NDT) [7,8], and concealed threat detection [9,10] to many new application scenarios and has attracted wide attention.However, several scenarios, such as unmanned aerial vehicle (UAV) SAR [11][12][13][14], automotive SAR [15][16][17], and freehand imaging [18][19][20] remain hampered by irregular scanning geometries that do not conform to the classic array geometries needed for high-resolution image reconstruction algorithms.
Previous work has proposed the use of traditional imaging methods for irregular SAR geometric arrays, such as the back-projection algorithm (BPA) [21][22][23] and the nonuniform fast Fourier-transform range-migration algorithm (NUFFT-RMA) [24][25][26][27][28].These two algorithms belong to the time-domain and frequency-domain imaging algorithms, respectively, which are essentially imaging-optimization algorithms within the framework of matched filtering (MF).The traditional MF method, as a simple linear filtering process, is limited by the system bandwidth and Nyquist sampling theorem, which make it hard to break the Rayleigh boundary without exploiting any prior information.Moreover, the primary lobes' expansion, higher sidelobes, and grating lobes inevitably occur in the case of sparse data cubes and a large amount of coherent speckle noise appears in the image, resulting in the degradation of SAR image quality.
In the past few decades, many methods have been developed to suppress the sidelobes/grating lobes and enhance the radar image, such as coherence factor (CF)-like methods [29,30] and sparsity recovery [22,[31][32][33].The CF algorithm can effectively suppress sidelobes/grating lobes and clutter by masking the target image with the coherence factor.However, it also weakens the weak targets around strong targets.Based on the CF algorithm, some scholars have also proposed related improvement methods, such as the exponential phase coherence factor (EPCF) method based on the phase difference of each channel echo at the primary lobes and the grating lobes, which makes full use of the phase information and offers a certain improvement on the gate-side lobe suppression effect [34].However, it does not completely solve the problem of mutual influence between the targets of the CF algorithm, and the computational complexity is significantly increased.
Sparse reconstruction algorithms such as the orthogonal matching pursuit (OMP) algorithm [31,35,36] have an advantage in many fields due to their sparse characteristic of requiring an extremely high signal-to-noise ratio and extremely accurate physical models to break through the limitation of the Nyquist sampling criterion and use a sub-Nyquist sampling rate to achieve high-fidelity imaging results [31,37,38].They can also effectively suppress clutter and sidelobes/grating lobes.However, they require conversion of the entire echo data matrix into vectors and optimization of iterations as well, which generates huge time-consuming iteration calculations.In addition, nonguaranteed stabilization and robustness through a mismatch of the sensing matrix with noise and interference prevents rapid implementation of 3D reconstruction.
In this article, taking into consideration the advantages and disadvantages of the approaches mentioned above, we propose a 3D sparse SAR direct imaging framework with irregular samples.In our framework, mmWave radar is used to detect the target and then perform 3D imaging.In contrast to existing imaging algorithms, we first perform initial imaging of the scene using the MF algorithm to determine the target location information, and then apply it as a priori information to reduce the matrix dimensionality and decrease the matrix storage via signal-processing techniques to speed up the computation.Our research will satisfy the demand for 3D direct imaging with irregular sampling, such as in UAV SAR, automotive SAR, concealed detection, nondestructive testing, etc.The main contributions of our paper are summarized as follows: 1.
A 3D sparse direct image reconstruction framework is established to solve the SAR imaging problem of irregular samples.

2.
We designed a NUFFT-RMA method to obtain the initial image position information.Subsequently, a threshold filtering method is put forward to reduce the dimensionality of the observation matrix.Finally, signal-processing techniques, such as storing the matrix in advance, are applied to achieve fine 3D sparse imaging.

3.
In contrast to previous compressed sensing methods, we propose replacing the interpolated data organized into vectors with raw data to avoid interpolation errors.

4.
This is the first attempt to develop a hybrid imaging algorithm of matched filtering combined with sparse reconstruction for 3D imaging in an irregular sampled scene.It can provide better performance in the suppression of sidelobes/grating lobes, the reduction of computation time and storage, and a significant improvement in the imaging quality of multiple targets.
The proposed algorithm was validated using simulation and empirical investigation.Our algorithm can perform 3D high-resolution SAR imaging under arbitrary array configurations.Within the paper, Section 2 introduces the relevant research theories, including the signal model, the BPA, the NUFFT-RMA, the OMP, and the novel sparse enhancement technique.Section 3 discusses the imaging results and parameter metrics of image quality presented to demonstrate the superiority of the proposed algorithm.Section 4 summarizes the whole article.

MF Model
The near-field 3D mmWave SAR imaging framework using irregular samples is shown in Figure 1.The computer controls the scanning of the radar antenna to form a randomly sampled synthetic aperture plane at z = 0, where the radar antenna transmits and receives the electromagnetic wave signal during the movement.Assuming that the transmitting signal of the system is a linear frequency-modulated continuous wave (FMCW) signal [21], the transmitted signal can be expressed as where f c is the instantaneous carrier frequency at time t = 0, K = B/T c is the chirp of the frequency slope, B is the sweep band of the signal, and T c is the duration of the fast time.
After dechirping, the complex echo signal is obtained as follows: where τ is the pulse round-trip time delay with σ representing the target scattering coefficient, f (t) = f c + Kt denotes the signal carrier frequency at each sampling point, and k(t) = 2π f (t)/c is the wave number.The last phase term is known as the residual video phase (RVP), which is ignored.Finally, the signal model can be represented as We model the target as a continuous set of point targets located in a volume V in the (x, y, z) space.Expanding the signal model in (3), we can obtain where f (x, y, z) represents the radar cross-section (RCS) coefficient of the target.The amplitude decay is considered here.Assuming the scattering points of the distributed target scene are closely located, |p − p | 2 = R 2 can be approximated as |p − p | [21].In this case, the raw data cube contains three-dimensional information in the range dimension, azimuthal direction, and altitude direction [39].The data cube in Equation ( 4) is obtained via range compression on a range direction with the Fourier transform.
Using the classic BPA, Equation ( 5) can be reformulated to recover the 3D complex image.
The Fourier-based algorithm in the subsequent analysis is known as the rangemigration algorithm (RMA); Equation ( 5) can be rewritten as Using the spherical wave as an approximate superposition of plane waves [40,41], Taking Equation ( 8) into (7), the signal can be written as where Leveraging the conjugate symmetry of the spherical wavefront, (10) can be rewritten in the following form to exploit the spatial Fourier transform on z: where (•) * is the complex conjugate operation.Using the 3D spatial Fourier transform and rearranging the terms in (10), Equation ( 10) becomes k z e j(k x x +k y y ) dk x dk y (11) The outer double integral above represents a 2D inverse Fourier transform.Hence, Equation (11) becomes The backscattered data spectrum s * (k x , k y , k) is assumed to be uniformly sampled in the k domain.Hence, the backscattered data are resampled to uniformly spaced positions in the k z domain.Therefore, the 3D RMA image reconstruction can be carried out as follows: where Stolt (k→k z ) means the interpolation to resample the data cube to uniformly spaced positions in k z domain.The 3D RMA image reconstruction for a 2D rectilinear synthetic array requires the signal to be uniformly sampled in the frequency and spatial domains.However, signals are randomly sampled in the spatial domain, which also prevents the RMA algorithm being used for 3D image reconstruction of general irregular sampled arrays.To overcome the above problem, the NUFFT-RMA is introduced to tackle the nonuniform spatial sampling in this paper [27].NUFFT is a fast algorithm for converting time-domain nonuniform to frequency-domain uniform sampling, changing the data directly via NUFFT without interpolating the data before imaging, which can reduce the errors and unnecessary operations caused by interpolation.Equation ( 13) can be reformulated as

The Proposed Algorithm
Compared with MF algorithms, sparse imaging algorithms can effectively suppress sidelobes and clutter, significantly improving the image quality.However, an observationmatrix-based sparse imaging algorithm demands conversion of the echo data matrix into vectors and reconstruction of the imaging scene, for which the computational complexity is high.Therefore, a direct 3D sparse reconstruction algorithm is proposed.Using the original data, instead of interpolating data, can effectively save time and reduce the error rates.Firstly, the SAR-based imaging model can be expressed as Equation ( 15) describes the decomposition of the measurement signal; it presents the signal formation process, where Y is the measurement vector (raw data), A represents the dictionary matrix, and X is the scattering coefficient matrix.The OMP algorithm mainly solves the minimization problem of the complex domain l 0 norm.
For our purposes, we use the l P norm to solve the above problem.By bringing the compressive sensing theory into Equation ( 15) and adding the regularization constraint [31], the sparse reconstruction of the observation scene can be expressed as where • p is the l p norm of a matrix, p(0 < p < 1) represents the shrinkage parameter, and µ is the regularization parameter [31].In actual imaging scenes, most areas of the scene are not essential, and the target areas of interest often account for a small percentage of the scene which is sufficient for 3D reconstruction.Therefore, we chose to perform dimension reduction on the observation matrix dictionary to improve the imaging speed and reconstruction accuracy.Imaging algorithms based on the matched filtering framework can all perform preliminary imaging of the target to initially select the target region.
Here, we use the NUFFT-RMA method to obtain the initial imaging results X and those characteristic of the target area with stronger scattering intensity to select the threshold τ.
The threshold τ is the value of the boundary between the initial screening target area and the clutter.Supposing the target area is C, Q loc is the local mean, Q all is the global mean, k is the constant of the scale factor, and β is the standard deviation of the image.Slide the 3 × 3 window of the image through the whole image in order, and we can obtain where X (x i , y i , z i ) is the N scattering objects which correspond to the target area C. A detailed explanation is given in Figure 2.While reducing the number of measurements due to the sparsity of the scene, a sensing matrix is introduced to multiply the dictionary matrix, where Θ is the sensing matrix [42,43].Assume that the original dimension of the observation matrix A is M × N, and the observation matrix A is now updated to M × N .We can rewrite Equation ( 17) as Taking the idea of Cetin's process [44], the cost function can be defined as Calculating the partial derivative of f with respect to X gives the following: where 2A H A + µpΛ can be assumed to be an approximate Hessian matrix and , where diag(•) indicates the diagonal matrix.We can use the Gaussian iterative method to obtain the 3D sparse reconstruction: The geometric explanation of the iterative operation is given in Figure 3, where ξ represents the iteration termination threshold and ∆ n+1 = (∆ n ) 0.8 is applied to accelerate the convergence of the iteration In fact, the 3D sparse direct imaging model without threshold selection presents large storage and computational challenges; assuming the size of the raw data in the 3D imaging dimension is 5 × 200 × 200, then M is approximately 2 × 10 5 units of storage.If the size of the reconstructed 3D scene is also 5 × 200 × 200, then N is approximately 2 × 10 5 units of storage.One storage cell is 8 bytes, meaning that the observation matrix A requires a storage space of about 298 GB.Generally, A cannot be stored, and it requires a huge amount of computation to compute A H A and A H Y in the iteration process.Because of the mentioned huge storage and computational challenges, it is impractical to apply 3D sparse reconstruction directly.Therefore, using threshold selection to reduce the observation matrix dimensions, we can obtain A H A = M × diag(ones(N )) and the dimension of and A H Y can be obtained after filtering Y with a threshold τ.By estimating A H A and A H Y in advance, the speed of the iterative operations can be enhanced and the amount of storage will be reduced, thus allowing implementation of 3D image reconstruction.Based on the above idea, the challenges of massive calculations related to the 3D sparse SAR direct reconstruction are resolved.The specific steps of the proposed sparse imaging algorithm are summarized in Table 1.Table 1.The steps of the proposed sparse imaging algorithm.

3D Sparse SAR Direct Imaging with Irregular Samples
1. Initial imaging by NUFFT-RMA to get the imaging result X .2. According to Equation ( 18), determine the target area C.
3. For the sparse SAR imaging a.According to the target area C and imaging result X ,take the idea of Cetin's process, and obtain the cost function f .b.Using the Gaussian iterative method to Simplify the derivative equation.

c. Obtain and store A
H A and A H Y in advance and Calculate X n+1 iteratively according to Equation ( 21), and get the final calculation results

Imaging Results and Evaluation
In this section, simulated and real measured data are shown to demonstrate the effectiveness and reliability of the proposed algorithm.The simulated and measured scanning sizes were 200 × 200 mm 2 and300 × 300 mm 2 .In the setup, all the computations were implemented using MATLAB in Windows 10 on an Intel Core i7-8250K CPU@1.8GHz and 16 GB memory.No parallel computing or GPU was adopted and the data acquisition time was not considered.The radar system parameters are shown in Table 2.

Imaging Results
Figure 4 indicates the geometric model of the nine ideal point targets distributed in space.Figure 5a-d    After obtaining the original data cube, these four algorithms were used to reconstruct the point target images.The image reconstruction by these algorithms used 35% random samples and added 60 dBW Gaussian white noise.By comparing these methods, it can be clearly seen that with the irregular sampling, the signal was severely attenuated with the MF algorithm and the images from the BPA and the NUFFT-RMA had greater ghosts and streaks on the contours of the object.However, the OMP and proposed algorithm effectively suppressed sidelobes and clutters.Although the OMP algorithm was able to almost completely reconstruct the image, its huge imaging computational complexity is not sufficient to satisfy the requirements of real-time imaging.Furthermore, compared with the image results based on these algorithms, the proposed algorithm improved the measurement accuracy of the target RCS and had better image-focusing performance than OMP.
In this paper, to verify the superiority of the proposed algorithm in real scenarios, the optical imaging of an experimental wrench model was targeted, as shown in Figure 6.The PC as the control center was connected to an IWR1843 mmWave radar, a highperformance DCA1000 raw data acquisition card, and a three-axis controllable stepper to capture the echo.The distance between the reference planar and the target center was 0.3 m.The platform operated by scanning a Z-shaped trajectory to collect data and carried out irregular sampling in the later signal processing.The mmWave radar imaged the wrench in the forward-and side-facing modes.The number of spatial sampling locations was chosen in an irregular array along x and y directions from −150 mm to 150 mm.The imaging dynamic range was set to 14 dB, and the operational frequency was set to 77-81 GHz with 256 points.
Using 30% random samples to reconstruct images, the 3D SAR imaging results of the various methods are presented in Figure 7.It can be clearly seen that the target was drowned out by clutter and sidelobes, and the image quality severely deteriorated when using the BPA and NUFFT-RMA.A comparison between Column 3 and Column 4 of Figure 7 indicates that a loss of detail and edge diffraction can be observed and the image quality of detail at the edges slightly deteriorated when using the OMP algorithm.In contrast, the measurement accuracy and image quality were improved and the image edge detail was recovered most effectively by the proposed algorithm.These results verify the practical applicability of the method.

Evaluation of Imaging Results' Characteristics
To further evaluate the performance gap between the proposed algorithm and the other algorithms, we discuss the characteristics of these methods below.
(1) Amplitude: In order to directly compare the effect of sidelobe suppression and detail retaining capability, we sliced the imaging amplitude results of the nine ideal point targets and the wrench along the azimuth and elevation directions, as shown in Figure 8.The results showed that the proposed method sharpened the principal lobes more than the traditional algorithms, which showed a narrower impulse response.
The mean differences of principal lobes in the results of the BPA, NUFFT-RMA, OMP, and proposed method were −6.24 dB, −2.23 dB, −18.96 dB, and −3.68 dB in the simulation, respectively.The mean differences in the results of the BPA, NUFFT-RMA, OMP, and proposed method were −3.18 dB, −7.62 dB, −9.26 dB, and −2.49dB in the real data, respectively, which shows that the measurement intensity and stability of the proposed method were improved compared to the other methods.In addition, the principal lobes based on the proposed algorithm were narrower than those of the other three algorithms, and the principal lobes were more focused.This also indicates the enhancement of the sidelobe clutter suppression effect and image details by the sparse method, as shown in Figure 8.The mean differences of principal lobes are shown in Table 3.  (2) Phase: Figure 9a-d shows the phases of BPA, NUFFT, OMP, and the proposed method, respectively.It can be clearly seen that the images underwent severe detail loss when using OMP, while there was a lot of clutter on the images when using MF-based approaches, which affected the image quality.Furthermore, these results indicate that the reconstruction of the image via the proposed algorithm maintained the phase information of the image, i.e., the complete complex scattered field was estimated instead of the field with only amplitude.(3) Image quality and time complexity: To quantitatively verify the image quality metrics and the computational load of these methods, the image entropy (IE), image contrast (IC), and computational time for each algorithm are shown in Table 4 [45].The image entropy shows the image-focus quality.The lower the entropy, the better the effectiveness of focusing.Image contrast is the difference in color within the image and indicates the texture characteristics of the image.The greater the contrast, the more visible the image details.Compared with the other three algorithms, the contrast of the images was improved and the image entropy was reduced by an order of magnitude when using the proposed method.These results show that the image quality of the proposed method was superior to that of the other algorithms and achieved a good focus on the target.Furthermore, it is obvious that the time complexity of the proposed method was only 1.5% that of the BPA, 0.3% that of the OMP, and close to that of the NUFFT-RMA.Its imaging efficiency can offer good performance in such scenarios as near-field real-time SAR imaging.

Imaging Evaluation with Multiple Targets
The imaging effect of multiple targets is determined by many factors, such as target scattering characteristics, data sparsity, etc.In this section, to further validate the strong robustness of the proposed algorithm for multiple targets, we randomly selected 50% of the 3DRIED [22] raw data and reconstructed the images from these imaging algorithms.Since the radar antenna spacing in 3DRIED is greater than half a wavelength, the SAR images based on conventional algorithms appeared as ghost images.Figure 10 shows the optical images of multiple targets.These include a pistol, knife, and stiletto.The 2D multi-target imaging results and the 3D imaging envelope are shown in Figure 11, and a brief analysis is provided.Observing the images for the multiple targets with different scattering coefficients in the space, their outlines on the reconstructed image were varied.The sidelobe/grating lobes clutter of the images obtained using different algorithms was diverse, which was reflected in the purity of the background.The sidelobe/grating lobes generated by the traditional MF method were relatively severe, and the operation using beam compensation efficiently decreased the sidelobe clutter [21]; this is outside the scope of this paper.The imaging results obtained using the OMP showed that a weak scattering target, such as a pistol, could not be imaged properly.Therefore, these results show that OMP is weak in robustness, requiring high SNR, and the computational efficiency of the algorithm is poor, which demands a large computational load.Compared with these methods, the proposed method effectively reduced the sidelobes/grating lobes and reconstructed the image of multiple targets with high efficiency by combining MF and sparse reconstruction algorithms.The image entropy and contrast of the imaging results obtained by the proposed method also reached an excellent level, as described in the previous paragraph, and the image quality was significantly improved.

Conclusions
In this paper, we describe the design of a 3D near-field millimeter-wave SAR imaging platform and propose a 3D near-field sparse direct imaging algorithm for irregular trajectories.Our proposed algorithm provided better performance in suppressing side lobes/grating lobes compared with traditional methods such as BPA in the time domain, NUFFT-RMA in the frequency domain, and OMP using compressed sensing theory; it required less computation time and storage space and provided significant improvement in the imaging quality of multiple targets.The algorithm has good robustness and obtained higher focusing performance.Use of the algorithm is feasible in typical applications such as freehand imaging, UAV SAR, and automotive imaging.We demonstrated through simulation and experimental studies that our algorithm achieved much higher SAR image quality compared with the traditional algorithms.Our goal is to contribute to the promotion of a complete and effective imaging process for industrial research.

Figure 2 .
Figure 2. The target scattering area selection map.

Figure 3 .
Figure 3.The geometric explanation of the iterative operation.
Figure 4 indicates the geometric model of the nine ideal point targets distributed in space.Figure 5a-d shows the simulation imaging results of BPA, NUFFT-RMA, OMP, and the proposed algorithm, respectively.The target center was 0.23 m away from the array aperture of the reference plane.The number of spatial sampling positions was selected along x− and y− orientations from −100 mm to 100 mm for an irregular array, and the operational frequency was 77-81 GHz with 201 points.

Figure 4 .
Figure 4.The geometric model of the point target.

Figure 5 .
Figure 5. Spatial target model of multiple points in an ideal space.(a) The 2D multiple-point imaging result of BPA.(b) The 2D multiple-point imaging result of NUFFT-RMA.(c) The 2D multiple-point imaging result of OMP.(d) The 2D multiple-point imaging result of the proposed method.

Figure 7 .
Figure 7. Imaging results.Row 1 (a-d) is the 2D wrench imaging results of the four imaging algorithms.Row 2 (e-h) is the 3D wrench imaging results of the four imaging algorithms.Column 1 is the BPA imaging results; Column 2 is the NUFFT-RMA imaging results; Column 3 is the OMP imaging results; Column 4 is the proposed method's imaging results.

Figure 8 .
Figure 8. Azimuthal and elevation profiles.(a) The azimuthal profiles of the multiple targets.(b) The elevation profiles of the multiple targets.(c) The azimuthal profile of the wrench.(d) The elevation profile of the wrench.

Figure 9 .
Figure 9. Imaging results of phase.(a) The imaging result of phase when using BPA.(b) The imaging result of phase when using NUFFT-RMA.(c) The imaging result of phase when using OMP.(d) The imaging result of phase when using the proposed method.

Figure 10 .
Figure 10.The optical images of multiple targets.

Figure 11 .
Figure 11.Imaging results of multiple targets using different algorithms.Rows 1 (a-d) and 3 (i-l) are the 2D imaging results; Rows 2 (e-h) and 4 (m-p) are the 3D imaging results; Rows 1 and 2 are the results for the pistol and knife; Rows 3 and 4 are the results for the pistol and stiletto; Column 1 is the BPA imaging results; Column 2 is the RMA imaging results; Column 3 is the OMP imaging results; Column 4 is the proposed method imaging results.

Table 2 .
The parameters of the SAR platform.

Table 3 .
The mean difference of principal lobes.

Table 4 .
Quality and time complexity.