A Long-Time Coherent Integration STAP for GEO Spaceborne-Airborne Bistatic SAR

Abstract

A geosynchronous spaceborne-airborne bistatic synthetic aperture radar (GEO SA-BSAR) system is an important technique to achieve long-time moving target monitoring over a wide area. However, due to special bistatic configuration of GEO SA-BSAR, two major challenges, i.e., severe range migration and space-variant Doppler parameters for moving targets, hinder the moving target indication (MTI) processing. Traditional SAR MTI methods, which do not take the challenges into consideration, will defocus the moving targets, leading to a loss of the signal-to-noise ratio (SNR). To focus moving targets and estimate motion parameters accurately, long-time coherent integration space-time adaptive processing (LTCI-STAP) is proposed for GEO SA-BSAR MTI in this paper. First, a modified adaptive spatial filtering based on the bistatic signal model is performed to suppress the clutter. Then, an LTCI filter bank is constructed to achieve range migration correction and moving target focusing, which yields the optimal output signal and filtering parameters. Finally, constant false alarm rate (CFAR) detection is carried out to determine the targets, and the space-variant Doppler parameters, solved from the filtering parameters, are used for estimating moving target positions and velocities. Simulations verify the effectiveness of our method.

1. Introduction

The geosynchronous spaceborne-airborne bistatic synthetic aperture radar (GEO SA-BSAR) system, consisting of a GEO transmitter and an airborne receiver, has attracted much attention in recent years. The GEO transmitter, covering thousands of kilometers and providing stable coverage for a long exposure time [1,2,3,4,5], is an ideal illuminator of opportunity. Furthermore, the airborne multichannel receiver can achieve different bistatic configurations for various applications without the constraints of a fixed orbit [6]. Most studies on GEO SA-BSAR are focused on the imaging of a stationary scene, including bistatic configuration design [7], resolution analysis [8,9,10], and imaging algorithms [11,12,13,14,15]. Various imaging modes, such as side-looking, squint-looking, and forwards-looking, were proved to be promising in recent research [8]. Various imaging modes can obtain the target’s scattering properties in different observation geometries, which are advantageous for moving target indication (MTI).

GEO SA-BSAR MTI has been studied for over 20 years. The idea of realizing MTI by GEO SA-BSAR was first proposed in 1997 [16], which only focused on the feasibility analysis of MTI. Subsequently, an MTI method for GEO SA-BSAR was investigated in [17], where two different airborne side-looking radars were employed to address azimuth ambiguity. In [18,19], the problems of azimuth Doppler ambiguity and azimuth spectrum aliasing of a moving target were solved. Then, the imaging of the moving target with unknown parameters was studied in [20]. However, in these methods, only a short synthetic aperture time (SAT) is considered, where the output signal-to-noise ratio (SNR) is not optimal, as most of the available pulses are not processed. In addition, most of the above studies focused on a simple bistatic configuration where the aircraft works in side-looking mode.

In this paper, we consider GEO SA-BSAR MTI with a long SAT and complex bistatic configurations, which can obtain higher SNR, higher azimuth resolution and abundant scattering characteristics. In fact, scenes can be illuminated by the GEO satellite for a long time, providing the airborne receiver with an opportunity to obtain a longer SAT (approximately 2~3 s [8]).

However, the long SAT and complex bistatic configurations bring about significant challenges to GEO SA-BSAR MTI processing. On the one hand, the long SAT will induce a high-order range migration. The classical method [18] did not consider such a long SAT, which will deteriorate the coherent integration. The methods in [21,22,23] performed SAR MTI with complex range migration. However, these methods are based on a monostatic circular SAR system, and its signal model is thus not suitable for the GEO SA-BSAR system. On the other hand, the complex bistatic geometry of GEO SA-BSAR brings about obvious spatial variation [24]. The optimal SAR MTI processor [25,26,27], i.e., imaging space-time adaptive processing (ISTAP), was based on the Doppler parameters of the scene center and cannot compensate for the space-variant range migration and Doppler parameters. In addition, ISTAP is designed for monostatic SAR systems, thus, its adaptive filter will mismatch the moving target signal of the GEO SA-BSAR, causing serious SNR losses and detection performance degradation. Therefore, a novel GEO SA-BSAR MTI method should be investigated to improve the MTI performance.

To address the above problems, this paper proposes a long-time coherent integration STAP (LTCI-STAP) method for the GEO SA-BSAR system. Based on the characteristics of GEO SA-BSAR, an LTCI method suitable for the GEO SA-BSAR system and the STAP method are combined to obtain the maximum output SNR. First, to adapt the long SAT, a high-order slant range model of moving targets, expressed as a function of position and velocity vectors, is presented. Second, based on the derived signal model, a modified spatial filter and an LTCI filter is proposed to correct range migration, focus moving targets, and adaptively optimize the output SNR along with corresponding filtering parameters. Finally, a constant false alarm rate (CFAR) detection is carried out to determine the targets, and the space-variant Doppler parameters, solved from the filtering parameters, are used for estimating moving target positions and velocities. Numerical computer simulation experiments verify our technique. The main contributions of this paper include the following two aspects:

• We extend the bistatic signal model and the spatial filter to the cases of long SAT and various bistatic configuration. Based on this, the moving targets’ spatial signals can be accumulated coherently and the radial velocity can be estimated accurately.
• We propose a novel LTCI-STAP method to deal with the problem of high-order range migration and space-variant Doppler parameters in the GEO SA-BSAR system. Via the method, the moving targets are well focused, and their positions and velocities are accurately estimated.

The rest of the paper is organized as follows. The multichannel signal model for the GEO SA-BSAR moving target is derived in Section 2. The LTCI-STAP method is proposed in Section 3. The experimental results are shown in Section 4, and Section 5 concludes this paper.

2. GEO SA-BSAR Signal Model for Moving Target

The observation geometry of a GEO SA-BSAR system for a moving target is shown in Figure 1. The origin of the coordinate system is at the center of the illuminated area by GEO SAR, with the z-axis vertical to the ground plane, the y-axis pointing to the aircraft’s flight direction, and the x-axis completing the right-handed system. The aircraft, whose position is $r_R$, flies along the y-axis with a constant velocity $v_R$. GEO SAR, whose position is $r_T$, has a velocity of $v_T$ and an acceleration of $a_T$ at the aperture center moment (ACM). $u_{RP}$ and $u_{TP}$ are the unit vectors of the target-to-receiver and target-to-transmitter lines of sight (LOSs) at the ACM, respectively. The target is at the position of $r_p$ and moves with constant velocity $v_p$. There are M receiving channels along the aircraft’s track. The distance between the first and m-th channels is $y_m$, where the uniformly spaced channels satisfy $y_m=\left(m-1\right)d$. The first channel is regarded as the reference channel.

2.1. Slant Range Model

In the traditional GEO SA-BSAR MTI method [18], the moving-target-to-receiver slant range and the moving-target-to-transmitter slant range are modeled as a quadratic polynomial and a linear polynomial, considering that the accumulation time is less than 1 second. This results in a low SNR and poor image resolution for a moving target. By appropriately increasing the SAT of the GEO SA-BSAR system, a higher SNR and improved image resolution can be obtained. However, the traditional GEO SA-BSAR slant range model has a large error for a long SAT, as shown in Figure 2a,b. The simulation parameters are presented in

Table 1. GEO SA-SAR system parameters and GEO SAR orbital elements.

System Parameters		GEO SAR Orbital Elements
Wavelength	0.0566 m	Eccentricity	0
Bandwidth	60 MHz	Inclination	16°
Receiver’s Height	10 km	Semimajor Axis	42,164 km
Receiver’s Velocity	200 m/s	Longitude Ascending Node	113°
Synthetic Aperture Time	3 s
Number of Channels	10	Argument of Perigee	-
Channel Spacing	0.3 m

.

The errors caused by the long SAT can be overcome by high-order expansion [28]. Thus, the higher-order polynomial term has remained in the slant range model. Figure 2c,d show the slant range model error of the second-order polynomial for the transmitter and third-order polynomial for the receiver. The errors are less than $\lambda /16$, where $\lambda$ is the wavelength, proving that the approximations meet the requirement of focusing. Thus, the slant range of the transmitter is expanded to second-order by Taylor expansion, written as

$$\begin{array}{l}{R}_T\left(t_a;r_p,v_p\right)=\Vert r_T+v_T{t}_a+a_T{t}_a^2/2-r_p-v_p{t}_a\Vert \\ \approx \Vert r_T-r_p\Vert +{\left(v_T-v_p\right)}^T{u}_{TP}{t}_a+\frac{1}{2}\left[a_T^T{u}_{TP}+\frac{{\Vert \left(v_T-v_p\right)\times u_{TP}\Vert }^2}{\Vert r_T-r_p\Vert }\right]t_a^2\end{array}$$

(1)

where $t_a$ is the slow time. Figure 3 shows the errors of the slant range model for different target velocities and target ranges. We can see that the errors are less than $\lambda /8$, even for the targets with the velocities of 50 m/s, which meets the focus requirement. The slant range model is also valid for the targets located at different positions (the target ranges are 11 km, 12 km and 13 km).

The receiver has equipped multiple channels along the track, and the slant range model of the m-th channel can be expressed as

$$\begin{array}{l}{R}_{R,m}\left(t_a;r_p,v_p\right)=\Vert r_R+y_m+v_R{t}_a-r_p-v_p{t}_a\Vert \\ \approx \Vert r_R-r_p\Vert +u_{RP}^T{v}_p\frac{y_m}{v_R}+{\left(v_R-v_p\right)}^T{u}_{RP}\left(t_a+\frac{y_m}{v_R}\right)\\ +\frac{{\Vert \left(v_R-v_p\right)\times u_{RP}\Vert }^2}{2\Vert r_R-r_p\Vert }{\left(t_a+\frac{y_m}{v_R}\right)}^2-\frac{{\left(v_R-v_p\right)}^T{u}_{RP}{\Vert \left(v_R-v_p\right)\times u_{RP}\Vert }^2}{2{\Vert r_R-r_p\Vert }^2}{\left(t_a+\frac{y_m}{v_R}\right)}^3\end{array}$$

(2)

Then, the slant range history of the m-th channel can be expressed as

$$\begin{array}{c}{R}_{bi,m}\left(t_a;r_p,v_p\right)=R_T\left(t_a-y_m/v_R;r_p,v_p\right)+\Delta R_{T,t}\left(t_a\right)\\ +R_{R,m}\left(t_a;r_p,v_p\right)\end{array}$$

(3)

where

$$\begin{array}{l}\Delta R_{T,t}\left(t_a\right)=-{\left(v_T-v_p\right)}^T{u}_{TP}\frac{y_m}{v_R}\\ \underset{\Delta R_{k_{T2}}}{\underbrace{-\frac{1}{2}\left[u_{TP}{a}_T^T+\frac{{\Vert \left(v_T-v_p\right)\times u_{TP}\Vert }^2}{\Vert r_T-r_p\Vert }\right]\frac{y_m^2}{v_R^2}-\left[u_{TP}{a}_T^T+\frac{{\Vert \left(v_T-v_p\right)\times u_{TP}\Vert }^2}{\Vert r_T-r_p\Vert }\right]\frac{y_m}{v_R}{t}_a}}\end{array}$$

(4)

Considering that the GEO SAR Doppler rate is very low, $\Delta R_{k_{T2}}$ in (4) is on the order of ${10}^{-4}$ m, far less than $\lambda /16$, and can be ignored. This yields

$$\begin{array}{l}{R}_{bi,m}\left(t_a;r_p,v_p\right)={\left(u_{RP}+u_{TP}\right)}^T{v}_p\frac{y_m}{v_R}\\ +R_{bi0}-v_T^T{u}_{TP}\frac{y_m}{v_R}+k_{bi1}\left(t_a+\frac{y_m}{v_R}\right)+k_{bi2}{\left(t_a+\frac{y_m}{v_R}\right)}^2+k_{bi3}{\left(t_a+\frac{y_m}{v_R}\right)}^3\end{array}$$

(5)

where

$$\begin{array}{l}{R}_{bi0}=\Vert r_T-r_p\Vert +\Vert r_R-r_p\Vert \\ k_{bi1}={\left(v_T-v_p\right)}^T{u}_{TP}+{\left(v_R-v_p\right)}^T{u}_{RP}\\ k_{bi2}=\frac{1}{2}\left[a_T^T{u}_{TP}+\frac{{\Vert \left(v_T-v_p\right)\times u_{TP}\Vert }^2}{\Vert r_T-r_p\Vert }+\frac{{\Vert \left(v_R-v_p\right)\times u_{RP}\Vert }^2}{\Vert r_R-r_p\Vert }\right]\\ k_{bi3}=-\frac{{\left(v_R-v_p\right)}^T{u}_{RP}{\Vert \left(v_R-v_p\right)\times u_{RP}\Vert }^2}{2{\Vert r_R-r_p\Vert }^2}\end{array}$$

(6)

The coefficient $k_{bi1}$, representing the Doppler centroid, is related to the target velocity’s projection at the vector $u_{TP}+u_{RP}$, which is regarded as the radial velocity. Because of the satellite’s very high orbit, the coefficient $k_{bi2}$, representing the Doppler rate, is mainly determined by the target’s cross-range velocity related to the airborne receiver.

2.2. Multichannel Signal Model

It is assumed that the linear frequency modulation signal is transmitted by GEO SAR and is expressed as

$$s_T\left(t\right)=\mathrm{rect}\left[\frac{t}{T_p}\right]\exp \left\{j\pi k_r{t}^2\right\}\exp \left\{j2\pi f_{c}t\right\}$$

(7)

where $t$ is the fast time and $T_p$, $k_r$ and $f_c$ denote the pulse width, frequency modulation rate and carrier frequency, respectively.

The moving target echoes of the GEO SA-BSAR system can be expressed as the time delay of the transmitted signal. The echo of the m-th channel after demodulation and range compression can be represented as

$$s_m\left(t,t_a\right)=\sigma \sin \mathrm{c}\left[B\left(t-R_{bi,m}\left(t_a;r_p,v_p\right)/c\right)\right]w_a\left(t_a\right)\exp \left\{-j2\pi R_{bi,m}\left(t_a;r_p,v_p\right)/\lambda \right\}$$

(8)

where $B$ is the signal bandwidth. The target’s scattering coefficient is defined as $\sigma$. $w_a(\cdot )$ is the azimuth envelope. The azimuth spectrum of the m-th channel is derived by Fourier transformation with respect to the slow time and is expressed as (see Appendix A)

$$\begin{array}{c}{s}_m\left(t,f_a\right)=\sigma \sin \mathrm{c}\left[B\left(t-{\tilde{R}}_{bi,m}\left(f_a;r_p,v_p\right)/c\right)\right]W_a\left(f_a\right)\\ \cdot \exp \left\{-j\psi \left(f_a\right)\right\}\exp \left\{j2\pi d_m\left(f_a\right)/\lambda \right\}\end{array}$$

(9)

where $f_a$ is the Doppler frequency. ${\tilde{R}}_{bi,m}\left(f_a;r_p,v_p\right)$ is the slant range in frequency domain. The phase $\psi \left(f_a\right)$ in (9) represents the common phase for the different channels and can be written as a high-order polynomial:

$$\psi \left(f_a\right)=2\pi R_{bi0}/\lambda -\pi \lambda /2k_{bi2}\cdot {\left(f_a+k_{bi1}/\lambda \right)}^2-\pi {\lambda }^2{k}_{bi3}/4k_{bi2}^3\cdot {\left(f_a+k_{bi1}/\lambda \right)}^3$$

(10)

The second exponential term in (9) represents the phase difference between the m-th channel and the reference channel, where $d_m\left(f_a\right)$ is related to the channel spacing and can be expressed as

$$d_m\left(f_a\right)=\left\{v_T^T{u}_{TP}-v_p^T\left(u_{TP}+u_{RP}\right)+\lambda f_a\right\}\frac{y_m}{v_R}$$

(11)

The phase differences between channels are determined by the radial velocity. The signal model of a stationary target can be obtained by setting the velocity to zero. Furthermore, the signal model of background clutter is the superposition of all scattered points’ echoes, which can be represented as

$$c_m^{rd}\left(t,f_a\right)={\displaystyle \sum _{r_p}\begin{array}{l}{\sigma }_c\sin \mathrm{c}\left[B\left(t-{\tilde{R}}_{bi,m}\left(f_a;r_p,0\right)/c\right)\right]W_a\left(f_a\right)\\ \cdot \exp \left\{-j{\psi }_c\left(f_a\right)\right\}\exp \left\{j2\pi d_{c,m}\left(f_a\right)/\lambda \right\}\end{array}}$$

(12)

where the expression of ${\psi }_c$ and $d_{c,m}\left(f_a\right)$ can be obtained by setting ${\mathbf{v}}_p=0$ in $\psi$ and $d_m\left(f_a\right)$, respectively.

To image the moving targets, the signal in (9) is transformed into a two-dimensional frequency domain, and the signal model can be expressed as

$$s_m\left(f_r,f_a\right)=\sigma W_r\left(f_r\right)W_a\left(f_a\right)\cdot \exp \left\{-j\phi \left(f_r,f_a\right)\right\}\exp \left\{j2\pi d_m\left(f_a\right)/\lambda \right\}$$

(13)

where

$$\begin{array}{c}\phi \left(f_r,f_a\right)=\frac{2\pi R_{bi0}\left(f_r+f_c\right)}{c}-\frac{\pi c}{2k_{bi2}\left(f_r+f_c\right)}{\left[f_a+\frac{k_{bi1}\left(f_r+f_c\right)}{c}\right]}^2\\ -\frac{\pi c^2{k}_{bi3}}{4k_{bi2}^3{\left(f_r+f_c\right)}^2}{\left[f_a+\frac{k_{bi1}\left(f_r+f_c\right)}{c}\right]}^3\end{array}$$

(14)

3. LTCI-STAP Method for GEO SA-BSAR MTI

This section proposes an LTCI-STAP algorithm for GEO SA-BSAR systems. The proposed method considers range-azimuth coupling because of the spatial variation. The processing steps are depicted in Figure 4. The LTCI-STAP method and parameter iteration processing will be detailed in the following.

3.1. Preprocessing Step

During the preprocessing step, the multichannel echoes are transformed into the range-Doppler domain via range compression and azimuth Fourier transformation. There are slant range differences between channels. However, it is less than a range resolution because of the small channel spacing. The multichannel range-Doppler signals are written in the form of vectors and can be expressed as

$$z\left(t,f_a\right)=s\left(t,f_a\right)+c\left(t,f_a\right)+n$$

(15)

where $n$ represents the noise. $s(\cdot )$ and $c(\cdot )$ denotes the multichannel range-Doppler signals of moving targets and background clutter with dimensions of $M\times 1$, respectively.

$$\begin{array}{l}s\left(t,f_a\right)={\left[\begin{array}{cccc}{s}_1^{rd}\left(t,f_a\right)& s_2^{rd}\left(t,f_a\right)& \cdots & s_M^{rd}\left(t,f_a\right)\end{array}\right]}^T\\ c\left(t,f_a\right)={\left[\begin{array}{cccc}{c}_1^{rd}\left(t,f_a\right)& c_2^{rd}\left(t,f_a\right)& \cdots & c_M^{rd}\left(t,f_a\right)\end{array}\right]}^T\end{array}$$

(16)

3.2. Clutter Suppression and Beamforming

To suppress clutter and accomplish beamforming, a modified adaptive spatial filter based on the derived signal model in (9) is constructed in the Doppler domain. The spatial filter consists of the clutter covariance matrix and the target spatial steering vector. Since the clutter is independent of the noise and different Doppler frequency cells are independent of each other, the clutter covariance matrix with the dimension of $M\times M$ can be estimated using data from several range bins as follows:

$$R_Q\left(f_a\right)={\sigma }_n^2{I}_M+\frac{1}{N_r}{\displaystyle \sum _{l=t_r}^{t_r+N_r-1}c\left(l,f_a\right)c{\left(l,f_a\right)}^H}$$

(17)

where $I_M$ is the M-dimensional identity matrix. ${\sigma }_n^2$ is the thermal noise variance. $N_r$ is the number of range bins used to estimate the covariance matrix. It should be greater than $2M$. $c\left(l,f_a\right)$ is the spatial clutter signal at fast time $l$.

The modified steering vector of the target is related to the moving target’s radial velocity. It can be expressed as

$$p_t\left(f_a\right)=\left[\begin{array}{c}\exp \left\{j\frac{2\pi y_1}{\lambda v_R}\left[-v_T^T{u}_{TP}+v_r-\lambda f_a\right]\right\}\\ ⋮\\ \exp \left\{j\frac{2\pi y_M}{\lambda v_R}\left[-v_T^T{u}_{TP}+v_r-\lambda f_a\right]\right\}\end{array}\right]$$

(18)

where $v_r$ is the target’s radial velocity and can be expressed as

$$v_r=v_p^T\left(u_{TP}+u_{RP}\right)$$

(19)

Then, an optimal adaptive processor is constructed in space according to the steering vector and the clutter’s covariance matrix. After applying the filter, we get

$$y\left(r,f_a\right)=p_t{\left(f_a\right)}^H{R}_Q^{-1}\left(f_a\right)z\left(r,f_a\right)$$

(20)

Next, we transform the output signal $y\left(r,f_a\right)$ into the range frequency domain to obtain $s_t\left(f_r,t_a\right)$.

3.3. LTCI Processing

The long slant range of the GEO SA-BSAR system introduces a low SNR in the range-Doppler domain. To improve the target detection performance, the moving target should be focused on before detection. However, because of space-variant range migration and Doppler modulation, the target is hard to focus on by the traditional azimuth filter, which is based on the Doppler parameters of the scene center. In this section, an LTCI filter is constructed in the two-dimensional frequency domain to achieve range-azimuth decoupling and moving targets focusing at different positions. In the LTCI filter, there are two parameters, which correspond to the position and velocities of the moving target. They are estimated iteratively by searching the maximum output SNR. Finally, the space-variant Doppler parameters of the moving target can be solved according to the LTCI parameters.

First, the Doppler centroid must be compensated for in the range frequency domain by

$$s_{t,com}=s_t\left(f_r,t_a\right)\exp \left\{j\frac{2\pi \left(f_r+f_0\right)}{c}{k}_{10}{t}_a\right\}$$

(21)

where $k_{10}$ is the first-order coefficient of the scene center’s slant range.

Next, we transform the output signal into the two-dimensional frequency domain by applying a range Fourier transformation. According to the signal model in (9), the LTCI filter is constructed to compensate for the high-order phase term with respect to $f_a$ and can be expressed as

$$h_{LTCI}=\exp \left\{j\left[\frac{c\cdot f_a^2}{\left(f_r+f_c\right)}\beta -\frac{c^2\cdot f_a^3}{{\left(f_r+f_c\right)}^2}\gamma \right]\right\}$$

(22)

where

$$\begin{array}{l}\beta =-\frac{\pi }{2k_{bi2}}-\frac{3\pi k_{bi1}{k}_{bi3}}{4k_{bi2}^3}\\ \gamma =-\frac{\pi k_{bi3}}{4k_{bi2}^3}\end{array}$$

(23)

Then, the two-dimensional IFFT is conducted to obtain the output results. According to the beamforming and LTCI processing, the parameters $v_r$, $\beta$ and $\gamma$ are determined by the moving target’s positions and velocities. However, the targets’ position and motion parameters are usually unknown in practice. Thus, the parameters $v_r$, $\beta$ and $\gamma$ are estimated by searching the maximum of the output SNR over the range of possible $v_r$, $\beta$ and $\gamma$. The ranges of $\beta$ and $\gamma$ can be obtained according to (23), where $k_{bi1}$, $k_{bi2}$ and $k_{bi3}$ are the boundary values considering all of the target’s positions and velocities. For each range-azimuth pixel, the maximum of the output SNR is also compared with a threshold. The threshold is obtained according to a cell averaging-constant false alarm rate (CA-CFAR) detector. Thus, moving targets detection can be achieved.

However, if all of the potential parameters $\beta$ and $\gamma$ are calculated to obtain the maximum SCNR, the LTCI processing is complex in computation. A suitable grid spacing of parameters $v_r$, $\beta$ and $\gamma$ is considered to ensure that the target SNR loss does not lead to detection failure (i.e., SNR loss less than 3 dB).

3.4. Parameter Estimation

After moving target detection after LTCI processing, the targets’ positions and velocities can be obtained according to the estimated parameters ${\widehat{v}}_r$, $\widehat{\beta }$ and $\widehat{\gamma }$. Because of the targets’ large range migrations caused by their motion, the targets have both azimuth and range displacement. Thus, the range and azimuth positions where the target is detected are not its real position. Although the targets’ radial velocities $v_r$ have been estimated by the above processing, the radial velocity directions are still unknown because of the targets’ false locations. The target’s real positions and motion parameters are estimated in the following step. After LTCI processing, the output signal can be expressed as

$$\begin{array}{c}{s}_{out}\left(r,t_a\right)=\sigma \sin \mathrm{c}\left[B\left(r-R_{bi0}+\frac{k_{bi1}^2}{4k_{bi2}}+\frac{k_{bi1}^3{k}_{bi3}}{8k_{bi2}^3}\right)\right]\sin \mathrm{c}\left[B_a\left(t_a+\frac{k_{bi1}}{2k_{bi2}}+\frac{3\pi k_{bi1}^2{k}_{bi3}}{8k_{bi2}^3}\right)\right]\\ \cdot \exp \left\{j\frac{2\pi R_{bi0}}{\lambda }\right\}\exp \left\{-j\frac{\pi k_{bi1}^2}{2\lambda k_{bi2}}\left[1+\frac{k_{bi1}{k}_{bi3}}{2k_{bi2}^2}\right]\right\}\end{array}$$

(24)

We assume that the detected target’s range and azimuth time are $\widehat{r}$ and ${\widehat{t}}_0$ at the output image, respectively. According to (23) and (24), the following equation set can be obtained:

$$\{\begin{array}{l}{R}_{bi0}-\frac{k_{bi1}^2}{4k_{bi2}}-\frac{k_{bi1}^3{k}_{bi3}}{8k_{bi2}^3}=\widehat{r}\\ -\frac{k_{bi1}}{2k_{bi2}}-\frac{3k_{bi1}^2{k}_{bi3}}{8k_{bi2}^3}+N\cdot T_s={\widehat{t}}_0\\ -\frac{\pi }{2k_{bi2}}-\frac{3\pi k_{bi1}^2{k}_{bi3}}{4k_{bi2}^3}=\widehat{\beta }\\ -\frac{\pi k_{bi3}}{4k_{bi2}^3}=\widehat{\gamma }\end{array}$$

(25)

where $T_s$ is the SAT and $N$ is the ambiguous number. It is worth noting that the displacement of the moving target caused by the radial velocity may be beyond the mapping area. It results in an ambiguous position for a moving target. Thus, we should estimate the ambiguous number $N$ first. Because the ambiguous position is totally caused by the radial velocity, the ambiguous number can be estimated by

$$N=\mathrm{ceil}\left(-\frac{{\tilde{k}}_{bi1}+{\widehat{v}}_r}{2T_s{\tilde{k}}_{bi2}}-\frac{1}{2}\right)$$

(26)

where ${\tilde{k}}_{bi1}$ and ${\stackrel{⌢}{k}}_{bi2}$ are approximated by

$${\tilde{k}}_{bi1}=\frac{\widehat{\beta }-\sqrt{{\widehat{\beta }}^2-6\pi \widehat{\gamma }{\widehat{t}}_0}}{3\widehat{\gamma }}$$

(27)

$${\tilde{k}}_{bi2}=-\frac{\pi }{2\sqrt{{\widehat{\beta }}^2-6\pi \widehat{\gamma }{\widehat{t}}_0}}$$

(28)

Then, the parameters $R_{bi0}$, $k_{bi1}$, $k_{bi2}$ and $k_{bi3}$ can be estimated by solving Equations (25). The results are

$${\widehat{R}}_{bi0}=\widehat{r}+\frac{1}{27\pi {\widehat{\gamma }}^2}\left[{\widehat{\beta }}^3-9\pi \widehat{\beta }\widehat{\gamma }\left({\widehat{t}}_0-N\cdot T_s\right)-{\left(\sqrt{{\widehat{\beta }}^2-6\pi \widehat{\gamma }\left({\widehat{t}}_0-N\cdot T_s\right)}\right)}^3\right]$$

(29)

$${\widehat{k}}_{bi1}=\frac{\widehat{\beta }-\sqrt{{\widehat{\beta }}^2-6\pi \widehat{\gamma }\left({\widehat{t}}_0-N\cdot T_s\right)}}{3\widehat{\gamma }}$$

(30)

$${\widehat{k}}_{bi2}=-\frac{\pi }{2\sqrt{{\widehat{\beta }}^2-6\pi \widehat{\gamma }\left({\widehat{t}}_0-N\cdot T_s\right)}}$$

(31)

$${\widehat{k}}_{bi3}=\frac{{\pi }^2\widehat{\gamma }}{2{\left(\sqrt{{\widehat{\beta }}^2-6\pi \widehat{\gamma }\left({\widehat{t}}_0-N\cdot T_s\right)}\right)}^3}$$

(32)

Because of the complex bistatic configuration, the target’s position coordinates $r_p$ have to be obtained by solving the following nonlinear equations:

$$\{\begin{array}{l}\Vert r_T-r_p\Vert +\Vert r_R-r_p\Vert ={\widehat{R}}_{bi0}\\ v_T^T{u}_{TP}+v_R^T{u}_{RP}={\widehat{k}}_{bi1}+{\widehat{v}}_r+k_{10}\end{array}$$

(33)

Thus, the target’s position coordinates are estimated as ${\widehat{r}}_p$ by Newton’s method. Once the target’s real position is obtained, its radial velocity direction can be determined, which is $\frac{r_T-{\widehat{r}}_p}{\Vert r_T-{\widehat{r}}_p\Vert }+\frac{r_R-{\widehat{r}}_p}{\Vert r_R-{\widehat{r}}_p\Vert }$.

The target velocities can also be estimated by solving the following nonlinear equations:

$$\{\begin{array}{l}{v}_p^T\left(\frac{r_T-{\widehat{r}}_p}{\Vert r_T-{\widehat{r}}_p\Vert }+\frac{r_R-{\widehat{r}}_p}{\Vert r_R-{\widehat{r}}_p\Vert }\right)={\tilde{v}}_r\\ \frac{1}{2}\left[a_T^T{u}_{TP}+\frac{{\Vert \left(v_T-v_p\right)\times \left(r_T-{\widehat{r}}_p\right)/\Vert r_T-{\widehat{r}}_p\Vert \Vert }^2}{\Vert r_T-{\widehat{r}}_p\Vert }+\frac{{\Vert \left(v_R-v_p\right)\times \left(r_R-{\widehat{r}}_p\right)/\Vert r_R-{\widehat{r}}_p\Vert \Vert }^2}{\Vert r_R-{\widehat{r}}_p\Vert }\right]={\widehat{k}}_{bi2}\end{array}$$

(34)

The target’s velocity can be estimated by Newton’s method. After processing, the target’s real position and velocity can be obtained.

4. Simulation Experiments

To verify the effectiveness of the proposed algorithm, we conducted some simulations. The system parameters are shown in

Table 2. GEO SA-BSAR MTI bistatic configuration.

Configuration	1 (Side Looking)	2 (Squint-Forwards Looking)	3 (Squint-Backwards Looking)
${\theta }_R$ (deg)	35	35	35
$\varphi$ (deg)	0	30	330
$\psi$ (deg)	0	20	340

. Because of GEO SAR’s large footprint, the airborne receiver can collect echoes at different bistatic configurations. Several bistatic configurations, where the aircraft works in side-looking mode, squint-forwards-looking mode and squint-backwards-looking mode, are considered in the simulations. The bistatic configuration parameters are shown in Table 2, where ${\theta }_R$ is the aircraft incident angle, $\varphi$ is the bistatic angle and $\psi$ is the angle between the velocities of the GEO satellite and the aircraft.

4.1. Simulations of Single Point

We set 36 stationary point targets, which are uniformly distributed with 200 m intervals in the x-direction and 200 m in the y-direction. Figure 5 shows the distribution of stationary targets. In addition, a moving target with a velocity of $\left(6,4,0\right)$ m/s at the position of $\left(200,200,0\right)$ m is set in the scene. It is assumed that the signal-to-clutter ratio (SCR) after range compression is −60 dB and that the signal-to-noise ratio (SNR) after range compression is −20 dB.

In the preprocessing step, the multichannel echoes are range compressed and transformed into the Doppler domain. Figure 6 shows the data in the Range-Doppler domain. It can be seen that the moving target is corrupted by clutter and cannot be detected. Figure 6a–c represent the signals of Conf. 1, Conf. 2 and Conf. 3, respectively. The clutter shows distinct shapes between the graphs because the different bistatic geometries lead to different range walk migrations in the GEO SA-BSAR system.

Then, the clutter is suppressed by clutter’s covariance matrix. The clutter covariance matrix is estimated using all of the range cells, and the results are shown in Figure 7. Clutter cannot be seen in the graphs, and the curves shown in the graphs represent the moving target. The moving target’s data face range-azimuth coupling, so the target cannot be focused only by azimuth filtering, as done in ISTAP.

In the next step, LTCI processing is carried out in the two-dimensional frequency domain to achieve range-azimuth decoupling and moving target focusing. The data are transformed into the image domain by 2D IFFT. Figure 8 gives the output results using the radial velocity and LTCI parameters to achieve the maximum SNR. The target has a high SNR after processing and is easily detected. An incorrect radial velocity and LTCI parameters are used in the processing, and the output results are shown in Figure 9. The incorrect LTCI parameters induce a defocusing for the moving target, leading to a low SNR and detection performance degradation. Thus, to detect the moving targets successfully, an LTCI filter bank is required to span a range of parameters $\beta$ and $\gamma$.

The target’s radial velocity and LTCI parameters can be estimated by searching the maximum output SNR. The estimated results are shown in

Table 3. Estimation results of the moving targets’ positions and velocities.

Estimated Parameters	${\mathit{v}}_{\mathit{r}}(\mathbf{m}/\mathbf{s})$	$\mathit{\beta }$ (s2/m)	$\mathit{\gamma }$ (s2/m2)	$\mathit{x}$ (m)	$\mathit{y}$ (m)	${\mathit{v}}_{\mathit{x}}(\mathbf{m}/\mathbf{s})$	${\mathit{v}}_{\mathit{y}}(\mathbf{m}/\mathbf{s})$
Conf. 1	4.2	−1.02	$3.58\times {10}^{-5}$	201.08	204.06	3.91	6.37
Conf. 2	4.8	−1.00	−0.001	199.61	190.14	4.16	5.89
Conf. 3	3.7	−1.02	$3.77\times {10}^{-4}$	200.05	196.93	4.05	5.98

. Note that the radial components are distinct in different configurations. According to the principle of error propagation, the estimated parameter accuracies are different. Details about the relationship between the bistatic configuration and parameter estimation accuracy can be seen in [29]. The moving target has been located to its real position, and its velocity has been estimated accurately.

4.2. Experiments with Scene Simulation

In this section, the effectiveness of LTCI-STAP processing is verified by experiments with a scene simulation. A stationary scene consisting of a road and its surrounding area is set, whose SAR imaging result is presented in Figure 10. Several moving targets are set along the road, whose velocities are either 8 m/s or 5 m/s in the direction parallel to the road. Their positions are presented in

Table 4. Set positions and velocities of the moving targets.

Target No.	T1	T2	T3	T4
x-axis (m)	−233	−75	64	188
y-axis (m)	−224	−116	−16	67
$v_x$ (m/s)	−4.11	4.11	6.58	−6.58
$v_y$ (m/s)	−2.84	2.84	4.55	−4.55

.

It is assumed that the SCR in range compression is −60 dB, and the SNR after range compression is −20 dB. The preprocessing results can be seen in Figure 11. The moving targets are submerged by noise and clutter in the range-Doppler domain and cannot be detected directly. Figure 12 shows the results after clutter suppression, where clutter cannot be seen, and the curves represent the moving targets.

After LTCI-STAP processing, different moving targets will be focused for their corresponding radial velocities and LTCI parameters. Figure 13 and Figure 14 show the output results of the processing for different LTCI parameters. It can be seen that high SNRs can be achieved only when the LTCI parameters match the set ones. Thus, the LTCI parameters should be determined as accurately as possible to improve the SNRs and thus detect all the moving targets, especially for the targets with low SNR.

For every detected target, its position and velocity parameters can be estimated by the radial velocity and LTCI parameters that maximize the target’s SNR. The detected targets are marked on the images in Figure 15. The red crosses represent the estimated locations of the moving targets, and the green circles are their real locations. It can be seen that all the targets’ positions are accurately estimated, and they are located near their real positions. The velocity estimation results and errors are shown in

Table 5. Position and motion parameters estimation results.

Parameter	Configuration	T1	T2	T3	T4
$\widehat{x}$ (m)	Conf. 1	−233.24	−74.65	65.05	187.41
	Conf. 2	−231.35	−75.85	63.65	187.19
	Conf. 3	−233.82	−74.23	62.69	188.61
$\widehat{y}$(m)	Conf. 1	−227.90	−109.59	−17.54	65.10
	Conf. 2	−232.41	−110.33	−15.35	66.63
	Conf. 3	−239.28	−109.42	−19.40	63.62
${\widehat{v}}_x$ (m/s)	Conf. 1	−4.06	4.01	6.59	−6.57
	Conf. 2	−4.01	4.06	6.58	−6.58
	Conf. 3	−3.89	4.02	6.64	−6.54
${\widehat{v}}_y$ (m/s)	Conf. 1	−2.46	3.43	4.82	−4.65
	Conf. 2	−2.47	2.61	4.59	−4.62
	Conf. 3	−3.17	2.90	4.60	−4.81

and

Table 6. The estimation errors for different configurations.

Parameter	Configuration	T1	T2	T3	T4
$\Delta x$ (m)	Conf. 1	0.24	0.35	1.05	0.59
	Conf. 2	1.65	0.85	0.35	0.81
	Conf. 3	0.82	0.77	1.31	0.61
$\Delta y$ (m)	Conf. 1	3.90	6.41	1.54	1.90
	Conf. 2	8.41	5.67	0.65	0.37
	Conf. 3	15.28	6.58	3.40	3.38
$\Delta v_x$ (m/s)	Conf. 1	0.06	0.10	0.01	0.01
	Conf. 2	0.10	0.06	0.01	0.02
	Conf. 3	0.22	0.09	0.06	0.04
$\Delta v_y$ (m/s)	Conf. 1	0.38	0.59	0.27	0.10
	Conf. 2	0.37	0.23	0.04	0.07
	Conf. 3	0.33	0.06	0.05	0.26

, respectively, indicating that all of the targets’ velocities and positions were estimated accurately.

Moreover, we have compared the proposed method with other SAR GMTI methods including the ISTAP method (the optimal processing method for spaceborne SAR MTI) and the VSAR method for GEO SA-BSAR MTI.

The simulation results are shown in Figure 16. The first column in Figure 15 presents the results of ISTAP processing. Since this algorithm is suitable for monostatic spaceborne SAR systems in a simple geometry (i.e., the trajectory is straight line and the sensor works in side-looking mode), the processor is seriously mismatched with the bistatic radar signal of the studied GEO SA-BSAR system. Thus, the moving targets in the GEO SA-BSAR system defocus, and the output SNR is very low.

The second column of Figure 16 shows the VSAR processing results. The VSAR method is designed for the GEO SA-BSAR system; thus, the performances are better than those of the classical ISTAP method. However, as the method only considers a short SAT and the aircraft is assumed to work in the side-looking mode, the targets are defocused to some degree, especially in the case of the aircraft in squint-forward-looking and squint-backward-looking mode.

The processing results of the proposed LTCI-STAP method are shown in the final column of Figure 16. It can be seen that, compared with the ISTAP and VSAR methods, the proposed method can output the highest SNR and improve the targets’ detection performance. The reason is that, due to the accurate signal model and the optimal processor, the proposed method can achieve coherent accumulation of moving targets in different bistatic configurations.

4.3. Computational Complexity Analysis

The computation of the proposed algorithm greatly depends on the size of raw data and the number of parameters to be searched, including the target’s radial velocity, LTCI parameters $\beta$ and $\gamma$. It is assumed that the size of the raw data is $N_a$ in azimuth and $N_r$ in range, and the number of channels is M. The search number of radial velocities, LTCI parameter $\beta$ and $\gamma$ is $N_{v_r}$, $N_{\beta }$ and $N_{\gamma }$, respectively. In the following text, we mainly analyze the time complexity of the clutter suppression, beamforming and LTCI processing.

The time complexity of clutter suppression mainly consists of the calculation of the covariance matrix. For a range-Doppler cell, if we estimate the covariance matrix (dimension is $M\times M$) by using the data of all range cells, the time complexity is $O\left(N_r\cdot M^2\right)$. Then, the covariance matrix is inverted and multiplied by the echo signal to complete clutter suppression, whose time complexity is $O\left(M^3\right)$. As all Doppler cells need to be processed, the time complexity of the whole clutter suppression processing is $O\left(N_a{M}^2\left(N_r+M\right)\right)$.

The time complexity of beamforming and LTCI processes depends on the search number of estimated parameters. For a given radial velocity, the beamforming is processed by multiplying the spatial steering vector and the results of clutter suppression; thus, the time complexity is $O\left(N_a{N}_{r}M\right)$. Then, for the given parameters $\beta$ and $\gamma$, the LTCI processing is conducted by multiplying the constructed LTCI filter with the echoes, whose time complexity is $O\left(N_a{N}_r\right)$. The time complexity of the whole processing is $O\left(N_{v_r}{N}_r\left(M+N_{\beta }{N}_{\gamma }{N}_a\right)\right)$.

In summary, the computational complexity of the proposed method is $O\left(N_a{M}^2\left(N_r+M\right)+N_{v_r}{N}_r\left(M+N_{\beta }{N}_{\gamma }{N}_a\right)\right)$. Using the Montage Jintide(R) C6248R 3.00 GHz processor for simulation, the whole processing takes 947.44 s.

4.4. Parameter Estimation Accuracy Analysis

We also evaluated the parameter estimation accuracy for different SNRs. The simulation parameters of the radar system and the moving target are the same as Section 4.1. The estimation errors are shown in

Table 7. The estimation errors for different SNRs.

SNR (dB)		−25	−20	−15	−10	−5	0
$\Delta x$ (m)	Conf. 1	1.01	1.20	0.99	1.35	0.99	0.89
	Conf. 2	0.83	0.39	0.0026	0.41	0.0026	0.041
	Conf. 3	1.61	1.66	1.68	1.62	1.62	1.65
$\Delta y$ (m)	Conf. 1	6.23	10.47	1.70	7.95	1.89	8.12
	Conf. 2	3.84	9.56	15.72	9.67	15.72	9.66
	Conf. 3	21.05	21.50	21.59	21.36	21.36	21.40
$\Delta v_x$ (m/s)	Conf. 1	0.087	0.07	0.17	0.0045	0.090	0.0016
	Conf. 2	0.074	0.16	0.25	0.16	0.25	0.16
	Conf. 3	0.33	0.32	0.33	0.32	0.32	0.32
$\Delta v_y$ (m/s)	Conf. 1	0.27	0.91	0.24	0.47	0.33	0.27
	Conf. 2	0.18	0.14	0.12	0.11	0.12	0.11
	Conf. 3	0.012	0.10	0.08	0.11	0.11	0.11

. $\Delta x$, $\Delta y$, $\Delta v_x$ and $\Delta v_y$ represent the accuracies of $x$, $y$, $v_x$ and $v_y$, respectively. Since the proposed method is able to estimate the motion parameters accurately, the signal with long synthetic aperture time can be coherently integrated to output a high SNR. The output SNR of the proposed method can reach 15 dB (ensuring that the moving targets can be detected [30]) when the SNR of the range-compressed signal is −25 dB. It can be seen that when the SNR is greater than −25 dB, the velocity measurement accuracy is better than 1 m/s, and the positioning accuracy is on the order of 10 m. Because the parameter estimation accuracy is mainly determined by the parameter search step, the SNR has little influence on the estimation accuracy.

5. Conclusions

This paper proposes an LTCI-STAP algorithm for the GEO SA-BSAR system. The GEO SA-BSAR moving target multichannel signals are modeled. The spatial filter and the LTCI filter are modified to overcome the problem of space-variant range migration and Doppler parameters and to obtain the maximum output SNR. The proposed MTI algorithm can detect and relocate moving targets in different GEO SA-BSAR bistatic configurations and achieve optimal SNRs without any prior information. Simulation experiments for different bistatic configurations are carried out to verify the effectiveness.