An Improved CZT Algorithm for High-Precision Frequency Estimation

Abstract

Estimating the frequencies of multiple superimposed exponentials in noise is an important problem due to its various applications in engineering. In order to obtain good inhibition of spectral leakage and improve the estimation accuracy, an improved Chirp-Z transform (CZT) algorithm is proposed for high-precision frequency estimation. Firstly, the proposed algorithm analyzes the characteristics of the CZT spectrum and utilizes the CZT spectrum to construct a bias correction factor for frequency bias estimation. Then, an expression between the bias correction factor and the frequency estimation error is derived to obtain a more accurate estimate of the frequency bias. Finally, the frequency estimate of the CZT is corrected by the estimated frequency bias so as to obtain a higher frequency estimation accuracy. Compared with the conventional CZT algorithm, the proposed improved CZT algorithm achieves a higher frequency estimation accuracy by correcting the frequency estimate of the CZT method using the estimated frequency bias. The proposed improved CZT algorithm is verified using simulation studies and experimental results, and the results show that it has a higher accuracy and better robustness than the existing methods.

1. Introduction

Frequency estimation is a problem frequently encountered in many fields of science and engineering, such as electrical power systems [1,2,3,4], radar systems [5], acoustic wave detection [6], and speech systems [7,8]. The utilization of frequency information has been extended to many fields, such as frequency-modulated continuous wave (FMCW) radar systems in which the distance between the radar and the target can be obtained by processing and analyzing the frequencies of the received beat signal [9] and in flow meters where fluid characteristics and flow states are obtained based on frequency information in the received signal [10].

There are many typical frequency estimation algorithms, such as the fast Fourier transform (FFT), the Chirp-Z transform (CZT) [11], and the phase method estimation [12]. Among these methods, the CZT method does not require filtering, and it facilitates a high-resolution narrowband analysis. Therefore, the CZT method is a commonly used algorithm for estimating frequency and has very good implementability in engineering applications. On the basis of these algorithms, many frequency estimation algorithms [13,14,15,16,17,18,19,20,21,22,23] have been proposed to improve the frequency estimation accuracy required in various applications. The zero-padding method [13] reduces the computational effort by using the number of non-integer periods in the discrete Fourier transform (DFT) kernel’s quadrature signal. From the perspective of correlation operations, the carrier frequency bias can be estimated from the auto-correlation of the received signal with flexible intervals in OFDM systems [14], and the cross-correlation of the received signal was utilized for frequency estimation by implementing frequency compensation [15]. The phase unwrapping method attempts to estimate the frequency by performing a linear regression on the phase of the received signal [16]. To cope with the problems of fence effects and spectral leakage in the classical FFT, the interpolated DFT (IpDFT) algorithm [17] is proposed to reduce these errors, and the transformed consecutive points of the DFT can be utilized to obtain a more accurate frequency estimation [18]. Candan [19] estimated the frequency in two stages. A coarse frequency estimate is obtained by utilizing the amplitude spectrum of the DFT in the first stage, and a bias correction term is estimated and utilized to improve the frequency estimation accuracy in the second stage. Aboutanios and Mulgrew proposed the A&M algorithm in [20], which shifted the DFT coefficients by half and resulted in an improved frequency estimation performance. Based on [20], Serbes proposed the q-Shift estimator (QSE) algorithm in [21], whose coefficients are equivalent to interpolating the DFT of the signal by $1/\left|q\right|$. The QSE algorithm can be considered a generalization of the A&M algorithm and has a better estimation performance than the A&M algorithm. Reshma [22] proposed combining the CZT method with the Candan algorithm [19] to improve the frequency resolution by deriving the bias correction term through the amplitude of the CZT spectral. In [23], Kayvan proposed the CZT-based algorithm for the detection of multiple targets for the FMCW radar.

As illustrated in [24], the CZT-based methods can improve the frequency resolution of the spectrum and reduce the measurement errors caused by spectral leakage by analyzing a narrow frequency band. Considering the advantages of the CZT method over other methods, we investigate the CZT method in depth in this paper. This paper proposes an improved CZT algorithm to further improve the frequency estimation accuracy of the CZT algorithm. The proposed algorithm first analyzes the characteristics of the CZT spectrum and utilizes the CZT spectrum to construct a bias correction factor for frequency bias estimation. Moreover, the relationship between the bias correction factor and the intrinsic frequency estimation error is derived, and from this, a more accurate estimate of the frequency bias is obtained. Then, the frequency bias estimate is used to modify the frequency estimate of the CZT algorithm to obtain a more accurate frequency estimate than the CZT algorithm. Finally, the proposed improved CZT algorithm is verified using simulation studies and experimental results, and the results show that it has a higher estimation accuracy and better robustness than the existing methods.

The remainder of this paper is organized as follows. In Section 2, the signal model is described, and the basic Chirp-Z transform is described. In Section 3, the proposed improved CZT algorithm is presented. In Section 4, the performance of the proposed algorithm is evaluated using simulation studies and experimental results. Finally, the conclusions are provided in Section 5.

2. Signal Model and Related Works

Suppose that there are P frequency components in the signal. It takes the following form:

$$x\left(t\right)=\sum _{i=1}^P{a}_{i}exp\left(j2\pi f_{ci}t+{\phi }_i\right)$$

(1)

where $a_i\left(i=1,\dots ,P\right)$ represents the amplitude, $f_{ci}\left(i=1,\dots ,P\right)$ represents the frequency, and ${\phi }_i\left(i=1,\dots ,P\right)$ represents the phase.

In order to facilitate the derivation of the proposed method, a single-tone signal model with $P=1$ in (1) is considered. It should be noted that the same derivation can be applied to the signal model of multiple frequency components. When $P=1$, the single-tone signal model takes the following form:

$$x\left(t\right)=a_{0}exp\left(j2\pi f_{c}t+{\phi }_0\right)$$

(2)

where $a_0$ is the signal amplitude, ${\phi }_0$ is the initial phase, and $f_c$ represents the frequency to be estimated. Then, the sampled signal $x\left(n\right)$ is

$$x\left(n\right)=a_{0}exp\left(j2\pi f_c\frac{n}{f_s}+{\phi }_0\right),n=0,1,\cdots ,N-1$$

(3)

where $f_s$ is the sampling rate, and N is the number of samples. The FFT algorithm is applied to $x\left(n\right)$ to obtain the spectrum $X\left(k\right)\left(k=0,1,\cdots ,N-1\right)$, which then obtains the estimate ${\widehat{f}}_p$, corresponding to the frequency index $k_p$ of the maximum of the FFT spectrum amplitude. The conventional CZT algorithm [11] is centered on the frequency estimate ${\widehat{f}}_p$, and the refined frequency range is chosen as $\left(f_1,f_2\right)$,

$$f=f_1+B_{\mathrm{CZT}}\frac{k}{M},k=0,\cdots ,M-1$$

(4)

where $B_{\mathrm{CZT}}$ is the bandwidth of $\left(f_1,f_2\right)$, and M represents the length of CZT. The CZT algorithm [11] is performed on $x\left(n\right)$ to obtain

$$X_{\mathrm{CZT}}\left(z_k\right)=\sum _{n=0}^{N-1}x\left(n\right)z_k^{-n},k=0,\cdots ,M-1$$

(5)

where $z_k=A{W}^{-k}$, $A=A_0{e}^{j{\theta }_0}$, $W=W_0{e}^{-j{\varphi }_0}$, ${\theta }_0$ is the starting sampling angle, and ${\varphi }_0$ is the angle between two adjacent sampling points. The frequency range $\left(f_1,f_2\right)$ of the CZT is determined based on the frequency index $k_p$ and the parameter q and is given by

$$f_1=\frac{f_s}{N}\left(k_p-q\right)$$

(6)

$$f_2=\frac{f_s}{N}\left(k_p+q\right)$$

(7)

Let $A_0=1$ and $W_0=1$ .Then, the CZT of $x\left(n\right)$ is given by

$$\begin{array}{c}{X}_{\mathrm{CZT}}\left(k\right)=\sum _{n=0}^{N-1}x\left(n\right)z_k^{-n}=\sum _{n=0}^{N-1}x\left(n\right)exp\left(-j{\theta }_0\right)exp\left(-j{\varphi }_{0}nk\right)\\ =a_{0}exp\left(j2\pi {\phi }_0\right)\times exp\left(j\pi \left(N-1\right)\left(\frac{f_c}{f_s}-\frac{{\theta }_0}{2\pi }-\frac{{\varphi }_{0}k}{2\pi }\right)\right)\times \frac{sin\left(\pi N\left(\frac{f_c}{f_s}-\frac{{\theta }_0}{2\pi }-\frac{{\varphi }_{0}k}{2\pi }\right)\right)}{sin\left(\pi \left(\frac{f_c}{f_s}-\frac{{\theta }_0}{2\pi }-\frac{{\varphi }_{0}k}{2\pi }\right)\right)},\\ k=0,\cdots ,M-1\end{array}$$

(8)

where ${\theta }_0=2\pi \left(k_p-q\right)/\phantom{\left(k_p-q\right)N}N$, and ${\varphi }_0=2\pi \left(2q\right)/\phantom{\left(2q\right)\left(MN\right)}\left(MN\right)$. In this case, the frequency resolution of the CZT is $2q{f}_s/\phantom{2q{f}_s\left(NM\right)}\left(NM\right)$. By substituting ${\theta }_0$ and ${\varphi }_0$ into (7), the CZT spectrum of $x\left(n\right)$ can be expressed as:

$$\begin{array}{c}{X}_{\mathrm{CZT}}\left(k\right)=a_{0}exp\left(j2\pi {\phi }_0\right)\times \hfill \\ exp\left(j\pi \left(N-1\right)\left(\frac{f_c}{f_s}-\frac{k_p-q}{N}-\frac{2qk}{MN}\right)\right)\times \frac{sin\left(\pi \left(\frac{f_{c}N}{f_s}-\left(k_p-q\right)-\frac{2qk}{M}\right)\right)}{sin\left(\pi \left(\frac{f_c}{f_s}-\frac{\left(k_p-q\right)}{N}-\frac{2qk}{NM}\right)\right)}\hfill \end{array}$$

(9)

Assuming that the frequency index at the maximum amplitude of the CZT spectrum is $k_m$, the frequency estimate is expressed as:

$${\widehat{f}}_m=\left(k_p-q\right)\frac{f_s}{N}+k_m\frac{2q{f}_s}{NM}$$

(10)

Because of the fence effect of the discrete signal, the maximum value ${\widehat{f}}_m$ of the CZT spectral in (9) is difficult to coincide with the true value $f_c$, and thus, there exists a frequency bias between them. In order to solve this problem, an improved CZT algorithm will be proposed in the next section for estimating this frequency bias to obtain a more accurate frequency estimate.

3. Improved CZT Algorithm

In this section, the design of the bias correction factor for frequency estimation is provided first. Then, the approximate solution of the bias correction factor $\mu$ is presented and utilized to estimate the frequency bias.

3.1. Design of Bias Correction Factor $\mu$ for Frequency Bias Estimation

Assuming that the true frequency to be estimated $f_c$ can be expressed as:

$$f_c=\left(k_p-q\right)\frac{f_s}{N}+\left(k_m+\delta \right)\frac{2q{f}_s}{NM}$$

(11)

where $2q{f}_s/\phantom{2q{f}_{s}N}N$ is the frequency bandwidth, $k_m$ is the frequency index corresponding to the maximum amplitude of the CZT spectrum, and $\delta \in \left[-0.5,0.5\right]$ is the bias between the true frequency and the frequency estimate obtained by the CZT algorithm. In order to accurately estimate the frequency bias $\delta$, the CZT spectrum is analyzed, and the bias correction factor $\mu$ is constructed based on the maximum point $X_{\mathrm{CZT}}\left(k_m\right)$ and its two adjacent spectral points, $X_{\mathrm{CZT}}\left(k_m-1\right)$ and $X_{\mathrm{CZT}}\left(k_m+1\right)$, as follows.

Following the principle of the CZT algorithm, the maximum point $X_{\mathrm{CZT}}\left(k_m\right)$ of the CZT spectrum can be expressed as:

$$X_{\mathrm{CZT}}\left(k_m\right)=a_{0}exp\left(j2\pi {\phi }_0\right)\times exp\left(j2\pi \left(\frac{q\delta \left(N-1\right)}{MN}\right)\right)\times \frac{sin\left(\pi \frac{2q\delta }{M}\right)}{sin\left(\pi \frac{2q\delta }{MN}\right)}$$

(12)

The left-hand point $X_{\mathrm{CZT}}\left(k_m-1\right)$ and the right-hand point $X_{\mathrm{CZT}}\left(k_m+1\right)$ of the maximum point $X_{\mathrm{CZT}}\left(k_m\right)$ can be expressed as:

$$X_{\mathrm{CZT}}\left(k_m-1\right)=a_{0}exp\left(j2\pi {\phi }_0\right)\times exp\left(j2\pi \left(\frac{q\left(\delta +1\right)\left(N-1\right)}{MN}\right)\right)\times \frac{sin\left(\pi \frac{2q\left(\delta +1\right)}{M}\right)}{sin\left(\pi \frac{2q\left(\delta +1\right)}{MN}\right)}$$

(13)

$$X_{\mathrm{CZT}}\left(k_m+1\right)=a_{0}exp\left(j2\pi {\phi }_0\right)\times exp\left(j2\pi \left(\frac{q\left(\delta -1\right)\left(N-1\right)}{MN}\right)\right)\times \frac{sin\left(\pi \frac{2q\left(\delta -1\right)}{M}\right)}{sin\left(\pi \frac{2q\left(\delta -1\right)}{MN}\right)}$$

(14)

By using $X_{\mathrm{CZT}}\left(k_m\right)$, $X_{\mathrm{CZT}}\left(k_m-1\right)$ and $X_{\mathrm{CZT}}\left(k_m+1\right)$, the bias correction factor $\mu$ is constructed as:

$$\mu =\frac{X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m+1\right)exp\left(j2\pi \frac{q\left(N-1\right)}{MN}\right)\frac{sin\left(\pi \frac{2q\left(\delta -1\right)}{MN}\right)}{sin\left(\pi \frac{2q\delta }{MN}\right)}}{X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m-1\right)exp\left(-j2\pi \frac{q\left(N-1\right)}{MN}\right)\frac{sin\left(\pi \frac{2q\left(\delta +1\right)}{MN}\right)}{sin\left(\pi \frac{2q\delta }{MN}\right)}}$$

(15)

By substituting (12), (13), and (14) into (15), $\mu$ is derived as:

$$\mu =\frac{sin\left(\pi \frac{2q\delta }{M}\right)-sin\left(\pi \frac{2q\left(\delta -1\right)}{M}\right)}{sin\left(\pi \frac{2q\delta }{M}\right)-sin\left(\pi \frac{2q\left(\delta +1\right)}{M}\right)}$$

(16)

Through simple operations, (16) can be rewritten as:

$$\mu =\frac{tan\left(\pi \frac{2q\delta }{M}\right)-tan\left(\pi \frac{2q\delta }{M}\right)cos\left(\pi \frac{2q}{M}\right)+sin\left(\pi \frac{2q}{M}\right)}{tan\left(\pi \frac{2q\delta }{M}\right)-tan\left(\pi \frac{2q\delta }{M}\right)cos\left(\pi \frac{2q}{M}\right)-sin\left(\pi \frac{2q}{M}\right)}$$

(17)

From (17), the relationship between $\mu$ and $\delta$ is given by:

$$\delta =arctan\left(\frac{\left(\mu +1\right)sin\left(\pi \frac{2q}{M}\right)}{\left(\mu -1\right)\left(1-cos\left(\pi \frac{2q}{M}\right)\right)}\right)\frac{M}{2q\pi }=arctan\left(\frac{\left(\mu +1\right)cos\left(\pi \frac{q}{M}\right)}{\left(\mu -1\right)sin\left(\pi \frac{q}{M}\right)}\right)\frac{M}{2q\pi }$$

(18)

Then, the frequency bias $\delta$ can be obtained from the bias correction factor $\mu$ in (18).

3.2. Approximate Solution of the Bias Correction Factor $\mu$

As $\theta$ approaches 0, there exists $sin\left(\theta \right)\approx \theta$. Therefore, when $MN$ is large, (15) can be approximated to:

$$\begin{array}{c}\mu \approx \frac{X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m+1\right)exp\left(j2\pi \frac{q\left(N-1\right)}{MN}\right)\frac{\left(\delta -1\right)}{\delta }}{X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m-1\right)exp\left(-j2\pi \frac{q\left(N-1\right)}{MN}\right)\frac{\left(\delta +1\right)}{\delta }}\hfill \\ =\frac{X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m+1\right)exp\left(j2\pi \frac{q\left(N-1\right)}{MN}\right)+X_{\mathrm{CZT}}\left(k_m+1\right)exp\left(j2\pi \frac{q\left(N-1\right)}{MN}\right)\frac{1}{\delta }}{X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m-1\right)exp\left(-j2\pi \frac{q\left(N-1\right)}{MN}\right)-X_{\mathrm{CZT}}\left(k_m-1\right)exp\left(-j2\pi \frac{q\left(N-1\right)}{MN}\right)\frac{1}{\delta }}\hfill \end{array}$$

(19)

From (12) and (14), the following can be obtained:

$$\begin{array}{c}{X}_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m+1\right)exp\left(j2\pi \frac{q\left(N-1\right)}{MN}\right)\hfill \\ =a_{0}exp\left(j2\pi {\phi }_0\right)\times exp\left(j2\pi \left(\frac{q\delta \left(N-1\right)}{MN}\right)\right)\times \left(\frac{sin\left(\pi \frac{2q\delta }{M}\right)}{sin\left(\pi \frac{2q\delta }{MN}\right)}-\frac{sin\left(\pi \frac{2q\left(\delta -1\right)}{M}\right)}{sin\left(\pi \frac{2q\left(\delta -1\right)}{MN}\right)}\right)\hfill \end{array}$$

(20)

From (12) and (13), the following can be obtained:

$$\begin{array}{c}{X}_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m-1\right)exp\left(-j2\pi \frac{q\left(N-1\right)}{MN}\right)\hfill \\ =a_{0}exp\left(j2\pi {\phi }_0\right)\times exp\left(j2\pi \left(\frac{q\delta \left(N-1\right)}{MN}\right)\right)\times \left(\frac{sin\left(\pi \frac{2q\delta }{M}\right)}{sin\left(\pi \frac{2q\delta }{MN}\right)}-\frac{sin\left(\pi \frac{2q\left(\delta +1\right)}{M}\right)}{sin\left(\pi \frac{2q\left(\delta +1\right)}{MN}\right)}\right)\hfill \end{array}$$

(21)

When $MN$ is large, (20) and (21) can be, respectively, approximated as follows:

$$X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m+1\right)exp\left(j2\pi \frac{q\left(N-1\right)}{MN}\right)\approx 0$$

(22)

$$X_{\mathrm{CZT}}\left(k_m\right)-X_{\mathrm{CZT}}\left(k_m-1\right)exp\left(-j2\pi \frac{q\left(N-1\right)}{MN}\right)\approx 0$$

(23)

Substituting (22) and (23) into (19), $\mu$ is approximated as follows:

$$\widehat{\mu }\approx -\frac{X_{\mathrm{CZT}}\left(k_m+1\right)exp\left(j2\pi \frac{q\left(N-1\right)}{MN}\right)}{X_{\mathrm{CZT}}\left(k_m-1\right)exp\left(-j2\pi \frac{q\left(N-1\right)}{MN}\right)}\approx -\frac{X_{\mathrm{CZT}}\left(k_m+1\right)}{X_{\mathrm{CZT}}\left(k_m-1\right)}cos\left(4\pi \frac{q\left(N-1\right)}{MN}\right)$$

(24)

Substituting (24) into (18), the approximation of $\delta$ can be derived as follows:

$$\widehat{\delta }=arctan\left(\frac{\left(\widehat{\mu }+1\right)cos\left(\pi \frac{q}{M}\right)}{\left(\widehat{\mu }-1\right)sin\left(\pi \frac{q}{M}\right)}\right)\frac{M}{2q\pi }$$

(25)

Substituting $\widehat{\delta }$ into (11), the frequency ${\widehat{f}}_c$ is estimated as follows:

$${\widehat{f}}_c=\left(k_p-q\right)\frac{f_s}{N}+\left(k_m+\widehat{\delta }\right)\frac{2q{f}_s}{NM}$$

(26)

It can be seen that the frequency estimate of the improved CZT method in (26) is obtained by correcting the frequency estimate of the CZT method in (10) using the frequency bias $\widehat{\delta }$. Therefore, the proposed improved CZT algorithm can obtain a higher frequency estimation accuracy than the basic CZT algorithm. The improved CZT algorithm is summarized in Algorithm 1. Moreover, the computational complexity of various methods is shown in

Table 1. The computational complexity of various methods.

Method	Computational Complexity
FFT	$O\left(NlogN\right)$
Candan	$O\left(NlogN\right)$
A&M algorithm	$O\left(NlogN\right)$
QSE algorithm	$O\left(NlogN\right)$
CZT-based algorithm	$O\left(PNlogN+PMlogN\right)$
CZT Candan	$O\left(NlogN+MlogN\right)$
Proposed improved CZT	$O\left(NlogN+MlogN\right)$

.

Algorithm 1: Improved CZT algorithm.

Input: the sampled signal $x\left(n\right)$, $n=0,1,\cdots ,N-1$;           Set the initial paraments q and M of CZT. Output: ${\widehat{f}}_c$ 1: N-point FFT algorithm is applied to $x\left(n\right)$ to obtain its spectrum $X\left(k\right)\left(k=0,1,\cdots ,N-1\right)$; Find its index $k_p$ corresponding to the maximum of $X\left(k\right)$ and obtain the frequency estimate ${\widehat{f}}_p$; 2: Apply the conventional CZT algorithm to $x\left(n\right)$ in the refined frequency range $\left(f_1,f_2\right)$; 3: Find the frequency index $k_m$ corresponding to the maximum amplitude of the CZT spectrum; 4: Calculate the approximation of $\mu$ using (24); 5: Calculate the approximation of $\delta$ using (25); 6: Obtain the frequency estimate ${\widehat{f}}_c$ using (26).

4. Experimental Analysis

In this section, the proposed improved CZT algorithm is verified using simulation studies and field test results for the frequency estimation in the FMCW radar. The performance of the improved CZT algorithm is compared with that of the FFT algorithm, the Candan algorithm [19], the A&M algorithm [20], the QSE algorithm [21], the CZT Candan algorithm [22], and the CZT-based algorithm [23]. In simulations, to ensure that the various algorithms are being compared under the same conditions, the initial values of the parameters of the various algorithms are kept consistent. The sampling rate and sampling points are, respectively, set to $f_s$ = 92,783 Hz and $N=1024$. The parameters q and M are the key parameters of the basic CZT algorithm and the proposed improved CZT algorithm. The parameter q determines the frequency band of the CZT algorithm, and M represents the number of points of the CZT algorithm. The smaller the q and the bigger the M, the more accurate the frequency estimation. The value of q is determined by the required frequency resolution. In order to verify the effect of q on the frequency resolution, the frequency estimates of the proposed method versus various values of q are provided in Figure 1. As can be seen from Figure 1, the smaller the q, the more accurate the frequency estimation.

In order to keep the same conditions for the methods based on the CZT algorithm in later simulations and the experimental analysis, the parameters q and M are set to $q=1$ and $M=32$. In addition, all results are obtained over 5000 independent Monte Carlo trials.

4.1. Simulation Studies

In the first simulation, we compared the root mean square error (RMSE) and variation ranges of the frequency estimates of the improved CZT algorithm versus the SNR. Figure 2 shows the means and variation ranges of the frequency estimates of the various methods versus the SNR when the frequency is $f_c=5070\mathrm{Hz}$; the compared algorithms include the conventional FFT, the Candan algorithm [19], the A&M algorithm [20], the QSE algorithm [21], the CZT Candan algorithm [22], the CZT-based algorithm [23], and the proposed improved CZT algorithm. Figure 2 shows the means and variation ranges of the frequency estimates of the various methods versus the SNR when $f_c=5070\mathrm{Hz}$, where (a) is the signal disturbed by the Gaussian noise, and (b) is the signal disturbed by the uniform noise. As can be seen from Figure 2, the mean of the frequency estimates of the proposed algorithm is closest to the true frequency of all of the methods, which indicates that the proposed improved CZT algorithm has a higher estimation accuracy than the other methods. Moreover, the proposed algorithm is affected by noise less than other methods, and the fluctuation range of the frequency estimates of the proposed algorithm is less than that of the other methods when the SNR is relatively high. In addition, the mean of the frequency estimates for the CZT-based algorithm, the CZT Candan, and the proposed algorithm is closer to the true frequency when the SNR is relatively high, and it is affected by noise. The variation ranges of frequency estimates become larger when the SNR is relatively low. Figure 3 shows the RMSE of the frequency estimates of various methods and the Cramer–Rao lower bound (CRLB) versus the SNR, where (a) is the signal disturbed by the Gaussian noise, and (b) is the signal disturbed by the uniform noise. In Figure 3, the formula of the theoretical CRLB has the following expression:

$$\mathrm{CRLB}=\frac{{\sigma }^2}{{a_0}^2{\displaystyle \sum _{i=0}^{N-1}}4{\pi }^2\frac{n^2}{{f_s}^2}{\left|exp\left(2\pi j\frac{n}{f_s}{f}_c\right)\right|}^2}$$

(27)

As can be observed from Figure 3, the larger the SNR, the smaller the RMSE of all algorithms. When $\mathrm{SNR}<-10\mathrm{dB}$ in the signal disturbed by the Gaussian noise, the RMSE of the proposed improved CZT algorithm is less than that of the CZT-based algorithm and is slightly greater than that of the FFT algorithm. When $\mathrm{SNR}\ge -10\mathrm{dB}$ in the signal disturbed by the Gaussian noise and in the signal disturbed by the uniform noise, the RMSE of the proposed algorithm is smaller than that of the existing algorithms, which indicates that the proposed improved CZT algorithm has a better estimation performance than other methods.

In the second simulation, the signal is disturbed by the Gaussian noise.

Table 2. Means of the frequency estimates versus the frequency when the SNR is 10 dB.

Preset	Mean of Frequency Estimates/Hz
Frequency/Hz	FFT	Candan	A&M Algorithm	CZT-based Algorithm	QSE Algorithm	CZT Candan	Proposed Algorithm
3000	2990.0933	2990.7082	2995.0410	3001.4171	3001.3411	3000.7109	3000.4748
3200	3171.3110	3179.0949	3185.6476	3199.6263	3199.6462	3199.8110	3199.8726
3400	3443.1377	3408.5396	3421.5705	3398.0031	3398.1102	3398.9988	3399.3300
3600	3624.3555	3619.3810	3612.1630	3601.7010	3601.6092	3600.8477	3600.5631
3800	3805.5732	3805.3871	3802.7203	3799.9102	3799.9149	3799.9543	3799.9689
4000	3986.7910	3987.9414	3993.3664	3998.1386	3998.2387	3999.0689	3999.3789
4200	4168.0088	4178.8213	4183.9950	4201.9565	4201.8515	4200.9798	4200.6542
4400	4439.8354	4415.9718	4419.9163	4400.1941	4400.1838	4400.0991	4400.0672
4600	4621.0532	4617.6298	4610.5368	4598.4021	4598.4880	4599.1998	4599.4658
4800	4802.2710	4802.2373	4800.9280	4801.9527	4801.8480	4800.9787	4800.6545
5000	4983.4888	4985.4033	4991.7410	5000.4779	5000.4520	5000.2367	5000.1557

shows the mean of the frequency estimates of the various methods versus the frequency $f_c$ when the SNR is 10 dB. As can be seen from Table 2, the mean of the frequency estimates of the proposed algorithm is closest to the true frequency of all of the methods, which indicates that the proposed algorithm has a better estimation performance than the other methods.

Table 3. RMSEs of the frequency estimates versus the frequency when the SNR is 10 dB.

Preset	RMSE of Frequency Estimates/Hz
Frequency/Hz	FFT	Candan	A&M Algorithm	CZT-Based Algorithm	QSE Algorithm	CZT Candan	Proposed Algorithm
3000	21.8910	20.5400	11.8729	3.1414	2.9734	1.6241	1.1728
3200	63.3943	46.2041	31.7640	0.8257	0.7829	0.5627	0.5766
3400	95.3217	19.1503	47.6832	4.9017	4.6463	2.5550	1.8097
3600	53.8185	42.8342	26.9534	3.7670	3.5642	1.9191	1.3526
3800	12.3152	11.9157	10.1294	0.1985	0.1925	0.3969	0.5168
4000	29.1880	26.6527	15.0581	4.1845	3.9612	2.1415	1.5042
4200	70.6913	46.8121	35.4033	4.4196	4.1842	2.2624	1.5841
4400	88.0247	35.3550	44.0314	0.4288	0.4084	0.4441	0.5366
4600	46.5215	38.9633	23.3930	3.5353	3.3457	1.8123	1.2881
4800	5.0182	4.9683	20.4068	5.1889	4.9222	2.7312	1.9435
5000	36.4850	32.2607	18.4582	1.0561	0.9995	0.6479	0.6154

shows the RMSE of the frequency estimates of the various methods versus the frequency when the SNR is 10 dB. As can be observed from Table 3, the proposed algorithm has the smallest RMSE for most of the frequencies, with the exception that it has a higher RMSE than the other methods at individual frequencies, such as $f_c=3800\mathrm{Hz}$ and $f_c=4400\mathrm{Hz}$, due to the fence effect of the CZT-based algorithm. From the results in Table 2 and Table 3, it can be seen that the frequency estimation accuracy of the proposed algorithm is better and more robust than for the other algorithms.

The third simulation is utilized to verify the performance of the proposed algorithm for the multiple-tone signal. Figure 4 shows the means of the frequency estimates of the multiple-tone signal versus the SNR when the frequencies are (a) $f_c=3050\mathrm{Hz}$, (b) $f_c$ =14,000 Hz, and (c) $f_c$ =18,000 Hz. The compared algorithms include the CZT Candan algorithm [22], the CZT-based algorithm [23], and the proposed improved CZT algorithm. As can be observed in Figure 4, three methods have the ability to recognize multiple frequencies. As the SNR increases, the estimated frequency results tend to be stable. The CZT-based algorithm and the CZT Candan algorithm perform better in (a), while the proposed algorithm performs better in (b) and (c). Therefore, the proposed algorithm outperforms the other algorithms for the multi-tone signal.

4.2. Field Test Results

In this section, a test system for the FMCW radar was built to verify the performance of the proposed algorithm in real-world engineering applications, as shown in Figure 5, where the resolution of the AD converter is 12-bit in the FMCW radar. The test system includes an FMCW radar, a laser rangefinder, a reflector plate mounted on the wall, and a rail for mounting the FMCW radar. In this test system, the rail is used to adjust the distance between the reflector plate on the wall and the FMCW radar, and the laser rangefinder is used to measure the reference distance between the FMCW radar and the reflector plate. According to the principle of the FMCW radar in [9], various distances between the reflector plate and the FMCW radar mean that the frequency value of the beat signal received by the radar varies, and the distance information can be obtained from the frequency estimates of the beat signal. In this test system, the accuracy of the frequency estimation determines the distance estimation accuracy, so this system is utilized to test the estimation performance of the proposed improved CZT algorithm in practical applications. The parameters of the test system are shown in

Table 4. Parameters of the FMCW test system.

System Parameter	Value
Initial frequency	23.6 GHz
Chirp bandwidth	999.807 MHz
Sampling points	1024
Sampling rate	123.2 KHz
Electromagnetic wave speed	299,709 Km/s

.

Utilizing the test system above, the distance between the FMCW radar and the reflector is adjusted between 1000 mm and 2000 mm.

Table 5. Frequency estimates versus various distances and frequencies of the various algorithms.

Distances	True	Frequency Estimates/Hz
/mm	Frequency/Hz	FFT	Candan	A&M Algorithm	CZT-Based Algorithm	QSE Algorithm	CZT Candan	Proposed Algorithm
1000.00	802.71	883.36	881.64	835.55	801.66	802.45	803.56	804.26
1050.00	842.84	883.36	886.15	888.73	846.79	846.52	844.78	844.04
1100.00	882.98	883.36	885.38	860.58	876.88	877.61	878.59	879.15
1200.00	963.25	883.36	968.65	906.90	955.11	955.95	957.40	958.05
1250.00	1003.38	1124.06	1092.32	1039.72	997.23	997.94	999.04	999.50
1300.00	1043.52	1124.06	1118.35	1064.30	1042.36	1042.75	1042.89	1042.97
1350.00	1083.65	1124.06	1124.30	1090.51	1087.49	1087.42	1086.06	1085.55
1450.00	1163.93	1124.06	1143.97	1128.64	1162.71	1163.14	1163.34	1163.58
1500.00	1204.06	1124.06	1215.90	1150.18	1207.85	1207.56	1205.25	1204.42
1550.00	1244.20	1364.77	1335.92	1285.75	1237.94	1238.65	1239.90	1240.59
1600.00	1284.33	1364.77	1360.45	1310.09	1283.07	1283.29	1282.75	1282.68
1650.00	1324.47	1364.77	1368.09	1374.77	1328.20	1327.82	1325.09	1324.10
1700.00	1364.60	1364.77	1369.62	1357.73	1358.29	1358.89	1359.78	1360.22
1750.00	1404.74	1364.77	1385.40	1371.85	1403.42	1403.56	1402.76	1402.48
1800.00	1444.87	1364.77	1459.35	1393.15	1448.55	1448.35	1446.23	1445.43
1900.00	1525.14	1605.47	1602.95	1557.62	1523.77	1524.17	1524.34	1524.51
1950.00	1565.28	1605.47	1606.32	1580.01	1568.91	1568.68	1566.45	1565.53
2000.00	1605.41	1605.47	1606.17	1575.14	1598.99	1599.65	1600.74	1601.40

gives the frequency estimates versus the various distances as well as the frequencies of the various algorithms, including the 1024-point FFT, the Candan algorithm [19], the A&M algorithm [20], the QSE algorithm [21], the CZT Candan algorithm [22], the CZT-based algorithm [23], and the proposed algorithm. Here, the means of the errors are defined as the difference between the frequency estimates and the true frequency. As can be seen in Table 5, the means of the proposed improved CZT algorithm are closer to the true values than those of other algorithms. Therefore, the proposed algorithm has a higher estimation accuracy than the other methods. Figure 6 shows the estimation errors of the frequency estimates under the different distances of the various algorithms. As can be seen in Figure 6, the maximum relative error of the proposed algorithm is 5.20 Hz, and the minimum error is 0.25 Hz while the maximum error of the CZT-based algorithm is 8.14 Hz, and the maximum error of the CZT Candan algorithm is 5.85 Hz. Moreover, the frequency estimates of the proposed algorithm are closer to the reference frequency than those of other algorithms. Therefore, the proposed improved CZT algorithm can effectively improve the frequency estimation capability of the FMCW radar and thus can improve the ranging accuracy of the FMCW radar system.

5. Conclusions

This paper proposes an improved CZT algorithm for high-precision frequency estimation, which can effectively suppress spectrum leakage. By theoretical analysis and derivation, the bias correction factor $\mu$ is constructed based on the CZT spectrum, and the approximate solution for the bias correction factor $\mu$ is derived according to the CZT spectrum. The relationship between the bias correction factor and the true frequency is used to make corrections to the frequency estimation. The simulation results show that the proposed algorithm has a lower estimated frequency mean error and lower RMSE than other algorithms and can effectively suppress the fence effect. In contrast to the CZT Candan algorithm, the proposed algorithm utilizes both the amplitude and phase information and shows a better performance. Moreover, the computation complexity of the proposed improved CZT algorithm is only slightly larger than the basic CZT algorithm, but it has a better frequency estimation performance than the basic CZT algorithm. In addition, the test system for the FMCW radar is applied to verify the effectiveness of the proposed algorithm in actual engineering applications.