Estimation of High-Frequency Vibration Parameters for Terahertz SAR Imaging Based on FrFT With Combination of QML and RANSAC

High-frequency vibration of motion platform leads to paired echo for synthetic aperture radar (SAR) imaging, especially in terahertz band due to its shorter wave length. Different from most existing parameters estimation methods only considering single component high-frequency vibration, in this paper a novel method considering multi-components vibration model is proposed based on fractional Fourier transform (FrFT) with combination of quasi-maximum likelihood (QML) and random sample consensus (RANSAC). Based on the model establishment of high-frequency vibration error in the echo, its instantaneous chirp rate (ICR) is firstly estimated by FrFT in sliding sub-aperture, followed which the vibration parameters are coarsely obtained through spectrum analysis and least square (LS) regression. To further refine the parameters estimates, QML is developed for compensating the deviation both caused by the frequency spectrum leakage and the error propagation effects by one-dimensional search over the vibration frequency. Meanwhile, RANSAC is adopted for avoiding the outlier of the ICR estimates in LS regression, especially at low signal-to-noise ratio (SNR). Thus, the refinement strategy based on the combination of QML and RANSAC is developed, whose utilization improves the estimation accuracy of vibration parameters. Finally, the paired echo in terahertz SAR (THz-SAR) imaging is effectively suppressed by the proposed method, and the high-quality THz-SAR imaging results are achieved. Both simulations of single component and multi-components high-frequency vibration are used to verify the validity of the proposed method. The simulation results show that the proposed method has higher estimation accuracy even at low SNR.


I. INTRODUCTION
Terahertz (THz) wave generally refers to the electromagnetic wave with frequency ranging from 0.1THz (100 GHz) to 10THz, which is the transition zone from electronics to The associate editor coordinating the review of this manuscript and approving it for publication was Zihuai Lin .
photonics [1], [2]. Due to its special characteristics different from microwave and infrared radiation, THz wave has been used in many applications [3]- [8]. Among the numerous applications of THz wave, THz radar imaging technique has been receiving more and more attention in the modern radar imaging fields. Compared to other radar imaging wavebands, THz wave is easier to achieve higher carrier frequency VOLUME 9, 2021 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ and larger absolute bandwidth, which provide a higher spatial resolution and more details of targets. Besides, THz wave has the ability of imaging at higher frame rates, which makes it be utilized in video synthetic aperture radar (ViSAR) imaging. These advantages make THz radar imaging technique have promising application prospects in maneuvering target surveillance and recognition in space and near space [9], [10], and ground moving target detection and tracking [11], [12]. Synthetic aperture radar (SAR) imaging systems acquire the target image through relative motion between the radar and the target. Thus, motion is the basis of SAR imaging and also the source of the problem. On the one hand, it means that SAR uses the relative motion of radar and target to generate Doppler width band signals and achieve high resolution image. On the other hand, the focusing quality of SAR image is affected by unstable motion of platform due to the influence of factors such as air flow and engine vibration. Moreover, the motion error should not exceed the sub-wavelength level to obtain a well-focused image. Therefore, the effect of high frequency vibration error negligible in microwave SAR due to small amplitude should be taken account into THz-SAR imaging. In the absence of high precise measurement for high-frequency vibration, one of the main problems for THz-SAR imaging is the estimation and compensation of high-frequency vibration.
Different from low-frequency motion error, highfrequency vibration error does not affect the track of THz-SAR due to its small vibration amplitude. However, it influences the phase of radar echo signal and produces paired echo in SAR imaging, which seriously affects the image quality [13], [14]. Currently, various paired echo suppression methods have been presented. Some methods are to directly estimate and compensate the phase errors induced by high-frequency vibration [15], [16]. Although the phase errors are compensated, the specific information about vibrations are still unknown. Thus, another kind of methods are proposed based on parameters estimation. Their main idea involves the estimation of high-frequency vibration parameters, such as vibration amplitude, frequency and initial phase, and the compensation of constructed phase errors with the above estimated parameters. In [17]- [19], the phases of vibration signal are directly utilized for parameter estimation. However, their performance is seriously affected by signalto-noise ratio (SNR) due to phase unwrapping step. And if the phase ambiguity is serious, the phase unwrapping process is even unable to accurately recover the phase without the prior information. In [20], an adaptive chirplet decomposition method for parameter estimation of high-frequency vibration is proposed. In [21], a novel THz-SAR high-frequency vibration estimation and compensation imaging algorithm based on local fractional Fourier transform (FrFT) is proposed. However, all these methods only consider single component high-frequency vibration. In [22], the vibration parameter estimation method based on discrete sinusoidal frequency modulation transformation is proposed, which can deal with multi-components vibration. However, it needs to search in multi-dimensional space, resulting in huge computational burden. Thus, the global optimization techniques such as the simulated annealing algorithm and the particle swarm optimization algorithm are used for reducing the computational load, but converging to a local minimum is inevitable when they are used in practice, which results in a large error in the parameter estimation results [23], [24].
In this paper, a novel high-frequency vibration parameters estimation method considering multi-components model is proposed. The instantaneous chirp rate (ICR) of vibration signal is estimated based on FrFT, which is robust to the noise influence. In the following, the discrete Fourier transform (DFT) is performed on the ICR estimates to obtain vibration frequency. And least square (LS) regression is used for the estimation of other parameters. Rough estimate of high-frequency vibration parameters is done. However, when the true frequency of high-frequency vibration is not equal to the integer multiples of frequency resolution, the estimated vibration frequency by DFT has deviation, which increases estimation error of other parameters due to the effects of error propagation. Thus, one-dimensional search over the vibration frequency is performed to estimate the optimal parameters, which is similar to the maximum likelihood (ML) and referred to as quasi-ML (QML). Furthermore, when performing LS regression, some of the ICR estimates are outliers causing inaccurate estimate, especially at low SNR. To avoid outliers in the LS regression, some alternative sampling tools should be considered when the positions of outliers are unknown. In this paper we achieve such goal with the random sample consensus (RANSAC), which is usually used for estimation of parameters in various signal processing fields [25], [26]. In our application, it is possible to find combination of samples without outlier for more accurate parameters estimate when performing multiple random selections of the ICR estimate samples. Thus, with the combination of QML and RANSAC, the refinement strategy can estimate the vibration parameters with high precision and strong robustness. Simulation results validate the effectiveness of the proposed method.
The remainder of this paper is organized as follows. The model of high-frequency vibration error in the echo is established in Section II. In Section III, the proposed high-frequency vibration parameters estimation method based on FrFT with combination of QML and RANSAC is described in detail. In Section IV, the simulation analysis is carried out to verify the effectiveness of the proposed method. Conclusions are drawn in Section V.

II. MODEL ESTABLISHMENT OF HIGH-FREQUENCY VIBRATION ERROR IN THE ECHO
The helicopter SAR imaging geometry model is shown in Fig. 1. Ideally, the radar platform flies horizontally along the X direction with a constant velocity of v and a reference height of H . For a point target P (vt 0 , y 0 , 0) in the imaging region, the ideal instantaneous slant range r (t) can be where t is the slow time, r 0 is the nearest slant range of point target, t 0 is the zero doppler time. However, in practical applications, affected by the rotation of rotor blades, it is inevitable to introduce periodically vibration error for the helicopter platform during the flight. Different from traditional motion error, the vibration error has high vibration frequency and small vibration amplitude. And the high-frequency vibration error of aircraft platform is close to a superposition of multiple normal sinusoidal signal components. The vibration error model is established as where I is the number of vibration components, e i is the ith vibration component, A i , f i and ϕ i are the vibration amplitude, frequency, and initial phase of the ith vibration component, respectively. Generally, A i is the centimeter-level or even millimeter-level, which is far smaller than the range cells and f i is usually within a few tens of Hz. The high-frequency vibration leads to variations of the instantaneous slant range from the radar to the scattering point targets in the imaging scene. The actual instantaneous slant range R (t) is expressed as the sum of r (t) and the projection d LOS (t) of the high-frequency vibration error e (t) on the line of sight of point targets. Thus, it is given as where c = e· r < 1 and is constant, and e and r represent the unit direction vectors of the high-frequency vibration error and of the radar to point target P, respectively. Supposed that the emitted signal by SAR antenna is linear frequency modulation pulse signal and is given as where τ is the fast time, w r (·) is the antenna pattern of range direction, f c is the center frequency, and K r is the chirp rate.
After mixing demodulation, the echo signal of point target received by SAR antenna is expressed as where c is the light speed, λ = c f c is the carrier wavelength, and w a (·) is the antenna pattern of azimuth direction.
Due to its small vibration amplitude, the high-frequency vibration error does not affect the track of THz-SAR. And the range cell migration caused by high-frequency vibration error can be ignored. Thus, after range compression and range cell migration correction, the echo signal is expressed as where B is the signal bandwidth. Then, the dechirp operation and the residual video phase compensation are performed [13], (6) is expressed as Finally, we can focus the azimuth signal in the Doppler domain through applying fast Fourier transform (FFT) to (7). However, due to the existence of high-frequency vibration, it introduces the additional phase modulation to the ideal echo, which exhibits a periodically sinusoidal characteristic and gives rise to paired echo in the imaging result. And the interval and intensity of paired echo depend on the vibration amplitude and frequency [13]. To eliminate the effects of paired echo, the periodically modulation phase must be compensated in the imaging procedure. Thus, the parameters of high-frequency vibration should be estimated first.

III. PROPOSED PARAMETERS ESTIMATION METHOD
In this section, we propose a novel method to estimate the parameters of high-frequency vibration signal. First, the FrFT is adopted for ICR estimation. The following step is to obtain the rough estimation by spectrum analysis and LS regression. Finally, to refine the rough estimation, the ML function is evaluated with the random selection of vibration frequency and ICR estimation instants, and the estimate results with the largest value of the ML function are adopted as final estimates.

A. FRFT-BASED ICR ESTIMATION
According to (7), the phase error signal caused by highfrequency vibration in each range bin is expressed as From (8), it is clear that the phase error signal is an example of nonstationary signals, which can be analyzed by means of using sliding small-time windows. In a short-time interval [t 0 , t 0 + t], d LOS (t) is approximated by second-order polynomials, namely is the initial vibration velocity and a (t 0 ) is the initial vibration acceleration. Therefore, in the short-time interval [t 0 , t 0 + t], (8) can be expressed as When t is much less than the vibration duration, the second-order approximation in (9) is fairly accurate. From (10), g (t) in a short-time window is approximately a chirp signal. According to the definition of FrFT, the FrFT of a chirp signal y (t) with the matched order p 0 is an impulse function at the fractional domain frequency u 0 . And the matched order p 0 and the fractional domain frequency u 0 are estimated by solving the optimal equation [27], i.e., where F p {y (t)} is the pth FrFT of the signal y (t).
Thus, the chirp rate k (t 0 ) is determined by the matched order p 0 , i.e.,k Performing FrFT on (8) in successive sliding short-time window, the ICRk (t) is estimated.

B. ROUGH ESTIMATE
In the above section, the ICRk (t) is estimated and contains information on all vibration parameters since the ICR of (8) is given by The LS regression is used for estimating the vibration parameters, and the details are shown as follows. First, using the sampling time t = 1 f PRF , where f PRF is the pulse repetition frequency, (10) is sampled as discrete. And the discrete ICR k (n) is denoted as where N is the number of samples.
The ICR in (14) is rewritten as Equation (13) shows that the frequency of ICR is identical to the vibration frequency. Both of them are obtained by applying DFT to the estimated ICR. That is to say, the vibration frequencyf i is obtained witĥ where DFT k (n) is the DFT of the signalk (n).
Then, with the estimatedf i , the other vibration parameters are estimated from the ICR estimates using the LS regression: And Q is the number of instants. Finally, the parameters Â i ,φ i are obtained aŝ On the one hand, vibration frequency is estimated mainly through direct observation of spectrum maximum position. Due to the existence of window truncation and fence effect, it makes the sinusoid signal frequency spectrum leakage and results in the problem of insufficient estimation precision. That is to say, the DFT frequency estimation accuracy depends on the sampling length [28], [29]. On the other hand, at low SNR environment some of ICR estimates are outlier causing inaccurate estimate in LS regression. Owing to these two factors, when compensating direct the phase error with the above estimated parameters, the SAR image quality will be influenced by residual phase error. Fortunately, we can adopt some technique for the refinement of the rough results. In this way, we are able to obtain results that are more robust to noise, while making the SAR image well-focused.

C. REFINEMENT STRATEGY
When the DFT is used for estimating the frequency of (14), the estimated valuef i has the following relationship with the true value f i , namely where f is the frequency resolution and f = f PRF N ; δ is the estimate error and |δ| ≤ 0.5. To summarize, the true value f i belongs to the set of Therefore, the improved ML approach is used for the determination of the optimal vibration frequency. For each vibration (18)- (20). Then, evaluate the ML function where g (n) is the discretization of the phase error signal g (t).
Here we use g (n) instead of k (n). It is due to the fact that we estimate the signal parameters to eliminate the phase of g (n) as much as possible.
The optimal value of vibration parameters maximizes the ML function In this way, instead of 3I times search over the parameter space in classical ML, we only need I times search over the set of vibration frequency. That is why the improved ML approach is referred to as the QML method. However, the error of ICR estimates has an effect on the QML. Large percentage of the ICR estimates samplesk (n) are outlier when the noise influence is emphatic. As noise level increases, percentage of outliers rapidly increases. However, we do not know where are the positions of outliers, so the outliers are inevitable when performing the LS regression on such data. Here, RANSAC is the good choice for avoiding outliers in the procedure of QML. Therefore, here we propose the QML approach with combination of RANSAC for the refinement of the rough estimates. In the refinement procedure, random selections of both vibration frequency and ICR estimates instants are employed in the LS regression. These random selections are performed multiple times and the final estimates are signal parameters that obtain the largest ML function value. The refinement algorithm is summarized as follows.
Step 1: Initialize the vibration frequency set F i and the iterations number M .
Step 2: Set the rough estimate values in Section B as initialization vibration parameters estimate and the current criterion function J is evaluated through (22).
Step 5. Calculate the ML function J (m) with the estimated parameters by (22). If J (m) > J , set J = J (m) and the corresponding parameters as the current estimates.
Step 6. If m < M , repeat Step3 to Step5. Otherwise, the current estimates are adopted as the signal parameters.
After perform the above refinement strategy, the effects of vibration frequency estimation precision and ICR estimates outliers are reduced as much as possible. The estimates with the largest ML function value are considered as the optimal results. According to the optimal estimates, the vibration error compensation function is constructed to compensate for the vibration phase error in (7). The paired echo in SAR images will be suppressed perfectly. In summary, the proposed parameters estimate method mainly involves three procedures: a) ICR estimation based on FrFT; b) Rough estimate with spectrum analysis and LS regression; c) Precise estimate with combination of QML and RANSAC. The flow chart of the proposed method in this paper is shown in Fig. 2. VOLUME 9, 2021

IV. SIMULATIONS RESULTS
In order to verify the validity of the proposed method based on FrFT with the combination of QML and RANSAC, numerical simulations of single component and multi-component high-frequency vibration parameters estimation with different SNR are given. The corresponding THz-SAR imaging results with the above vibration error are also presented. The simulation parameters of THz-SAR system are given in Table 1.

A. SIMULATION ANALYSIS OF SINGLE COMPONENT VIBRATION
Considering the single component high-frequency vibration signal with the form of (8) and supposing that f 1 = 8.3 Hz, A 1 = 2.5 × 10 −3 m and ϕ 1 = π 4. In real THz-SAR imaging applications, the dominant scatter point is first chosen by the intensity of range bins. Then, a coarsely focused SAR image is achieved by azimuth FFT, and the strongest response is circularly shifted to the SAR image center. After that, the azimuth windowing is implemented to cut out all the energy of dominant point, and the intercepted image segment is processed by the inverse FFT. And the data sequences used for vibration error estimation is obtained. However, the data sequences extracted from dominant scatter point are affected by other scatter points and noises in the same range bin. Here, the clutter and noise are both modeled as Gaussian white noise. Therefore, we assume that the signal is corrupted with Gaussian white noise of 20dB. In order to demonstrate the excellent performance of the proposed method, the estimate results based on the LS method and the QML method are also given here, respectively. First, the signal ICR estimated by FrFT and its spectrum are shown in Fig. 3. From the  comparison between the estimated ICR and the true ICR, it can be found that the ICR estimation errors are larger at the extreme value points. It is due to the fact that the approximation error is larger when the second order approximation is adopted. It is obvious from the frequency spectrum diagram that the ICR frequency is 8.0645 Hz. Compared to the true vibration frequency value, there is a considerable estimate error. With the estimated vibration frequency, the other vibration parameters are obtained, and all parameters estimates and ML function values are listed in Table 2. Fig. 4 shows the estimated vibration displacement constructed with the estimated parameters. It is clear from the Table 2 that the estimated values by the LS method deviate from the true values due to the error propagation effects from the vibration frequency towards vibration amplitude and initial phase, which makes the estimated vibration displacement deviate from the true vibration displacement, as shown in the blue line of Fig. 4. The ML value also shows that the LS method is poor. In the QML method and the proposed method, the error propagation effects are both suppressed, which makes the parameters estimate values be almost consistent with the theoretical values. Fig. 4 also does confirm this. Meanwhile, from Table 2, we find that the proposed method is slightly better than the QML method.
The THz-SAR system parameters in Table 1 and the single component high-frequency vibration error in the above simulation are used to generate the echo signal for multi-point targets. Firstly, the range-Doppler (RD) algorithm is used for focusing on imaging, as shown in Fig. 5(a). It can be seen that the imaging results have many paired echoes and the SAR image quality seriously degrades. Then, with the estimated vibration parameters based on the above methods, the vibration error compensation functions are constructed to compensate for the vibration phase error, respectively. After phase error compensation, the imaging results are shown in Fig. 5(b)-(d). Fig. 6 shows the phase estimate error between the true phase and the estimated phase. Compared to Fig. 5(a),   the image quality in Fig. 5(b) is improved, but there are still some paired echoes in the image. It is due to the fact that the residual phase errors are much larger than π/4, as shown in the blue line of Fig. 6. However, the residual phase errors estimated by the QML method and the proposed method are smaller than π/4, as shown in the green line and red line of Fig. 6. Thus, we see from Fig. 5(c)-(d) that the paired echoes are completely suppressed.
In order to quantitatively evaluate the focus quality of the point target, a rectangular window is used for intercepting the center point target of Fig. 5(b)-(d). After up-sampling eight times, the resulting azimuth-dimensional amplitude profiles are shown in Fig. 8. And the imaging quality indices of  point target, such as impulse response width (IRW), peak sidelobe ratio (PSLR) and integral sidelobe ratios (ISLR), for different parameters estimation methods are listed in Table 3. As can be seen from the azimuth profiles and the imaging quality indices, the LS method is too poor to suppress the paired echoes well. In the contrary, the QML method and the proposed method are both better and the paired echoes are suppressed well. The PSLR of the QML method is slightly higher than that of the proposed method. It is because that the residual phase error of the QML method is slightly larger than that of the proposed method, as shown in Fig. 6. This suggests that the performance of the proposed method is slightly better than that of the QML method at high SNR.
In the following, we consider the single component highfrequency vibration with the same parameter values, but it is corrupted with Gaussian white noise of 5dB. First, the signal ICR estimated by FrFT and its spectrum are shown in Fig. 8. From the comparison between the estimated ICR and the true ICR, most points are almost consistent with the theoretical values, but some points deviate from the true values. Compared to Fig. 3 at high SNR, there are some outliers affecting the parameters estimate performance. Then we obtain that the ICR frequency is 8.0645 Hz, which is the same as the estimated value at SNR = 20dB. However, the spectrum of the estimated ICR has been affected by noise. With the estimated vibration frequency, the other vibration parameters are obtained, and all parameters estimates and ML function values are listed in Table 4. Fig. 9 shows the estimated vibration displacement. It is clear from the Table 4 that the estimated values by the LS method and the QML method both deviate from the true values. On the one hand, it is due to the error propagation effects mentioned in the above simulation. On the other hand, the outliers of the estimated ICR have aggravated the estimate error. Although one dimensional search for optimal vibration frequency is adopted in  the QML method, the estimates are still no improving. The ML values also show that these two methods are poor. Fortunately, the proposed method not only suppresses the error propagation effects, but also eliminates the effects of outliers, making the parameters estimate values be consistent with the theoretical values, as shown in Table 4. The estimated vibration displacement in the red line of Fig. 9 also does confirm this. Similarly, we use the THz-SAR system parameters in Table 1 and the single component high-frequency vibration error to generate the echo signal for multi-point targets. Fig. 10 gives the imaging results of point targets. Fig. 11 shows that the phase estimate error between the true phase and the estimated phase. Because of no compensation for vibration error, Fig. 10(a) shows that the SAR image quality is affected by paired echoes and seriously degrades. After phase error compensation by the LS method and the QML method, the image quality in Fig. 10(b)-(c) are no improving. It is due to the fact that the residual phase errors are much larger than π 4, as shown in the blue line and green line of Fig. 11. The exciting part is that the paired echoes are completely suppressed after phase error compensation by the proposed method. It is obvious from Fig. 11 that the phase estimate error by the proposed method is much smaller than that by the LS method and the QML method. The resulting azimuth-dimensional amplitude profiles after up-sampling eight times are shown in Fig. 12. And the imaging quality indices of point target for different methods are listed in Table 5. As can be seen from the azimuth profiles and the imaging quality indices, only the proposed method can greatly suppress the paired echo and obtain the image quality equal to that in the ideal case.

B. SIMULATION ANALYSIS OF MULTI-COMPONENTS VIBRATION
In this section, we consider the multi-component highfrequency vibration with the form of (8) and I = 2. Supposing that the first component vibration parameters are the  same as that of simulation in section A and the second component vibration parameters are that f 2 = 15.0 Hz,A 2 = 0.3 × 10 −3 m and ϕ 2 = π 4. In addition, we assume that the signal is corrupted with Gaussian white noise of 20dB. The signal ICR estimated by FrFT and its spectrum are shown in Fig. 13, respectively. It can be seen from Fig. 3 and 13, similar to the ICR estimation results of single component vibration, the ICR estimations are also consistent with the true values. Then we obtain that the signal ICR frequencies are 8.0645 Hz and 16.1290 Hz. Compared to the true vibration frequency, there are considerable estimate errors. With the estimated vibration frequency, the other vibration parameters are obtained, and all parameters estimates and ML function values are listed in Table 6. Fig. 14 gives the estimated vibration displacement. It is clear from the Table 6 that the estimated values based on the LS method deviate from the true values because of the error propagation effects, which makes the estimated vibration displacement deviate from the true vibration displacement, as shown in the blue line of Fig. 14. The ML value also shows that the LS method is not ideal. However, in the QML method and the proposed method, the error propagation effects are both suppressed, which makes the parameters estimate values are almost consistent with the theoretical values, as shown in Table 6.  The THz-SAR system parameters in Table 1 and the two components high-frequency vibration error are used to generate the echo signal for multi-point targets. Firstly, the RD algorithm is used for focusing on imaging, as shown in Fig. 15(a). It can be seen that the SAR image have many paired echoes and the image quality seriously degrades. Then, with the estimated vibration parameters by the different methods, the vibration phase errors are compensated respectively, and the imaging results are given in Fig. 15(b)-(d). Fig. 16 shows the phase estimate error between the true phase and the estimated phase. Compared to Fig. 15(a), the image quality in Fig. 15(b) is improved, but there are still some paired echoes in the image. It is due to the fact that the residual phase errors are much larger than π 4, as shown in the blue line of Fig. 16. However, since the residual phase errors compensated by the QML method and the proposed method are smaller than π 4, it can be seen from Fig. 15(c)-(d) that the paired echoes are completely suppressed. VOLUME 9, 2021   Similarly, the resulting azimuth-dimensional amplitude profiles after up-sampling eight times are shown in Fig. 17.
And the imaging quality indices of point target for different parameters estimation methods are listed in Table 7. As can be seen from the azimuth profiles and the imaging quality indices, the LS method is too poor to suppress the paired echoes well. In the contrary, the QML method and the proposed method are both better and the paired echoes are suppressed well. The PSLR of the QML method is slightly higher than that of the proposed method. It is because that the residual phase error of the QML method is slightly larger than that of the proposed method, as shown in Fig. 16. It suggests that the proposed method is slightly better than the QML method at high SNR. In the following, we consider the two components highfrequency vibration with the same parameter values, but it is corrupted with Gaussian white noise of 5dB. First, the signal ICR estimated by FrFT and its spectrum are given in Fig. 18, respectively. From the comparison between the estimated ICR and the true ICR, most points are almost consistent with the theoretical values, but some points deviate from the true values. Compared to Fig. 13(a) at high SNR, there are some outliers, which affect the parameters estimate. Although the frequency spectrum of estimated ICR is affected by noise, we still obtain that the signal ICR frequencies are 8.0645 Hz  Table 8. Fig. 19 shows the vibration displacement constructed with the estimated parameters. Similar to the single component vibration at low SNR, the estimated values based on the LS method and the QML method both deviate from the true values due to the two factors explained in the section of the single component vibration at low SNR. The ML values also do confirm this. However, as shown in Table 8, the parameters estimate values by the proposed method are almost consistent with the theoretical values, which suggests that the proposed method can suppress both the effects of the two factors. The vibration displacement in the red line of Fig. 19 also does confirm this. Similarly, Fig. 20 gives the imaging results of point targets with two components vibration error. Fig. 21 shows the phase estimate error between the true phase and the estimated phase. When direct imaging without vibration error compensation, the SAR image quality is affected by paired echoes and seriously degrades, as shown in Fig. 20(a). After vibration error compensation based on the LS method and the QML method, the image quality in Fig. 20(b)-(c) are still no improving. It is due to the fact that the residual phase errors are much larger than π 4, as shown in the blue line and green line of Fig. 21. The exciting part is that the paired echoes are completely suppressed after phase compensation by the proposed method, as shown in Fig. 20(d).
Finally, the resulting azimuth-dimensional amplitude profiles after up-sampling eight times are shown in Fig. 22.
And the imaging quality indices of point target with two components vibration for different parameters estimation    methods are listed in Table 9. As can be seen from the azimuth profiles and the imaging quality indices, the LS method and the QML method are both too poor to suppress the paired echoes well. In the contrary, the imaging quality indices of point target are almost consistent with the theoretical values, which shows that the proposed method is good at low SNR.
Remark: In the estimation of single component and two components high-frequency vibration, no matter the SNR is high or low, the LS method is not ideal, but the proposed method is good. Only when the SNR is high, the QML method can obtain ideal estimation results, but its estimation performance is still worse than that of the proposed method.

V. CONCLUSION
In this paper, a high-frequency vibration parameter estimation method based on FrFT with combination of QML and RANSAC was proposed to effectively suppress the paired echo in THz-SAR imaging. The proposed method not only can deal with the multi-components vibration error but also has high estimation accuracy. No matter the SNR is high or low, based on the proposed method, the parameters estimation accuracy can fulfill image well-focused requirements. Both simulation results of single component and multicomponents high-frequency vibration for different SNR are presented to verify the validity of the proposed method. The method also can be extended to solve the problem of high-frequency vibration error for other wavebands, such as the microwave band. Until now, the premise of the proposed method is that there are dominant scatter points in the imaging scene. Otherwise, the proposed method will fail. Thus, the case where there is no dominant scatter point is worth studying in future work. Furthermore, the existing paired echo suppression algorithms assume that the low-frequency motion error has been compensated. The next step is to study the simultaneous compensation of low-frequency motion error and high-frequency vibration error.
YINWEI LI (Member, IEEE) received the bachelor's degree in electronic information engineering from the University of Electronic Science and Technology of China, Chengdu, China, in 2009, and the Ph.D. degree in information and communication engineering from the University of Chinese Academy of Sciences, Beijing, China, in 2014. He was with the Shanghai Radio Equipment Research Institute, Shanghai, China, from 2014 to 2019. He is currently with the Terahertz Technology Innovation Research Institute, University of Shanghai for Science and Technology, Shanghai, and with the Shanghai Institute of Intelligent Science and Technology, Tongji University, Shanghai. He is also an Associate Professor with the University of Shanghai for Science and Technology. His research interests include (inverse) synthetic aperture radar (SAR/ISAR) imaging, interferometric SAR system design and signal processing, and terahertz radar imaging.
QI WU received the bachelor's degree in communication engineering from Chizhou University, Anhui, China, in 2019. He is currently pursuing the master's degree with the University of Shanghai for Science and Technology, Shanghai, China. His research interests include terahertz SAR imaging and motion compensation.
JIAWEI WU received the bachelor's degree in electronic information engineering from Shanghai Dianji University, Shanghai, China, in 2020. He is currently pursuing the master's degree with the University of Shanghai for Science and Technology, Shanghai. His research interests include terahertz ISAR imaging and motion compensation.
PING LI received the Ph.D. degree in engineering from Northwestern Polytechnical University, in 2005. She is currently presiding over one project of the National Key Research and Development plan and one project of the Shanghai Science and Technology Commission. Her research interests include millimeter-wave and terahertz-wave technology, especially with the applications to target detection and imaging.
QIBIN ZHENG received the bachelor's degree from the Northwest University, Xi'an, China, in 2011, and the Ph.D. degree from the University of Science and Technology of China, Hefei, China, in 2016. He is currently a Lecturer with the University of Shanghai for Science and Technology, Shanghai, China. His current research interest includes high speed electronics. His research interests include developing optimization readout electronics and signal processing for physical experiments and industrial application.