一種快速高精度正弦信號頻率估計新方法

楊建

【摘? ?要】? ?針對受高斯白噪聲干擾的正弦波信號，提出了一種基于DFT 插值的信號頻率估計的新方法。在總結經典插值算法的基礎上，充分利用DFT后頻譜峰值及其左右譜線的實部和虛部值進行高精度頻率估計。經理論分析和仿真實驗證明，當信噪比（SNR）在-5 dB以上時，頻率估計的方差趨近于克拉-美羅下界（CRLB）。與經典Rife法對比，相同條件下該算法綜合性能更好、更穩定，并且計算復雜度更低，便于實現和應用。

【關鍵詞】? ?頻率估計；快速傅里葉變換；正弦波信號；高斯白噪聲；信噪比

A Fast and Accurate Frequency Estimation Algorithm of Sinusoid Signal

Yang Jian

（National Key Laboratory of Electromagnetic Space Security， Chengdu 610036， China）

【Abstract】? ? For the sinusoidal signal disturbed by Gaussian white noise，a frequency estimation algorithm based on DFT interpolation algorithm is proposed in this paper. Based on the classical interpolation algorithm， the algorithm of this paper takes full use of the Peak Spectral Frequency and its neighbor spectral lines to estimate the frequency of the signal. The analysis and simulation results indicate that the performance of the algorithm approaching to the Cramer Rao lower bound （CRLB） when the signal to noise ratio （SNR） is higher than -5 dB. Compared with the classical Rife algorithm， the proposed algorithm is better and more stable in performance under the same conditions， and the algorithm has low computational complexity and is easy to implement and apply.

【Key words】? ? ?frequency estimation; fast Fourier transform （FFT）; single sine signal; Gauss white noise; SNR

〔中圖分類號〕? TN971? ? ? ? ? ? ? ?〔文獻標識碼〕? A ? ? ? ? ? ? ?〔文章編號〕 1674 - 3229（2024）01- 0042 - 05

0? ? ?引言

目前應用于通信領域的數字接收機的信號處理要求具有良好的帶通采樣和高效的下變頻功能，這就需要測量信號頻率以便于實施后續的信號和信息處理。接收機需要對待接收的電磁信號進行采樣并下變頻后精確估計其頻率。對接收的信號（如通信信號等）載頻精確估計，是幾乎所有數字接收機的基本而關鍵的技術。測頻算法包括時域法、頻域法以及時頻分析法等等。這些算法中，多信號分離能力強且信噪比要求較低的頻域算法，非常適合在多種通信信號混疊環境以及復雜通信信號調制的情況使用。頻譜分析為傅里葉變換的最重要功能，頻域算法大多采用離散傅里葉變換（DFT，discrete fourier transform）進行頻譜分析。DFT算法成熟且所得頻譜的含義明確，通過快速傅里葉變換（FFT，fast fourier transform）可以使計算量降低到可接受范圍，因此頻域分析法估計頻率得到了廣泛應用。

頻域法估計信號頻譜時，在信號功率未被噪聲功率覆蓋（信噪比高于某一閾值）情況下，定位幅度譜峰值位置可以計算出信號載頻。離散化的兩譜線距離為[fsN]（[fs]是信號采樣率，N是FFT點數），介于譜線間的頻率都會劃歸到最近的譜線上，故FFT后直接估計載頻的誤差范圍是[0～[fs2N]]。應用于通信等領域的數字接收機需要更高的頻率測量精度，因此需要采用更好的頻率估計算法來提高測頻的精度，以保證信息傳輸質量和接收效果。

針對頻域算法進行信號頻率估計，學者們已經提出很多方法。文獻［1］和文獻［2］提出了頻率的最大似然估計（ML）算法，該算法對受高斯白噪聲干擾的正弦信號的估計誤差能達到克拉-美羅下限（CRLB，Cramér-Rao low bound），是一種最優估計，該算法缺點是需要采用搜索運算，計算量很大。文獻［3］提取出Rife算法，利用最大的兩根譜線進行插值運算，然后對頻率進行估計，該算法優點是在兩離散量化頻率間的中心位置附近時，頻率估計非常接近真實值，并且估計方差值趨近于CRLB，缺點是在信號頻率位于離散的量化頻率點附近時頻率估計誤差很大，甚至大于最簡單的FFT法。文獻［4］和文獻［5］提出修正Rife算法，優點是在最大譜線之間所有區域上估計方差都接近CRLB，缺點是需要進行迭代以提高精度，使運算量大大增加。

上述算法各有優點，但應用于數字接收機卻受其缺點限制，很難兼顧高穩定性、計算量小和精度高這三方面。實際應用中需要的是一種在兩量化頻率之間整個區域都有很高估計精度、計算量小的頻率估計算法。對此本文提出了一種三譜線頻率估計插值算法，充分利用信號FFT計算結果的實部與虛部值，進行內插計算便可獲得很高的估計精度，計算量比Rife法稍大，但頻率估計穩定性優于Rife算法，且兼顧了高估計精度，便于實現和應用。

1? ? ?經典頻域插值測頻算法

應用DFT進行譜估計，硬件實現中可用成熟的FFT降低算法復雜度。目前的頻率估計插值法主要是DFT變換并結合幅度譜進行插值以提高測頻準確度。目前公認高效的頻率估計插值算法包括：經典的雙線幅度法（Rife法）、二次頻率內插法、雙線性對稱頻率內插法、單線相位法等［6］。

單載頻實正弦信號數字化為

[x（n）=asin（2πnf0/fs）]? ? ? ? ? ? ? ?（1）

其中[a]是信號幅度，[fs]是采樣頻率，[f0]是信號頻率。該信號N點的DFT為

[X0（k）=n=0N-1x（n）e-j2πkn/N=a2j1-ej2π（k0-k）1-ej2π（k0-k）/N-1-ej2π（k0+k）1-ej2π（k0+k）/N? ? ? ? ? ? ? ? ? ? ? ? ?（2）]? ? ? ? ? ? 其中，[k0=Nf0/fs]為信號真實頻率位置。由信號的DFT變換得到其幅度譜

[S（k）=X0（k）]? ? ? ? ? ? ? ? ?（3）

然后用幅度譜進行各種插值運算，使測頻精度在直接FFT法基礎上有了不同程度提高。

1.1? ?雙線幅度法

雙線幅度法即Rife方法是一種經典的頻域插值頻率估計算法，它的修正因子為

[Δk=-S（km-1）S（km）+S（km-1）? ? ? ? ? ? ? ? 0S（km+1）S（km）+S（km+1）當S（km-1）>S（km+1）當S（km-1）=S（km+1）當S（km-1）

[f0=（km+Δk）*fs/N]? ? ? ? ? ? ? ? ? ?（5）

再根據式（5）計算估計頻率[f0]。該算法的修正因子算法決定，當[f0]離幅度值最大譜線[km?fs/N]稍遠時，由Rife插值算法估計的[f0]性能較好；缺點是當信號頻率[f0]接近最大譜線[km?fs/N]時，估計誤差很大甚至大于直接FFT后的頻率估計誤差。

1.2? ?二次頻率內插法

該算法假設二次方程[S（k）=ak2+bk+c]是幅度譜的逼近，取[（km-1，S（km-1））][（km，S（km））][（km+1，S（km+1））]三點組成一個方程組求解可以得到

[Δk=-b2a=12*S（km+1）-S（km-1）2S（km）-S（km-1）-S（km+1）]? ?（6）

二次逼近頻率內插法算法復雜度低、運算量小且估計精度較高有一定應用價值。

1.3? ?雙線性對稱頻率內插法

該算法利用了信號FFT后的離散化頻譜之間的時間間隔一般很小的特點，在此基礎上最大譜線兩端次大值譜線間隔時間也接近。由此可以推導出本算法的頻率修正因子為

[Δk=-S（km-1）-S（km+1）S（km）-S（km+1）? ? ? ? ? ? ? ?0S（km+1）-S（km-1）S（km）-S（km-1）當S（km-1）>S（km+1）當S（km-1）=S（km+1）當S（km-1）

再根據式（5）計算估計頻率[f0]。雙線性頻率內插法主要基于近似認為次大譜線關于最大譜線對稱，然后進行頻率修正因子計算。這種插值法在二次頻率內插法優點的基礎上進一步提高了測頻的精度。

1.4? ?單線相位法

該算法適合已知初相的情況，主要是利用最大譜線和提取其DFT后的相位來估計頻率。去掉初相的方法是取兩不同點數的來自原信號序列，然后對這兩個序列做DFT運算，將所得的頻譜最大值相位相減即可去掉相位影響。設[N、M]為兩序列的點數，[km0]、[km1]為[N、M]點FFT后所得頻譜的峰值位置，[β]是兩個序列的相位差。綜合本算法的頻率估計為

[f0=fsN-MN-1Nkm0-M-1Mkm1-βπ]? ?（8）

2? ? ?頻率估計新方法

Rife算法只利用了FFT后頻譜的幅度信息進行內插計算修正因子，并沒有充分利用FFT后頻譜實部與虛部兩部分的信息。先驗信息并未得到充分利用也是Rife算法精度不穩定的原因。本文提出了一種充分利用FFT后頻譜實部和虛部先驗信息的內插算法，計算過程如下：

（1） 對（1）式信號進行N點FFT計算，得到頻譜的復數序列[Xk]。

（2） 對所得的復數頻譜序列[Xk]求絕對值并取其前N/2點，然后找到絕對值最大譜線所對應的序號[ko]，其中0≤[ko]≤N/2。

（3） 計算得到[ko]的修正因子

[δk=-realXko+1-Xko-12Xko-Xko+1-Xko-1]? ? ?（9）

（4） 將[δk]帶入式（5）得到頻率的估計值

[f0=（ko+δk）*fs/N]? ? ? ? ? ? ? ? （10）

式（9）是新算法核心，新算法充分利用了FFT后頻譜峰值的實部與虛部值來進行修正因子的計算，可以充分利用實部和虛部信息不斷修正FFT計算所得的頻率，由于真實頻率在最大譜線和次大譜線之間所有區域對應了不同的頻譜峰值復數值，所以修正因子計算更穩定，保證了算法的穩定性。

3? ? ?算法建模仿真與性能分析

用MATLAB對新算法進行建模和仿真，為了保證仿真結果的穩定性和可靠性，應用了蒙特卡羅方法。仿真中觀測到信號序列的點數為N，信號的采樣率為[fs]，信號頻率[f0=（ko+Δk）*fs/N]，其中[1≤ko≤N2-1]且為整數，[-0.5≤Δk≤0.5]。信號經過DFT變換后對應的最大譜線位置為[ko+Δk]，由于[Δk]在-0.5～0范圍與其在0～0.5范圍對應頻率估計性能的變化曲線關于[Δk=0]對稱，所以只研究后者對應曲線即可。仿真中可以通過固定[ko]并使[Δk]在0～0.5范圍，此時經過算法前兩步得到最大譜線位置為[ko]。取合適的增量使[Δk]從0遞增到0.5，繼續本文算法的后兩步，便可以在仿真中得到信號真實頻率偏離FFT所得頻譜最大譜線范圍為0～[fs/2N]時，新算法得到的修正因子[δk]以及頻率的估計值[f0]。為了量化測頻精度，用歸一化后頻率估計誤差的方差來衡量新算法的估計精度，以下簡稱為歸一化頻偏方差，它的定義式為

[var=Ef0-fo2]? ? ? ? ? ? ? ? （11）

文中以歸一化頻偏方差的大小作為對算法的估計性能優劣的衡量標準。 此外本文要將新算法與經典的Rife算法進行性能比較以論證其是否具有應用價值。用受高斯白噪聲干擾的單載頻正弦波信號作為仿真實驗待估計頻率信號

[x=asin2π（ko+Δk）*fs/N+Noise]? ?（12）

其中[a]為信號幅度，[Noise]為所加高斯白噪聲。

3.1? ?估計性能

由采樣定律可知載波頻偏[fo≤12fs]時，本算法適用；當采樣速率[fs=1T]時， 新算法的估計范圍為[fo≤fs2]。圖1是[Δk]在0～0.5范圍變化時， 歸一化頻偏方差隨[Δk]變化的曲線，4條曲線分別是信噪比為-5 dB、0 dB、5 dB、10 dB時歸一化頻偏方差變換規律。圖1表明在幾種信噪比下，新算法的歸一化頻偏方差隨[Δk]取值不同可基本保持平穩。

為進一步評價新估計方法的估計精度，將新算法歸一化頻偏方差與克拉-美羅界進行比較。載波頻率估計的修正克拉-美羅下界是估計能達到的最小值，對于相位和幅度沒有給出時，正弦信號的歸一化頻率估計誤差方差的下界是：

[CRB（fo）=6f2s2π2N（N2-1）SNR]? ? ? ? ? （13）

式中SNR為信噪比，且滿足[SNR=A2/2δ2]。

圖2表明信噪比低于某個門限值時， 歸一化頻偏方差值急劇增大，估計精度大幅降低。原因是信噪比低于門限時，噪聲淹沒信號頻譜，最大幅值位置隨機變化，使估計誤差非常大。綜合圖2各條曲線隨信噪比變化的規律，可以得出新算法的信噪比門限值為-5 dB，使得新算法可以對低信噪比含噪信號的頻率估計更準確。

3.2? ?與Rife算法性能比較

首先將新算法與Rife算法進行性能穩定性比較。信噪比為5 dB時，新算法和Rife算法的歸一化頻率估計方差隨[Δk]變化曲線如圖3所示。新算法在[0≤Δk≤0.25]范圍內性能有了明顯的改善，解決了Rife算法在這個范圍內估計精度較差的問題。這是因為新算法充分利用最大幅度及其左右的幅度進行修正因子的計算（式（9），不需要對左右幅度值比較，也不必再次用大幅度值對修正因子計算，很大程度上避免了Rife算法當真實頻率在量化頻率附近時最大幅度的誤選，故新算法穩定性比Rife算法好。

其次將新算法與Rife算法進行高精度估計時的信噪比門限比較。由圖3仿真結果及前述分析可知，Rife法在[Δk]取值在0.25～0.5范圍時，估計精度接近CRLB性能，因此有更好的信噪比門限特性。如圖4所示，仿真中取Rife算法估計精度非常高的[Δk]=0.40和估計精度低的[Δk]=0.15進行實驗。

圖4說明Rife算法對噪聲更為敏感：[Δk]=0.40時，Rife算法性能較好，其噪聲門限達到-8 dB，而新算法信噪比門限約為-8.2 dB；[Δk]=0.15時Rife算法性能較差，信噪比門限約為-7 dB，而新算法信噪比門限約為-9 dB，盡管這種條件下Rife算法信噪比門限接近本文算法，但是其估計性能遠低于新算法。綜合比較上述兩種條件下仿真結果，新算法的工作信噪比門限比Rife算法好。

最后將新算法與Rife算法進行運算量比較。通過表1得出對信號進行FFT運算由位數N決定，所以新算法的運算量僅比Rife算法多1次復乘和3次復加運算。而M-Rife算法相比前兩者增加了大量復乘以及復加運算，其高精度的頻率估計建立在很高的計算復雜度上。綜合來看本文算法性能優于經典Rife算法并且幾乎沒有增加運算量，與M-Rife算法相比運算復雜度低，在一般情況下本文算法更容易硬件實現。

4? ? ?結語

本文提出了一種高精度頻率插值估計新算法，仿真實驗表明，估計精度非常高近修正克拉-美羅限，并且穩定性比Rife算法好。此外，新算法信噪比門限比Rife算法低3 dB左右，當信噪比不理想時能更精確測得頻率值。與此同時，本文算法的計算量與Rife算法計算量相當，較當前性能研究很多的M-Rife算法，新算法的運算量大大降低。新算法的關鍵是快速傅里葉變換求頻譜，在硬件上非常容易實現，因此新方法可以在通信領域數字接收機的信號處理中獲得應用并提高其解析和處理信號的能力。

［參考文獻］

［1］ Rife D C′ ，Boorstyn R R．Single-tone parameter estimation from discrete-time observation［J］. IEEE Trans Inform Theory，1974，20（5）：591-598．

［2］ Rife D C′，Boorstyn P R．Multiple tone parameter estimation from discrete-time observation［J］.Bell Syst Tech J，1976，55（9） ： 1389-1410．

［3］ Rife D C′ ，Vincent G A.Use of the discrete Fourier transform in the measurement of frequencies and levels of tones［J］.Bell Syst Tech J，1970，49（2） ： 197-228．

［4］ 鄧振淼，劉渝，王志忠.正弦波頻率估計的修正Rife算法［J］.數據采集與處理，2006，21（4）：473-477．

［5］ 王旭東，劉渝，鄧振淼.基于修正Rife 算法的正弦波頻率估計及FPGA 實現［J］. 系統工程與電子技術，2008，30（4）： 621-624．

［6］? 張昌菊，唐斌.單頻信號快速頻率估計算法比較及改進［J］.電訊技術，2005，1（5）：72-76．

［7］ 孫宏軍，王小威.基于幅值-相角判據的修正Rife正弦波頻率估計算法［J］.天津大學學報（自然科學與工程技術版），2018，51（8）：810-816.

［8］ 周龍健，羅景青，房明星.基于IIN算法和Rife算法的正弦波頻率估計算法［J］.數據采集與處理，2013，28（6）：839-842.

［9］ 趙強，侯孝民，廉昕.基于Zoom-FFT的改進Rife正弦波頻率估計算法［J］.數據采集與處理，2017，32（4）：731-736.

［10］ 王曉峰，周瑞，田潤瀾，等.用于正弦波頻率估計的修正Kay算法［J］.電訊技術，2017，57（5）：518-522.

［11］ 崔誦祺，王愛華，湯樂，等.基于雙門限判決修正Rife算法的正弦信號頻率估計［J］.北京理工大學學報，2017，37（3）：292-296+303.

［12］ 陳役濤，劉渝，鄧振淼.基于相關累加的正弦波頻率估計算法［J］.數據采集與處理，2008，23（6）：729-733.

［13］ 陳毅濱.一種低信噪比下正弦波頻率估計方法［J］.通信技術，2021，54（10）：2316-2320.

［14］ 王宏偉，趙國慶，齊飛林.一種實時精確的正弦波頻率估計算法［J］.數據采集與處理，2009，24（2）：208-211.

［15］ 鄧振淼，劉渝.正弦波頻率估計的牛頓迭代方法初始值研究［J］.電子學報，2007（1）：104-107.

［16］ 周喜慶，趙國慶，王偉.實時準確正弦波頻率估計綜合算法［J］.西安電子科技大學學報，2004（5）：657-660+681.