**1. Introduction**

Historically, the application of inverse synthetic aperture radar (ISAR) imaging in target recognition is limited, owing to the fact that the technique only captures the projected two-dimensional (2-D) characteristics of the target. To overcome this drawback, interferometric ISAR (InISAR) that provides three-dimensional (3-D) information of the target has been developed [1–11]. InISAR systems are generally composed of several fixed channels, with one as both a transmitter and receiver and the rest as full-time receivers. The 3-D geometry of the target can be reconstructed based on the phase difference of different ISAR images. Until now, nearly all researches regarding InISAR imaging have been carried out in the microwave band, while limited research about InISAR imaging with terahertz (THz) radars is available in the current literature. Compared to InISAR imaging with microwave radars, THz radars can more easily achieve a higher carrier frequency and wider absolute bandwidth, which provide higher spatial resolution and more detailed information of the target. Furthermore, THz InISAR imaging holds large potential in maneuvering target surveillance and recognition in space and near space, and it is of great significance to study and advance this technique.

With the noticeable progress of THz sources and detectors over the past few decades, it has become possible to achieve imaging and recognition with a THz radar system. For now, we can summarize the imaging radar system in the THz band into three main categories: raster-scanning radar system, mixed-scanning radar system, and an SAR/ISAR system. The raster-scanning radar system focuses the beam to a fixed area using several lenses, and the 3-D imaging result is obtained by recording the data of each scanning area. Representatives of such systems include Jet Propulsion Laboratory [12–17] and Pacific Northwest National Laboratory [18,19]. The imaging frame rate is determined by the number of scanning pixels and oscillation frequency of the scanning mirrors, which makes it time-consuming to obtain a target image. Besides, the system structure is complex and can be easily damaged. Compared to the raster-scanning radar system, the mixed-scanning radar system substitutes one-dimensional (1-D) raster scanning with 1-D electrical scanning or 1-D mechanical scanning, and the imaging speed can be greatly improved. Representatives of mixed-scanning radar systems include Chinese Academy of Sciences [20–22] and China Academy of Engineering Physics [23,24]. The mixed-scanning radar system can achieve 3-D imaging near real time, but some defects still exist (e.g., the imaging field of view is limited, and the target needs to be stationary). The SAR/ISAR system acquires the target image through relative motion between radar and target. The resolution depends on the bandwidth of the sweep signal and the length or relative rotation angle of the synthetic aperture. Representatives include FGAN-FHR [25–27], the US Defense Advanced Research Projects Agency (DARPA) [28], and so on [29–33]. Compared to the former two THz radar systems, the SAR/ISAR system has a relatively simple system structure and low hardware cost, and it has no limitation of target distance. However, the SAR/ISAR system can only achieve 2-D imaging. The THz InSAR/InISAR system has the advantages of the SAR/ISAR system and also has the ability to achieve 3-D imaging. Thus, it is meaningful to establish an InSAR/InISAR system in the THz band.

As we know, the translational motion of a target will lead to image misregistration in InISAR imaging, and an image registration process is necessary to be carried out first. In recent years, some research regarding image registration in InISAR systems have been presented. These image registration methods can be summarized into three categories. The first is based on the correlation coefficient of ISAR images such as the time domain correlation method and the frequency domain searching method [4]. However, in order to guarantee registration precision, the searching step of motion parameters should be controlled within a very limited range in the THz band, which significantly increases the computational complexity. The second is based on the selection of a reference distance, such as the respective reference distance compensation method and the reference distance deviation compensation method [5,6]. In these methods, the reference distance is chosen as the total distance from the transmitter to reference center and the reference center to corresponding receiver. There is no need to estimate the motion parameters, and the method is suitable for both three-antenna and multiantenna configurations. However, in a real InISAR imaging scene, it is highly challenging to acquire the accurate distances from the reference center to the other receivers (except the one as both a transmitter and receiver). The third is based on estimation of the motion parameters of target. Reference [7] presents a procedure based on a multiple antenna-pair InISAR imaging system with nine antennas to estimate angular motion parameters. The angular motion parameters measurement is based on the spatial spectrums of the whole target, and the phase wrapping of motion measurement curves is avoided by using a pair of antennas with a short baseline. However, there are some shortcomings in this method. Firstly, the multiple antenna-pair configuration increases the system complexity and hardware cost, especially in high-frequency applications. Secondly, the spatial spectrums in a fixed range cell usually contain information of several scattering centers, which would introduce the "angle glint" phenomenon to the motion measurement curves. Thirdly, the method does not consider the influence of noise, but the motion parameter measurement precision of the weak scattering centers is sensitive to noise. In addition to the limitations of the above methods, most research just gives the simulation results and are short on experimental validation.

Unlike SAR images, high-frequency ISAR images of moving targets usually consist of several dominating reflectors such as corner reflectors formed by the tail, fuselage, and wings of the aircraft. Taking the computing efficiency, practicability, robustness, and precision into consideration, a strong scattering centers fusion (SSCF) technique is proposed in this paper to estimate translational motion parameters using these dominating reflectors, which overcomes the defects in the aforementioned image registration methods. With a one-input three-output THz InISAR system, translational motion parameters both in the azimuth and height directions are accurately estimated. Strong scattering centers (SSCs) are extracted in the image domain with a rectangular filter operation, which eliminates most of the noise. Motion measurement curves are derived from the spatial spectrums of mutually independent SSCs so that the "angle scintillation" phenomenon can be effectively suppressed. Then, the translational motion parameters are obtained by fitting the motion measurement curves with phase unwrapping and intensity-weighted fusion processing. Finally, image registration is achieved by compensating the estimated motion parameters to the radar echoes, and the InISAR imaging results are obtained with a simple interference operation.

This paper is organized as follows: in Section 2, the signal model is established. The SSCF technique is described in detail in Section 3. In Section 4, the simulations of the point targets model under different signal-to-noise ratio (SNR) and equivalent verification experiments with a multichannel THz radar system are carried out. A discussion is given in Section 5. Conclusions are presented in in Section 6.

#### **2. Signal Model of Interferometric Inverse Synthetic Aperture Radar (InISAR) Imaging**

The configuration of the InISAR system is demonstrated in Figure 1. Antenna *O* acts as both the transmitter and receiver, while antennas *A* and *B* operate in the receiving mode only. *L*<sup>1</sup> and *L*<sup>2</sup> denote the lengths of the baselines *OA* and *OB*, respectively. A target coordinate system *xoz* is built referencing the right-angle layout of the three antennas. *P*(*xp*, *yp*, *zp*) is an arbitrary scattering center located on the target whose projections on planes *xoy* and *yoz* are *P1* and *P2*, respectively. *y*<sup>0</sup> denotes the initial distance from antenna *O* to the reference center *o*, and *RAP*, *RBP*, and *ROP* denote the initial distances from *P* to three antennas, respectively. Suppose the transmitting linear frequency modulated (LFM) signal from antenna *O* is

$$\mathbf{s}(\hat{t}, t\_m) = \text{rect}(\frac{\hat{t}}{T\_p}) \exp[j2\pi(f\_c t + \gamma \hat{t}^2/2)],\tag{1}$$

where

$$\text{rect}\left(u\right) = \begin{cases} 1 & |u| \le 0.5\\ 0 & |u| > 0.5 \end{cases} \tag{2}$$

then, *Tp* is the pulse width, *fc* is the carrier frequency, γ is the chirp rate, *tm* is the slow time, ˆ*t* is the fast time, and *t* = *tm* + ˆ*t* is the full time.

**Figure 1.** Schematic of the three-receiver interferometric inverse synthetic aperture radar (InISAR) imaging configuration.

Ignoring the signal envelope, the received signal at receiver *i* (*i* = *O*, *A*, *B*) from *P* is

$$s\_{\hat{\mathbf{t}}}(\hat{\mathbf{t}},t\_{m}) = \sigma\_{P} \exp\{j2\pi[f\_{\hat{\mathbf{c}}}(t-\frac{R\_{\mathrm{OP}}(\mathbf{f},t\_{m})+R\_{\mathrm{iP}}(\mathbf{f},t\_{m})}{c}) + \frac{\mathcal{V}}{2}(\hat{\mathbf{t}} - \frac{R\_{\mathrm{OP}}(\mathbf{f},t\_{m}) + R\_{\mathrm{iP}}(\mathbf{f},t\_{m})}{c})^{2}\}),\tag{3}$$

where *c* is the wave propagation velocity, σ*<sup>P</sup>* is the reflection coefficient of *P*, and *RiP*(ˆ*t*, *tm*) represents the distance from receiver *i* to *P* at time *t*. These distances are:

$$R\_{\rm OP}(t) = \sqrt{\left(\mathbf{x}\_{\rm P} + \boldsymbol{\Delta R}\_{\rm x}(t)\right)^2 + \left(y\_{\rm P} + y\_0 + \boldsymbol{\Delta R}\_{\rm y}(t)\right)^2 + \left(z\_{\rm P} + \boldsymbol{\Delta R}\_{\rm z}(t)\right)^2},\tag{4}$$

$$R\_{AP}(t) \approx R\_{OP}(t) + \frac{L\_1^2 - 2L\_1(\mathbf{x}\_P + \Delta R\_x(t))}{2y\_0},\tag{5}$$

$$R\_{\rm BP}(t) \approx R\_{\rm OP}(t) + \frac{L\_2^2 - 2L\_2(z\_P + \Delta R\_z(t))}{2y\_0},\tag{6}$$

where Δ*Rx*(*t*), Δ*Ry*(*t*), and Δ*Rz*(*t*) represent the displacement of *P* from time instant 0 to time *t* in the *x*, *y*, and *z* directions, respectively. The approximations of *RAP*(*t*) and *RBP*(*t*) are based on the assumption that the imaging scene satisfies the far-field condition.

Based on the 'stop–go' model approximation and de-chirping processing, the received signal after range compression is

$$\begin{array}{rcl} s\_i(\boldsymbol{y}, t\_m) = & \sigma\_P \text{sinc} \mathbf{c} [\frac{2\gamma T\_p}{c} \left[ \mathbf{y} - \frac{(\mathcal{R}\_{\text{OP}}(t\_m) + \mathcal{R}\_{\text{iP}}(t\_m) - 2\mathcal{R}\_{\text{Ov}}(t\_m))}{2} \right] \\ & \times \exp\{-j\frac{2\pi}{\lambda} \left[ \mathcal{R}\_{\text{OP}}(t\_m) + \mathcal{R}\_{\text{iP}}(t\_m) - 2\mathcal{R}\_{\text{Ov}}(t\_m) \right] \end{array} \tag{7}$$

where *ROo*(*tm*) is the reference distance

$$R\_{\rm Oo}(t\_{\rm m}) = \sqrt{\Delta R\_{\rm x}^{\cdot 2}(t\_{\rm m}) + \left(y\_0 + \Delta R\_{\rm y}(t\_{\rm m})\right)^2 + \Delta R\_{\rm z}^{\cdot 2}(t\_{\rm m})}.\tag{8}$$

Within the short time of data acquisition, the translational velocities of the aircraft are assumed as constant (i.e., Δ*Rx*(*tm*) = *Vxtm*, Δ*Ry*(*tm*) = *Vytm*, Δ*Rz*(*tm*) = *Vztm*). After migration through range cell (MTRC) correction and azimuth compression, the ISAR images of three channels can be obtained as

$$\cos(y\_\prime f\_a) = A\_P \sin \mathbf{c} [\frac{2\gamma T\_p}{c}(y - y\_P)] \sin \mathbf{c} [T\_a [f\_a + \frac{2\mathbf{P} \mathbf{V}}{\lambda}] \exp(j\rho\_0)] \tag{9}$$

$$s\_A(y, f\_a) = A\_P \sin \mathbf{c} |\frac{2\gamma T\_p}{\mathbf{c}}(y - y\_P)| \sin \mathbf{c} |T\_a[f\_a + \frac{2\mathbf{P}\mathbf{V}}{\lambda y\_0} - \frac{L\_1 V\_x}{\lambda y\_0}]| \exp[j(\varphi\_0 + \varphi\_1)],\tag{10}$$

$$s\_B(y, f\_a) = A\_P \text{sinc} \{ \frac{2\gamma T\_p}{c} (y - y\_P) \} \text{sinc} \{ T\_a [f\_a + \frac{2\text{PV}}{\lambda y\_0} - \frac{L\_2 V\_z}{\lambda y\_0}] \} \exp[j(\varphi\_0 + \varphi\_2)],\tag{11}$$

where *AP* is the scattering intensity of *P* in the ISAR images, *Ta* is the total data acquisition time, *P* = (*xP*, *yP*, *zP*) is the coordinate vector, and *V* = (*Vx*, *Vy*, *Vz*) is the translational motion vector. ϕ<sup>0</sup> and ϕ*<sup>i</sup>* (*i* = 1, 2) denote the constant phase and interferometric phases, respectively, and they are expressed as

$$\varphi\_0 = -\frac{2\pi}{\lambda} (y\_P + \frac{x\_P^2 + y\_P^2 + z\_P^2}{y\_0}),\tag{12}$$

$$
\varphi\_1 = \frac{2\pi}{\lambda} \frac{2L\_1 \mathbf{x}\_P - L\_1^2}{2y\_0},
\tag{13}
$$

$$
\rho\_2 = \frac{2\pi}{\lambda} \frac{2L\_2 z\_P - L\_2^2}{2y\_0}.\tag{14}
$$

From Equations (9) to (11), it can be seen that the Doppler shifts of ISAR images are *L*1*VxTa*/λ*y*<sup>0</sup> and *L*2*VzTa*/λ*y*0, respectively. Obviously, the Doppler shifts will lead to pixels misregistration between ISAR images. Therefore, image registration must be accomplished before interferometric imaging. It is clear that image registration is essentially the compensation of the translational motion parameters *Vx* and *Vz*. Once these motion parameters are achieved, the Doppler shifts can be eliminated by compensating the radar echoes of antenna *A* and antenna *B* with them.

### **3. The Strong Scattering Centers Fusion (SSCF) Technique**

In this subsection, the proposed SSCF technique is described in detail. This method can settle problems such as the "angle glint", sensitivity to noise, and phase wrapping in motion parameter estimation.

The details of the proposed SSCF technique are described as follows.

Step 1) After data preprocessing (i.e., range alignment, autofocus [34,35], and MTRC correction), the spatial spectrums of radar echoes are obtained as *sO*(*m*, *k*), *sA*(*m*, *k*), and *sB*(*m*, *k*), 1 ≤ *m* ≤ *M* and 1 ≤ *k* ≤ *N*, where *M* and *N* denote the number of samples and the number of pulses, respectively. Then, the ISAR images can be obtained as *sO*(*m*, *n*), *sA*(*m*, *n*), and *sB*(*m*, *n*) after azimuth compression, where 1 ≤ *n* ≤ *N*.

Step 2) Based on a fixed threshold, strong scattering areas on the object can be extracted from the three ISAR images. Firstly, the strongest scattering center in the first strong scattering area of antenna *O* can be easily found. Subsequently, the corresponding scattering centers in ISAR images of antenna *A* and antenna *B* can be confirmed, since they are distributed in the same range cell. This scattering center is then extracted in the image domain with a rectangular filter, whose length is determined by the main lobe width (3 dB) of the extracted scattering center. The noise is filtered in this step. This process is iterated to search for the localized strongest scattering center in all the strong scattering areas until the strongest scattering center in the last strong scattering area has been extracted.

Step 3) The extracted data of all SSCs in the image domain are rearranged to form new image matrices *sO*(*m* , *n*), *sA*(*m* , *n*), and *sB*(*m* , *n*), where 1 ≤ *m* ≤ *M* , and *M* denotes the number of SSCs. By performing an inverse Fourier transformation on variable *n*, the new spatial spectrums are obtained as *sO*(*m* , *k*), *sA*(*m* , *k*), and *sB*(*m* , *k*), respectively. Here, each row of the spatial spectrums contains information of only one scattering center. The phase difference curves *s*∗ *O*(*m* , *k*)*sA*(*m* , *k*) and *s*∗ *O*(*m* , *k*)*sB*(*m* , *k*) can be calculated as

$$q\_{xm'}(k) = -\frac{2\pi}{\lambda} \frac{L\_1^2 - 2L\_1(x\_{m'} + V\_{xm'}k\Delta t)}{2y\_0},\tag{15}$$

$$q\_{zm'}(k) = -\frac{2\pi}{\lambda} \frac{L\_2^2 - 2L\_2(z\_{m'} + V\_{zm'}k\Delta t)}{2y\_0},$$

where Δ*t* is the pulse repetition interval, and *Vxm* and *Vzm* are the estimated translation velocity of the *m* th strong scattering center.

Step 4) In the THz band, the wavelength is very short. On the other hand, in order to guarantee the precision of altitude measurement, a relative long baseline is required [7]. Therefore, the values of ϕ*xm*(*k*) and ϕ*zm*(*k*) usually exceed the range [−π, π], and the phase unwrapping operation is necessary. Here, the one-dimensional path integral method is adapted to achieve phase unwrapping of the motion measurement curves, and the theory is

$$\Lambda q(k) = \begin{cases} q(k+1) - q(k) & \left| q(k+1) - q(k) \right| \le \pi \\ q(k+1) - q(k) + 2\pi & q(k+1) - q(k) < -\pi \\ q(k+1) - q(k) - 2\pi & q(k+1) - q(k) > \pi \\ q\_{\text{new}}(k+1) = q\_{\text{new}}(k) + \Delta q(k) \end{cases} \tag{17}$$

where ϕ(*k*) and ϕnew(*k*) are phases before and after phase unwrapping, respectively.

Step 5) The motion measurement curves along the *x* axis and *z* axis are calculated as

$$R\_{\rm xm'}(k) = q\_{\rm xm' \rm new}(k) \cdot \frac{\lambda y\_0}{2\pi L\_1} + \frac{L\_1}{2} = \mathbf{x}\_{\rm m'} + V\_{\rm xm'} k \Delta t\_\prime \tag{18}$$

$$R\_{zm'}(k) = \varphi\_{zm' \text{new}}(k) \cdot \frac{\lambda y\_0}{2\pi L\_2} + \frac{L\_2}{2} = z\_{m'} + V\_{zm'} k \Delta t. \tag{19}$$

The estimated velocities of each SSC can be obtained by fitting these time-dependent curves. Finally, the estimation values of translational motion parameters are acquired by intensity-weighted fusion of all SSCs:

$$V\_x = \frac{\sum\_{m'=1}^{M} V\_{xm'} A\_{m'}}{\sum\_{m'=1}^{M} A\_{m'}} \; \; V\_z = \frac{\sum\_{m'=1}^{M} V\_{zm'} A\_{m'}}{\sum\_{m'=1}^{M} A\_{m'}} ; \tag{20}$$

where *Am* is the mean scattering intensity of the *m* th SSC in the ISAR images.

By compensating the radar echoes of antenna *A* and antenna *B* with the estimated motion parameters, the coordinates of the target are finally obtained from *s*∗ *<sup>O</sup>*(*y*, *fa*) *sA*(*y*, *fa*) and *s*∗ *<sup>O</sup>*(*y*, *fa*)*sB*(*y*, *fa*),

$$\exp = \frac{\lambda y\_0}{2\pi L\_1} q \nu\_1 + \frac{L\_1}{2} \nu\_2 = \frac{\lambda y\_0}{2\pi L\_2} q \nu\_2 + \frac{L\_2}{2} . \tag{21}$$

The derivation of the value of *yP* is omitted here, as it can be directly obtained from the range scale in the ISAR images. To summarize the above analysis, a block scheme of the InISAR imaging procedure based on the SSCF technique is presented in Figure 2.

**Figure 2.** The flowchart of InISAR imaging based on the strong scattering centers fusion (SSCF) technique.

Regarding the InISAR imaging theory, the Doppler shift between ISAR images should be less than one-eighth of the Doppler cell after image registration to guarantee a relatively high InISAR imaging precision. Thus, the velocity estimation error Δ*V* should satisfy

$$|\Delta V| < \frac{y\_0 \lambda}{8LT\_a}.\tag{22}$$

Equation (22) is applied as the criteria to evaluate the performance of the SSCF technique in this paper.

#### **4. Results and Analysis**

#### *4.1. The Point Target Simulation Results*

In order to verify the effectiveness of the proposed SSCF technique, an InISAR imaging simulation of a moving airplane model is presented. The parameters in the simulation are shown in Table 1. The airplane model contained 64 scattering centers, and the size of the airplane model was 21, 24, and 7.5 m in length, width, and height, respectively. For a vivid visualization, the model was depicted in four different views from four visual angles in Figure 3, with (a), (b), (c), and (d) corresponding to the 3-D view and projections on the *xoy*, *xoz*, and *yoz* planes, respectively. Three SSCs were assigned at the tail, wing, and fuselage, respectively, as highlighted with red circles in Figure 3. The ratio of scattering intensity between these SSCs and the others was 3:1. As a result, several ordinary scattering centers were distributed in the same range cell with the SSCs. The target was assigned to move at a constant velocity along the *x* direction, which led to image misregistration between ISAR images of antenna *O* and antenna *A*.

**Table 1.** Parameter settings of the radar system and target.

**Figure 3.** Four different views of the airplane model. (**a**) 3-D view; (**b**) *xoy* view projection; (**c**) *xoz* view projection; and (**d**) *yoz* view projection.

After range alignment, autofocus, and MTRC correction processing, three ISAR images were obtained. Based on the system parameters in Table 1, we calculated the Doppler shift between ISAR images of antenna *A* and antenna *O* at 5.5 cells. The range profiles of ISAR images at *y* = 0 of antenna *A* and antenna *O* are shown in Figure 4. It was obvious that there was a Doppler shift of six cells on the peak positions. To compare performances between the conventional method in [7] and our method, we added a pair of antennas with 0.04 m length of baseline to the conventional method because it used short baselines to achieve phase unwrapping of motion measurement curves.

**Figure 4.** ISAR images at *y* = 0 of antenna *O* and antenna *A* before image registration.

Figure 5 shows the spatial spectrums of radar echoes, and Figure 6 shows the phase difference curves of the spatial spectrums, both with (a) corresponding to the conventional method and (b) corresponding to our method. The SNR was 0 dB in this simulation. It can be seen from Figure 5a that the multiple scattering centers in a fixed range cell introduced a serious "angle glint" phenomenon. The "angle glint" phenomenon was eliminated by extracting the SSCs, as shown in Figure 5b. The influence of the multiscattering centers in motion parameter estimations is visualized in Figure 6a. Firstly, the "angle glint" phenomenon and noise introduced serious nonlinearity to the phase difference curves, which destroyed the real phase difference relationship and deteriorated the estimation precision of motion parameters. Secondly, without the short baseline, different scattering centers had different phase wrapping positions, and the one-dimensional path integral method was not applicable in this condition. In contrast, the curves in Figure 6b were linear, and the phase wrapping positions of each curve were constant. Noise was effectively filtered with the filter operation in the image domain, which ensured an excellent linearity of the motion trajectory even under a low SNR.

**Figure 5.** The spatial spectrums of radar echoes. (**a**) Conventional method. (**b**) Our method.

**Figure 6.** Phase difference curves of spatial spectrums. (**a**) Conventional method. (**b**) Our method.

Based on the one-dimensional path integral method, the phase difference curves of the three SSCs after phase unwrapping were acquired, as shown in Figure 7. All curves were absolutely linear with no fluctuation. Subsequently, translational motion parameters were easily estimated based on polynomial curve fitting and intensity-weighted fusion operation. In this simulation, the estimated velocities of the conventional method and our method were 68.678 and 300.057 m/s, respectively, which meant the conventional method was invalid under this low-SNR simulation environment. The estimated velocity of our method was nearly the same as the real velocity, which verified the effectiveness of the proposed SSCF technique. By compensating the estimated velocity of our method to the echo signals of antenna *A*, the final imaging results are shown in Figure 8, with (a) corresponding to the range profiles at *y* = 0 of antenna *A* and antenna *O* and (b) corresponding to the final InISAR imaging results. It can be seen from Figure 8a that the Doppler shift was eliminated, and image registration was achieved. In Figure 8b, the red circles are the real positions of the target, and the blue dots are the InISAR imaging results. It was clear that the InISAR image and the target model were precisely overlapped.

**Figure 7.** Phase difference curves after phase unwrapping based on the SSCF technique.

**Figure 8.** Imaging results after image registration. (**a**) ISAR images at *y* = 0 of antenna *O* and antenna *A*. (**b**) InISAR imaging results.

At the end of this simulation, the antinoise abilities of the two methods were evaluated. With our simulation parameters, the maximum allowed mean absolute error (MAE) of the estimated velocity in this simulation was derived as 6.8 m/s on the basis of Equation (22). The MAE of the estimated velocity is defined as

$$\text{MAE} = \sum\_{i=1}^{N} |V\_i - V| / N\_\prime \tag{23}$$

where *Vi* is the *i*th estimated velocity, *V* is the real velocity, and *N* is the number of simulation times. Here, the MAE of the estimated velocity was provided based on 100 Monte Carlo simulations under different SNR environments, as shown in Figure 9. The MAE of the conventional method was larger than 6.8 m/s when the SNR was under 30 dB, which meant the antinoise performance of the conventional method was very poor. This method was only suitable for extremely high SNR environments. The MAE of our method was under 6.8 m/s when the SNR was larger than −30 dB. The low MEA in estimated velocity proved the excellent antinoise ability of the proposed SSCF technique.

**Figure 9.** The mean absolute error (MAE) under different signal-to-noise ratio.

In order to show the superiority of the proposed SSCF technique over the past methods, performance comparisons of the frequency domain searching (FDS) method in [4], respective reference distance compensation (RRDC) method in [5] and [6], angular motion parameters estimation (AMPE) method in [7], and proposed SSCF technique in this paper are listed in Table 2. The algorithm reconstruction times (excluding the echo generation process) were obtained on a desktop computer with Intel core i7-7820X 3.60GHz CPU and 32GB RAM using Matlab codes. The reconstruction time of the FDS method depended on the initial value and searching step of motion parameters. In this

experiment, the initial value and searching step were set as 200 and 3 m/s, respectively, and a typical single-parameter optimization algorithm was adapted. From Table 2, it was seen that the FDS method was too time-consuming to fulfill the requirement of real time imaging. The RRDC method had the fastest imaging speed, but the radar needed to have ranging function. Besides, it was highly challenging to acquire accurate distances from the reference center to the other receivers (except the one as both a transmitter and receiver). The AMPE method needed a complex system structure, and the antinoise ability was poor. The SSCF technique proposed in this paper was advantageous in its high computing efficiency, excellent antinoise performance, and simple system structure. Therefore, the SSCF technique was more suitable for image registration in InISAR imaging.


**Table 2.** Comparison of algorithm performances.
