A fast and gridless STAP algorithm based on mixed‐norm minimisation and the alternating direction method of multipliers

National Natural Science Foundation of China, Grant/Award Number: 61372133 Abstract Sparse recovery (SR) based space‐time adaptive processing (STAP) methods have received much attention recently due to their dramatically reduced requirements of training samples. However, most of the existing SR‐STAP algorithms suffer from the off‐ grid effect induced by the discretisation of angle‐Doppler plane. To eliminate this effect and improve the performance of clutter suppression, a novel gridless mixed‐norm minimisation (MNM) based STAP algorithm is proposed. Unlike previous grid‐based SR‐STAP methods, the clutter covariance matrix (CCM) is estimated in continuous domain by solving a semidefinite programme (SDP), which is established by utilising the properties of the CCM and the compact formulation of MNM. Furthermore, to reduce the computational burden of solving this SDP, a fast iterative algorithm via the framework of the alternating direction method of multipliers (ADMM) is also derived, where the unknown parameters are iteratively updated with closed‐form expressions. Simulation results show that our proposed algorithm achieves better performance of clutter suppression and target detection with lower computational complexity.


| INTRODUCTION
Space-time adaptive processing (STAP) [1][2][3][4], which can separate the target and clutter from a joint angle-Doppler domain, plays an important role in airborne radar moving target detection. The key challenge of most STAP approaches is how to accurately estimate the clutter covariance matrix (CCM) of the cell under test (CUT). According to the Reed-Mallett-Brennan (RMB) rule [2], the number of independent and identically distributed (IID) training samples must be greater than twice the system degrees of freedom (DOFs, which is defined as the product of the number of array elements and the number of coherent pulses). Only then can the obtained signal-to-interference-plus-noise-ratio (SINR) loss be within 3 dB. However, due to the terrain variations and array geometry configurations, sufficient training samples are difficult to be obtained in practical applications.
To reduce the requirement of training samples, different types of approaches have been proposed in the last few decades, such as non-adaptive reduced-dimension (RD) [5][6][7] and adaptive reduced-dimension approaches [8][9][10][11][12][13][14], parametric adaptive matched filter (PAMF) approaches [15], direct data domain (DDD) [16] methods, and knowledge-aided (KA) algorithms [17][18][19]. However, although the number of training samples required by these approaches is further reduced, they all have their limitations. The requirement of these RD approaches is still large, especially for systems with large DOFs, the order for PAMF is hard to decide, the benefit of DDD is achieved at the cost of reduced DOFs, and prior environmental knowledge of the illuminated region is hard to obtain for KA-based approaches.
Recently, based on the intrinsic sparsity of clutter data, many sparse recovery based STAP (SR-STAP) approaches [20][21][22][23][24][25][26][27][28] have been proposed. In these kinds of methods, the continuous angle-Doppler plane is first discretised into a finite set of grid points. Under the assumption that the clutter components can be precisely represented by the steering vectors corresponding to the predefined grid points, the SR-STAP methods [20][21][22][23][24][25][26][27][28] can achieve near-optimal performance with a few training samples. In practical applications, however, the clutter ridge cannot locate on the predefined grid points exactly, which results in the so-called off-grid/grid mismatch, and thereby severely reduces the performance of clutter suppression [29][30][31][32][33][34]. A natural idea to mitigate the off-grid effect is to decrease the grid interval. But too dense a grid will enhance the correlation between adjacent atoms of the spacetime dictionary matrix, which might result in the power leakage of strong clutter components and therefore degrade the target detection performance. Besides, discretisation with a dense grid will lead to unaffordable computational workloads, which might hinder the real-time application of the algorithms.
To mitigate this off-grid effect, several remedies [29][30][31][32][33][34] including grid-based and gridless STAP algorithms have been proposed most recently. Generally, the grid-based remedies [29][30][31][32][33] attempt to generate new grid points that match the location of the clutter ridge by dictionary learning [29], local searching [30][31][32], or prior platform information [33]. Nevertheless, discretisation is still assumed in [29][30][31][32][33], and hence the off-grid effect cannot be avoided completely. By utilising the structure and low-rank properties of CCM, a gridless remedy has been proposed in [34], which attempts to reconstruct the clutter subspace in a continuous two-dimensional plane by solving an atomic norm minimisation-based optimisation problem. However, the performance of the atomic norm minimisation algorithm is sensitive to the low-rank property, and the computational complexity of this approach is very high due to the interior point method-based solving process [34].
In this work, to efficiently eliminate the off-grid effect, we investigate a new gridless STAP approach inspired by the mixed-norm minimisation (MNM) [35] criterion (i.e. MNM-STAP). The proposed MNM-STAP approach estimates the CCM in the continuous domain by solving a semidefinite programme (SDP), which is constructed by the block-Toeplitz and positive semidefinite (PSD) properties of the clutter covariance matrix (CCM), and the compact formulation of MNM. Moreover, to improve the computational efficiency, we also derive a fast iterative algorithm for solving this SDP problem via the alternating direction method of multipliers (ADMM) [36], in which the unknown parameters are iteratively updated with closed-form expressions. Numerical experiments with simulated data show that the proposed algorithm has better clutter suppression performance and higher computational efficiency. The main contributions of this paper are summarised as follows. The rest of this paper is organised as follows: Section 2 introduces the signal model and the grid-based SR-STAP algorithm for airborne radar. Section 3 details the proposed gridless MNM-STAP algorithm and its fast implementation based on ADMM. Numerical experiments are performed in section 4. Conclusions are drawn in section 5.
Notation: ℝ and ℂ denote the sets of real and complex values respectively. D and D þ denote the sets of diagonal and non-negative diagonal matrices respectively. Boldface uppercase letters denote the matrices, and Boldface lowercase letters denote the vectors. ⊗ denotes the kronecker product. For a matrix A, A T and A H denote the transpose and conjugate transpose of A respectively. trðAÞ denotes the trace of A. A ≥ 0means A is a positive semidefinite matrix.
, and Frobenius norms respectively. For a vector x, diagðx Þ denotes a diagonal matrix with the elements in x on its main diagonal. For a scalar x, |x| denotes the absolute value of x.

| Signal model
Consider an airborne phased-array radar system with a uniform linear array (ULA) consisting of N elements. The interelement distance of this array is d ¼ λ=2, in which λ is the wavelength. The geometry of the airborne radar system is shown in Figure 1. The altitude and velocity of the platform are H and v, respectively. The symbols θ, φ, and φ 0 stand for the elevation angle, azimuth angle, and crab angle, respectively. M pulses are transmitted in a coherent processing interval (CPI) at a constant pulse repetition frequency (PRF) f r .
It is well known that the radar targets detection is a binary hypothesis problem where hypothesis H 0 denotes the absence of the target and H 1 denotes the presence of the target.
where x t denotes the target signal; x l;c ¼ ∑ N c k¼1 α k sðf d;k ; f s;k Þ stands for the clutter data of the lth range gate [4]; N c denotes the number of clutter patches evenly distributed in azimuth; α k , f d,k and f s,k are the complex amplitude, normalised Doppler frequency and normalised spatial frequency of the kth clutter patch, respectively; n ∼ CN ð0; σ 2 n I MN Þ is the complex white Gaussian noise with the covariance matrix σ 2 n I MN ; sðf d;k ; f s;k Þ is a MN-dimensional space-time steering vector, which can be represented as where s t ðf d;k Þ and s s ðf s;k Þ are the corresponding temporal steering vector and spatial steering vector, respectively, and where 8 > > > < > > > : Suppose the components in x l are mutually uncorrelated. Then the ideal clutter-plus-noise covariance matrix (CNCM) can be represented as where R c represents the CCM and R n represents the noise covariance matrix. Under the criterion of maximum SINR, the weighting vector of the optimal STAP filter is given by where sðf d;0 ; f s;0 Þ ¼ s t ðf d;0 Þ ⊗ s s ðf s;0 Þ is the assumed target steering vector, f d,0 and f s,0 are the normalised Doppler frequency and normalised spatial frequency of the target, respectively. However, the ideal covariance matrix R is unknown in practical applications, and it needs to be estimated from the IID secondary data, that is, where L denotes the number of available IID training samples.

| Grid-based SR-STAP algorithms
In the framework of grid-based SR-STAP algorithms, we first need to discretise the continuous angle-Doppler plane into sufficiently small grid points. Here, the whole angle-Doppler plane is uniformly discretised into Q ¼ N s N d grid points, where N s ¼ ρ s N is the number of grids in the spatial domain, is the number of grids in the temporal domain, and ρ s ; ρ d > 1 denote the resolution scales. Let s q ðq ¼ 1; 2; ⋯; QÞ denote the space-time steering vector corresponding to the qth grid point, then the STAP dictionary matrix Φ ∈ ℂ MN�Q , which is defined as the set of space-time steering vectors corresponding to all discrete grid points, can be represented as Let X ¼ ½x 1 ; x 2 ; ⋯; x L � ∈ ℂ MN�L denote the received data matrix composed of L IID training samples, and then the SR-STAP model in the case of multiple samples can be formulated as: where N ¼ ½n 1 ; n 2 ; ⋯; n L � ∈ ℂ MN�L is the noise matrix, and B ¼ ½β 1 ; β 2 ; ⋯; β Q � T ∈ ℂ Q�L is the unknown angle-Doppler profile matrix with each row representing a possible clutter component [27].
Assuming L IID training samples share the same intrinsic sparsity of clutter, which means the support index set for every column of B is identical. Then the unknown angle-Doppler profile matrix can be recovered by solving the following ℓ 2,1 mixed-norm minimisation problem [37,38].
where the superscript sr denotes the CNCM estimated using SR-technologies, and Γ ¼ diagðγ 1 ; γ 2 ; ⋯; γ Q Þ ∈ ℝ Q�Q is a diagonal matrix with diagonal elements representing the power values of corresponding grid points, in which

| PROPOSED ALGORITHMS
In this section, we will derive a new gridless MNM-based STAP algorithm (i.e. MNM-STAP) to eliminate the off-grid effect as well as a fast implementation method (i.e. FMNM-STAP) by utilising the framework of ADMM to improve its computational efficiency.

| MNM-STAP algorithm
Before deriving the MNM-STAP algorithm, we first look at Equation (11). It can be observed from Equation (11) that the parameters we are really interested in are fγ q g Q q¼1 which represents the power values of corresponding clutter components and can be regarded as functions of the row-norms of B, that is, where ω q ∈ ℂ and c q ∈ ℂ L�1 , similar to β q . According to the inequality of arithmetic and geometric means, the ℓ 2 norm of β q must satisfy [35].
with equality holding if and only if |ω q | 2 ¼ ‖c q ‖ 2 2 . By exploiting (12), the row-norm ‖β q ‖ 2 [35] can be rewritten as And the optimal solution of Equation (13) holds that Based on Equation (13), the ℓ 2;1 mixed-norm of the angle-Doppler profile matrix B ¼ ½β 1 ; β 2 ; ⋯; β Q � T can be represented as [35]. (15) into (10) yields the following equivalent optimisation problem For a fixed matrix Ω, the minimiser C has the following closed-form expression Substituting Equation (17) into Equation (16) and carrying out some basic reformulations, we can obtain the following representation of Equation (16) as By exploiting Equation (7) and the expression Equation (18) can be reformulated as the following equivalent form Based on Equations (14) and (19), we can observe that the matrix Λ ¼ diagðλ 1 ; λ 2 ; ⋯; λ Q Þ in Equation (20) represents the row-norms of B on its diagonal with (20) is also known as the compact formulation of Equation (10) in [35]. The derivation process of Equation (20) in this article is inspired by [35], where we extend the one-dimensional (1D) frequency estimation method in [35] to two-dimensional (2D) STAP model. Once the optimal solution of Equation (20) is obtained, we can reconstruct the CNCM aŝ However, it should be noted that, up to now, the reconstruction of R in Equation (21) is still based on the discretisation, which just changes the parameters of interest from B to the row-norms of B.
To avoid the discretisation and eliminate the off-grid effect, let's reconsider the optimisation problem in Equation (20). Since the radar system considered in this article is equipped with a ULA, and the PRF is constant, the atom in Φ, that is, s q ; q ¼ 1; ⋯; Q is thus a 2D-Vandermonde vector. Based on the Λ ∈ D þ and 2D-Vandermonde structure of s q , the product ΦΛΦ H in Equation (20) is therefore a matrix with PSD and block-Toeplitz properties. Different from the grid-based SR-STAP approaches, we here parameterise the matrix [39]. Then the M � M block-Toeplitz matrix ΦΛΦ H can be parameterised by T ðU Þ from U as where the block matrix U m is a N � N Toeplitz matrix defined as Substituting Equation (22) into Equation (20), we get the following gridless minimisation problem where trðΛÞ ¼ 1 MN trðT ðU ÞÞ. Compared to the grid-based optimisation problem Equation (20), the derived gridless minimisation problem Equation (24) transforms the interest parameters from the row-norms of B into a parameter matrix U by exploiting the PSD and block-Toeplitz properties of ΦΛΦ H , which avoids the discretisation and the corresponding off-grid effect.
In addition, according to the Schur complement [40,41], the convex optimisation problem in Equation (24) Then, an SDP solver, such as SDPT3 [42], can be used here to calculate the solution of Equation (25). When the parametric matrix U is obtained, we can directly reconstruct T ðU Þ using Equations (22) and (23). However, since the di- (25) is therefore not an accurate estimation of CCM. In other words, the clutter subspace is well estimated from T ðU Þ, but the corresponding power values require some form of correction. Consequently, we here re-estimate the CCM with the received data X and eigenvectors of T ðU Þ by [34].
where the superscript gl denotes the CCM is estimated by the proposed gridless algorithm, and Ρ is given by T ðU Þ with eigendecomposition, that is, Then, the final CNCM can be calculated as and the weighting vector of the proposed algorithm can be calculated as For clarity, the procedure of the proposed MNM-STAP approach is summarised in Table 1.

| FMNM-STAP algorithm
In the previous subsection, the MNM-STAP approach is derived to eliminate the off-grid effect. In the proposed approach, the optimisation problem Equation (25) is solved by the off-the-shelf solver SDPT3, where the interior-point method is performed. Note that, however, this solver tends to be slow, especially for large-scale problems. To reduce the computational time of solving Equation (25), we next derive a fast iterative method (i.e. FMNM-STAP) by utilising the framework of ADMM [36], which can be viewed as a version of the method of multipliers where a single Gauss-Seidel pass [43][44][45] over decomposed variants is used instead of the usual joint minimisation. Before proceeding to FMNM-STAP, we rewrite Equation (25) in the following equivalent form.
Then the augmented Lagrangian expression of Equation (31) is defined as where ϒ ∈ ℂ ð2MNþLÞ�ð2MNþLÞ is the Lagrangian multiplier and ρ > 0 is a penalty parameter, 〈ϒ; Θ 〉 ¼Reðtrðϒ H ΘÞÞ. In the (t+1)-th iteration, the implementation of ADMM includes the following updated steps: In order to implement Equation (33), we first partition the matrices Θ and ϒ as follows: Computing the partial derivatives of Equation (33) with respect to Z and T ðU Þ respectively, we obtain , and let Equations (38) and (39) be 0, we then have By exploiting Equation (41) and the prior block-Toeplitz property of ðT ðU ÞÞ tþ1 , we can obtain where E M;N is a ð2M − 1Þ � ð2N − 1Þ matrix, of which the (M,N)-th element is 1 and the other elements are 0. Besides, the symbol T * ð·Þ in Equation (42)  Obviously, the optimal solution of Equation (44) without the constraint Θ ≥ 0 is Let Θ tþ1 un ¼ Gdiagðfδ g gÞG −1 ; g ¼ 1; 2; ⋯; 2MN þ L, then the updated step for Θ is given as where fδ g g þ is to set all negative eigenvalues to zero. Note that Θ tþ1 un in Equation (45) is a block diagonal matrix, thus a blockmatrix eigendecomposition method [46] can be used here to further reduce the computational complexity, where the computational complexity of updating Θ can be reduced from oðð2MN þ LÞ 3 Þ to oððMN þ LÞ 3 þ ðMNÞ 3 Þ.
Finally, we can update ϒ according to Equation (35). The iterative steps will be terminated if either of the following conditions is reached: (a) The number of iterations t reaches the predefined maximum iterative number, that is, t ≥ T max . (b) The relative change in U is smaller than a predetermined threshold ε s , that is, For clarity, we summarise the updated steps in Table 2.

| NUMERICAL EXPERIMENTS
In this section, numerical experiments are conducted to assess the clutter suppression and target detection performance of the proposed FMNM-STAP approach. For comparison purposes, we also present the performance of conventional OS-SBL-STAP [32], IAA-STAP [21], DDD-STAP [16], the diagonal T A B L E 2 The updated steps of the proposed FMNM-STAP approach Input: Received data X ¼ ½x 1 ; x 2 ; ⋯; x L � and system parameters

Until a stopping criterion is satisfied
Output: The weighting vector w gl via Equation (29) and the obtained U.
The main parameters of the simulated scenarios are listed in Table 3, the number of clutter patches is N c ¼ 360, the amplitude for each clutter patch follows a complex Gaussian distribution, the thermal noise is modelled as a zero-mean Gaussian process with σ 2 n ¼ 1, and the crab angle is assumed to be 0°. Additionally, in the following experiments, for the proposed FMNM-STAP approach, η ¼

| Convergence of the proposed FMNM-STAP
We first assess the convergence of the proposed FMNM-STAP algorithm from the perspective of numerical experiments via the metric of average SINR loss, which is defined as the mean of the SINR loss for f d ∈ ð−0:5; −0:1Þ ∪ ð0:1; 0:5Þ. Here the value of SINR loss is calculated as [4].
where w is the estimated STAP weighting vector, and R is the known CNCM. The simulated result is shown in Figure 2, which plots the average SINR loss versus the number of iterations. It can be observed from Figure 2 that the performance of FMNM-STAP increases with the number of iterations, and it can converge to a near-optimal state with about 50 iterations. Moreover, [36] also gives a rigorous mathematical convergence proof for convex-based ADMM methods. As the constructed optimisation problem Equation (25) of FMNM-STAP is well known to be convex, thereby the proposed algorithm can also converge to the optimal solution theoretically. The interested reader is referred to [36] for more detailed analyses.

| Clutter suppression performance
In this experiment, we compare the clutter suppression performance of the proposed algorithms with that of OS-SBL, IAA, LSMI, and OPT approaches via the metric of SINR loss. The results are shown in Figure 3, where the curve of MNM-STAP is omitted as it is visually identical to FMNM-STAP. Figure 3(a) shows the SINR loss versus the normalised Doppler frequency, where we can observe that, in the side-lobe region, the three SR-STAP algorithms (i.e. FMNM, OS-SBL, and IAA) achieve the near-optimal performance, which is much better than that of conventional DDD and LSMI. In the main-lobe region, the proposed FMNM-STAP approach still obtains the near-optimal performance, however, the performance of grid-based SR-STAP methods degrades significantly due to the off-grid effect. Especially, the performance of IAA is even poorer than that of DDD. By iterative searching scheme, the OS-SBL approach can mitigate such off-grid effect and achieve better SINR loss performance compared with IAA. Nevertheless, discretisation is still assumed in OS-SBL; therefore, its performance is lower than that of the proposed gridless FMNM-STAP method. Figure 3(b) depicts the average SINR loss versus the number of training samples, namely, the convergence speed. As shown in Figure 3(b), FMNM-STAP and OS-SBL-STAP exhibit the fastest converge speed and achieve a suboptimal performance with only 2 training samples. In addition, with the increase of the number of training samples, the proposed FMNM-STAP method can achieve optimal performance, while the OS-SBL-STAP method cannot due to the off-grid effect, which is consistent with the previous analysis.  To better illustrate the advantage of the proposed FMNM-STAP method, Figure 4(a)-(d) plot the Capon spectra of corresponding SR-STAP and OPT-STAP approaches. From Figure 4(a)-(d), we can see that the spectrum of the proposed FMNM-STAP is the closest to the optimal without any clutter power leakage, but the spectra of two grid-based STAP have much expansion in the main-lobe area, especially the IAA-STAP method, which correspond to the simulation results in

| Target detection performance
In this subsection, we access the target detection performance via the probability of detection (PD) versus signal-tonoise (SNR) curves, which are achieved by utilising the adaptive matched filter (AMF) detector [47]. In this simulation, the probability of false alarm rate (PFA) is set as 10 −3 , the threshold and probability of detection estimates  are based on 1000 samples. The target is assumed in the boresight with the normalised Doppler frequency 0.35 in Figure 5(a) and 0.08 in Figure 5(b), respectively. As depicted in Figure 5(a), the detection performance of FMNM, OS-SBL, and IAA is close to the optimal performance, and is much better than that of DDD and LSMI in the case of f d;0 ¼ 0:35 (i.e. the target speed is relatively fast). In the case of f d;0 ¼ 0:08 (i.e. the target speed is relatively slow), as shown in Figure 5(b), the proposed FMNM still obtains the near-optimal detection performance; however, the detection performance of OS-SBL and IAA degrades to a great extent, especially the IAA, which is even much lower than that of DDD. The reason is that discretisation of the angle-Doppler plane results in the off-grid effect, which expends the clutter spectra of grid-based OS-SBL and IAA approaches and yields wider clutter notches as well as weak detection performance for slow-moving targets.
To further investigate the detection performance of the proposed FMNM-STAP approach, we next assess the receive operating characteristic (ROC) (i.e. PD vs. PFA) of each method with different SNR levels, that is, SNR = 0dB, SNR = −4dB, and SNR = −10dB. The detection It can be seen that for fast-moving targets (see Figure 6 (a)-(c)), the proposed FMNM-STAP method shows the best detection performance in all three cases, although the detection performance decreases in varying degrees with the reduction of the target SNR. Besides, the detection performance of OS-SBL-STAP and IAA-STAP is very close to the proposed method, and is much better than the conventional STAP methods, since the grid-based SR-STAP can estimate CCM more accurately than the conventional ones with limited training samples.
For slow-moving targets (see Figure 7 (a)-(c)), the proposed algorithm still shows the best performance, however, the detection performance of OS-SBL-STAP and IAA-STAP degrades significantly. This is because the two gridbased approaches are not good at estimating the clutter components around the main-lobe region. As previously presented (see Figure 4 (c) and (d)), in the main-lobe region, the estimated spectra of OS-SBL-STAP and IAA-STAP are wider than that of the optimal and proposed approaches due to the off-grid effect, which might result in the selfcancelation of slow-moving targets and therefore degrade its detection performance. In addition, we want to mention that the proposed FMNM-STAP method also exhibits a robust performance in the case of non-Gaussian distributed clutter.

| Computational complexity analysis
Finally, we analyse the computational complexity of the proposed FMNM-STAP method and compare it with the other two SR-STAP methods, that is, OS-SBL-STAP and IAA-STAP. As the calculation of STAP filters for the three algorithms share the same complexity, we only consider the computational complexity of CNCM estimation in terms of the number of complex multiplications. The results are shown in Table 4, where N FMNM , N OS−SBL , and N IAA denote the number of iterations of FMNM-STAP, OS-SBL-STAP and IAA-STAP.
To make the computational complexity comparison more illustrative, we examine the computational time of each algorithm via tic and toc instruction embedded in MATLAB. Note that all simulations are implemented by MATLAB 2017b on a computer equipped with Intel Xeon E5-2620 CPU @2.10 GHz�2 and 64GByte RAM. Figure 8(a) presents the computational time of different methods versus the system DOFs, where we set ρ s ¼ ρ d ¼ 4. Figure 8(b) plots the computational time versus the resolution scale ρ s , where we set M ¼ N ¼ 8 and ρ d is always equal to ρ s . It can be observed from Figure 8 that the computational time of the proposed FMNM-STAP approach is much lower than that of OS-SBL-STAP and IAA-STAP with the increase of DOFs, and it is independent of the resolution scale, which illustrates that the proposed FMNM-STAP method has an advantage in practical real-time applications.
Besides, we would like to briefly analyse the memory requirement of different methods. As MATLAB stores all the other two SR-STAP methods since it needs to estimate the mean value (the memory requirement of which is 0.1611 MB) and the covariance matrix (the memory requirement of which is 17 MB) of sparse coefficients B in each iteration.

| CONCLUSIONS
We studied the problem of clutter covariance matrix estimation with limited training samples. To avoid the discretisation of the angle-Doppler plane, we first propose a gridless MNM-STAP approach by exploiting the PSD and block-Toeplitz properties of CCM, and the compact formulation of MNM. To improve the computational efficiency, we also derive a fast iterative method, that is, FMNM-STAP, by exploiting the framework of ADMM. Numerical results demonstrate that the proposed FMNM-STAP method can achieve better performance of clutter suppression and target detection with lower computational complexity.