Three-dimensional SAR imaging with sparse linear array using tensor completion in embedded space

Due to the huge data storage and transmission pressure, sparse data collection strategy has provided opportunities and challenges for 3D SAR imaging. However, sparse data brought by the sparse linear array will produce high-level side-lobes, as well as the aliasing and the false-alarm targets. Simultaneously, the vectorizing or matrixing of 3D data makes high computational complexity and huge memory usage, which is not practicable in real applications. To deal with these problems, tensor completion (TC), as a convex optimization problem, is used to solve the 3D sparse imaging problem efficiently. Unfortunately, the traditional TC methods are invalid to the incomplete tensor data with missing slices brought by sparse linear arrays. In this paper, a novel 3D imaging algorithm using TC in embedded space is proposed to produce 3D images with efficient side-lobes suppression. With the help of sparsity and low-rank property hidden in the 3D radar signal, the incomplete tensor data is taken as the input and converted into a higher order incomplete Hankel tensor by multiway delay embedding transform (MDT). Then, the tucker decomposition with incremental rank has been applied for completion. Subsequently, any traditional 3D imaging methods can be employed to obtain excellent imaging performance for the completed tensor. The proposed method achieves high resolution and low-level side-lobes compared with the traditional TC-based methods. It is verified by several numerical simulations and multiple comparative studies on real data. Results clearly demonstrate that the proposed method can generate 3D images with small reconstruction error even when the sparse sampling rate or signal to noise ratio is low, which confirms the validity and advantage of the proposed method.

perfectly by a tucker decomposition in the embedded space and subsequently the 3D imagery can be focused by the various 3D imaging approaches. Finally, the effectiveness and the accuracy of the proposed algorithm have been evaluated by various experiment results of simulated and measured data.
The paper sections are shown below. Section 2 briefly introduces the notations and formulas of tensor as the preliminaries and backgrounds. In Sect. 3, the construction of 3D sparse array SAR signal model in the tensor space is given. In Sect. 4, we discuss the low rank and sparse property and then present a 3D imaging algorithm using TC in embedded space for 3D sparse array SAR. Section 5 shows experimental findings that demonstrate the validity of the presented method. Finally, conclusions are drawn in Sect. 6.

Notation and preliminaries
This section first introduces the notations and formulas related to tensor or tensor operation as a preparation for the following algorithm derivation, as shown in Table 1. Figure 1 illustrates the imaging geometric relationship of 3D sparse array SAR. Suppose that the aircraft flies at the altitude H with a speed V a . Under the wings of the aircraft, a number of antenna arrays suspended at unequal intervals made up a sparse linear array with the length of L y . In the 3D space coordinate system, three orthogonal axes represent the azimuth, cross-tack and altitude, corresponding to the aircraft flight direction, linear array direction and radar irradiation direction, respectively.

Definitions Notations and formulas Notes
Vector b Matrix B Tensor A A (i) -the ith matrix in a sequence Unfolding/matricization The i-mode unfolding/matricization: A (i) Unfolding a tensor to a matrix Frobenius norm �A� F = m1 m2 · · · m I a 2

Tensor multiplication
The i-mode product of tensor A and matrix B : (2) · · · × I F (I) C is core tensor Multi-linear tensor product Multi-linear tensor product: is a set of the factor matrices Multi-linear tensor product with the i-th mode excluded:

Rank-one tensor
If A can be expressed as an outer product of vectors, i.e., A = x (1) • x (2) • · · · • x (I) , it satisfies rank-one tensor • denotes the outer product of vectors Tensor rank rank(A) The minimum quantity of rank-one tensors which make up A by their sum From the theory of scattering center model [29], the large volume target consists of a finite number of scattering points. The activate antennas transmit a step frequency signal. This modulated signal consists of a set of K pulses with a reference frequency f 0 . The step frequency increases sequentially in a constant frequency increment f . The signal reaches target B and then is reflected back to be received by the activate antennas. After a double-way propagation distance, the receiving antenna can obtain the echo signal E as where A is the amplitude. γ denotes the backscattering coefficient of the target. τ is the time variable. The one-way distance of signal propagation R d is indicated as the distance between the transmitting/receiving antenna P t = (V a t m , y n , H) and the target where R 0 = (V a t m ) 2 + y 2 n + H 2 indicates Zero-Doppler distance. The first approximately equal sign is established according to Fresnel approximation. It can be seen that the third part can be neglected due to the far-field condition, i.e., R 0 ≫ x B , y B , z B .
Since the sample on a 3D signal grid is discrete, the received data by these scatterers can be expressed as Further simplification, the part of constant exp(−4πjR 0 ) is ignorable in the below content. Let three variables δ x , δ y , δ z are introduced as Therefore, the echo signal can be represented as follow format. It can be view as a generic signal model for 3D SAR imaging.
In order to obtain the 3D imagery in tensor space, the above signal model can be represented into tensor form. The imaging scene can be divided into grids at equal intervals with the size of P × Q × L . That is to say, if there is a scattering point in the grid, γ is not equal to 0. Thus, Eq. (5) can be re-described as The tensor representation of the signal model can be described as where E is tensor format of the echo E with the order of M × N × K , G is tensor format of the backscattering coefficient γ with the order of P × Q × L , × i , i = 1, 2, 3 represents tensor multiplication. x is the azimuth steering vector matrix, y is the cross-track steering vector matrix, z is the range steering vector matrix.

Methods
Before reconstruction the target, we consider recovering the lost data elements first. Then, a traditional 3D imaging algorithm is exploited to obtain a more excellence 3D imagery. However, the echo data acquired by the sparse linear array can be viewed as a 3-order incomplete tensor with randomly missing slices. In this case, it is difficult for the traditional TC approach to recover the lost data elements, which cause the unacceptable imaging performance.
In this section, we present a novel 3D imaging algorithm for sparse array SAR using TC in embedded space to solve the above problem. Figure 2 depicts the flow diagram of this 3D sparse array SAR imaging method based on the echo tensor with missing slices. The algorithm includes two parts: tensor completion and 3D imaging, where the part of tensor completion contains MDT, low-rank tensor approximation, and inverse MDT steps.

Property analysis for signal tensor
First, we discuss the sparse property based on the signal model in (6), which helps to perform 3D image reconstruction from a sparse signal tensor. For 3D sparse array SAR, massive quantity of non-target areas exists in the 3D scene. That means the sum of the backscattered responses of few prominent scatters constitutes the 3D radar image G , which means that G is expected to be sparse. Then, the prerequisite of tensor completion is the low-rank property. The signal data for one scattering point B is constructed as Thus, the tensor representation of signal E is expressed as The operator of outer product indicates as • . Since E is expressed as an outer product of vectors, it satisfies a rank-one tensor. It also shows that E is a low-rank tensor.
Due to the superposition of signals from several scatterers, the final signal E can be characterized as a linear combination of the corresponding 3-mode rank-one tensor.
The minimum quantity of rank-one tensors which make up E by their sum (see Sect. 2) denotes the rank of E . That means the rank of tensor E cannot exceed B. In addition, a few number of strong scattering points form the image scene, thus B ≪ PQL , and further derived that rank(E) ≤ B ≪ PQL . Hence, if the target is sparse, then the tensor E has the property of low CP-rank.

MDT and inverse MDT
Given that the 3-order echo tensor E is transformed by a delay embedding with parameters ξ = {ξ 1 , ξ 2 , ξ 3 } ∈ N 3 and = {M, N , K } . This processing includes two steps: duplication step and folding step [28], as shown in Fig. 3.
First, the MDT produces the low order tensor E into a duplicated high order tensor, which is called as "Hankelization" [30]. Consider the duplication matrices are satisfied as Then the MDT can be obtained by denotes a folding operator from a low order tensor to a high order one. Here a 6-th order tensor with In contrast, the inverse MDT transform can be decomposed into two steps: a matricization operation (also known as unfolding) and the Moore-Penrose pseudoinverse D † = (D T D) −1 D T . Thus, the Hankel tensor E H after the inverse MDT processing is presented as

Low-rank tensor approximation
According to (13), the incomplete tensor E ∈ C M×N ×K and its mask tensor M ∈ {0, 1} M×N ×K can be transformed by MDT, which are given by where I = 6 because the Hankel tensor E H is a 6-th order tensor here. The zero elements in M H correspond to the lost entry; otherwise, the one element corresponds to the available entry.
Here, we can solve the low-rank tensor approximation by tucker decomposition. Hence, the solution of the optimization problem can be converted into the following form.
where ⊛ denotes the element wise Hadamard product.
The above equation is not a convex issue and the solution is not unique [20]. For the tensor with the complete elements, its stationary point can be efficiently obtained by (13)  the alternating least squares (ALS) [20] for tucker decomposition, whereas for the tensor with lost elements, we can address this optimization problem effectively by gradient descent method [31] and manifold optimization [32]. It is obvious that the step-size parameter affects the efficiency of both algorithms. Therefore, an auxiliary function [28] is introduced as where a parameter set α = {C, F (1) , F (2) . . . , F (I) } , T α = C × {F} denotes a tucker decomposition, and M H as a complement set of M H is equal to 1 − M H . According to [28], we can transform the auxiliary function into Clearly, there is a two-step processing composing the auxiliary function minimization solution. First, the auxiliary tensor X can be calculated by Then, the factor matrices {F} and the core tensor C are updated by using the ALS to optimize In order to the non-uniqueness of the tensor T solution, the rank increment strategy is integrated into tucker-based completion, which has been discussed in [28]. Specifically, a very low-rank tucker decomposition is first set up and used as initialization to obtain a higher-rank decomposition. Following that the rank is updated iteratively until the noise condition is less than a threshold. In short, the algorithm of low-rank tensor approximation can be summarized step by step as follows.
1: Set a low rank sequence R i = 1 where i = 1, . . . , 6 as the initial value. 2: According to Eqs. (19) and (20), compute C and F as a noise condition is not larger than a noise threshold parameter η . If it is satisfied, the algorithm is terminated; otherwise, the algorithm continues to execute. 4: update the parameter F and increment R i ′ , and then go back to step 2.

3D image reconstruction
After the completely sampled 3D data has been recovered by TC mentioned above, 3D imagery can be focused exactly by the Fourier transform-based technologies. Further, some super-resolution imaging algorithms such as the spectrum estimation strategies [33] also can be employed. To summarize, the proposed method for 3D sparse array SAR imaging using TC in embedded space is shown in Algorithm 1. It is clear that the sparse tensor is transformed to an incomplete high order Hankel tensor by MDT first. Next, lowrank tensor approximation is leveraged to complete the higher order tensor and in the next step converted to the full-sampled data tensor by inverse MDT. Last, the 3D imagery can be focused by applying the 3D Range Doppler (RD) algorithm. More details of 3D RD algorithm can be found in [34].
Algorithm 1 3D SAR sparse imaging using tensor completion in embedded space

Results and discussion
In this section, the validity and accuracy of the proposed 3D imaging algorithm is evaluated by the simulated and real datasets. The running environment carried out in our experiments is listed in Table 2.

Verification of the proposed algorithm using simulated dataset
First, the simulations of the imaging distributed scene for 3D sparse array SAR are presented to verify the validity of the proposed algorithm. An X-band SAR and a single scatterer located at (3 m, 0 m, − 1 m) are utilized. Table 3 lists the parameters used in the simulations. A full-sample 3D data in this simulation is collected by a uniform virtual linear array. The virtual array of 120 elements is built by using 4 transmit and 30 receive elements with MIMO technique. SNR = 10 dB. The first line in Fig. 4 shows the 3D imaging results by using 3D RD method with full-sampled data. The sparse sampling rate (SSR) is denoted by the number ratio of the observed samples to the all samples. Here, we use a sparse linear array by randomly selecting 60 elements to produce the sparse data, i.e., SSR = 50%. This sparse data can be seen as a 3D tensor data with the missing slice. We pad the lost samples with zeros and then apply the different imaging algorithms to acquire 3D images. Figure 4 shows the 3D imaging results with different methods. Images are obtained by the conventional 3D RD method and the proposed method in the second and third line, respectively. The proposed method can recover the missing slices exactly. We set ξ = (32, 1, 1) and a (120, 200, 120) echo data was converted into a (32, 89, 1, 200, 1, 120) tensor. This Hankel tensor was re-expressed as a fourth-order tensor with (32,89,120,200). Thus, we set the rank sequences as L 1 =   From Fig. 4, the scatterer is recovered with high-level side-lobes along the crosstrack direction by the 3D RD algorithm. For comparison, the proposed method can focus the scatterer precisely because of TC processing with MDT. Noted that there is no significant change along both azimuth and height directions because the missing slice is only appeared along the cross-tack direction. The peak sidelobe ratio (PSLR) and integrated sidelobe ratio (ISLR) as evaluated indexes are leveraged, which are shown in Table 4 for the point target. Here the best values are emphasized in bold font. As we can see from the cross-track direction, PSLR and ISLR values of 3D RD with sparse data deteriorate significantly, whereas the proposed method obtains the best performances, which is much closer to the corresponding indexes with 100% data. Moreover, PSLR and ISLR values in the other two directions obtain the similar results among these algorithms.

Performance comparison with different algorithms
Furthermore, the performance of the proposed method is compared with those of other different algorithms including the conventional TC algorithms: high accuracy low rank tensor completion (HaLRTC) and tucker decomposition (TDC). Consider a target   consisting of 10 scattering points with the same amplitude. 3D RD method is chosen as a baseline method for comparison and the imaging results are illustrated in Fig. 5 with full-sampled data. Figure 6 illustrates the results obtained by various methods. From top to bottom, images are acquired by the 3D RD method, the HaLRTC-based method, the TDC-based method and the proposed method, respectively. The 3D tensor data with 80% random missing slices, i.e., SSR = 80%. From left to right, the 3D result and its 3D views on crosstrack-azimuth slice, cross-track-height slice and azimuth-height slice are displayed, respectively. Comparing Fig. 6 with referenced Fig. 5, it can be observed evidently that although the scatter points can be resolvable using the conventional 3D RD method with sparse data, higher side-lobes is obvious along the cross-track direction because of the  inadequate sampling. As shown in the second and third lines, HaLRTC and TDC fail to reconstruct the lost elements, which lead to the invalid suppression of side-lobes in 3D imaging result. By contrast, the proposed method can effectively suppress the sidelobes and can get similar 3D imaging performance with 100% data in Fig. 5. Similar to the results from Fig. 4, the imaging performance along the azimuth and height direction changes between different methods are not obvious because the remaining samples on the azimuth-height slices are full sample.

Performance with different SSRs and SNRs
Different SSRs and SNRs would bring the changes in the imaging performance. To quantitatively evaluate these changes, we compute the mean squared error between the 3D image reconstructed with sparse data T sparse and the full-sampled data T full .
Apparently, we prefer a smaller MSE value, which illustrates the reconstructed signal is much closer to the original signal.
We set up multiple sets of experiments with different SSRs from 10 to 90% under Monte Carlo simulation. The SNR is chosen to be 10 dB. The Monte-Carlo trials is set as 50 to evaluate average MSE. Figure 7 compares average MSEs between four different methods based on sparse data with varying SSRs. It can be expected clearly that the average MSEs of the 3D RD method, the HaLRTC-based method and the TDC-based method overlap totally because the HaLRTC-based and the TDC-based methods are invalid completely. When the SSR is improved, the overall trend for average MSE of all methods is declining, which means the reconstruction accuracy becomes higher. Among these curves, when the SSR is more than nearly 30%, the error is extremely tiny (almost less than 0.1) for the proposed method, whereas the average MSE increases drastically when the SSR is less than 30%. It demonstrates that, if SSR is not excessively low (e.g., 30%), the proposed method gains a reliable imaging performance.
The SNRs are examined ranging from − 20 to 20 dB under Monte Carlo simulation. The number of Monte Carlo trials is 50. Figure 8 draws the trend of average MSE on Fig. 7 Trend of average MSE between different methods based on sparse data with various SSRs. The SNR is set to be 10 dB each SNR levels. Although we found that the fluctuation of the curves with 3D RD method, HaLRTC-based method and TDC-based method are not prominent, the average MSEs are slightly reduced with the increasing of SNR, as shown in the top right corner of Fig. 8. Moreover, with the increased SNR, the MSE goes down remarkably. When the SNR reaches − 20 dB, the value of average MSE for the proposed method is still less than 0.1, which implies that the imaging performance can be satisfied with a relatively low SNR. Even so, we can observe from the detail view in the lower right corner that there appears a transition point at the curve of − 15 dB, which indicates when the SNR is very low (less than − 15 dB), the increment of error becomes larger.

Verification of the proposed algorithm using real dataset
The validation and evaluation of imaging performance for the proposed algorithm are carried out on the real dataset. The experimental principle of real data acquisition is shown in Fig. 9. A pair of antennas moving along the 2-D track can transmit or receive signal, which can synthesize a virtual 2-D array. Figure 10 displays the signal transmitter, antenna track, foamy plate and eight targets covered with tinfoil. The experiment parameters on real data are shown in Table 5.
In the experiment, the size of the tensor is (50, 161, 1601). The proposed method is applied with ξ = (32, 1, 1) . The rank parameters are configured as L 1 =  Figure 11 demonstrates the 3D images reconstructed by various methods based on sparse data with missing slices. Similar to the above conclusions, the targets cannot be satisfactorily reconstructed by ordinary low-rank model. However, the proposed method gets very clear results even though its accuracy is not high. At the same time, the comparison results by four methods with different SSRs are also shown in Fig. 11. With the decrease number of samples, all of the methods obtain more blurred images. When the SSR reaches 20%, the targets cannot distinguish at all by 3D RD and HaLRTC-based method and TDCbased method produced similar results. By contrast, the proposed method obtains significant improvements compared with the other methods. Table 6 shows the average MSE results with three different SSRs. It is demonstrated that the proposed method performs better than the other methods and it was also very robust in terms of the SSR.

Conclusion
In this manuscript, we have provided a new idea to realize the 3D imaging of targets for 3D sparse array SAR in embedded space. On the basis of the sparsity and low-rank, the unsampled elements are reconstructed in the multiway delay embedded space by utilizing the tucker decomposition. After tensor complementation, the satisfactory 3D images can be easily achieved by any conventional algorithm. In comparison with other matched filter-based methods, the proposed algorithm has the advantage of acquiring target images  with high-resolution and low side-lobes. Through extensive simulations and experiments on real datasets, the results clearly show that the proposed algorithm can significantly curb the adverse effects resulting from sparse data. In contrast to 3D RD and state-of-the-art tensor completion algorithms, the proposed algorithm can reach similar imaging performance with complete data. Experimental results at different SSR and SNR values also illustrate that the proposed method is robust in the field of SSRs and SNRs. Note that the theoretical derivation of signal model in this paper is based on the first-order Taylor expansion approximation. It is only suitable for the research of point targets in 3D scenes. In the future, 3D imaging will be further deeply studied for area targets with the distributed characteristics.