一種新的正弦波頻率估計算法

張剛兵，錢顯毅

（常州工學(xué)院 電子信息與電氣工程學(xué)院，江蘇 常州 231002）

一種新的正弦波頻率估計算法

張剛兵＊，錢顯毅

（常州工學(xué)院 電子信息與電氣工程學(xué)院，江蘇 常州 231002）

研究了高斯白噪聲條件下大樣本點單一正弦波信號的頻率估計方法.首先利用離散傅里葉變換確定頻率粗估計，然后以該值為參考頻率構(gòu)造本地信號，將原信號下變頻至基帶，對基帶信號分段求和，能得到一個新的正弦波信號，該信號的頻率為真實頻率與參考頻率之差，最后利用最小二乘法估計新信號的載頻，修正粗估計值就能得到原信號頻率的最優(yōu)估計.推導(dǎo)了算法的漸近方差與克拉美－羅限之間的關(guān)系.仿真結(jié)果表明，本算法能適用于整個頻段范圍，頻率估計的精度接近正弦波頻率估計的克拉美－羅限.

大樣本；離散傅里葉變換；相關(guān)積累；頻率估計

正弦波頻率估計算法廣泛應(yīng)用于雷達、通信以及電子對抗等信號處理領(lǐng)域.Rife［1］最先提出了被加性復(fù)高斯白噪聲污染的正弦波信號頻率估計算法——最大似然（Maximum Likelihood，ML）估計，雖然其性能接近正弦波頻率估計的克拉美－羅限（Cramer－Rao Lower Bound，CRLB），但需要進行一維頻率搜索，計算量太大，不利于工程采納.多位學(xué)者對正弦波信號頻率估計問題作了進一步研究，相繼提出了多種頻率估計算法.這些算法大致可分為兩類，一類是基于譜線插值，另一類是基于相位信息.

文獻［2］利用信號的最大兩根譜線插值進行頻率估計，即Rife算法.當信號頻率位于量化頻率附近時，由于插值方向錯誤，會導(dǎo)致頻率估計性能下降.針對這一問題，文獻［3］先對信號進行頻移，使新信號的頻率位于兩個相鄰量化頻率點的中心區(qū)域，提出了修正Rife算法，再以該算法的估計值為初始值進行一次牛頓迭代（Sinusoid Frequency Estimation Based on Newton's Method，SFENM）.當初始值位于收斂區(qū)域時，迭代收斂，其性能穩(wěn)定，否則會導(dǎo)致頻率估計精度下降.文獻［4］提出了基于DFT相位的正弦波頻率和初相估計方法（Based on Phase of DFT，BPDFT），利用分段DFT頻譜的相位差消除了初相對頻率估計的影響且避免了相位測量模糊問題.文獻［5］指出，當信號真實頻率與DFT量化頻率差為某一范圍時Rife算法精度不高，并研究了噪聲背景中插值FFT估計正弦信號頻率估計的問題.文獻［6］利用與最大譜線對應(yīng)的量化頻率點相差半個量化頻率的兩根譜線進行插值，提出了一種迭代插值（（Iterative Frequency Estimation by Interpolation on Fourier Coefficients，IFEIFC））頻率估計算法，性能接近克拉美－羅限.

Tretter［7］將加性復(fù)測量噪聲等效為加性實相位噪聲，利用最小二乘法對展開的相位估計信號的頻率.高斯噪聲條件下的最小二乘估計等價于最大似然估計，因此，文獻［7］的算法在高斯噪聲條件下的性能接近正弦信號頻率估計的克拉美－羅限，但該算法模型成立的信噪比門限約為10dB，降低了其在低信噪比條件下的性能.為了避免相位展開，文獻［8］利用相位平均和相位加權(quán)平均提出了基于相位差分的頻率估計算法，在高信噪比條件下的性能與最小二乘法相當，也接近克拉美－羅限，但該算法也不適用于信噪比低于10dB的場合.文獻［9－11］相繼提出了3種基于信號自相關(guān)函數(shù)的頻率估計算法，雖然改善了信噪比門限附近的性能，使之不會急劇惡化，但高信噪比條件下的性能卻無法接近克拉美－羅限，甚至還會縮小頻率估計范圍［11］.文獻［12］反復(fù)利用低通濾波、抽取，線性預(yù)測以及數(shù)字差拍變頻進行頻率估計（Iterative Lin－ear prediction，ILP），大大降低了算法的信噪比門限，信噪比門限以上的性能最多超過克拉美－羅限0.7dB.文獻［13］多次利用自相關(guān)函數(shù)進行頻率估計（Autocorrelation－based algorithm，AA），當信號樣本長度為1 024時，其頻率估計的均方根誤差比克拉美－羅限大0.14dB.

當信號樣本點數(shù)較大時，受硬件條件的限制，難以一次性對全部接收信號進行處理.如果先將信號分成等長的L段，再對每段信號分別進行頻率估計，最后對各頻率估計值取算術(shù)平均，那么只能得到非相干頻率估計值，其估計方差僅為各段方差的L分之一，與整段信號頻率估計的克拉美－羅限（為各段頻率估計方差的L3分之一）相差甚遠，此時的估計精度可能難以滿足系統(tǒng)的性能指標.為了能估計大樣本點信號的頻率參數(shù)，同時保證參數(shù)估計的精度滿足系統(tǒng)設(shè)計要求，必須研究新的頻率估計方法.

受文獻［3－6］的啟發(fā)，本文利用離散傅里葉變換（Discrete Fourier Transform，DFT）的信噪比增益作用降低信噪比門限，結(jié)合高斯噪聲條件下最小二乘等價最大似然的特性，提出了一種正弦波頻率估計算法.首先利用離散傅里葉變換得到頻率粗估計值，然后以該值作為參考頻率將原信號下變頻至基帶，分段求和后能得到一個新的正弦波信號.以該正弦波的頻率修正參考頻率就能得到原信號的頻率估計值.與基于相位信息的算法相比，本算法的信噪比門限更低、均方誤差更接近克拉美－羅限，與基于譜線插值的算法相比，計算效率更高.

1 算法介紹

設(shè)信號模型為

式中，A為信號的幅度，φ0為初始相位，f為信號的頻率，ts為采樣間隔，N為信號樣本點數(shù).ε是均值為0、方差為σ2的復(fù)高斯白噪聲，定義信號的信噪比（signal－to－noise ratio，SNR）為SNR＝A2／σ2.

對式（1）定義的正弦波采樣序列x（n），假定已得到信號頻率f的粗估計值＾f0，現(xiàn)構(gòu)造序列

將（1）式和（2）式共軛相乘

式中，ε′是均值為0、方差為σ2的復(fù)高斯白噪聲，Δf0＝f－為粗估計的誤差.

將式（3）表示的信號z（n）按M點分為一組并求和，能得到一個點數(shù)為L＝N／M的新序列，假設(shè)N／M為整數(shù)，有

式（4）中的S（m）是一個頻率為Δf0、采樣間隔為Mts、樣本點數(shù)為L點的正弦波序列，對S（m）進行頻率估計可得到頻偏Δf0的估計值Δ＾f0，修正粗估計能得到信號頻率f的精確估計.無模糊估計頻率要求頻偏Δf0滿足奈奎斯特采樣定理，即

序列S（m）的信噪比SNRout?MA2／σ2＝M×SNR，與正弦波采樣序列x（n）相比，此時的信噪比增益提高了10lg MdB.累加點數(shù)M越多，信噪比增益越大，頻偏Δf0估計的信噪比門限越低.Tretter在文獻［7］中指出，當復(fù)正弦序列的輸入信噪比遠遠大于1時，加性復(fù)測量噪聲可近似為等效的加性實相位噪聲，且相位噪聲的方差為正弦波輸入信噪比倒數(shù)的一半.當原測量信號信噪比不太小、信號樣本點數(shù)較大時，經(jīng)過求和積累后，式（4）的信噪比SNRout有可能遠遠大于1.因此，式（4）可近似為

式中，ζm是均值為0、方差為1／（2SNRout）的實高斯白噪聲.

高斯白噪聲條件下的最小二乘估計等價于最大似然估計，如果能得到序列S（m）無模糊的相位值φm，就能用最小二乘估計得到頻偏Δf0的最優(yōu)估計.設(shè)

式（8）中X的最小二乘估計為

利用頻偏的估計值修正頻率粗估計就能得到信號頻率的精確估計，綜合以上分析，本算法的實現(xiàn)步驟為

2）按照公式（2）、（3）、（4）得到序列S（m），m＝0，1，…，L－1；

3）利用文獻［7］的算法對S（m）進行載頻估計，得到頻偏Δf0的估計值Δ

4）利用頻偏Δf0的估計值得到頻率估計值

2 性能分析

現(xiàn)從漸近方差（Asymptotic Error Variance，AEV）以及計算量對本算法進行定量分析.

2.1 漸近方差（AEV）

利用最大譜線對應(yīng)的量化頻率點作為頻率粗估計，在不出現(xiàn)頻率模糊的條件下，粗估計誤差滿足

對式（3）表示的信號z（n）進行M 點累加以后，新序列的信噪比增加了10lg（M）dB，M 越大，信噪比增益越大，影響頻偏估計的信噪比門限就越低.累加之后要求信號S（m）有兩個以上的樣本點，M 的取值必須滿足M ≤N／2，利用式（10），有

因此，累加之后能滿足式（5）對無模糊頻率估計的要求.

頻率估計的精度由頻偏Δf0的估計精度決定，而頻偏估計是無偏估計，估計的均方誤差為［7］

式中，E（·）表示取數(shù)學(xué)期望.

正弦波信號頻率估計＾f的克拉美 －羅限（CRLB）為［1］

由式（12）和式（13）有

從式可以看出，累加點數(shù)M越大，性能損失越大，當M＝N／2時，性能下降最嚴重，此時的性能與克拉美－羅限相差1.25dB.但當M分別為N／4和N／8時，累加求和后的信號樣本分別為4點和8點，性能相對克拉美－羅限分別下降0.28dB和0.07dB.在M?N的條件下，本算法的均方誤差與克拉美－羅限相等.

2.2 計算量分析

現(xiàn)分析本算法需要的計算量，假設(shè)N是2的整數(shù)次冪.在頻率粗估計時利用快速傅里葉變換（Fast Fourier Transform，F(xiàn)FT）確定信號的頻譜，需要N／2·log2N次復(fù)數(shù)乘法、Nlog2N次復(fù)數(shù)加法，確定最大譜線的位置需要N次復(fù)數(shù)乘法，式（4）需要N次復(fù)數(shù)乘法，式需要N（M－1）／M復(fù)數(shù)加法，獲得相位需要N／M次反正切運算、N／M次相位展開，估計頻偏則需要N／M次實數(shù)乘法和N／M－1次實數(shù)加法，作一次迭代需要N次復(fù)數(shù)乘法、N（M－1）／M復(fù)數(shù)加法、N／M 次反正切運算、N／M 次相位展開、N／M次實數(shù)乘法和N／M－1次實數(shù)加法.

3 仿真結(jié)果分析

為了驗證本算法的測頻性能，對其進行計算機仿真，信號樣本點數(shù)為1 024，采樣頻率為100MHz，信噪比為0dB，累加點數(shù)M為128，信號頻率從25MHz開始，到25MHz＋1／（2 NΔt）結(jié)束，將頻率變化范圍等分20份，每個離散頻率點上各仿真2 000次.其均方根誤差與克拉美 －羅限之比隨頻率變化的性能曲線如圖1所示.

圖1的仿真結(jié)果表明，基于自相關(guān)函數(shù)相位信息的頻率估計算法的均方根誤差最大，基于譜線迭代插值和以修正Rife算法估計值為初值進行一次牛頓迭代的頻率估計算法性能相當，都能在整個頻段范圍內(nèi)進行頻率估計.

下面驗證本算法在不同信噪比條件下的性能.采樣頻率為100MHz，信號的頻率為f＝21.123 42MHz，信號樣本點數(shù)N為1 024，累加點數(shù)M為128.對本算法的性能進行仿真，給出了各算法在相同條件下的性能對比，每種條件下各進行2 000次獨立的仿真.均方根誤差（RMSE）與信噪比（SNR）之間的關(guān)系如圖2所示.圖中給出了各算法的均方根誤差及克拉美－羅限，圖中橫軸為線性坐標，縱軸為對數(shù)坐標.

圖1 測頻范圍性能對比Fig.1 Comparison of frequency range

圖2 正弦波頻率估計的性能曲線Fig.2 Performance comparison of frequency estimation

圖2的仿真結(jié)果表明，本算法與文獻［6］算法一樣，具有相同的信噪比門限，且都低于文獻［3－4－12－13］算法.因此，本算法對信噪比要求更低，能估計較低信噪比條件下信號的頻率.在較低信噪比條件下，牛頓迭代的初始值位于收斂區(qū)域之外，迭代不收斂，此時的估計性能不如本算法.隨著信噪比的增加，初始值進入收斂區(qū)域，迭代后的性能接近克拉美－羅限.文獻［12－13］都是基于迭代自相關(guān)函數(shù)估計信號的頻率，其信噪比門限均高于基于DFT的頻率估計算法.在整個信噪比變化區(qū)間內(nèi)，本算法的均方根誤差都接近克拉美－羅限，表明本算法性能穩(wěn)定，對信噪比變化不敏感.

4 結(jié)束語

本文提出了一種適用于單一信號的頻率估計算法，本算法具有較低的信噪比門限，其性能均勻分布在整個測頻范圍內(nèi).當累加求和后的信號樣本點數(shù)為8時，頻率估計的均方誤差相對克拉美－羅限下降0.07dB.當累加后的信號樣本點數(shù)遠小于原測量信號樣本點數(shù)時，本算法的均方誤差與克拉美－羅限相同.考慮信噪比門限、測頻范圍以及估計性能，本算法都優(yōu)于基于相位和自相關(guān)函數(shù)的頻率估計器，因此，本算法具有一定的工程應(yīng)用價值.

［1］Rife D C，Boorstyn R R.Single－tone parameter estimation from discrete－time observations［J］.IEEE Trans on Information Theory，1974，20（5）：591－598.

［2］Rife D C，Bowstyn R R.Multiple tone parameter estimation from discrete rime observation［J］.Bell Syst Tech J，1976，55（9）：1389－1410.

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

［4］齊國清，賈欣樂.基于DFT相位的正弦波頻率和初相的高精度估計方法［J］.電子學(xué)報，2001，29（9）：1164－1167.

［5］齊國清，賈欣樂.插值FFF估計正弦信號頻率的精度分析［J］.電子學(xué)報，2004，32（4）：625－629.

［6］Aboutanios E，Mulgrew B.Iterative frequency estimation by interpolation on Fourier coefficients［J］.IEEE Trans on Signal Processing，2005，53（4）：1237－1242.

［7］Tretter S.Estimating the frequency of a noisy sinusoid by linear regression［J］.IEEE Trans on Information Theory，1985，31（6）：832－835.

［8］Kay S M.A fast and accurate single frequency estimator［J］.IEEE Trans on Acoustics，Speech，and Signal Processing，1989，37（12）：1987－1990.

［9］Fitz M P.Further results in the fast estimation of a single frequency［J］.IEEE Trans Communications，1994，42：862－864.

［10］Luise M，Reggiannini R.Carrier frequency recovery in alldigital modems for burst－mode transmissions［J］.IEEE Trans Communications，1995，43：1169－1178.

［11］Fowler M L，Johnson J A.Extending the threshold and frequency range for phase－based frequency estimation［J］.IEEE Trans on Signal Processing，1999，47（10）：2857－2863.

［12］Brown T，Wang M M.An iterative algorithm for singlefrequency estimation［J］.IEEE Trans on Signal Processing，2002，50（11）：2671－2682.

［13］Xiao Yang－Can，Wei Ping，Tai Heng－Ming.Autocorrelation－based algorithm for single－frequency estimation［J］.Signal Processing，2007：1224－1233.

Novel frequency estimation method for sinusoidal signal

ZHANG Gangbing，QIAN Xianyi
（Changzhou Institute of Technology，School of Electronic Information and Electric Engineering，Changzhou，Jiangsu 231002）

A novel algorithm for frequency estimation of sinusoid with many samples in the complex white Gaussian noise was proposed.A coarse frequency estimate was obtained through discrete Fourier transform and a sequence was constructed with the coarse frequency estimate being a reference frequency.Then the original sequence was converted into baseband by down－conversion.A new sinusoid was acquired after correlation accumulation，of which the frequency was the difference between the original frequency and reference one.Linear regression made it possible to get the optimal frequency difference estimate.The relationship between the asymptotic error variance（AEV）and the Cramer－Rao lower bound（CRLB）was derived.Simulation results show that the performance of the proposed algorithm approaches the CRLB of the sinusoidal signal when the signal－to－noise ratio（SNR）is higher than the SNR threshold.

large sample；DFT；correlation accumulation；frequency estimation

TN911

A

1000－1190（2012）01－0040－05

2011－07－10.

江蘇省高校自然科學(xué)基金項目（10KJD480003）；常州工學(xué)院校級基金項目（YN1101）.

＊E－mail：caeppub＠sohu.com.