A Novel Non-Stationary Clutter Suppression Approach for Space-Based Early Warning Radar Using an Interpulse Multi-Frequency Mode

Abstract

The non-stationary clutter of space-based early warning radar (SBEWR) is more serious than that of airborne early warning radar. This phenomenon is primarily attributed to the Earth’s rotation and range ambiguity. The increase in clutter degrees of freedom (DOFs) and the significant widening of the clutter suppression notch are not conducive to moving target detection near main lobe clutter. This paper proposes an effective approach to suppress non-stationary clutter based on an interpulse multi-frequency mode for SBEWR. Using the orthogonality of the uniform stepping frequency signal, partial range ambiguity can be effectively suppressed, and the clutter DOFs will be reduced. Subsequently, joint pitch-azimuth-Doppler three-dimensional spacetime adaptive processing and slant range preprocessing are used to perform clutter suppression. This combination not only curtails the estimation error associated with the clutter covariance matrix but also enhances the overall detection capabilities of the system. The simulation results verify the effectiveness of the proposed approach.

1. Introduction

Compared with airborne early warning radar (AEWR) and ground-based early warning radar, the most significant difference of space-based early warning radar (SBEWR) is that its orbit height is located in space. This feature results in unique advantages for SBEWR which are different from those of AEWR and ground-based radar. Firstly, SBEWR employs a satellite as carrying platform. It will not be restricted by national boundaries and climate, and its security is higher. Secondly, the high orbit of SBEWR can obtain a wider detection area and a longer warning time. In addition, SBEWR adopts a downward-looking mode to detect ground and air targets, and thus the influence of ground shelters is smaller [1,2,3,4]. While SBEWR boasts various advantages, it also faces the challenge of clutter suppression. The vast detection area leads to more serious range ambiguity. And influenced by the Earth’s rotation, the non-stationary clutter will be more serious. The original overlapping range ambiguity will be separated, and the main lobe clutter will widen. Therefore, the clutter degrees of freedom (DOFs) will increase. Additionally, the radar needs to increase the pulse repetition frequency (PRF) to reduce clutter Doppler ambiguity. This will also lead to more serious range ambiguity for the clutter, and slow-moving targets located near the main lobe clutter will be submerged in the clutter. Based on these problems, the suppression notch of the main lobe clutter will be widened. If the target is located near the notch, then the energy will be suppressed, and the target signal cannot be detected correctly. Therefore, non-stationary clutter suppression is a key challenge for SBEWR [5,6,7,8,9].

AEWR often uses compensation, interpolation and prediction methods to solve the non-stationary clutter of a non-side-looking array [10,11,12,13,14,15]. But in the case of range ambiguity, it will instead aggravate the non-stationary clutter. These methods are not applicable to SBEWR, whose range ambiguity is more severe. The authors of [2,16] elaborated upon the effects of Earth’s rotation and range ambiguity on SBEWR based on the geometric model of space-based radar (SBR). The crab angle and crab amplitude were introduced to measure the influence of Earth’s rotation. They found that the effects of Earth’s rotation on clutter are related to the orbital inclination and the latitude of the nadir point of the radar platform. The essence of range ambiguity is to increase the clutter DOFs, and the original controllable system DOFs are not enough to suppress range ambiguity clutter [17,18,19]. There are two main ideas for solving non-stationary clutter: one is to reduce the clutter DOFs, and the other is to increase the system DOFs. The increase in system DOFs will increase the number of training samples required for clutter covariance matrix estimation, which can be solved with dimension reduction STAP [8,9]. Therefore, finding methods to reduce clutter DOFs is the first task. The authors of [16] proposed the orthogonal pulse method, which transmits a string of mutually orthogonal pulse signals and uses the orthogonality of the echoes of the range ambiguity for matched filtering, thus separating the range ambiguity clutter from the echo to be measured. This method requires a high level of waveform orthogonality, and the number of orthogonal waveforms is limited.

The proposed new system of a waveform diversity array changes the way radar obtains information. The signals of each transmitting antenna differ in space, time, frequency, polarization and modulation mode, increasing the DOFs of the radar system. Radar can better obtain information and optimize target detection performance. It is mainly divided into frequency division, time division and code division in terms of the technical implementation [20,21,22,23]. By introducing the DOFs of the transmitting dimension, the transmitting dimension, jointly with the receiving dimension, distance dimension, Doppler dimension and polarization dimension, can improve the flexibility of radar system resources, and the new method of multi-dimensional domain signal processing is the main development direction of future radar to obtain target information. Domestic and foreign scholars have proposed the concept of three-dimensional spacetime adaptive processing (3D-STAP), which introduces a coding dimension, waveform dimension, pitch dimension and carrier frequency dimension through its system design in addition to two dimensions of space and time. The controllable DOFs of the new dimensions are used to suppress the range ambiguity of SBEWR [24,25,26,27]. A multi-domain cascaded adaptive processing method was proposed in [28]. It uses a digital antenna array to suppress some of the clutter energy by introducing the transmit pulse coding domain. In addition, the rest of the influential clutter energy is suppressed by introducing the elevation receiving channel, and then it carries out two-dimensional adaptive processing. Due to the complexity of orthogonal coding and orthogonal waveform design and the limited number of orthogonal signals, the range ambiguity suppression ability is limited.

Compared with code division and time division, frequency division quadrature is easier to implement. Based on this, the concept of frequency diversity is now a popular research area. A frequency diverse array (FDA) forms an angle-distance-time-dependent emission direction map by introducing a small carrier frequency offset between the array elements. It uses this distance-angle-dependent direction map to improve the STAP performance of range ambiguity. The range ambiguity clutter at different distances can be separated by the DOFs of the distance dimension. Forming a notch in the emission direction map is used to suppress the range ambiguity clutter echoes [15,29,30,31]. To make better use of the DOFs of the transmitter, many studies have combined multi-input multi-output (MIMO) with an FDA. Element pulse coding-MIMO (EPC-MIMO) radar is a technique for phase modulating each transmitting array element between pulses. It can achieve waveform diversity using a transmitting beam that varies with slow time [32]. EPC-MIMO introduces new DOFs in the frequency domain of the transmit dimension. The matched filtering of the original data can recover the transmitting DOFs of the FDA due to the use of orthogonal waveforms. Beam forming at the receiver side and distance-dependent compensation are performed so that the echo signals can be distinguished in the space–frequency domain. By performing beam forming at the transmitter side on the compensated data, the suppression of range ambiguous signals is achieved, and thus the echo signal in each range gate can be extracted. In [33], FDA-phase MIMO was proposed on the basis of the traditional STAP. The multi-domain joint reduced STAP not only achieves a better clutter suppression performance, but it also reduces the number of required auxiliary channels. In addition to setting the carrier frequency offset between transmitting arrays, intra-pulse multi-frequency mode is also a frequency diversity mode which was first proposed in the field of synthetic aperture radar (SAR). The authors of [34] used multiple sub-pulses whose frequency bands were different to improve the performance of range ambiguity suppression, called multi-frequency sub-pulses (MFSPs). Since the range ambiguity is assigned in non-overlapping frequency bands, range ambiguity suppression is achieved by frequency domain bandpass filtering rather than spatial filtering. However, beam modulation between sub-pulses requires fast phase switching between sub-pulses. What is more, the designs of beam pointing and the sub-pulse bandwidth are more complicated. Instead of discrete beam modulation from sub-pulse to sub-pulse, the intra-pulse frequency scanning array technique allows the transmission beam of the pitch dimension to continuously scan a single transmit pulse [35,36,37]. SSF-STAP was proposed in [38], which uses an intra-pulse multi-frequency mode to introduce a dimension related to the slant range. It reduces the influence of range ambiguity clutter through a two-step adaptive clutter suppression method. The sub-pulse number of the intra-pulse multi-frequency mode is limited, only being suitable for airborne radar with a few range ambiguities, and the clutter suppression performance of SBEWR will be degraded.

The above methods introduce the change in carrier frequency and bring about our new idea. This paper proposes an effective approach to suppressing non-stationary clutter based on the inter-pulse multi-frequency mode for SBEWR, which can suppress some of the range ambiguity. The other parts of this paper are arranged as follows. Section 2 establishes a transceiver model for inter-pulse multi-frequency signals based on the signal model of the single-carrier frequency. Section 3 utilizes the orthogonality of inter-pulse multi-frequency signals and matched filtering at the receiving platform to suppress some of the range ambiguity, which reduces the clutter DOFs. Section 4 discusses the clutter suppression process of the inter-pulse multi-frequency mode and analyzes 3D-STAP and covariance matrix estimation in detail. Combining slant range preprocessing with 3D joint domain localization (JDL), STAP is proposed to perform clutter suppression, and it further reduces the clutter DOFs. Section 5 shows simulations to illustrate the effectiveness of the proposed approach in this paper. Section 6 concludes the paper.

2. Signal Model

The SBEWR orbit is located in space, and it employs beam scanning techniques in both the pitch and azimuth directions for downward detection. The signal is modeled as a one-dimensional uniform line array in side-looking mode, as shown in Figure 1. The azimuth beam scanning range of the phased array is about ${120}^{\circ }$. Limited by the angle range of the antenna beam scanning, its detection area is a partial sector of concentric circles with the nadir point as the center. The satellite orbit is approximated as a standard circular orbit. Its orbital height is H. The pitch angle and azimuth angle between the radar platform and the ground scattering unit are ${\theta }_{EL}$ and ${\theta }_{AZ}$, respectively. The echo signal of SBEWR is two-dimensional, where the time domain is composed of the phase generated by the pulse interval and the space domain is composed of the phase generated by the array spacing of the phased array antenna at the receiving end.

SBEWR transmits the narrow-band signal

$$x\left(t\right)=AW\left(\frac{t}{T_p}\right)exp\left\{j2\pi f_{c}t\right\}$$

(1)

where A is the amplitude of the transmitted signal, $W\left(·\right)$ is the envelope of the transmitted signal, $T_p$ is the bandwidth of the transmitted signal, $f_c$ is the carrier frequency of the transmitted signal and t is the fast time. Seventy percent of the Earth’s surface is ocean, and thus the clutter signals of satellite radar are usually sea clutter. In the subsequent study of this paper, the fluctuation of sea clutter within a CPI is ignored.

When the transmitted signal reaches the ground scattering unit q, the echo returns to the receiving end. The echo signal from the kth pulse received by the nth array element of the receiving end is

$$x_{n,k,q}\left(t\right)=\sum _{m=1}^N\begin{array}{c}{A}_{n,k,q}W\left(\frac{t-{\Delta }_{n,k,m,q}}{T_p}\right)·exp\left\{j2\pi f_c\left(t-{\Delta }_{n,k,m,q}\right)\right\}\hfill \end{array}$$

(2)

where m represents the sequence number of the array elements of the transmitting end and ${\Delta }_{n,k,m,q}$ is the time delay of the signal from transmitting end to the receiving end, which is specifically expressed as

$${\Delta }_{n,k,m,q}=\frac{2R_{n,k,m,q}}{C}$$

(3)

The slant distance of the echo which is received by the nth array element of the receiving end and from the kth pulse transmitted by the mth array element of the transmitting end is

$$\begin{array}{c}\begin{array}{c}{R}_{n,k,m,q}=\frac{R_{T,m,q}\left(t_k\right)+R_{R,n,q}\left(t_k\right)}{2}\hfill \\ \approx R_{0,q}-V{t}_k{\rho }_{c}cos{\varphi }_q^{{}^{\prime }}-\frac{d_m}{2}cos{\varphi }_q-\frac{d_n}{2}cos{\varphi }_q+\frac{{\left(V{t}_k+d_n\right)}^2}{4R_{0,q}}+\frac{{\left(V{t}_k+d_m\right)}^2}{4R_{0,q}}\hfill \end{array}\hfill \end{array}$$

(4)

where $d_m$ is the relative distance between the mth element of the transmitting end and the reference element $d_n$ is the relative distance between the nth element of the receiving end and the reference element. If the first element is the reference element, then $d_n$ is $d_n=(n-1)d$. Here, d is element spacing, $R_{0,q}$ is the initial slant distance between the reference array element of the radar platform and the scattering unit q, $t_k=(k-1)T_r$ is the slow time, $f_r=1/T_r$ is the PRF, ${\varphi }_q$ is the incident cone angle, which is $cos{\varphi }_q=sin{\theta }_{EL,q}cos{\theta }_{AZ,q}$, ${\rho }_c$ is the crab amplitude, $cos{\varphi }_q^{{}^{\prime }}$ is the cosine of the incident cone angle, which is corrected by crab angle, $cos{\varphi }_q^{{}^{\prime }}=sin{\theta }_{EL}cos\left({\theta }_{AZ}+{\varphi }_c\right)$ [2], and V is the radar velocity. The direction of the velocity is tangent to the orbit at that time. During the coherent accumulation time, the moving distance of the platform is much less than the slant distance (i.e., $V{t}_k+d_n\ll R_{0,q}$ and $V{t}_k+d_m\ll R_{0,q}$). Therefore, the slant distance can be simplified as follows:

$$R_{n,k,m,q}\approx R_{0,q}-V{t}_k{\rho }_{c}cos{\varphi }_q^{{}^{\prime }}-\frac{d_m}{2}cos{\varphi }_q-\frac{d_n}{2}cos{\varphi }_q$$

(5)

By substituting the above equation into Equation (3), Equation (6) can be obtained:

$$x_{n,k,q}\left(t\right)=\sum _{m=1}^N\begin{array}{c}{A}_{n,k,q}W\left(\frac{t-{\Delta }_{n,k,m,q}}{T_p}\right)exp\left\{j2\pi f_{c}t\right\}exp\left\{-j2\pi f_c\frac{2R_{0,q}}{C}\right\}exp\left\{j2\pi f_c\frac{d_{m}cos{\varphi }_q}{C}\right\}\hfill \\ exp\left\{j2\pi f_c\frac{d_{n}cos{\varphi }_q}{C}\right\}exp\left\{j2\pi f_c\frac{2V{t}_k{\rho }_{c}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}\hfill \end{array}$$

(6)

The echo signals are processed. Firstly, the phase information of the transmitting end is represented by the transmitting direction map ${\omega }_{T,k,q}$ as shown in (7):

$$\begin{array}{c}{x}_{n,k,q}\left(t\right)=A_{n,k,q}W\left(\frac{t-{\Delta }_{n,k,m,q}}{T_p}\right)w_{T,k,q}exp\left\{j2\pi f_{c}t\right\}exp\left\{-j2\pi f_c\frac{2R_{0,q}}{C}\right\}\hfill \\ \begin{array}{ccc}& & \end{array} exp\left\{j2\pi f_c\frac{d_{n}cos{\varphi }_q}{C}\right\}exp\left\{j2\pi f_c\frac{2V{t}_k{\rho }_{c}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}\hfill \end{array}$$

(7)

Then, after pulse compression and Fourier transform of the fast time dimension, the echo signal can be obtained as shown in Equation (8), where the amplitude information of the echo is uniformly represented by $A$:

$$\begin{array}{c}{x}_{n,k,q}=A_{n,k,q}{w}_{T,k,q}exp\left\{-j2\pi f_c\frac{2R_{0,q}}{C}\right\}exp\left\{j2\pi f_c\frac{d_{n}cos{\varphi }_q}{C}\right\}exp\left\{j2\pi f_c\frac{2V{t}_k{\rho }_{c}cos{\varphi }_q^{{}^{\prime }}}{C}\right\}\hfill \\ \begin{array}{cc}& \end{array} =A_{n,k,q}{w}_{T,k,q}exp\left\{j2\pi f_c\frac{d_{n}cos{\varphi }_q}{C}\right\}exp\left\{j2\pi f_c\frac{2V{t}_k{\rho }_{c}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}\hfill \end{array}$$

(8)

The complete echo signal vector of the scattering unit q is

$$\begin{array}{c}{\mathbf{x}}_q=\left[\begin{array}{cccccc}{x}_{11,q}& x_{21,q}& \cdots & x_{N1,q}& \cdots & x_{1k,q}\end{array}\right.{\left.\begin{array}{cccc}\cdots & x_{Nk,q}& \cdots & x_{NK,q}\end{array}\right]}^T\hfill \\ ={\mathrm{A}}_q{\mathbf{s}}_t\otimes {\mathbf{s}}_s\hfill \end{array}$$

(9)

where ${\mathbf{s}}_t$ and ${\mathbf{s}}_s$ are the time steering vector and space steering vector in Equation (10), respectively. The phase information of the steering vector includes the Doppler and angle of the echo signal:

$$\left\{\begin{array}{c}{\mathbf{s}}_t={\left[\begin{array}{ccccc}1& \cdots & exp\left\{j\pi \frac{\left(k-1\right)f_c·4V{T}_r{\rho }_{c}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}& \cdots & exp\left\{j\pi \frac{\left(K-1\right)f_c·4V{T}_r{\rho }_{c}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}\end{array}\right]}^T\\ {\mathbf{s}}_s={\left[\begin{array}{cccccc}1& exp\left\{j2\pi f_c\frac{dcos{\varphi }_q}{C}\right\}& \cdots & exp\left\{j2\pi f_c\frac{\left(n-1\right)dcos{\varphi }_q}{C}\right\}& \cdots & exp\left\{j2\pi f_c\frac{\left(N-1\right)dcos{\varphi }_q}{C}\right\}\end{array}\right]}^T\end{array}\right.$$

(10)

According to Equations (9) and (10), the echo signal of the traditional single-carrier frequency SBEWR can be expressed as the Kronecker product of the space steering vector and the time steering vector. When SBEWR adopts an interpulse multi-frequency mode, the carrier frequency of different pulses is different, as shown in Figure 2. The pulse width is $\tau$. The pulse repetition interval is $T_r$.

The carrier frequency of the kth transmitted signal is $f_k$, which is substituted into Equation (8), and the echo signal is

$$\begin{array}{c}{x}_{n,k,q}={\mathrm{A}}_{n,k,q}{w}_{T,k,q}exp\left\{-j\pi f_k\frac{4R_{0,q}}{C}\right\}exp\left\{j\pi f_k\frac{2d_{n}cos{\varphi }_q}{C}\right\}exp\left\{j\pi {\rho }_c{f}_k\frac{4V{t}_{k}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}\hfill \end{array}$$

(11)

The echo signal vector of the ground scattering unit q is

$${\mathbf{x}}_q={\left[\begin{array}{cccccccccc}{x}_{11,q}& x_{21,q}& \cdots & x_{N1,q}& \cdots & x_{1k,q}& \cdots & x_{Nk,q}& \cdots & x_{NK,q}\end{array}\right]}^T$$

(12)

Since the carrier frequency of interpulse multi-frequency SBEWR changes, and the spacetime coupling relationship changes, the signal vector can be expressed as

$${\mathbf{x}}_q=kron\left({\mathbf{b}}_r\odot {\mathbf{b}}_t,{\mathrm{I}}_{N\times 1}\right)\odot \mathbf{A}$$

(13)

where ⊙ is the dot product, ${\mathrm{I}}_{N\times 1}$ means that all elements of the vector are one and its size is $N\times 1$. The slant steering vector ${\mathbf{b}}_r$ and the time steering vector ${\mathbf{b}}_t$ are shown in Equation (14):

$$\left\{\begin{array}{c}{\mathbf{b}}_r={\left[\begin{array}{cccccc}1& exp\left\{-j\pi f_1\frac{4R_0}{C}\right\}& \cdots & exp\left\{-j\pi f_2\frac{4R_0}{C}\right\}& \cdots & exp\left\{-j\pi f_k\frac{4R_0}{C}\right\}\end{array}\right]}^T\\ {\mathbf{b}}_t={\left[\begin{array}{ccccc}1& \cdots & exp\left\{j\pi \frac{\left(k-1\right)f_k·4V{T}_r{\rho }_{c}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}& \cdots & exp\left\{j\pi \frac{\left(K-1\right)f_K·4V{T}_r{\rho }_{c}cos{{\varphi }_q}^{{}^{\prime }}}{C}\right\}\end{array}\right]}^T\end{array}\right.$$

(14)

In the following, $\mathbf{A}$ is the space phase information:

$$\mathbf{A}={\left[\begin{array}{ccccc}{\mathbf{a}}_1& \cdots & {\mathbf{a}}_k& \cdots & {\mathbf{a}}_K\end{array}\right]}^T$$

(15)

where ${\mathbf{a}}_k$ is the space steering vector of the kth carrier frequency, as shown in Equation (16):

$${\mathbf{a}}_k={\left[\begin{array}{cccccc}1& exp\left\{j2\pi f_k\frac{dcos{\varphi }_q}{C}\right\}& \cdots & exp\left\{j2\pi f_k\frac{\left(n-1\right)dcos{\varphi }_q}{C}\right\}& \cdots & exp\left\{j2\pi f_k\frac{\left(N-1\right)dcos{\varphi }_q}{C}\right\}\end{array}\right]}^T$$

(16)

The echo signal of a range gate received by radar is the superposition of the echo of all ground scattering units on this range ring; in other words, we have

$$\mathbf{x}=\sum _q{\mathbf{x}}_q$$

(17)

It can be seen from the above echo signal model that the multi-frequency mode introduces a new dimension to the echo, which is the slant distance dimension. Different from the traditional single-carrier frequency SBEWR, the relationship between the three dimensions of space, time and slant distance is no longer a Kronecker product. Different pulses correspond to different phases of the slant distance and space. In the case of Earth’s rotation, the three dimensions of the ground scattering unit q are in Equation (18):

$$\begin{array}{c}\left\{\begin{array}{c}{\vartheta }_{r,q}=exp\left\{-j2\pi f_k\frac{2R_{0,q}}{C}\right\}=exp\left\{-j2\pi f_0\frac{2R_{0,q}}{C}\right\}·exp\left\{-j2\pi \left(k-1\right){\delta }_f\frac{2R_{0,q}}{C}\right\}\\ {\vartheta }_{t,q}=exp\left\{j\pi {\rho }_c{f}_k\frac{\left(k-1\right)4V{T}_{r}cos{\varphi }_q^{{}^{\prime }}}{C}\right\}=exp\left\{j\pi {\rho }_c{f}_0\frac{\left(k-1\right)4V{T}_{r}cos{\varphi }_q^{{}^{\prime }}}{C}\right\}·exp\left\{j\pi {\rho }_c{\delta }_f\frac{{\left(k-1\right)}^{2}4V{T}_{r}cos{\varphi }_q^{{}^{\prime }}}{C}\right\}\\ {\vartheta }_{s,q}=exp\left\{j\pi f_k\frac{2\left(n-1\right)dcos{\varphi }_q}{C}\right\}=exp\left\{j\pi f_0\frac{2\left(n-1\right)dcos{\varphi }_q}{C}\right\}·exp\left\{j\pi {\delta }_f\left(k-1\right)\frac{2\left(n-1\right)dcos{\varphi }_q}{C}\right\}\end{array}\right.\hfill \end{array}$$

(18)

In this paper, the uniform step frequency is adopted. The carrier frequency of the kth pulse is

$$f_k=f_0+\left(k-1\right){\delta }_f$$

(19)

where $f_0$ is the initial carrier frequency and ${\delta }_f$ the step frequency.

3. Range Ambiguity Suppression Based on Interpulse Multi-Frequency Signals

The radar periodically transmits signals. When the echoes of previously transmitted signals from a long distance are received simultaneously with the echoes from the current scattering unit, range ambiguity is generated, and the real slant distance cannot be discerned. For interpulse multi-frequency signals, the carrier frequencies of the echo signals from different range ambiguities and the target range gate echo signals are different. Due to the orthogonality of the multi-carrier frequency signals, the echoes of a partial range ambiguity can be suppressed by matching filters. The specific principle is shown in Figure 3.

The beginning of the first pulse is taken as the initial moment. Here, $t=0$ represents the start time of the first pulse emitted to the range ring $R_1$. At the initial moment $t=0$, the antenna transmits the pulse whose carrier frequency is $f_1$. After a delay, the signal reaches $R_1$. At the same time, the echo reflected from $R_2,R_3\dots R_K$ also just reaches the echo wavefront of $R_1$, which represent the previously transmitted pulses with the carrier frequency $f_2,f_3\dots f_K$. At this time, when the echo of the first transmitted pulse returns to the receiving end, the echo signals of multiple carrier frequencies from different range gates also arrive at the receiving end at the same time. The received signal is $c_1{s}_1\left(t\right)+c_2{s}_K\left(t\right)+\dots +c_{K-1}{s}_3\left(t\right)+c_K{s}_2\left(t\right)+\dots$. The carrier frequency of $s_i\left(t\right)$ is $f_i$. The receiving end adopts $h_i\left(t\right)=W\left(t/\phantom{t{T}_p}{T}_p\right)exp\left\{j2\pi f_{i}t\right\}$ to perform matched filtering. The preprocessed signal is

$$\begin{array}{c}{\mathbf{x}}_1=c_1{\rho }_{11}{\mathbf{a}}_{R_1}\left(f_1\right)+c_2{\rho }_{K1}{\mathbf{a}}_{R_2}\left(f_K\right)+\dots +c_K{\rho }_{21}{\mathbf{a}}_{R_K}\left(f_2\right)+\dots \hfill \end{array}$$

(20)

where ${\rho }_{ij}={\int }_0^{T_r}{s}_i\left(t\right)h_j^H\left(t\right)\mathrm{d}t$ is the response of the matching filter, ${\mathbf{a}}_{R_i}$ represents the echo from a distance $R_i$ received by an array of N elements and its carrier frequency is $f_j$, ${\mathbf{a}}_{R_i}\left(f_j\right)={\left[\begin{array}{cccccc}{s}_{i1}\left(t\right)& s_{i2}\left(t\right)& \dots & s_{i,n}\left(t\right)& \dots & s_{i,N}\left(t\right)\end{array}\right]}^T$, as shown in Equation (11). Similarly, the echo of the second pulse is $c_1{s}_2\left(t\right)+c_2{s}_1\left(t\right)+c_3{s}_K\left(t\right)+\dots +c_K{s}_3\left(t\right)$. It adopts $h_2\left(t\right)$ to perform matched filtering. The preprocessed signal is

$$\begin{array}{c}{\mathbf{x}}_2=c_1{\rho }_{22}{\mathbf{a}}_{R_1}\left(f_2\right)e^{j\pi {\omega }_{d,12}}+c_2{\rho }_{12}{\mathbf{a}}_{R_2}\left(f_1\right)·e^{j\pi {\omega }_{d,21}}+\dots +c_K{\rho }_{32}{\mathbf{a}}_{R_K}\left(f_3\right)e^{j\pi {\omega }_{d,K3}}\begin{array}{c}+\cdots \end{array}\hfill \end{array}$$

(21)

where ${\omega }_{d,ij}$ is the Doppler frequency of the echo from $R_i$ and its carrier frequency is $f_j$:

$${\omega }_{d,ij}={\rho }_c{f}_j\frac{4V{T}_{r}cos{\varphi }_i^{{}^{\prime }}}{C}$$

(22)

Thus, the echo of the Kth pulse is $c_1{s}_K\left(t\right)+c_2{s}_{K-1}\left(t\right)+c_3{s}_{K-2}\left(t\right)+\dots +c_K{s}_1\left(t\right)$, which adopts matched filtering, and the array output of the received signal is

$$\begin{array}{c}{\mathbf{x}}_K=c_1{\rho }_{KK}{\mathbf{a}}_{R_1}\left(f_K\right)e^{j\pi \left(K-1\right){\omega }_{d,1K}}+c_2{\rho }_{K-1,K}{\mathbf{a}}_{R_2}\left(f_{K-1}\right)e^{j\pi \left(K-1\right){\omega }_{d,2,K-1}}+\cdots +c_K{\rho }_{1K}{\mathbf{a}}_{R_K}\left(f_1\right)e^{j\pi \left(K-1\right){\omega }_{d,K1}}+\cdots \hfill \end{array}$$

(23)

After performing matched filtering of the echoes of K pulses, the echo vector can be written as shown in Equation (24). In Equation (24), ${\mathbf{b}}_i\left({\omega }_d\right)$ represents the time steering vector of the echo from $R_i$, ${\mathbf{b}}_i\left({\omega }_d\right)={\left[\begin{array}{cccc}0& e^{j\pi {\omega }_{d,i1}}& \cdots & e^{j\pi \left(K-1\right){\omega }_{d,i,K}}\end{array}\right]}^T$. ${{\mathbf{A}}_{R_i}}^{NK\times 1}={\left[\begin{array}{cccc}{\mathbf{a}}_{R_i}\left(f_1\right)& {\mathbf{a}}_{R_i}\left(f_2\right)& \dots & {\mathbf{a}}_{R_i}\left(f_K\right)\end{array}\right]}^T$, and ${\mathbf{a}}_{R_i}\left(f_j\right)$ is the space steering vector of the echo from $R_i$, whose carrier frequency is $f_j$. Its size is $N\times 1$:

$$\begin{array}{c}\mathbf{x}=\left[\begin{array}{c}{\mathbf{x}}_1\\ {\mathbf{x}}_2\\ ⋮\\ {\mathbf{x}}_K\end{array}\right]=c_1·kron\left(\left[\begin{array}{c}{\rho }_{11}\\ {\rho }_{22}\\ ⋮\\ {\rho }_{KK}\end{array}\right]\odot {\mathbf{b}}_1\left({\omega }_d\right),{\mathbf{I}}_{N\times 1}\right)\odot {\mathbf{A}}_{R_1}+c_2·kron\left(\left[\begin{array}{c}{\rho }_{K1}\\ {\rho }_{12}\\ ⋮\\ {\rho }_{K-1,K}\end{array}\right]\odot {\mathbf{b}}_2\left({\omega }_d\right),{\mathbf{I}}_{N\times 1}\right)\odot {\mathbf{A}}_{R_2}\hfill \\ \begin{array}{ccccccc}& & & & & & \end{array}+\cdots +c_K·kron\left(\left[\begin{array}{c}{\rho }_{21}\\ {\rho }_{32}\\ ⋮\\ {\rho }_{1K}\end{array}\right]\odot {\mathbf{b}}_K\left({\omega }_d\right),{\mathbf{I}}_{N\times 1}\right)\odot {\mathbf{A}}_{R_K}+\cdots \hfill \end{array}$$

(24)

Figure 4 shows the orthogonality of linear frequency-modulated (LFM) signals at different carrier frequencies. The step frequency is 0.8 MHz, the initial carrier frequency is 0.55 GHz, and the signal bandwidth is 0.8 MHz. When $i=j$, the correlation coefficient is ${\rho }_{ij}\approx 1$. For orthogonal signals with multiple frequencies, when $i\ne j$, the correlation coefficient is ${\rho }_{ij}\approx 0$ in the ideal case. This is equivalent to the amplitude modulation of the echo, suppressing the echo from other range gates such as $R_2$ and $R_3$.

Equation (24) can be simplified to Equation (25). Because there are K pulses in one CPI, $f_{K+1}=f_1$ is the first pulse of another CPI:

$$\begin{array}{c}\mathbf{x}=\left[\begin{array}{c}{\mathbf{x}}_1\\ {\mathbf{x}}_2\\ ⋮\\ {\mathbf{x}}_K\end{array}\right]=c_1·kron\left(\left[\begin{array}{c}{\rho }_{11}\\ {\rho }_{22}\\ ⋮\\ {\rho }_{KK}\end{array}\right]\odot {\mathbf{b}}_1\left({\omega }_d\right),{\mathbf{I}}_{N\times 1}\right)\odot {\mathbf{A}}_{R_1}\hfill \\ \begin{array}{ccc}\begin{array}{cc}& \end{array}& & \end{array}+c_{K+1}·kron\left(\left[\begin{array}{c}{\rho }_{11}\\ {\rho }_{22}\\ ⋮\\ {\rho }_{KK}\end{array}\right]\odot {\mathbf{b}}_{K+1}\left({\omega }_d\right),{\mathbf{I}}_{N\times 1}\right)\odot {\mathbf{A}}_{R_{K+1}}+\cdots \hfill \end{array}$$

(25)

According to the above equation, the echo of other range gates with a different carrier frequency from the echo of $R_1$ will be suppressed, while the echo of other range gates with the same carrier frequency will be retained.

The interval between adjacent range ambiguity units is one pulse interval $r_u=\frac{C}{2f_r}$. Combined with the above principle of orthogonality of interpulse multi-frequency signals, the carrier frequencies of the signals of range ambiguity $R_2\dots R_k$ are different from that of the target range gate $R_1$. The received data is shown in Figure 5.

In addition to the range gate $R_0$, the data of this range gate also includes its range ambiguities. For a single frequency, the echoes of these range gates cannot be separated. Under the influence of Earth’s rotation, the clutter spectrum widens, and the clutter suppression notch widens, which is rather unfavorable to the detection of targets near the main lobe clutter. For a multi-frequency mode, multiple range ambiguities are superimposed together. Due to the orthogonality of interpulse multi-frequency signals, the carrier frequencies of echoes from other range ambiguities are different from the echoes from range gate $R_0$. Therefore, the echoes from some range ambiguities can be suppressed. To sum up, the interpulse multi-frequency mode reduces the influence of range ambiguity and makes the clutter suppression notch narrower.

There are L range ambiguities for the qth range gate, and thus the echo is as shown in Equation (26):

$$x_{m,n,k,q}=\sum _{l=-L_1}^{L_2}\begin{array}{c}{\mathrm{A}}_{n,k,q,l}{w}_{T,k,q,l}exp\left\{-j\pi f_k\frac{4\left(R_{0,q}+l·r_u\right)}{C}\right\}exp\left\{j\pi {\rho }_c{f}_k\frac{\left(k-1\right)4V{T}_{r}cos{\varphi }_{q,l}^{{}^{\prime }}}{C}\right\}\hfill \\ exp\left\{j\pi f_k\frac{\left(n-1\right)2dcos{\varphi }_{q,l}}{C}\right\}\hfill \end{array}$$

(26)

If there are K pulses and K carrier frequencies, then the initial numbers of the forward range ambiguity and the backward range ambiguity are

$$\left\{\begin{array}{c}{L}_1=\frac{R_0-R_{min}}{r_u}\\ L_2=\frac{R_{max}-R_0}{r_u}\end{array}\right.$$

(27)

The total number for the range ambiguity is $L=L_1+L_2$. When the multi-frequency signal is orthogonal, the minimum range ambiguity number is

$$\left\{\begin{array}{c}{L}_{f1}=\frac{R_0-R_{min}}{K·r_u}\\ L_{f2}=\frac{R_{max}-R_0}{K·r_u}\end{array}\right.$$

(28)

The total number for the range ambiguity is $L_f=L_{f1}+L_{f2}$ after range ambiguity suppression, and the size of the range non-ambiguity resolution unit becomes $r_{u,new}=K·r_u$. The signal after range ambiguity suppression is as shown in Equation (29). The range ambiguity is greatly reduced to $L_f\approx \frac{L}{K}$:

$${\tilde{x}}_{n,k,q}=\sum _{l=-L_{f1}}^{L_{f2}}\begin{array}{c}{\mathrm{A}}_{n,k,q,l}{w}_{T,k,q,l}exp\left\{-j\pi f_k\frac{4\left(R_{0,q}+l·r_{u,new}\right)}{C}\right\}exp\left\{j\pi {\rho }_c{f}_k\frac{\left(k-1\right)4V{T}_{r}cos{\varphi }_{q,l}^{{}^{\prime }}}{C}\right\}\hfill \\ exp\left\{j\pi f_k\frac{\left(n-1\right)2dcos{\varphi }_{q,l}}{C}\right\}\hfill \end{array}$$

(29)

4. Residual Clutter Suppression

4.1. Clutter Suppression Process

STAP adopts the linear constraint minimum variance (LCMV) principle:

$$\left\{\begin{array}{c}min{\mathbf{w}}^H{\widehat{\mathbf{R}}}_c\mathbf{w}\\ s.t.\begin{array}{c}\end{array}{\mathbf{w}}^H\mathbf{S}=1\end{array}\right.$$

(30)

The adaptive weight vector of clutter suppression is

$$\mathbf{w}=\frac{{\widehat{\mathbf{R}}}_c^{-1}\mathbf{S}}{{\mathbf{S}}^H{\widehat{\mathbf{R}}}_c^{-1}\mathbf{S}}$$

(31)

where $\mathbf{S}$ is the preset ideal steering vector and ${\widehat{\mathbf{R}}}_c$ is the clutter covariance matrix obtained by the maximum likelihood estimation of the training samples near the sample to be tested:

$${\widehat{\mathbf{R}}}_c=\frac{1}{P}\sum _{i=1}^P{\tilde{\mathbf{x}}}_i{\tilde{\mathbf{x}}}_i^H$$

(32)

Based on the above formulas, the clutter suppression process is shown in Figure 6:

$\mathbf{Step}\mathbf{1}:$ Based on the orthogonality of multi-frequency signals, matched filtering is used to suppress the partial range ambiguity clutter.

$\mathbf{Step}\mathbf{2}:$ After pulse compression and a Fourier transform, the echo is converted to dimensions located in space, time and range, where the time dimension is the Doppler frequency. For 2D-STAP, the spatial dimension is the azimuth beam dimension, and for 3D-STAP, the spatial dimension is the azimuth beam dimension and pitch beam dimension.

$\mathbf{Step}\mathbf{3}:$ A range gate is selected as the sample to be tested, and a channel for testing is selected.

$\mathbf{Step}\mathbf{4}:$ According to the principle of reduction STAP, several auxiliary channels are selected to construct a reduction matrix, which will be used to reduce the dimension of the echo data.

$\mathbf{Step}\mathbf{5}:$ The reduction matrix is used to reduce the dimension of the echo, and then the covariance matrix is estimated.

$\mathbf{Step}\mathbf{6}:$ Solve the adaptive weight vector and carry out clutter suppression.

4.2. Discussion of 3D JDL-STAP

After interpulse multi-frequency matching filtering, part of the range ambiguity can be suppressed, and the remaining amount of range ambiguity is as shown in Equation (28). The interpulse multi-frequency mode makes the echo signal introduce the phase of the slant distance dimension, but it is different from traditional multi-dimensional adaptive processing. According to Equation (18), the slant distance dimension, Doppler domain and azimuth domain are not independent of each other. The change in carrier frequencies couples these three dimensions together, and thus it can be called “pseudo three dimensions”, which cannot be used directly as three dimensions to perform carrier frequency-azimuth beam-Doppler 3D-STAP.

When 2D-STAP is used for clutter suppression, the phase of the slant distance can assist the filter to form a clutter suppression notch at the range gate to be tested. After the preprocessing of slant distance phase compensation, the preset ideal steering vector is $\mathbf{S}={\left[S_{11}\dots S_{nk}\dots S_{NK}\right]}^T$:

$$\begin{array}{c}{S}_{nk}=exp\left\{j\pi \left(k-1\right)f_k\frac{{\omega }_d}{f_0}\right\}·exp\left\{j\pi \left(n-1\right)f_k\frac{{\omega }_s}{f_0}\right\}\hfill \end{array}$$

(33)

The following analysis takes one forward range ambiguity and three backward range ambiguities as an example. The slant distance of the range gate to be tested is r, and its normalized Doppler frequency is ${\omega }_{d,0}$. The slant distances of the range ambiguities are $r-r_{u,new}$, $r+r_{u,new}$, $r+2r_{u,new}$ and $r+3r_{u,new}$. Their normalized Doppler frequencies are ${\omega }_{d,-1}$, ${\omega }_{d,1}$, ${\omega }_{d,2}$ and ${\omega }_{d,3}$. The phase information of the slant distance and Doppler frequency of the range ambiguity are shown in Equation (34), where $\Delta {\phi }_l={\delta }_f\frac{4l·r_{u,new}}{C}\left(l=-1,1,2,3\right)$, which will be compensated by the normalized Doppler frequency ${\omega }_d$ shift as shown in Equation (35):

$$\left\{\begin{array}{c}exp\left\{-j\pi f_k\frac{4l·r_{u,new}}{C}\right\}=exp\left\{-j\pi f_0\frac{4l·r_{u,new}}{C}\right\}·exp\left\{-j\pi \left(k-1\right){\delta }_f\frac{4l·r_{u,new}}{C}\right\}\left(l=-1,1,2,3\right)\\ exp\left\{j\pi \left(k-1\right)f_k\frac{{\omega }_d}{f_0}\right\}=exp\left\{j\pi \left(k-1\right){\omega }_d\right\}·exp\left\{j\pi {\left(k-1\right)}^2{\delta }_f\frac{{\omega }_d}{f_0}\right\}\end{array}\right.$$

(34)

$$\begin{array}{ccc}r-r_{u,new}& \to & {\omega }_{d,-1}\\ r& \to & {\omega }_{d,0}\\ r+r_{u,new}& \to & {\omega }_{d,1}\\ r+2r_{u,new}& \to & {\omega }_{d,2}\\ r+3r_{u,new}& \to & {\omega }_{d,3}\end{array}\begin{array}{cc}\to & {\omega }_{d,-1}+\Delta {\omega }_{d,-1}\Rightarrow r\\ \to & {\omega }_{d,0}\\ \to & {\omega }_{d,1}+\Delta {\omega }_{d,1}\Rightarrow r\\ \to & {\omega }_{d,2}+\Delta {\omega }_{d,2}\Rightarrow r\\ \to & {\omega }_{d,3}+\Delta {\omega }_{d,3}\Rightarrow r\end{array}$$

(35)

where $\Delta {\omega }_{d,l}\approx \Delta {\phi }_l$. According to the above analysis, although the interpulse multi-frequency mode can suppress most of the range ambiguity and form a clutter suppression notch in the range gate to be tested, due to the particularity of its structure and the coupling relationship of the three dimensions, the normalized Doppler frequency will compensate the phase difference between the range ambiguities and the range gate to be tested through a frequency shift, which will form several additional clutter suppression notches.

There are still the quadratic terms $exp\left\{j\pi {\left(k-1\right)}^2{\delta }_f\frac{\Delta {\omega }_{d,l}}{f_0}\right\}$ in the compensated data. The spacetime phase of range ambiguity cannot be completely obtained with a frequency shift. Because the phase of this quadratic term is small, the clutter suppression notches can still be formed. The extra clutter suppression notches are more shallow than the main lobe clutter suppression notch.

The non-stationary clutter caused by residual range ambiguity makes the clutter DOFs greater than that of the system. Generally, there are two solutions: one is to further reduce the clutter DOFs, and the other is to increase the system DOFs to suppress the residual ambiguity clutter. After the previous matching filtering of the interpulse multi-frequency mode, the clutter DOFs can no longer be reduced. Therefore, when considering increasing the system DOFs, this paper introduces the pitch dimension to improve the system DOFs.

In addition to traditional 2D-STAP, the pitch array elements can also be used for clutter suppression. The snapshot data of clutter is the two-dimensional planar of the azimuth-Doppler STAP. When considering the pitch dimension, the snapshot data of clutter is the three-dimensional data of the azimuth-pitch-Doppler STAP. The signal of the pitch-azimuth-Doppler STAP is shown in Equation (36):

$${\tilde{x}}_{m,n,k,q}=\sum _{l=-L_{f1}}^{L_{f2}}\begin{array}{c}{\mathrm{A}}_{\mathrm{m},\mathrm{n},\mathrm{k},\mathrm{q},\mathrm{l}}{w}_{T,k,q,l}exp\left\{-j\pi f_k\frac{4\left(R_{0,q}+l·K·r_u\right)}{C}\right\}exp\left\{j\pi {\rho }_c{f}_k\frac{\left(k-1\right)4V{T}_{r}cos{\varphi }_{q,l}^{{}^{\prime }}}{C}\right\}\hfill \\ exp\left\{j\pi f_k\frac{\left(n-1\right)2dcos{\varphi }_{q,l}}{C}\right\}·exp\left\{j\pi f_k\frac{\left(m-1\right)2dcos{\theta }_{EL,q,l}}{C}\right\}\hfill \end{array}$$

(36)

According to Equation (13), the signal vector can be written as

$${\tilde{\mathbf{x}}}_q=kron\left({\mathbf{b}}_r\odot {\mathbf{b}}_t,{\mathrm{I}}_{MN\times 1}\right)\odot \Theta$$

(37)

where $\Theta$ represents the pitch-azimuth phase information:

$$\Theta ={\left[\begin{array}{ccccc}{\Theta }_1& \cdots & {\Theta }_k& \cdots & {\Theta }_K\end{array}\right]}^T$$

(38)

where ${\Theta }_k=kron\left({\mathbf{a}}_{EL,k},{\mathbf{a}}_{AZ,k}\right)$ represents the space steering vector of the echo with the carrier frequency $f_k$ in Equation (39):

$$\left\{\begin{array}{c}{\mathbf{a}}_{EL,k}={\left[\begin{array}{cccccc}1& exp\left\{j2\pi f_k\frac{dcos{\theta }_{EL}}{C}\right\}& \cdots & exp\left\{j2\pi f_k\frac{\left(m-1\right)dcos{\theta }_{EL}}{C}\right\}& \cdots & exp\left\{j2\pi f_k\frac{\left(M-1\right)dcos{\theta }_{EL}}{C}\right\}\end{array}\right]}^T\\ {\mathbf{a}}_{AZ,k}={\left[\begin{array}{cccccc}1& exp\left\{j2\pi f_k\frac{dcos{\varphi }_q}{C}\right\}& \cdots & exp\left\{j2\pi f_k\frac{\left(n-1\right)dcos{\varphi }_q}{C}\right\}& \cdots & exp\left\{j2\pi f_k\frac{\left(N-1\right)dcos{\varphi }_q}{C}\right\}\end{array}\right]}^T\end{array}\right.$$

(39)

The ideal spacetime steering vector of clutter suppression is

$$\begin{array}{c}\mathbf{S}={\left[\begin{array}{ccccc}{S}_{111}& \cdots & S_{mnk}& \cdots & S_{MNK}\end{array}\right]}^T\hfill \\ S_{mnk}=exp\left\{j\pi \left(k-1\right)f_k\frac{{\omega }_d}{f_0}\right\}·exp\left\{j\pi \left(n-1\right)f_k\frac{{\omega }_{s,az}}{f_0}\right\}·exp\left\{j\pi \left(m-1\right)f_k\frac{{\omega }_{s,el}}{f_0}\right\}\hfill \end{array}$$

(40)

The diagram of 3D-STAP clutter suppression according to Equations (30)–(32) and the channel selection of 3D-JDL is shown in Figure 7. The red cube represents the selected channel to be tested, with the channels around it serving as auxiliary channels. For 3D-JDL, the selected auxiliary channels are the cubes in the grid-like cube shown below. The introduction of the pitch dimension can suppress non-stationary clutter by using additional controllable system DOFs and also can form the clutter suppression notch located at the pitch angle corresponding to the range gate to be tested, which can better solve the problem of clutter suppression notch frequency shift.

4.3. Discussion of Clutter Covariance Matrix Estimation

The covariance matrix is estimated by the training samples near the sample to be tested, which is one of the key steps to achieving effective clutter suppression. In addition to the spacetime phase information, the echo signal also has the phase information of the slant distance dimension as shown in Equation (18). One can take the signal vectors in Equation (41) as an example to analyze the error of covariance matrix estimation. Here, ${\mathbf{s}}_1$ and ${\mathbf{s}}_{-1}$ are training samples around the sample to be tested ${\mathbf{s}}_0$. The training samples ${\mathbf{s}}_1$ and ${\mathbf{s}}_{-1}$ are used to estimate the clutter covariance matrix of the sample to be tested ${\mathbf{s}}_0$:

$$\left\{\begin{array}{c}{\mathbf{s}}_0={\left[\begin{array}{cccc}{s}_{111}& s_{211}& \cdots s_{mnk}& \cdots s_{MNK}\end{array}\right]}^T\\ {\mathbf{s}}_1={\left[\begin{array}{ccc}{s}_{111,1}& \cdots s_{mnk,1}& \cdots s_{MNK,1}\end{array}\right]}^T\\ {\mathbf{s}}_{-1}={\left[\begin{array}{ccc}{s}_{111,-1}& \cdots s_{mnk,-1}& \cdots s_{MNK,-1}\end{array}\right]}^T\end{array}\right.$$

(41)

If the slant distance phase is arranged diagonally as a diagonal matrix, then the samples that introduce the slant phase are ${\mathbf{G}}_0{\mathbf{s}}_0$, ${\mathbf{G}}_1{\mathbf{s}}_1$ and ${\mathbf{G}}_{-1}{\mathbf{s}}_{-1}$, where

$$\begin{array}{c}{s}_{mnk,i}=\mathrm{A}exp\left\{j2\pi f_k\frac{d_{n}cos{\varphi }_i}{C}\right\}·exp\left\{j2\pi f_k\frac{d_{m}cos{\theta }_{EL}}{C}\right\}·exp\left\{j2\pi f_k\frac{2{\rho }_{c}V{t}_{k}cos{{\varphi }_i}^{{}^{\prime }}}{C}\right\}\hfill \end{array}$$

(42)

The one element of the slant distance matrix ${\mathbf{G}}_i$ is

$$g_{mnk,i}=exp\left\{-j2\pi f_k\frac{2\left(R_0+i·{\Delta }_R\right)}{C}\right\}$$

(43)

The ideal covariance matrix is

$$\mathbf{R}=\left({\mathbf{G}}_0{\mathbf{s}}_0\right){\left({\mathbf{G}}_0{\mathbf{s}}_0\right)}^H={\mathbf{G}}_0{\mathbf{s}}_0{\mathbf{s}}_0^H{\mathbf{G}}_0^H={\mathbf{R}}_{G0}\odot {\mathbf{R}}_0$$

(44)

The estimated covariance matrix is as shown in Equation (45):

$$\begin{array}{c}\widehat{\mathbf{R}}=\frac{1}{2}\left({\mathbf{G}}_1{\mathbf{s}}_1+{\mathbf{G}}_{-1}{\mathbf{s}}_{-1}\right){\left({\mathbf{G}}_1{\mathbf{s}}_1+{\mathbf{G}}_{-1}{\mathbf{s}}_{-1}\right)}^H\hfill \\ =\frac{1}{2}\left({\mathbf{G}}_1{\mathbf{s}}_1{\mathbf{s}}_1^H{\mathbf{G}}_1^H+{\mathbf{G}}_{-1}{\mathbf{s}}_{-1}{\mathbf{s}}_{-1}^H{\mathbf{G}}_{-1}^H+{\mathbf{G}}_1{\mathbf{s}}_1{\mathbf{s}}_{-1}^H{\mathbf{G}}_{-1}^H+{\mathbf{G}}_{-1}{\mathbf{s}}_{-1}{\mathbf{s}}_1^H{\mathbf{G}}_1^H\right)\hfill \\ =\frac{1}{2}\left({\mathbf{R}}_{G1}\odot {\mathbf{s}}_1{\mathbf{s}}_1^H+{\mathbf{R}}_{G-1}\odot {\mathbf{s}}_{-1}{\mathbf{s}}_{-1}^H+{\mathbf{R}}_{G1,-1}\odot {\mathbf{s}}_1{\mathbf{s}}_{-1}^H+{\mathbf{R}}_{G-1,1}\odot {\mathbf{s}}_{-1}{\mathbf{s}}_1^H\right)\hfill \\ =\frac{1}{2}\left({\mathbf{R}}_{G1}\odot {\mathbf{R}}_1+{\mathbf{R}}_{G-1}\odot {\mathbf{R}}_{-1}+{\mathbf{R}}_{G1,-1}\odot {\mathbf{R}}_{1,-1}+{\mathbf{R}}_{G-1,1}\odot {\mathbf{R}}_{-1,1}\right)\hfill \end{array}$$

(45)

where ${\mathbf{R}}_{Gi}={\mathbf{g}}_i{\mathbf{g}}_i^H$, ${\mathbf{R}}_{G,i,j}={\mathbf{g}}_i{\mathbf{g}}_j^H$, $\mathbf{g}=vec\left(\mathbf{G}\right)$, ${\mathbf{R}}_i={\mathbf{s}}_i{\mathbf{s}}_i^H$ and ${\mathbf{R}}_{i,j}={\mathbf{s}}_i{\mathbf{s}}_j^H$.

The phase information of one element of the cross-correlation matrix ${\mathbf{R}}_{i,j}={\mathbf{s}}_i{\mathbf{s}}_j^H$ consists of the three parts in Equation (46), where i and j are the sequence numbers of the training samples and $n_1$, $n_2$, $k_1$ and $k_2$ are the sequence numbers of pulses and array elements corresponding to the position of the element in the matrix, respectively:

$$\left\{\begin{array}{c}\begin{array}{c}{\Gamma }_1=exp\left\{j\pi {\rho }_c\frac{4V{T}_r}{C}cos\left({\theta }_{AZ}+{\varphi }_c\right)\left[f_{k_1}\left(k_1-1\right)sin{\theta }_{EL,i}-f_{k_2}\left(k_2-1\right)sin{\theta }_{EL,j}\right]\right\}\hfill \\ \begin{array}{cc}& \end{array}\approx exp\left\{j\pi {\rho }_c\frac{4V{T}_r}{C}cos\left({\theta }_{AZ}+{\varphi }_c\right)\left[f_{k_1}\left(k_1-1\right)-f_{k_2}\left(k_2-1\right)\right]sin{\theta }_{EL,0}\right\}\hfill \end{array}\\ \begin{array}{c}{\Gamma }_2=exp\left\{j\pi \frac{2d_n}{C}cos{\theta }_{AZ}\left[f_{k_1}\left(k_1-1\right)sin{\theta }_{EL,i}-f_{k_2}\left(k_2-1\right)sin{\theta }_{EL,j}\right]\right\}\hfill \\ \begin{array}{cc}& \end{array}\approx exp\left\{j\pi \frac{2d_n}{C}cos{\theta }_{AZ}\left[f_{k_1}\left(k_1-1\right)-f_{k_2}\left(k_2-1\right)\right]sin{\theta }_{EL,0}\right\}\hfill \end{array}\\ \begin{array}{c}{\Gamma }_3=exp\left\{j\pi \frac{2d_m}{C}\left[f_{k_1}\left(k_1-1\right)cos{\theta }_{EL,i}-f_{k_2}\left(k_2-1\right)cos{\theta }_{EL,j}\right]\right\}\hfill \\ \begin{array}{cc}& \end{array}\approx exp\left\{j\pi \frac{2d_m}{C}\left[f_{k_1}\left(k_1-1\right)-f_{k_2}\left(k_2-1\right)\right]cos{\theta }_{EL,0}\right\}\hfill \end{array}\end{array}\right.$$

(46)

According to the Brennan criterion and maximum likelihood estimation, in order to ensure that the loss of the output of the signal-to-clutter-and-noise ratio (SCNR) is no more than 3 dB, the number of training samples is required to be more than twice the system DOFs; that is, the estimated covariance matrix obtained by a sufficient number of training samples can accurately estimate the clutter covariance matrix of the samples to be tested. When the number of independent identically distributed (i.i.d.) training samples is sufficient, the ≈ symbol is valid, and the phase error introduced by the cross-correlation terms can be negligible.

An element of the autocorrelation matrix ${\mathbf{R}}_{Gi}={\mathbf{g}}_i{\mathbf{g}}_i^H$ of the slant distance phase is in Equation (47). One of the elements of the cross-correlation matrix ${\mathbf{R}}_{G,i,j}={\mathbf{g}}_i{\mathbf{g}}_j^{*}$ of the slant distance phase is shown in Equation (48):

$$\begin{array}{c}{g}_{m_1{n}_1{k}_1,i}·g_{m_2{n}_2{k}_2,i}^{*}=exp\left\{-j2\pi f_{k_1}\frac{2\left(R_0+i·{\Delta }_R\right)}{C}\right\}exp\left\{j2\pi f_{k_2}\frac{2\left(R_0+i·{\Delta }_R\right)}{C}\right\}\hfill \\ \begin{array}{ccccc}& & & & \end{array} =exp\left\{-j2\pi \left(f_{k_1}\frac{2\left(R_0+i·{\Delta }_R\right)}{C}-f_{k_2}\frac{2\left(R_0+i·{\Delta }_R\right)}{C}\right)\right\}\hfill \\ \begin{array}{ccccc}& & & & \end{array} =exp\left\{-j2\pi {\delta }_f\frac{2R_0}{C}\left(k_1-k_2\right)\right\}·exp\left\{-j2\pi {\delta }_f\frac{2i·{\Delta }_R}{C}\left(k_1-k_2\right)\right\}\hfill \end{array}$$

(47)

$$\begin{array}{c}{g}_{m_1{n}_1{k}_1,i}·g_{m_2{n}_2{k}_2,j}^{*}=exp\left\{-j2\pi f_{k_1}\frac{2\left(R_0+i·{\Delta }_R\right)}{C}\right\}exp\left\{j2\pi f_{k_2}\frac{2\left(R_0+j·{\Delta }_R\right)}{C}\right\}\hfill \\ \begin{array}{cccc}& & & \end{array} =exp\left\{-j2\pi \left(f_{k_1}\frac{2\left(R_0+i·{\Delta }_R\right)}{C}-f_{k_2}\frac{2\left(R_0+j·{\Delta }_R\right)}{C}\right)\right\}\hfill \\ \begin{array}{cccc}& & & \end{array} =exp\left\{-j2\pi {\delta }_f\frac{2R_0}{C}\left(k_1-k_2\right)\right\}·exp\left\{-j2\pi f_0\frac{2\left(i-j\right)·{\Delta }_R}{C}\right\}·exp\left\{-j2\pi {\delta }_f\frac{2{\Delta }_R}{C}\left[\left(k_1-1\right)i-\left(k_2-1\right)j\right]\right\}\hfill \end{array}$$

(48)

It can be seen that the estimation error of the slant distance phase information is large and cannot be ignored, as shown in Equation (49). Therefore, due to the existence of slant distance phase information, the cross-correlation terms will be caused by the echoes of different training samples and affect the estimation accuracy of the clutter covariance matrix:

$$exp\left\{-j2\pi f_0\frac{2\left(i-j\right)·{\Delta }_R}{C}\right\}·exp\left\{-j2\pi {\delta }_f\frac{2{\Delta }_R}{C}\left[\left(k_1-1\right)i-\left(k_2-1\right)j\right]\right\}\ne exp\left\{-j2\pi {\delta }_f\frac{2i·{\Delta }_R}{C}\left(k_1-k_2\right)\right\}$$

(49)

According to Equations (47)–(49), errors are quantified to more intuitively analyze the influence of the slant distance phase on covariance matrix estimation. In the case of $k_1=k_2$, the following three values are used to measure the degree of error caused by cross-correlation terms in the covariance estimation in the three dimensions of time, space and slant distance:

$$\left\{\begin{array}{c}{\varphi }_1=\frac{4{\rho }_{c}V{T}_r}{\lambda }cos\left({\theta }_{AZ}+{\varphi }_c\right){\Delta }_{EL,ij}\\ {\varphi }_2=\frac{2d}{\lambda }\left(cos{\theta }_{AZ}{\Delta }_{EL,ij}+{\delta }_{EL,ij}\right)\\ {\varphi }_3=\frac{4\left(i-j\right)·{\Delta }_R}{C}{\delta }_f\end{array}\right.$$

(50)

where ${\Delta }_{EL,ij}=\left(sin{\theta }_{EL,i}-sin{\theta }_{EL,j}\right)$ and ${\delta }_{EL,ij}=\left(cos{\theta }_{EL,i}-cos{\theta }_{EL,j}\right)$.

A range gate is taken as a reference range gate, while 20 range gates are taken from its left and right sides each to analyze the phase difference generated by the cross-correlation terms. The slant distances of the reference range gate are 600 km, 1000 km and 1800 km respectively, representing a short distance, medium distance and far distance, respectively. However, ${\varphi }_1$ and ${\varphi }_2$ are not complete expressions of the cross-correlation terms. They are expressed here only to show that the variation in pitch angle of different range gates has little influence on the estimation results, and the error is minimal. Here, ${\varphi }_3$ is the phase of this term in which the step frequency is used rather than the initial carrier frequency, because the value of the term containing the initial carrier frequency is a constant and has no effect on the matrix.

The reason why only $k_1=k_2$ is considered is that from Equations (47)–(49), it can be seen that when $k_1=k_2$, the error is larger. Here, if the error of $k_1=k_2$ cannot be ignored, then it also means that the error will be more serious when $k_1\ne k_2$. In Figure 8, it can be seen that the error term of the time domain and the space domain was almost zero, and the error of the slant distance domain was quite large. If we want to reduce the estimation error, then we need to reduce the number of training samples because $abs(i-j)$ of ${\varphi }_3$ is smaller. However, the number of training samples is related to the system DOFs. If it is too low, then the estimation accuracy of the time domain and the space domain will be affected.

In order to eliminate the influence of the slant distance phase on the estimation of the covariance matrix, a phase compensation matrix was constructed to compensate for the slant distance:

$$\mathbf{\Phi }=diag\left(kron\left(\mathbf{r},{\mathrm{I}}_{MN\times 1}\right)\right)$$

(51)

where $\mathbf{r}=exp\left\{-j\pi \left(f_0+{\delta }_f·{\left[0:K-1\right]}^T\right)\frac{4r}{C}\right\}$. The echo data after the slant distance phase compensation are ${{\tilde{\mathbf{x}}}^{\prime }}_q={\mathbf{\Phi }}^{\mathbf{H}}{\tilde{\mathbf{x}}}_q$. One of the elements of the signal vector will be modified to Equation (52):

$${\tilde{x}}_{m,n,k,q}^{\prime }=\sum _{l=-L_{f1}}^{L_{f2}}\begin{array}{c}{\mathrm{A}}_{\mathrm{m},\mathrm{n},\mathrm{k},\mathrm{q},\mathrm{l}}{w}_{T,k,q,l}·exp\left\{-j\pi f_k\frac{4l·r_{u,new}}{C}\right\}·exp\left\{j\pi {\rho }_c{f}_k\frac{\left(k-1\right)4V{T}_{r}cos{\varphi }_{q,l}^{{}^{\prime }}}{C}\right\}\hfill \\ ·exp\left\{j\pi f_k\frac{\left(n-1\right)2dcos{\varphi }_{q,l}}{C}\right\}·exp\left\{j\pi f_k\frac{\left(m-1\right)2dcos{\theta }_{EL,l}}{C}\right\}\hfill \end{array}$$

(52)

Then, the compensated echo is used for clutter suppression according to Equations (30)–(32). The detailed process is shown in the Figure 9.

5. Simulation

The parameters in

Table 1. Simulation parameters.

Parameter	Value
Orbit height H	500 km
Pulse number K	16
The azimuth dimension of the antenna $L_{AZ}$	50 m
The pitch dimension of the antenna $L_{EL}$	3 m
Element spacing d	0.25 m
Detection range	800–810.5 km
Range gate to be tested	809 km
Orbit inclination $\eta$	45${}^{\circ }$
Bandwidth B	0.8 MHz
Pulse repetition frequency $f_r$	6000 Hz
Number of channels in azimuth dimension $N_{ch}$	10
Initial carrier frequency $f_0$	0.55 GHz
Step frequency ${\delta }_f$	0.8 MHz
Noise coefficient $F_n$	2 dB
Clutter background	Sea and ground

were used for the clutter simulation. In the azimuth dimension, the overlap mode was used for channel division, which was divided into 10 channels. At the orbit altitude of 500 km, the range ambiguity corresponding to the single carrier frequency radar was about 80. When there were 16 carrier frequencies in the interpulse multi-frequency mode, the range ambiguity number was reduced to about five. The slant distance of the RD spectrum with 800–810.5 km was selected for the simulation, and the range gate of 809 km was selected for the simulation in the output SCNR. Sea clutter and ground clutter were analyzed in the following simulations.

In the case of a single frequency, the energy of the moving target located near the suppression notch of the main lobe clutter will be suppressed, and the target cannot be detected correctly. The echo signals with different carrier frequencies from different range ambiguities were filtered by using the interpulse multi-frequency mode and the orthogonality of different carrier frequency signals. The clutter subspace of a range ring was spanned by the sampling vectors corresponding to all clutter units on this range ring, and the total clutter subspace was the union of the clutter subspaces of each range ambiguity. Therefore, in the case of range ambiguity, the dimension of clutter subspace increased; that is, the clutter DOFs increased. According to the above parameters, the amount of range ambiguity was changed from 80 to about 5 after matched filtering, and thus the non-stationary characteristics of the clutter were improved. In the eigenvalue spectrum, the clutter eigenvalues are greater than the noise eigenvalues (i.e., the clutter DOFs). As shown in Figure 10, the clutter DOFs were reduced from about 130 to 40.

In Figure 11 and Figure 12, the main lobe clutter in the RD spectrum is preprocessed by a matching filter and becomes narrower. Because the slant distance information led to the covariance estimation error, the clutter subspace did not match the clutter subspace of the ideal covariance matrix, and multiple suppression notches appeared. According to Equation (51), slant range preprocessing was carried out. After slant range compensation, the slant distance phase was eliminated, and the main clutter ridge was clearer. However, according to Equation (35), a Doppler shift caused the clutter notches of the range ambiguities to move, and there were still several suppression notches on both sides of the main lobe clutter. If the target is located in these notches, then the target energy will be suppressed. Therefore, the residual range ambiguities also need to be further processed.

In Figure 13 and Figure 14, 3D-STAP of the pitch-azimuth-Doppler were adopted for clutter suppression. The DOF of the pitch dimension was introduced, and the remaining range ambiguities were suppressed by additional controllable system DOFs, which effectively solved the frequency shift of the clutter suppression notch caused by range ambiguity echo. The main lobe clutter ridge was clear, and the clutter suppression notch was narrower after the preprocessing with slant range compensation. This is more advantageous for moving target detection near the main lobe clutter.

In the 20th range gate, it was assumed that there was a moving target with a normalized Doppler frequency of −0.2 near the main lobe clutter. The Doppler shift of the moving target was small, and its Doppler frequency was close to that of the main lobe clutter. Under the influence of the crab angle and range ambiguity, the clutter suppression notch was wider, which is unfavorable for the detection of targets near the notch. As shown in Figure 15, clutter suppression was carried out by single frequency + 2D-STAP, multi-frequency + slant range preprocessing + 2D-STAP and multi-frequency + slant range preprocessing + 3D-STAP. The matching filtering of the multi-frequency signal suppressed the partial range ambiguity and reduced the clutter DOFs. After slant range preprocessing and 3D-STAP, it became evident that the clutter energy was significantly diminished, and the target became prominently discernible.

6. Conclusions

Aiming at the serious problem of non-stationary clutter in SBEWR, an efficient approach was proposed for suppressing non-stationary clutter in this paper which can effectively reduce clutter DOFs. Firstly, an interpulse multi-frequency signal transmission mode was proposed. The carrier frequencies of echoes from different range ambiguities are different. By using the orthogonality of multi-carrier frequency signals and adopting the matching filter which is consistent with the carrier frequency of the target range gate, the range ambiguity echoes with different carrier frequencies can be suppressed so as to achieve the effect of suppressing partial range ambiguities. Then, the residual clutter was suppressed by slant range preprocessing together with pitch-azimuth-Doppler 3D-STAP. Slant range preprocessing can compensate the slant distance phase to avoid the error caused by cross-correlation term of covariance matrix estimation. After preprocessing, more controllable system DOFs can be obtained by 3D-STAP, and additional controllable pitch dimension DOFs are used to suppress residual clutter. The simulation results verify the effectiveness of the proposed method. The clutter ridge is clear, and the notch of clutter suppression is obviously narrowed, which is more favorable for moving target detection.

While non-stationary clutter suppression has advantages, new signal transmission modes will also bring new research problems. For example, the total bandwidth of a signal with multiple carrier frequencies is wider, and its impact on the narrow-band hypothesis needs to be considered. In addition, the carrier frequencies of different pulses are different, which puts forward new requirements for the method of target energy accumulation. Signals of different frequency bands need to adopt appropriate non-coherent accumulation methods for target detection. These problems need to be solved for the interpulse multi-frequency mode in future studies.