*3.1. RCM, Doppler Centre Shift and QDFM Compensation*

In accordance with Equation (5), the discrete form of *s*1(*f*,*tn*) omitting the azimuth window function is expressed as

$$\begin{array}{lcl} s\_1(m, n) &= \text{rect}\left(\frac{m\Delta f}{\mathcal{B}}\right) \exp\left[-\frac{j4\pi}{c}(m\Delta f + f\_c) \\ &\cdot \left(R\_0 - v\_c n \Delta t\_n + \frac{(v - v\_a)^2 - R\_0 a\_c}{2R\_0} n^2 \Delta t\_n^2 + \frac{v\_c(v - v\_a)^2 + R\_0 a\_c (v\_b - v)}{2R\_0^2} n^3 \Delta t\_n^3\right)\right] \end{array} \tag{7}$$

where *m*(*m* = −*M*/2,−*M*/2 + 1, ··· ,−*M*/2 − 1, *M*/2) denotes the discrete range frequency number index related to continuous range frequency *f*, *M* is assumed to be an even integer, Δ*f* denotes the range frequency sampling interval, *n*(*n* = −*N*/2,−*N*/2+1, ··· , *N*/2 − 1, *N*/2) represents the discrete slow time number index related to continuous slow time *tn*, *N* is assumed to be an even integer and Δ*tn* indicates the pulse repetition interval.

Given the symmetrical property of the discrete slow time sequence, a new signal after reversing the discrete slow time sequence for each range frequency yields [36,37]

$$\begin{aligned} s\_1(m, \overset{\leftarrow}{n}) &= s\_1(m, -n) = \text{rect}\left(\frac{m\Delta f}{\text{B}}\right) \exp\left[-\frac{\text{i}4\pi}{\text{c}}(m\Delta f + f\_\text{c}) \\ &\cdot \left(R\_0 + v\_\text{c}n\Delta t\_n + \frac{\left(v - v\_\text{s}\right)^2 - R\_0\mu\_\text{c}}{2R\_0}n^2\Delta t\_n^2 - \frac{v\_\text{c}\left(v - v\_\text{s}\right)^2 + R\_0\mu\_\text{s}\left(v\_\text{s} - v\right)}{2R\_0^2}n^3\Delta t\_n^3\right)\right] \end{aligned} \tag{8}$$

where "←" denotes the sequence reversing operation. As shown in Equations (7) and (8), the effects of *tn*-term and *t* 3 *<sup>n</sup>*-term can be removed by multiplying Equation (7) by Equation (8). Thus, the corresponding result yields

$$\begin{split} s\_2(m, n) &= s\_1(m, n) \cdot s\_1(m, \overleftarrow{n}) \\ &= \text{rect}\left(\frac{m\Delta f}{\mathcal{B}}\right) \text{exp}\left[-\frac{j8\pi}{c}(m\Delta f + f\_c) \left(R\_0 + \frac{\left(v - v\_x\right)^2 - R\_0 a\_c}{2R\_0} n^2 \Delta t\_n^2\right)\right]. \end{split} \tag{9}$$

As described in Equation (9), the FRCM, TRCM and QDFM are effectively compensated. In the meantime, the Doppler centre blur and spectrum distributions belonging to cases II and IV are avoided given that the Doppler centre shift is effectively removed. Nevertheless, the effect induced by *t* 2 *<sup>n</sup>*-term still exists. Thus, the SRCM and FDFM should be effectively eliminated.

As verified in [16,33,38,39], the SRCM of the common metre-level range resolution SAR system depends on the SAR platform velocity given that the target along-track velocity and cross-track acceleration are considerably smaller than the velocity of SAR platform. Therefore, the SRCM can be removed by constructing the correction function on the basis of the SAR platform velocity *v* as long as the residual SRCM correction error is smaller than one range resolution bin. In accordance with Equation (9), an improved SRCM correction function is constructed as follows:

$$H\_{RCM}(m,n) = \exp\left(\frac{j4\pi m\Delta f v^2 n^2 \Delta t\_n^2}{cR\_0}\right). \tag{10}$$

The SRCM correction error by applying the correction function in Equation (10) is obtained as

$$
\Delta R\_{RCM} = \left| \left[ \frac{\left(\upsilon\_a^2 - 2\upsilon\upsilon\_a\right) - R\_0 a\_c}{4R\_0} \right] \left(T\_a\right)^2 \right|. \tag{11}
$$

We suppose that a maneuvering target, which is denoted by Target A, is considered. The target parameters are *va* = 30 m/s, *aa* = <sup>−</sup>1 m/s2, *vc* = 32 m/s and *ac* = <sup>−</sup>1 m/s2. The main radar parameters are *fc* = 10 GHZ, *B* = 80 MHZ, *PRF* = 1400 HZ, *v* = 250 *m*/*s*, *R*<sup>0</sup> = 6000 m and *Ta* = 1.2 s. The SRCM correction error Δ*RRCM* is calculated as 0.486 m, which is smaller than one range resolution bin. Example A without noise is also presented. Figure 3a shows the target trajectory after range compression. The target suffers from severe RCM effect. As illustrated in Figure 3b, only the SRCM remains after FRCM and TRCM correction. Figure 3c exhibits the result of SRCM compensation by using the correction function in Equation (10). The SRCM is effectively removed, and the trajectory of the moving target is located on the same range bin. The residual SRCM correction error by using the SAR platform velocity can be ignored and the correction function in Equation (10) is considered effective for the common metre-level range resolution SAR system.

**Figure 3.** Results of Example A: (**a**) trajectory after range compression; (**b**) trajectory after first-order range cell migration and third-order range cell migration correction; (**c**) result of second-order range cell migration correction by using the correction function in Equation (10).

After multiplying Equation (9) by Equation (10), we have

$$\begin{split} s\_3(m, n) &= s\_2(m, n) \cdot H\_{\text{RCM}}(m, n) \\ &\approx \text{rect}\left(\frac{m\Delta f}{\overline{\mathcal{B}}}\right) \exp\left(-\frac{j8\pi}{c} \left(m\Delta f + f\_c\right)R\_0\right) \exp\left[-\frac{j8\pi}{\lambda} \frac{(v - v\_a)^2 - R\_0 a\_c}{2R\_0} n^2 \Delta t\_n^2\right]. \end{split} \tag{12}$$

After the range IFFT is applied to Equation (12), the continuous time expression is obtained as

$$s\_3(t, t\_n) = \text{sinc}\left[B\left(t - \frac{4R\_0}{c}\right)\right] \exp\left[-\frac{j8\pi}{\lambda} \frac{\left(\upsilon - \upsilon\_a\right)^2 - R\_0 a\_0}{2R\_0} t\_n^2\right].\tag{13}$$

According to Equation (13), the moving target is focused in the same range bin, and the target signal along azimuth slow time dimension can be modelled as a second-order phase signal (i.e., an LFM signal). The chirp rate (i.e., second-order phase coefficient) of the target signal in Equation (13) should be effectively estimated and compensated to focus the moving target.

#### *3.2. Estimation of Chirp Rate by 1D SCFT*

Several typical methods, such as, discrete chirp-Fourier transform (DCFT) [34,40], fractional Fourier transform (FRFT) [29,41], Lv's distribution (LVD) [26,42], and CICPF [17,31], have been proposed to estimate the chirp rate of the target signal in Equation (13). As for the parameter searching-based methods, DCFT and FRFT are computationally prohibitive due to the brute-force grid searching operation. Subsequently, the LVD and CICPF have been proposed without the searching process. However, these methods must transform the 1D LFM signal into a 2D parameter space to obtain the final estimation result. These methods still have a large computational burden due to the 2D data processing operation. Considering that the chirp rate of signal in Equation (13) has been doubled, the constant factor ε is introduced to expand the application scope. The corresponding 1D SCFT with constant factor ε yields

$$\begin{split} s\_3(\frac{4k\_0}{\varepsilon}, f\_{\frac{t^2}{n}}) &= \int s\_3(\frac{4k\_0}{\varepsilon}, t\_n) \exp\{-j2\pi f\_{t\_n^2} \varepsilon t\_n^2\} dt\_n^2 \\ &= \delta(f\_{t\_n^2} + \frac{2(v - v\_a)^2 - 2R\_0 a\_\varepsilon}{\varepsilon \Lambda R\_0}) \end{split} \tag{14}$$

where *<sup>s</sup>*3( <sup>4</sup>*R*<sup>0</sup> *<sup>c</sup>* , *tn*) = exp# <sup>−</sup> *<sup>j</sup>*8<sup>π</sup> λ (*v*−*va*) 2 −*R*0*ac* <sup>2</sup>*R*<sup>0</sup> *t* 2 *n* \$ , *f t* 2 *<sup>n</sup>* is the azimuth scaled frequency variable corresponding to *t* 2 *<sup>n</sup>* and δ(·) denotes the Dirac delta function. The selection criterion of constant factor ε is discussed in Appendix A.

In accordance with Equation (14), the chirp rate is estimated as follows:

$$\frac{2\left(v-v\_a\right)^2 - 2R\_0a\_c}{\lambda R\_0} = -\varepsilon f\_{t\_n'} \tag{15}$$

where ˆ *f t* 2 *<sup>n</sup>* denotes the peak position in the azimuth scaled frequency domain.

Example B without noise is provided to validate the constant factor ε. We consider a target that is denoted by Target B. The parameters of Target B are as follows: *va*<sup>2</sup> = <sup>−</sup>29.5 m/s and *ac*<sup>2</sup> = 0.8 m/s2. The radar parameters are the same as those of Example A. Figure 4a,b show the result using the 1D SCFT with the constant factor ε = 1. The defocusing result is obtained in Figure 4a,b because the chirp rate of the signal for Target B exceeds the scope of parameter estimation. Figure 4c,d depict the result using the 1D SCFT with the constant factor ε = 2. Given that the constant factor ε = 2 expands the parameter estimation scope, a clear peak appears in Figure 4c,d. Therefore, the results of Example B demonstrate the validity of the constant factor ε.

**Figure 4.** Results of Example B: (**a**) result using constant factor ε = 1; (**b**) dB version of amplitude for Figure 4a; (**c**) result using constant factor ε = 2; (**d**) dB version of amplitude for Figure 4c.

#### *3.3. Azimuth Focusing by Matched Filter Based on Estimated Chirp Rate*

The signal in Equation (9) is transformed into the range and azimuth frequency domain yields

$$s\_2(f, f\_{l\_n}) = \text{rect}\left(\frac{f}{B}\right) \exp\left[-\frac{j8\pi}{c}(f+f\_c)R\_0\right] \exp\left\{j\pi \frac{cR\_0}{4(f+f\_c)\left[\left(\upsilon-\upsilon\_a\right)^2 - R\_0 a\_c\right]}f\_{l\_n}^2\right\}.\tag{16}$$

By substituting Equation (15) into Equation (16), we have

$$s\_2(f\_\prime, f\_{t\_n}) = \text{rect}\left(\frac{f}{B}\right) \exp\left[-\frac{j8\pi}{c}(f+f\_c)R\_0\right] \exp\left\{-j\pi \frac{c}{2(f+f\_c)\varepsilon\lambda\hat{f}\_{t\_n}^2}f\_{t\_n}^2\right\}.\tag{17}$$

In accordance with Equation (17), the matched filter based on the estimated chirp rate is constructed as follows:

$$H\_1(f, f\_{t\_n}, \hat{f}\_{t\_n^2}) = \exp\left\{j\pi \frac{c}{2(f + f\_c)\varepsilon\lambda \hat{f}\_{t\_n^2}} f\_{t\_n}^2\right\}.\tag{18}$$

After Equation (17) is multiplied by Equation (18) and 2D IFFT is performed, we have

$$s\_4(t, t\_n) = \text{sinc}\left[B\left(t - \frac{4R\_0}{c}\right)\right] \text{sinc}(f\_D t\_n),\tag{19}$$

where *fD* denotes the bandwidth of the signal in Equation (9).

As shown in Equation (19), the RCM, Doppler centre shift and DFM of the maneuvering target are effectively compensated. After 2D IFFT is performed, the moving target is refocused in the range–azimuth time domain; accordingly, the subsequent moving target processing operations are facilitated [22].

Figure 5 provides the flowchart of the presented algorithm. As shown in Figure 5, the proposed algorithm after range compression mainly includes three steps. The steps are summarised as follows:

	- (1.1) Apply the range FFT to the range compressed signal *s*1(*t*, *tn*), and obtain *s*1(*f*, *tn*).
	- (1.2) Perform the sequence reversing operation to *s*1(*f*, *tn*), and calculate Equation (9) to obtain *s*2(*f*, *tn*).
	- (1.3) Construct the SRCM correction function *HRCM*(*f*, *tn*), and calculate Equation (12) to obtain *s*3(*t*, *tn*).
	- (3.1) Construct the azimuth match filter *H*1(*f*, *ftn* , ˆ *f t* 2 *n* ), and multiply Equation (16) by *H*1(*f*, *ftn* , ˆ *f t* 2 *n* ) to obtain *s*4(*f*, *ftn* ).
	- (3.2) Perform 2D IFFT to *s*4(*f*, *ftn* ), and obtain the final focused result *s*4(*t*, *tn*),

where step 1.1 is used to transform the range compressed signal into the range frequency domain, step 1.2 is applied to correct FRCM and TRCM, step 1.3 is performed to remove SRCM, step 3.1 is utilised to compensate LDFM and step 3.2 is used to obtain the final refocused result.

**Figure 5.** Flowchart of the presented algorithm.

#### **4. Analysis Related to the Proposed Algorithm**

#### *4.1. Analysis of the SRCM Correction Function Mismatch*

Considering that the SRCM correction error Δ*RRCM* may be larger than one range resolution cell in some high-range resolution SAR systems, the mismatch of SRCM correction function in Equation (10) will appear. As for this case, we can utilise the SOKT to correct the SRCM. However, the Doppler bandwidth of Equation (9) is easily larger than one PRF band because the chirp rate has been doubled after the previous processing steps. Then, the target trajectory will split into multiple parts after directly using SOKT due to the Doppler spectrum ambiguity (i.e., the azimuth spectrum distribution belonging to case V). Therefore, the Doppler bandwidth of Equation (9) should be effectively compressed. In [25,34], the deramp functions are constructed to compress the azimuth spectrum. However, these deramp functions ignore the effects of the target unknown motion parameters; this condition leads to the performance degradation, especially for fast-maneuvering targets. In accordance with Equation (9), the Doppler bandwidth is expressed as follows:

$$f\_D = \left| \frac{4}{\lambda} \left[ \frac{\left(\upsilon - \upsilon\_a\right)^2 - R\_0 a\_c}{R\_0} \right] T\_a \right|. \tag{20}$$

A new deramp function with constant factor ϕ is created to compress the Doppler bandwidth as follows:

$$H\_{drmump}(f, t\_l, q) = \exp\left[\frac{j\pi}{c}(f + f\_c)\frac{qPRF\lambda}{T\_a}t\_n^2\right].\tag{21}$$

By multiplying Equation (9) by Equation (21), we have

$$\begin{split} s\_{\mathfrak{F}}(f, t\_{\mathfrak{n}}) &= s\_{\mathfrak{2}}(f, t\_{\mathfrak{n}}) \cdot \mathcal{H}\_{\mathrm{dorm}}(f, t\_{\mathfrak{n}}, q) \\ &= \mathrm{rect}\Big(\frac{f}{\mathfrak{F}}\Big) \exp\Big[-\frac{j8\pi}{c}(f + f\_{\mathfrak{c}})\mathcal{R}\_{\mathbb{0}}\Big] \exp\Big[-\frac{j8\pi}{c}\Big(f + f\_{\mathfrak{c}}\Big)\Big(\frac{(p - v\_{\mathfrak{a}})^2 - \mathcal{R}\_{\mathbb{0}}a\_{\mathfrak{c}}}{2\mathcal{R}\_{\mathbb{0}}} - \frac{q\mathcal{R}\mathcal{R}\mathcal{F}\mathcal{A}}{8\Gamma\_{\mathfrak{a}}}\Big)\Big], \end{split} \tag{22}$$

In accordance with Equation (22), the residual Doppler bandwidth is written as

$$f\_D f\_{D-pre} = \left| \frac{4}{\lambda} \left[ \frac{\left(\upsilon - \upsilon\_d\right)^2 - R\_0 a\_c}{R\_0} - \frac{\varrho PRF\lambda}{4T\_a} \right] T\_a \right| = \left| f\mathbf{D} - \varrho PRF \right|. \tag{23}$$

As shown in Equation (23), the main part of Doppler bandwidth is removed. The Doppler bandwidth *fD* is extremely compressed. An appropriate constant factor ϕ is chosen to make the remaining Doppler bandwidth *fD*−*pre* less than one PRF. The value of constant factor ϕ depends on the Doppler bandwidth of signal in Equation (9). If ω*PRF* < *fD* < (ω + 1)*PRF*, ω = 1, 2, 3, ··· , then the constant factor ϕ will be chosen as ω. After the pre-processing step in Equation (23), the SOKT (i.e., (*f* + *fc*)*t* 2 *<sup>n</sup>* = *fc*ξ<sup>2</sup> [22,23]) can be applied to remove the SRCM. As a result, we have

$$s\_{\mathbb{S}}(f,\xi) = \text{rect}\left(\frac{f}{B}\right) \exp\left[-\frac{j8\pi}{c}(f+f\_c)R\_0\right] \exp\left[-\frac{j8\pi}{\lambda}\left(\frac{\left(\upsilon-\upsilon\_a\right)^2 - R\_0a\_c}{2R\_0} - \frac{q\alpha PRFA}{8T\_a}\right)\xi^2\right].\tag{24}$$

As exhibited in Equation (24), the coupling influence between the *t* 2 *<sup>n</sup>*-term and range frequency variable *f* is eliminated. The SRCM is accurately removed after the SOKT is used, which is beneficial to the subsequent refocusing processing of a moving target. The target trajectory split is avoided because the Doppler bandwidth is effectively compressed. The Doppler spectrum ambiguity is further eliminated by performing a new deramp function. In the meantime, the constant factor ϕ is introduced to change the Doppler bandwidth of the deramp function. Thus, the new deramp function has a wide applicability.

Some noise-free simulation results are presented to validate the above-mentioned processing step. The range bandwidth of the simulated radar is set to 300 MHz. Other main radar parameters are the

same as those of Example A. In Example C, a target, which is denoted by Target C, is considered. The parameters of the target are as follows: *va* = <sup>−</sup>6.6 m/s, *aa* = <sup>1</sup> m/s2, *vc* = 24 m/s and *ac* = <sup>−</sup>2.5 m/s2.

Figure 6 shows the result of Example C. An evident SRCM appears in Figure 6a. Figure 6b displays the result of SRCM compensation by using the correction function in Equation (10). The residual SRCM correction error cannot be ignored because it is larger than one range bin. Therefore, the mismatch of the SRCM correction function in Equation (10) appears. The Doppler spectrum distribution of Target C after FRCM and TRCM correction is presented in Figure 6c. The Doppler spectrum bandwidth is larger than one PRF and the Doppler spectrum occupies several PRF bands because the chirp rate has been doubled after previous processing steps. The trajectory splits into several parts after directly performing the SOKT, as shown in Figure 6d. Figure 6e illustrates the result of applying the new deramp function with constant factor ϕ = 1. The target Doppler spectrum bandwidth is extremely compressed, and the target Doppler spectrum is smaller than one PRF band. Figure 6f exhibits the result of performing the deramp and SOKT operations. Given that the target spectrum is located on one PRF band, the trajectory is focused in the same range bin after applying SOKT. Thus, the simulation results verify the new deramp function and SOKT operations.

**Figure 6.** Results of Example C: (**a**) trajectory of Target C after FRCM and TRCM compensation; (**b**) result using the correction function in (10); (**c**) trajectory of Target C in the range and Doppler domain; (**d**) result by directly using the second-order keystone transform (SOKT); (**e**) result after applying the new deramp function in Equation (21); (**f**) result after performing new deramp function and SOKT.

#### *4.2. Analysis of Multiple Target Focusing*

In the previous section, the case of one moving target is considered. However, multiple targets may be present in the observation scene. In the case of multiple targets, the effect of cross term induced by nonlinear operation in Equation (9) should be discussed. We assume that the number of targets is *D*. Then, the range compressed signal in Equation (5) is expressed as

$$s\_{mul,1}(f, t\_n) = \sum\_{i=1}^{D} \text{rect}\left(\frac{f}{B}\right) \exp\left[-j\frac{4\pi}{c}(f+f\_c)\left(\mathbb{R}\_{i,0} + \beta\_{i,1}t\_n + \beta\_{i,2}t\_n^2 + \beta\_{i,3}t\_n^3\right)\right] \tag{25}$$

where *Ri*,0 is the nearest slant range of the *<sup>i</sup>*th target. <sup>β</sup>*i*,1 <sup>=</sup> <sup>−</sup>*vi*,*c*, <sup>β</sup>*i*,2 <sup>=</sup> (*v* − *vi*,*a*) <sup>2</sup> <sup>−</sup> *Ri*,0*ai*,*<sup>c</sup>* /(2*Ri*,0) and β*i*,3 = *vi*,*c*(*v* − *vi*,*a*) <sup>2</sup> <sup>+</sup> *Ri*,0*ai*,*a*(*vi*,*<sup>a</sup>* <sup>−</sup> *<sup>v</sup>*) / 2*R*<sup>2</sup> *i*,0 represent the first-, second- and third-order phases of the *i*th target, respectively.

The target signal after RCM correction (RCMC) is written as follows:

$$s\_{\text{small},2}(f, t\_{\text{fl}}) = \sum\_{i=1}^{D} \text{rect}\left(\frac{f}{B}\right) \exp\left(-j\frac{8\pi}{\lambda}\beta\_{i,2}t\_{\text{n}}^{2}\right) \exp\left[-j\frac{8\pi}{c}(f+f\_{\text{c}})R\_{i,0}\right] + s\_{\text{cross},2}(f, t\_{\text{n}}),\tag{26}$$

where

$$\begin{array}{lcl} s\_{\text{cross},2}(f, t\_n) &= \sum\_{i=1}^{D} \sum\_{j=1, i \neq j}^{D} \text{rect}\left(\frac{f}{\theta}\right) \exp\left[-j\frac{4\pi}{c}(f + f\_c)\left(R\_{i,0} + R\_{j,0}\right)\right] \\ & \cdot \exp\left[-j\frac{4\pi}{c}(f + f\_c)\left(\beta\_{i,1} - \beta\_{j,1}\right)t\_n\right] \exp\left[-j\frac{4\pi f\_c}{c}\left(\beta\_{i,2} + \beta\_{j,2}\right)t\_n^2\right] \\ & \cdot \exp\left[-j\frac{4\pi}{c}(f + f\_c)\left(\beta\_{i,3} - \beta\_{j,3}\right)t\_n^3\right] \end{array} \tag{27}$$

and *scross*,2(*f*, *tn*) denotes the cross terms of the signal in Equation (26). Given that the proposed method contains the nonlinear operation, the signal in Equation (26) includes auto and cross terms. According to Equations (26) and (27), the RCMs of auto terms are effectively eliminated, and only the second-order phase remains. Considering that the 1D SCFT is a linear transform, the second-order phases of the auto terms can be effectively estimated by the 1D SCFT. The matched filtering operation in Equation (18) is also a linear processing step. Thus, the auto terms can be refocused using the corresponding matched filters. With regard to the cross terms, β*i*,1 β*j*,1 and β*i*,3 β*j*,3 generally hold because *vi*,*c*, *ai*,*c*, *vi*,*<sup>a</sup>* and *ai*,*<sup>a</sup>* usually differ. Not only the second-order phases, but also the first- and third-order phases of the cross terms remain. Therefore, the RCMs of the cross terms still exist, and the energy of cross terms spreads along the range dimension. The cross terms are typically defocused in the 1D SCFT domain. In summary, the cross terms may not affect the refocusing of auto terms in this case.

However, moving targets may have the same first- and third-order phases, namely, β*i*,1 = β*j*,1 and β*i*,3 = β*j*,3, respectively, in a particular case. Therefore, the signal in Equation (26) is simplified and transformed into range and azimuth slow time domain as follows:

$$\begin{array}{rcl} s\_{\text{mul},2}(t, t\_n) & = \sum\_{i=1}^{D} \text{sinc} \Big[ B \Big( t - \frac{4\mathcal{R}\_{i,0}}{c} \Big) \Big] \exp \Big( -j \frac{8\pi}{\lambda} \beta\_{i,2} t\_n^2 \Big) \\ & + \sum\_{i=1}^{D} \sum\_{j=1, j \neq j}^{D} \text{sinc} \Big( B \Big[ t - \frac{2\left( \mathcal{R}\_{i,0} + \mathcal{R}\_{j,0} \right)}{c} \Big] \Big) \exp \Big[ -j \frac{4\pi}{\lambda} \Big( \beta\_{i,2} + \beta\_{j,2} \big) t\_n^2 \Big]. \end{array} \tag{28}$$

In accordance with Equation (28), the RCMs of auto and cross terms are all removed. The auto terms are focused at the positions of *t* = 4*Ri*,0/*c*, and the cross terms are focused at the positions of *t* = 2 *Ri*,0 + *Rj*,0 /*c*. Then, the 1D SCFT is performed to the auto and cross terms. Accordingly, we have

$$(s\_{\text{nuto,mol,2}}(\frac{4R\_{i,0}}{\varepsilon}, f\_{t\_n^2}) = \sum\_{i=1}^{D} \delta(f\_{t\_n^2} + \frac{4\beta\_{i,2}}{\varepsilon\lambda}),\tag{29}$$

$$\mathbb{E}\_{\text{cross},\text{null},2} \left[ \frac{2\left(R\_{i,0} + R\_{j,0}\right)}{c}, f\_{l\_n^2} \right] = \sum\_{i=1}^{D} \sum\_{j=1, j \neq j}^{D} \delta \left[ f\_{l\_n^2} + \frac{2\left(\beta\_{i,2} + \beta\_{j,2}\right)}{\varepsilon \lambda} \right]. \tag{30}$$

According to Equations (29) and (30), the peak positions of the auto and cross terms in the 1D SCFT domain are *f t* 2 *<sup>n</sup>* = <sup>−</sup>4β*i*,2/ελ and *<sup>f</sup> t* 2 *<sup>n</sup>* = <sup>−</sup><sup>2</sup> β*i*,2 + β*j*,2 /ελ, respectively.

The peaks of the cross terms may not affect the determination of the auto term peaks according to the previous analysis, but they may lead to spurious peaks. As shown in Equations (29) and (30), the peak positions of the auto terms, namely, *t* = 4*Ri*,0/*c*, *f t* 2 *<sup>n</sup>* = <sup>−</sup>4β*i*,2/ελ and *t* = 4*Rj*,0/*c*, *f t* 2 *<sup>n</sup>* = <sup>−</sup>4β*j*,2/ελ, are symmetric with respect to the peak positions of the cross terms, that is, *t* = 2 *Ri*,0 + *Rj*,0 /*c*, *f t* 2 *<sup>n</sup>* = <sup>−</sup><sup>2</sup> β*i*,2 + β*j*,2 /ελ, in the range time and 1D SCFT domain.

In summary, Figure 7 shows the symmetric characteristics between the auto and the cross terms. As shown in Figure 7a, trajectories of auto terms A and B are represented by a straight purple line, and that of the corresponding cross term C is denoted by a straight red line. After RCMC, auto terms A and B and cross term C are focused in the range time domain. The positions of the auto terms in the range time domain, namely, *t* = 4*Ri*,0/*c* and *t* = 4*Rj*,0/*c*, are symmetric with respect to the position of corresponding cross term *t* = 2 *Ri*,0 + *Rj*,0 /*c*. Figure 7b–d depict the 1D SCFT result of auto term A, cross term C and auto term B, respectively. The positions of the auto terms in the 1D SCFT frequency domain, namely, *f t* 2 *<sup>n</sup>* = <sup>−</sup>4β*i*,2/ελ and *<sup>f</sup> t* 2 *<sup>n</sup>* = <sup>−</sup>4β*j*,2/ελ, are symmetric with respect to the position of the cross term, that is, *f t* 2 *<sup>n</sup>* = <sup>−</sup><sup>2</sup> β*i*,2 + β*j*,2 /ελ. Thus, the symmetric properties of the auto and cross terms in the range time and 1D SCFT domain can be used to preliminarily identify the potential spurious peak induced by the cross term.

**Figure 7.** Schematic of the symmetric characteristics between the auto and cross terms: (**a**) result of the auto and cross terms after range cell migration correction (RCMC); (**b**) 1D scaled Fourier transform (SCFT) result of auto term A; (**c**) 1D SCFT result of cross term C; (**d**) 1D SCFT result of auto term B.

A spurious peak recognition procedure based on the 1D SCFT is presented to confirm the potential spurious peak caused by the cross term. Firstly, a recognition function is proposed in the range frequency and azimuth slow time domain as follows:

$$\begin{array}{lcl} s\_{\text{tr}-\text{null}}(f,t\_{\text{n}}) &= s\_{1-\text{null}}(f,t\_{\text{n}}) \cdot s\_{1-\text{null}}^{\*}(f,t\_{\text{n}}) \\ &= \sum\_{i=1}^{D} \text{rect}\left(\frac{f}{\mathcal{B}}\right) + \sum\_{\substack{i,j=1, i \neq j}}^{D} \text{rect}\left(\frac{f}{\mathcal{B}}\right) \exp\left[-j\frac{4\pi}{c}(f+f\_{c})\{\mathcal{R}\_{i,0} - \mathcal{R}\_{j,0}\}\right] \\ & \cdot \exp\left[-j\frac{4\pi}{c}(f+f\_{c})\{\beta\_{i,1} - \beta\_{j,1}\}t\_{\text{n}}\right] \exp\left[-j\frac{4\pi}{c}(f+f\_{c})\{\beta\_{i,2} - \beta\_{j,2}\}t\_{\text{n}}^{2}\right] \\ & \cdot \exp\left[-j\frac{4\pi}{c}(f+f\_{c})\{\beta\_{i,3} - \beta\_{j,3}\}t\_{\text{n}}^{3}\right]. \end{array} \tag{31}$$

After the SOKT is applied to the signal in Equation (31), we have

$$s\_{re-mul}(f, \xi) = \sum\_{i=1}^{D} \text{rect}\left(\frac{f}{B}\right) + s\_{re-cross}(f, \xi),\tag{32}$$

where

$$\begin{array}{lcl} \mathbf{s}\_{n\leftarrow\text{cross}}(f,\xi) &=& \sum\_{i,j=1,i\neq j}^{D} \text{rect}\left(\frac{f}{\mathfrak{k}}\right) \exp\left[-j\frac{4\pi}{c}(f+f\_{\varepsilon})\left(\mathbf{R}\_{i,0}-\mathbf{R}\_{j,0}\right)\right] \\ & \cdot \exp\left[-j\frac{4\pi}{c}\left(\boldsymbol{\beta}\_{i,1}-\boldsymbol{\beta}\_{j,1}\right)(f\_{\varepsilon})^{\frac{1}{2}}(f+f\_{\varepsilon})^{\frac{1}{2}}\xi\right] \exp\left[-j\frac{4\pi f\_{\varepsilon}}{c}(\boldsymbol{\beta}\_{i,2}-\boldsymbol{\beta}\_{j,2})\xi^{2}\right] \\ & \cdot \exp\left[-j\frac{4\pi}{c}\left(\boldsymbol{\beta}\_{i,3}-\boldsymbol{\beta}\_{j,3}\right)(f+f\_{\varepsilon})^{-\frac{1}{2}}(f\_{\varepsilon})^{\frac{2}{2}}\xi^{3}\right]. \end{array} \tag{33}$$

The signal in Equation (32) contains *D* auto terms and *D*(*D* − 1) cross terms. After the range IFFT is performed, we have

$$s\_{n \leftarrow mul}(t, \xi) = \sum\_{t=1}^{D} \text{sinc}(Bt) + s\_{n \leftarrow \text{cross}}(t, \xi), \tag{34}$$

where *sre*−*cross*(*t*, ξ) represents the result of transforming the cross terms in Equation (33) into the range and azimuth slow time domain. If the first- and third-order phases do not satisfy β*i*,1 = β*j*,1 and β*i*,3 = β*j*,3, as shown in Equation (33), then not only the second-order phases, but also the first and third-order phases of cross terms remain for the recognition function. Therefore, the RCMs and DFMs of the cross terms exist, and the energy of cross terms is still defocusing. The first-, second- and third-order phases of the auto terms are removed. The auto terms are focused at the position of *t* = 0 in the range time domain.

As for the special case, namely, β*i*,1 = β*j*,1 and β*i*,3 = β*j*,3, Equation (32) is rewritten as follows after the range IFFT is applied:

$$s\_{\rm re-mul}(t,\xi) = \sum\_{i=1}^{D} \text{sinc}(Bt) + \sum\_{i,j=1, i\neq j}^{D} \text{sinc}\left[B\left(t - \frac{2\{R\_{i,0} - R\_{j,0}\}}{c}\right)\right] \exp\left[-j\frac{4\pi f\_c}{c}(\beta\_{i,2} - \beta\_{j,2})\xi^2\right].\tag{35}$$

According to Equation (35), the auto terms of the recognition function in Equation (35) are focused at *t* = 0 and the cross terms of the recognition function in Equation (35) are focused at *t* = 2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c* and *t* = −2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c* in the range domain. Then, the 1D SCFT is applied to the cross terms of the recognition function in Equation (35). As a result, we have

$$s\_{\text{cross,rc}-\text{null}} \left[ \frac{2\left(R\_{i,0} - R\_{j,0}\right)}{c}, f\_{\xi^2} \right] = \sum\_{i=1}^{D} \sum\_{j=1, j\neq j}^{D} \delta \left[ f\_{\xi^2} + \frac{2\left(\beta\_{i,2} - \beta\_{j,2}\right)}{\varepsilon\lambda} \right]. \tag{36}$$

As shown in Equation (36), the peak positions of cross terms are *f* <sup>ξ</sup><sup>2</sup> = −2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ) and *f* <sup>ξ</sup><sup>2</sup> = 2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ), respectively, in the 1D SCFT domain.

Therefore, the peak positions of the cross terms of the recognition function in Equation (35), namely, *t* = 2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*, *f* <sup>ξ</sup><sup>2</sup> = −2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ) and *t* = −2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*, *f* <sup>ξ</sup><sup>2</sup> = 2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ), are symmetric with respect to the position of *t* = 0, *f* <sup>ξ</sup><sup>2</sup> = 0.

In summary, Figure 8 shows that, if the peak is the spurious peak, then the recognition function will have evident symmetrical peaks with the same amplitudes at positions of *t* = 2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*, *f* <sup>ξ</sup><sup>2</sup> = −2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ) and *t* = −2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*, *f* <sup>ξ</sup><sup>2</sup> = 2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ), which are symmetric with respect to the position of *t* = 0, *f* <sup>ξ</sup><sup>2</sup> = 0. Otherwise, Figure 9 depicts that, if the peak is the target peak, then the evident symmetrical peaks at positions of *t* = 2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*, *f* <sup>ξ</sup><sup>2</sup> = −2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ) and *t* = −2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*, *f* <sup>ξ</sup><sup>2</sup> = 2 <sup>β</sup>*i*,2 <sup>−</sup> <sup>β</sup>*j*,2 /(ελ) are absent in the output result of the recognition function. *Remote Sens.* **2019**, *11*, 2214

**Figure 8.** Result of spurious peak recognition function for spurious peak: (**a**) RCMC result; (**b**) 1D SCFT result for signal in *t* = 2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*; (**c**) 1D SCFT result for signal in *t* = −2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*.

**Figure 9.** Result of spurious peak recognition function for target peak: (**a**) RCMC result; (**b**) 1D SCFT result for signal in *t* = 2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*; (**c**) 1D SCFT result for signal in *t* = −2 *Ri*,0 <sup>−</sup> *Rj*,0 /*c*.

The identification procedure of the potential spurious peak is summarised as follows:

Step 1) Extract all peak positions, which are denoted by *tz*, *f z*,*t* 2 *n* , *z* = 1 ··· *Z*, where *Z* is the number of peaks.

Step 2) Determine whether *tw* = *tu*+*tx* <sup>2</sup> and *f w*,*t* 2 *<sup>n</sup>* <sup>=</sup> *<sup>f</sup> u*,*t* 2 *n* + *f x*,*t* 2 *n* <sup>2</sup> , *u* = 1 ··· *Z*, *x* = 1 ··· *Z*, *w* = 1 ··· *Z*, *w u x* is satisfied; if so, then the peak at position of *tw*, *f w*,*t* 2 *n* may be the spurious peak of the cross term. Apply the recognition function (i.e., Steps 3–5) to identify the potential spurious peak; if not, then all peaks are the target peaks.

Step 3) Calculate Equation (31), and obtain *sre*−*mul*(*f*, *tn*).

Step 4) Perform the SOKT and range IFFT to *sre*−*mul*(*f*, *tn*), and obtain *sre*−*mul*(*t*, ξ).

Step 5) Apply azimuth 1D SCTF to *sre*−*mul*(*tw* − *tu*, ξ) and *sre*−*mul*(*tw* − *tx*, ξ). If evident symmetrical peaks with the same amplitudes described in Figure 8 are present, then the peak at *tw*, *f w*,*t* 2 *n* mentioned 

in Step 2 will be a spurious peak. Otherwise, the peak at *tw*, *f w*,*t* 2 *n* is the target peak.

The flowchart of potential spurious peak recognition procedure is given in Figure 10.

**Figure 10.** Flowchart of the potential spurious peak identification procedure.

Two examples are provided to validate the analysis of multiple target focusing and spurious peak recognition procedure. The radar parameters are the same as those of Example A. The signal-to-noise (SNR) is 7 dB. In Example D, we consider three targets with different motion parameters, and that are represented by Targets D, E and F. The phase parameters of these target signals are

 $\beta\_{D,1} = -19.8 \text{ m/s}$   $\beta\_{D,2} = 1.2 \text{ m/s}^2$ , and  $\beta\_{D,3} = 0.5 \text{ m/s}^3$ ,  $\beta\_{E,1} = 15.6 \text{ m/s}$ ,  $\beta\_{E,2} = 2.4 \text{ m/s}^2$ , and  $\beta\_{E,3} = -0.6 \text{ m/s}^3$ ,  $\beta\_{F,1} = 30.5 \text{ m/s}$ ,  $\beta\_{F,2} = 3.6 \text{ m/s}^2$ , and  $\beta\_{F,3} = 1.2 \text{ m/s}^3$ .

Figure 11 displays the processing results of Example D. Figure 11a depicts the RCMC result. The background noise is removed from the obtained result to show the target trajectories clearly. Three straight trajectories related to Targets D, E and F appear because the RCM is effectively corrected. However, the cross terms are still defocused due to the serious effect of the RCM, thereby helping in suppressing the interference of the cross term. Figure 11b–d exhibit the 1D SCFT results for Targets D, E and F in the 147th, 177th and 207th range sample bins, respectively. Evident peaks with respect to Targets D, E and F appear in the corresponding figure. Figure 11e,f illustrate the 1D SCFT results of data in the 157th and 187th range sample bins, respectively, to indicate the defocusing of cross terms without loss of generality. The evident peaks are absent in the processing results. The cross terms still suffer from the effect of defocusing induced by the residual first- and third-order phase errors. The positions of peaks in Figure 11b–d satisfy the symmetric properties described in Figure 7. The potential spurious peak recognition procedure is utilised to identify whether the peak in Figure 11c is a spurious peak. Figure 12a,b depict the corresponding recognition result. Considering that clear peaks, which satisfy the symmetric feature shown in Figure 8, do not emerge in the 30th and −30th range time sample bins of the recognition function, the peak in Figure 11c is determined as the target peak.

**Figure 11.** Results of Example D: (**a**) result after RCMC; (**b**) 1D SCFT result for Target D in the 147th range sample bin of Figure 11a; (**c**) 1D SCFT result for Target E in the 177th range sample bin of Figure 11a; (**d**) 1D SCFT result for Target F in the 207th range sample bin of Figure 11a; (**e**) 1D SCFT result of data in the 157th range sample bin of Figure 11a; (**f**) 1D SCFT result of data in the 187th range sample bin of Figure 11a.

**Figure 12.** Potential spurious peak recognition results of Example D: (**a**) 1D SCFT result of the −30th range time sample bin of the recognition function; (**b**) 1D SCFT result of the 30th range time sample bin of the recognition function.

In Example E, we consider two targets that are denoted by Targets G and H. The phase parameters of these target signals are

 $\beta\_{\rm G,1} = 32.6 \text{ m/s}, \beta\_{\rm G,2} = 1.2 \text{ m/s}^2, \text{and}$  $\beta\_{\rm G,3} = 0.8 \text{ m/s}^3, \beta\_{\rm H,1} = 32.6 \text{ m/s}, \beta\_{\rm H,2} = 3.6 \text{ m/s}^2, \text{and}$  $\beta\_{\rm H,3} = 0.8 \text{ m/s}^3.$ 

Figure 13 shows the obtained results of Example E. The RCMC result without noise is exhibited in Figure 13a to display the target trajectories clearly. Three straight trajectories are observed in Figure 13a because the first- and third-order phases of Target G equal those of Target H. Figure 13b–d illustrate the 1D SCFT result of the data in the 147th, 177th, and 207th range sample bins, respectively. Considering that peak positions in Figure 13b–d satisfy the symmetric characteristics described in Figure 7, the peak in Figure 13c may be a spurious peak. Figure 14 depicts the potential spurious peak recognition results. Given that evident symmetric peaks with the same amplitudes appear in the 30th and −30th range time sample bins of the recognition function, as shown in Figure 14a,b, the peak in Figure 13c is confirmed as the spurious peak.

**Figure 13.** Results of Example E: (**a**) result after RCMC; (**b**) 1D SCFT result of Target G in the 147th range sample bin of Figure 13a; (**c**) 1D SCFT result of cross term in the 177th range sample bin of Figure 13a; (**d**) 1D-SCFT result of Target H in the 207th range sample bin of Figure 13a.

**Figure 14.** Potential spurious peak recognition results of Example E: (**a**) 1D SCFT result of the −30th range time sample bin of the recognition function; (**b**) 1D SCFT result of the 30th range time sample bin of the recognition function.

#### **5. Experimental Results**

In this section, several simulation experimental results in the presence of Gaussian background and real data processing results are analysed to verify the proposed method.

#### *5.1. Simulation Experimental Result Analysis*

The simulation radar parameters are listed in Table 1. Two maneuvering targets, which are denoted by T1 and T2, are considered. The motion parameters of the two targets are summarised in Table 2. The SNR after range compression is 8 dB. T1 is a Doppler spectrum ambiguity target, and its azimuth Doppler spectrum bandwidth is larger than PRF/2. The azimuth Doppler spectrum of T1 distributes into two PRF bands. The Doppler ambiguity number of T2 is −1, and its azimuth Doppler spectrum is still located on one PRF band. The constant factor is set to ε = 2 following the constant factor selection strategy in Appendix A. The FOKT-based [25], stationary phase-based [28] and SOKT-GHHAF methods [32] are used for comparison.


**Table 1.** Basic simulation radar parameters.



Figure 15a depicts the result of range compression. Two curved trajectories are observed in the figure. Target energy also distributes into several range sample bins, thereby leading to severe RCM. The result by directly applying azimuth FFT is illustrated in Figure 15b. Notably, the energy of targets also spreads along the azimuth Doppler dimension, which induces serious DFM. The RCM and DFM result in severe defocusing effects. In addition, the azimuth Doppler spectrum of T1 spans over two PRF bands. The azimuth Doppler spectrum of T2 occupies one PRF. Figure 15c shows the result after RCMC, and the background is removed to indicate the trajectories of two targets clearly. The RCM is effectively eliminated, and the trajectory split is avoided. The energy of the target is focused in the corresponding range bin. In the meantime, the RCMs of cross terms remain, and cross terms still suffer from the effect of defocusing, thereby helping in suppressing the cross terms. Figure 15d exhibits the 1D SCFT result of T1. An evident peak with respect to T1 appears in Figure 15d. From the peak position, a well-focused result is obtained in the range–azimuth time domain by using the corresponding matched filter, as shown in Figure 15e–h. For the same reason, T2 is also accumulated as a peak in Figure 15i. With the peak position, T2 is well focused in the range–azimuth time domain by the matched filtering, as exhibited in Figure 15j–m.

**Figure 15.** *Cont.*

(**m**)

**Figure 15.** Results of simulation experiment: (**a**) range compression result; (**b**) azimuth Doppler spectrum distributions of two targets; (**c**) RCMC result; (**d**) 1D SCFT result of T1; (**e**) focusing result of T1 by using the proposed method; (**f**) profile for focusing result of T1; (**g**) dB version of amplitude for Figure 15f; (**h**) stereogram of Figure 15e; (**i**) 1D SCFT result of T2; (**j**) focusing result of T2 by performing the proposed method; (**k**) profile for focusing result of T2; (**l**) dB version of amplitude for Figure 15k; (**m**) stereogram of Figure 15j.

Figure 16a,b plot the SRCM correction result using the SOKT-GHHAF method before and after Doppler centre shifting by PRF/2 operation. The background is removed from the result to illustrate the target trajectory clearly. As presented in Figure 16a, the trajectory of T1 splits into two parts due to the azimuth Doppler spectrum of T1 spanning over two PRF bands. Given that the Doppler spectrum bandwidth is larger than PRF/2, the operation of Doppler centre shifting by PRF/2 in SOKT-GHHAF method is invalid. Figure 16b exhibits that the trajectory still splits into two parts, and this condition affects the performance of RCMC and induces the coherent integration loss. Figure 16c–e show the focusing result of T1 by applying the FOKT-based method. A defocused result appears in Figure 16c–e because the along-track velocity and acceleration motions are ignored. As presented in Figure 16f–h, the focusing performance of the stationary phase-based method deteriorates significantly given that the third-order phase is neglected.

**Figure 16.** Processing results of compared methods for T1: (**a**) SRCM compensation result using second-order keystone transform-based generalised Hough-high-order ambiguity function (SOKT-GHHAF) method before Doppler centre shifting by PRF/2; (**b**) SRCM correction result using SOKT-GHHAF method after Doppler centre shifting by PRF/2; (**c**) focusing result by performing the first-order keystone transform-based method; (**d**) profile for focusing result of FOKT-based method; (**e**) dB version of amplitude for Figure 16d; (**f**) focusing result by applying the stationary phase-based method; (**g**) profile for focusing result of the stationary phase-based method; (**h**) dB version of amplitude for Figure 16f.

In summary, the result of simulation experiment validates that the proposed method can effectively compensate the RCM and DFM of the maneuvering target, and a well-focused result can be achieved regardless of Doppler ambiguity including Doppler centre blur and spectrum ambiguity. The performance of presented method is better than that of the FOKT-based and stationary phase-based methods. This result is attributed to the fact that the proposed method considers the along-track velocity and acceleration motions and can deal with the high-order RCM and DFM induced by the third-order

phase. The proposed method is also more robust to Doppler ambiguity than the SOKT-GHHAF method because it can effectively solve the problems of Doppler centre shift and Doppler spectrum broadening.
