A Knowledge-Based Auxiliary Channel STAP for Target Detection in Shipborne HFSWR

: The broadened ﬁrst-order sea clutter in shipborne high frequency surface wave radar (HFSWR), which will mask the targets with low radial velocity, is a kind of classical space–time coupled clutter. Space–time adaptive processing (STAP) has been proven to be an effective clutter suppression algorithm for space-time coupled clutter. To further improve the efﬁciency of clutter suppression, a STAP method based on a generalized sidelobe canceller (GSC) structure, named as the auxiliary channel STAP, was introduced into shipborne HFSWR. To obtain precise clutter information for the clutter covariance matrix (CCM) estimation, an approach based on the prior knowledge to auxiliary channel selection is proposed. Auxiliary channels are selected along the clutter ridge of the ﬁrst-order sea clutter, whose distribution can be determined by the system parameters and regarded as pre-knowledge. To deal with the heterogeneity of the spreading ﬁrst-order sea clutter, an innovative training samples selection approach according to the Riemannian distance is presented. The range cells that had shorter Riemannian distances to the cell under test (CUT) were chosen as training samples. Experimental results with measured data veriﬁed the effectiveness of the proposed algorithm, and the comparison with the existing clutter suppression algorithms showed the superiority of the algorithm.


Introduction
In recent years, with the increasing emphasis on marine resources, long-distance detecting, tracking, and early warning of maritime targets have become the top priority of the current maritime rights protection work. High frequency surface wave radar (HFSWR), which achieves over-the-horizon detection of sea-surface targets and sea-skimming targets at low altitude, is an efficient way for ocean surveillance. Compared with shore-based HF-SWR, shipborne HFSWR has better survivability, flexibility, and maneuverability. However, the movement of the platform causes the spread of the first-order sea clutter, which will obscure targets with low radial velocity like ships. How to detect the targets submerged in the clutter has become the focus of research on signal processing in shipborne HFSWR.
The premise of clutter suppression is that the clutter and the target are distinguishable in a certain dimension. The characteristics of the clutter and target should be studied before designing a clutter suppression approach. The properties of the spread spectrum of first-order sea clutter has been analyzed and a technique for ship target detection by using the sum-and-difference beamformer was proposed [1]. The space-time distribution of the first-order sea clutter for shipborne HFSWR has been demonstrated and the linear expectation operator is denoted by E [·]. The rounding operator, which indicates rounding to the nearest integer, is denoted by int{·}. The absolute value of a scalar or the determinant of a square matrix is denoted by |·|. The 2 norm and Frobenius norm of a matrix are denoted by · 2 and · F , respectively. The trace of a square matrix is represented by tr(·).

Shipborne High Frequency Surface Wave Radar (HFSWR) Signal Model
Assume that the receiving subsystem of the radar consists of N receiving antennas that form a uniform linear array (ULA). The transmitting subsystem transmits a burst of M pulses during the coherent processing interval (CPI). For each channel and each pulse, the received array data are sampled and preprocessed to obtain range samples. The kth range cell takes the form where x m = [x 1,m x 2,m . . . x N,m ] T , m = 1, 2, . . . M represents the data collected from the mth pulse and [·] T denotes transposition. The received data snapshot x(k) ∈ C MN×1 is the sum of three parts, target, clutter, and noise, denoted s(k), c(k), and n(k), respectively. It is mathematically written as x(k) = s(k) + c(k) + n(k). ( The target snapshot, denoted s(k) ∈ C MN×1 , can be expressed as the product of complex amplitude σ tgt and space-time steering vector v tgt , and can be written as where the space-time steering vector of the target, denoted v tgt ∈ C MN×1 , is the Kronecker product of the spatial steering vector a θ tgt ∈ C N×1 and the temporal steering vector b f tgt ∈ C M×1 , and is written as where the mathematical symbol ⊗ is the Kronecker product. The spatial steering vector and the temporal steering vector take the form, where f s = (d/λ)sinθ tgt and f t = f tgt · t R are the normalized spatial and Doppler frequency, respectively. The symbols λ and d represent the wavelength and the interelement spacing of the ULA, respectively. The azimuth of the target and the Doppler frequency of the target are denoted by θ tgt and f tgt , respectively. The pulse repetition interval (PRI) is denoted by t R . Thus, the covariance matrix of the target signal can be derived as where E[·] indicates the expectation operator. The clutter snapshot c(k) ∈ C MN×1 in HFSWR is mainly composed of the first-order and second-order sea clutter, and the energy of the second-order sea clutter is usually 20-45 dB smaller than that of the first-order [33]. Therefore, the second-order sea clutter can be ignored. The clutter in a certain range cell can be regarded as the superposition of a large number of uncorrelated clutter patches [34,35], which has the expression as follows, (8) where q = 1, 2, . . . , N c denotes the qth clutter patch in the kth range bin; the complex amplitude of the receding (caused by ocean waves going away from the radar) and approaching (caused by ocean waves coming toward the radar) first-order ocean clutter are denoted by σ 1,q and σ 2,q , respectively; their corresponding Doppler frequency are denoted by f 1,q and f 2,q , respectively; and the azimuth of the qth clutter patch in the kth range bin is denoted by θ q . Under ideal conditions, there is a linear relationship between the Doppler frequency f i,q , i = 1, 2 and the azimuth cosine cos θ q of the clutter patch, as described in [2], where v p denotes the platform velocity. The symbol f B ≈ 0.102 f 0 is the first-order Bragg frequency, and is calculated using the system carrier frequency f 0 in MHz. The the clutter covariance matrix (CCM) was analyzed in [34,35] and has the following expression For a fully-developed sea where the energy supplied by the wind is equal to the energy lost in breaking waves, its surface is a random surface. According to [34,36,37], the random amplitudes of the harmonic components of electromagnetic waves scattered off a random surface are uncorrelated. Thus, the first-order sea clutter from different ocean patches can be assumed as uncorrelated. The complex amplitudes of the preceding and approaching Bragg components from the same clutter patch are also assumed to be uncorrelated. The CCM is simplified to The rank of the CCM for shipborne HFSWR is derived in a precise expression in [35], where the factors α and β are calculated as α = 2d λ and β = 2v p t R d ; the PRI t R in β is calculated as t R = 1 f R ; and the operator int{·} indicates rounding to the nearest integer.
The noise snapshot n(k) ∈ C MN×1 can be assumed as a complex Gaussian distribution with zero-mean and the power σ 2 n . The noise covariance matrix is

Knowledge-Based Auxiliary Channel Space-Time Adaptive Processing (STAP)
Auxiliary channel STAP is one of the GSC-structure STAP methods. The core of auxiliary STAP is to estimate the clutter in the main channel using the clutter in auxiliary channels. In this section, a knowledge-based auxiliary channel STAP algorithm is presented for spreading sea clutter suppression in shipborne HFSWR.

Auxiliary Channel Space-Time Adaptive Processing (STAP)
For the kth range cell data x(k), assume that the azimuth of the target and the Doppler frequency of the target are denoted by θ tgt and f tgt , respectively.
The main channel, denoted by x main (k) ∈ C 1×1 , is calculated as where T main ∈ C MN×1 is the transformation matrix, which has the same representation with the space-time steering vector of the target, Assume that D auxiliary channels are selected, the auxiliary channels, denoted by x aux (k) ∈ C D×1 , are calculated as, where T aux ∈ C MN×D is the transformation matrix for auxiliary channels.
After transforming into the angle-Doppler domain, the main channel data x main (k) is the sum of the target, the first-order sea clutter, and noise. As the following equation denotes, x main (k) = s main (k) + c main (k) + n main (k). (17) If the clutter in the auxiliary channels is related to that in the main channel, it can be estimated from the auxiliary channels. This estimation, denoted by y aux (k) ∈ C 1×1 , is named as the output of auxiliary channels, and can be denoted by the following equation, Subtracting the estimated clutter plus noiseˆ c main (k) +ˆ n main (k) from the main channel and the output of the auxiliary channel STAP processor, y(k) ∈ C 1×1 is where ε(k) = c main (k) + n main (k) −ĉ main (k) −n main (k) is the residual term of the clutter and noise. The weights vector w aux can be determined by the solution of the following least squares problem The solution is given by where R aux is the self-correlation matrix of the auxiliary channel, and r aux−main is the cross-correlation vector between the auxiliary channels and the main channel. It is worth noting that R aux and r aux−main are unknown in practice and the maximum likelihood estimate is commonly used where x aux (k), k = 1, 2, . . . , K, x main (k), k = 1, 2, . . . , K are the secondary data. The framework of the auxiliary channel STAP algorithm is illustrated in Figure 1. As shown in Figure 1, for the input range cell data x(k), after a transformation (Doppler processing and beam forming), the data are converted into the beam-Doppler space, and the auxiliary channels are selected in it. After that, the adaptive weight vector of the auxiliary channels is calculated according to the minimum mean square error (MMSE) criterion and the output of the auxiliary channels is calculated. Finally, the output of the auxiliary channels is subtracted from the main channel output to suppress the clutter.
The performance of the auxiliary channel STAP method is mainly affected by how the auxiliary channels are selected and the estimation accuracy of the weights vector w aux . The selection of auxiliary channels is determined as the transformation matrix T aux and the estimation accuracy of w aux is determined by the training samples used when calculating the self-correlation matrix R aux and cross-correlation vector r aux−main . Therefore, for these two directions, a knowledge-based auxiliary channel selection approach and a training sample selection approach are proposed in Sections 3.2 and 3.3 to improve the performance of the auxiliary channel STAP method.

Knowledge-Based Auxiliary Channels Selection
For shipborne HFSWR, the distribution of the first-order sea clutter is two straight lines on the angle-Doppler plane. Setting f i,q as abscissa and cos(θ q ) as the vertical coordinate, Equation (9) is rewritten as The wavelength λ, the platform velocity v p , and the Bragg frequency f B are constant. Hence, two single straight lines (clutter ridge) are expected in the angle-Doppler domain. One is the negative part y = λ 2v p x + λ 2v p f B , and the other is the positive part y = λ 2v p x − λ 2v p f B . These two straight lines are determined by system parameters and can be seen as prior knowledge.
A data acquisition experiment was conducted in the Yellow Sea of China in September 1998. The data acquisition system consists of a transmitting subsystem, a receiving subsystem, and a data acquisition subsystem. The transmitting subsystem is made of an omnidirectional antenna, and the signal used is frequency modulated interrupted continuous wave (FMICW). The receiving subsystem includes seven antennas that form a uniform linear broadside array. The geometry of the receiving array is illustrated in Figure 2. The whole system was mounted on a barge, which was tugged by a tugboat. As described in [1,3,7,12], the system parameters are shown in Table 1.  Figure 3 illustrates the range-Doppler map and the angle-Doppler map of the measured data. The dashed line in Figure 3a represents the spreading region of the first-order sea clutter. The dashed line in Figure 3b denotes the clutter ridge. For auxiliary channel STAP, auxiliary channels should contain clutter with the same characteristics as the main channel. Due to the wider beam width in the spatial domain, the main channel will contain some clutter components from the clutter ridge. Thus, we can utilize this prior knowledge and select auxiliary channels along the clutter ridge of the first-order sea clutter. Figure 4 shows the sketch map of the main channel and auxiliary channels. If the targets fall in the spreading range of the first-order sea clutter, the clutter plus noise component c main (k) + n main (k) in the main channel can be approximated from the clutter ridge.
According to the prior knowledge, the auxiliary channels selecting matrix T aux were chosen as a set of the space-time steering vectors of the first-order sea clutter.
where f c i , i = 1, 2, . . . , D and θ c i , i = 1, 2, . . . , D are the Doppler frequency and the azimuth of the clutter patches, respectively. They satisfy Equation (9). The quantity of auxiliary channels, denoted by D, is discussed as follows. The total number of auxiliary channels available is determined by how the clutter ridge is defined. Figure 5 illustrates the theoretical clutter ridge and two different extracted clutter ridges from the measured data.  Taking the theoretical clutter ridge in Figure 5a as an example, the quantity of auxiliary channels available when only one beam channel is selected for a determined Doppler frequency, can be calculated as where d res represents the Doppler resolution. In this situation, the number of auxiliary channels is equal to the number of Doppler channels. According to Equation (9) λd res in the positive part or the negative part. Hence, the quantity of auxiliary channels is int 4v p λd res . This usually numbers between tens and hundreds, which means that there are plenty of auxiliary channels for selecting.
It is worth noting that T aux also works as the blocking matrix and should satisfy T H aux T main = T H aux v tgt = 0. Here, we chose one column of T aux to calculate the inner product v H c i v tgt . Figure 6 illustrates the value of v H c i v tgt when f c i and θ c i changes. M is 1024 and N equals 7. The x-axis represents the difference of the normalized temporal frequency between the target and clutter; y-axis ∆ f s = (d/λ)sinθ tgt − (d/λ)sinθ c i is the difference of the normalized spatial frequency between the target and clutter. It can be seen from the result that as long as ∆ f t is not zero, v H c i v tgt is close to zero. Therefore, if the auxiliary channels and the main channel are not in the same Doppler bin, the transforming matrix used for the auxiliary channels can block the target signal.
Based on the analysis, the principle of auxiliary channel selection is to select auxiliary channels from the clutter ridge and the Doppler frequency of the auxiliary channels should be different from the main channel.

Training Sample Selection
When the auxiliary channels are determined, the performance of the auxiliary channel STAP is determined by the estimation accuracy of the auto-correlation matrix R aux and the cross-correlation vector r aux−main . Assuming that the noise is not related to the signal and clutter, and the blocking matrix T aux can block the target signal well, the estimation accuracy is discussed below.
For the cross-correlation vector r aux−main , it can be decomposed as follows: The first part 1 K ∑ K k=1 c aux (k) s * main (k) may not be equal to a zero-vector due to limited training samples, especially when the target appears in multiple range bins or the target power is relatively strong. Under these circumstances, the target signal of the main channel will leak into the auxiliary channels and lead to target self-cancellation. This phenomenon can be reduced by setting guard range cells around the cell under test (CUT) and selecting training samples without targets.
For the auto-correlation matrix R aux , it can be decomposed aŝ Considering the decomposition results, the estimation accuracy is mainly affected by the training samples used. For shipborne HFSWR, the heterogeneity of clutter in range has been analyzed based on measured data in our previous work [12]. The qualified independent and identically distributed (IID) training samples are limited. How to select the most homogeneous clutter sample in the finite training samples becomes important. The clutter in the training samples should be as similarly distributed as possible to the CUT. To measure the similarity, the Riemannian distance is introduced and a training samples selecting method is proposed here.
First, the explanation of the Riemannian distance is given. Different from the Euclidean distance used in [22,38], the Riemannian distance puts the covariance matrix on the Riemannian manifold, and the length of the shortest curve connecting two matrix R 1 and R 2 is defined as the Riemannian distance. The calculation of the Riemannian distance is given in [39], and the expression is as follows, where tr(·) is the trace of the matrix; · 2 F represents the matrix's Frobenius norm and is calculated as R 2 F =< R, R >= tr RR T ; and the eigenvalue of R Then, calculate the Riemannian distances between the clutter's self-correlation matrix of CUT and other range cells using the following equations After that, sort the distances η k = [η k,1 , η k,2 , · · · , η k,K ] in ascending order and choose L range cells with lower geometric distances. For the convenience of explaining the sample selection approach, a selecting matrix κ(k), which satisfies x(k) = Xκ(k), is defined. It has the following expression Hence, the training sample selecting matrix Λ ∈ C K×L can be written as (2)) · · · κ(index(L))], (35) [ [?] η k , index] = sort(η k , ascend ), (36) where index is the same size as η k and describes the arrangement of the elements of η k into [?] η k along the sorted dimension and [?] η k = η k (index). The selected training samples X training ∈ C MN×L can be denoted as X training = XΛ. (37) As for the number of the training samples, it has a lower limit. The rank of the self-correlation matrixR aux should satisfy where R c represents the local covariance matrix of the clutter (after transforming into the angle-Doppler domain). By Equation (12), the local DOF of clutter can be approximated as where K b and K d are the beams and Doppler bins in the local processing region (LPR), respectively. In this paper, there was only one beam and one Doppler bin in the LPR. Thus, the local DOF of clutter was approximated to 2. As stated by the Reed-Mallett-Brennan (RMB) criterion [40], the secondary data used to estimate the self-correlation matrix and cross-correlation vector should be twice the DOF of the clutter, which means that at least 2rank( R c ) = 4 IID secondary data are needed.
In conclusion, the procedural of the training sample selection approach can be summarized as follows:
Select L (L ≥ 2rank( R c )) range cells, which has the lowest values of η k,l to formulate the selecting matrix Λ using Equation (35), and the training samples are obtained as Equation (37).

Experimental Results with Measured Data
To evaluate the performance of the proposed algorithm, both simulated targets and real target were tested using the measured data. Comparisons with existing algorithms are also demonstrated.

Measured Data with Simulated Target
Four simulated targets with the detailed parameters shown in Table 2 were injected to the measured data. It is worth noting that target 2 was alongside the theoretical position of the negative sea clutter patch where the Doppler frequency was 0.2039 Hz, while target 4 was next to the position of the positive sea clutter patch with the Doppler frequency of 0.2650 Hz, when the azimuth was set to 10 • . These two targets were injected to test the performance when the targets and the first-order sea clutter were relatively close. The expression of signal-to-clutter-plus-noise ratio (SCNR) is as follows, SCNR = 10 log 10 ( P s P c+n ), where P s denotes the power of the signal, and P c+n represents the average power of clutter plus noise. For injected targets, the SCNR is calculated on the target's range cell data x(k).
The injected target signal is s(k) = σ tgt v tgt , the target's power P s is calculated as P s = σ 2 tgt , and the clutter plus noise power P c+n is calculated as P c+n = 1 MN ∑ N n ∑ M m |x nm | 2 .  The Riemannian distances between the CUT and other range cells are illustrated in Figure 7. As described in Table 2, the injected targets were distributed in the positive and negative broadened areas separately. For target 1 and target 2, the clutter ridge with the negative Doppler frequency was used as the auxiliary channels, while the clutter ridge with the positive Doppler frequency was used for target 3 and target 4. The targets were in the 17th range cell and this range cell was chosen as the reference. The first three range cells were in the blind zone and could be ignored. The 16th and 18th range cells were set as guard cells. In order to facilitate the comparison with the conventional JDL algorithm, in which the local process region consists of three doppler bins and three beam bins, L was set was 18. Range cells that had a lower Riemannian distance were selected as training samples, as the points marked with square signs denote in Figure 7. The angle-Doppler and range-Doppler maps are demonstrated in Figure 8. Figure 8b,d denotes the angle-Doppler map and range-Doppler map after clutter suppression, respectively. Compared with Figure 8a,c, processed using 2-dimensional fast Fourier transform and digital beamforming (FFT-DBF) and set as the reference, it is obvious that the broadened first-order sea clutter has been suppressed and targets submerged in the clutter appear.
SCNR Improve = SCNR after − SCNR before (dB).  Figure 9 illustrates the comparison of Doppler profiles with simulated targets between the existing clutter-suppression algorithms (conventional JDL [7], IOW [8], IOP [9]) and the proposed algorithm. The FFT-DBF method was set as the reference. The dashed vertical lines represent the Doppler frequency of the injected targets. For all four injected targets, the proposed algorithm worked well. The conventional JDL failed with target 4, the IOW failed with target 2 and target 4, and the IOP algorithm failed with target 3. The average attenuation of the spreading sea clutter was about 25 dB, while the targets almost maintained power for the proposed algorithm. In addition, the proposed algorithm does not introduce frequency estimation deviation, while the other three algorithms will encounter frequency estimation deviation. The SCNR improvement, which indicates the improvement of the SCNR after clutter suppression, is introduced to quantitatively evaluate the performance of different algorithms. The expression of SCNR improvement is as follows.
For HFSWR, it has a high Doppler frequency resolution and the range resolution is relatively poor, so unlike airborne radar, the Doppler profile is often used to show the performance of clutter suppression. In this paper, the SCNR improvement was calculated based on the Doppler profile. As Equation (40) denotes, P s is calculated as the power of the target's Doppler bin. P c+n is calculated as the average power of several Doppler bins around the target.
The results are illustrated in Table 3. For the simulated targets, the proposed algorithm had the highest average SCNR improvement. The average SCNR improvement reached about 26 dB. When the target and the first-order sea clutter became close to each other (target 2 and target 4), the IOW algorithm failed due to target gain loss. The IOP algorithm also suffered from performance degradation. The conventional JDL algorithm failed for target 4. Fortunately, the proposed algorithm was robust and had good clutter suppression performance while retaining the target gain.

Measured Data with Real Target
In the recorded measured data, there was a non-cooperative target of a passenger liner. The parameters of this target were: range 78 km, azimuth −5 • , and radial velocity 5.68 m/s with the corresponding Doppler frequency of −0.2 Hz [3]. The proposed algorithm was used to suppress the clutter and the results are shown in Figure 10. Figure 10b,d illustrates the angle-Doppler map and range-Doppler map after clutter suppression, respectively. Figure 10a shows the initial angle-Doppler map and Figure 10c illustrates the initial range-Doppler map. The training samples used were the same as the selected samples in Figure 7a.  Figure 11 shows the comparison of Doppler profiles between the existing cluttersuppression algorithms. All four algorithms could suppress the broadening first-order sea clutter, but all of them suffered target energy attenuation. The average attenuation of the spreading sea clutter was about 25 dB, while the real target energy attenuation was about 5 dB. The other three algorithms had more obvious attenuation of clutter energy, but at the same time, the energy attenuation of the target became more serious. The target's SCNR improvement was about 11.54 dB (conventional JDL), 8.51 dB (IOW), 13.91 dB (IOP), and 15.75 dB (proposed algorithm), separately.

Discussion
The proposed auxiliary channel STAP algorithm utilizes prior knowledge of the first-order sea clutter and presents an auxiliary channel selection approach. Unlike the conventional JDL algorithm, which has no special requirements for the distribution of clutter in the spatial and frequency domain, the GSC-structure STAP method is suitable for scenario where clutter has a specific spatial or temporal distribution, particularly for mobile platforms where the distribution of ground clutter or sea clutter can be obtained by the system parameters. The proposed training sample selection method can also be extended to other STAP algorithms where the clutter is not homogenous and training samples need to be selected.

Conclusions
In this paper, an auxiliary channel STAP algorithm was presented for suppressing the spreading first-order sea clutter in shipborne HFSWR. To obtain better clutter suppression results, this paper improved the STAP algorithm from two aspects. First, a method for selecting auxiliary channels was presented, which employs the knowledge of the spacetime distribution of the first-order sea clutter; second, to deal with the heterogeneity of the clutter in range dimension and obtain an accurate estimation of the CCM, a training sample selection method was proposed based on the Riemannian distance between the CUT and other range cells. The range cells that were more similar to the CUT were selected. The performance of the proposed algorithm was evaluated with the measured data. Results showed that for the simulated targets, there was average of 26 dB SCNR improvement after conducting the proposed clutter suppression algorithm, and for the real target, the SCNR improvement reached 16 dB. Compared with FFT-DBF, IOW, IOP, and the conventional JDL algorithm, the proposed algorithm had a better performance when evaluating by the SCNR improvement for both the simulated and real targets, which indicates that the proposed algorithm is effective and superior.