基于改進(jìn)的Rife測頻算法及其FPGA實(shí)現(xiàn)*

司偉建，郝 鑫，趙忠凱，陳 濤

（哈爾濱工程大學(xué)信息與通信工程學(xué)院，哈爾濱 150001）

0 引言

正弦波頻率估計(jì)是信號參數(shù)估計(jì)中的經(jīng)典問題，也是電子對抗領(lǐng)域研究的重點(diǎn)課題之一。文獻(xiàn)[1]提出了利用DFT進(jìn)行頻率估計(jì)的同時在時域加類sinc窗函數(shù)的頻譜校正方法，這種方法加寬了窗譜寬度，使得兩條相鄰DFT譜線產(chǎn)生測頻模糊，因此需要DFT點(diǎn)數(shù)較大才可以保證測頻精度。文獻(xiàn)[2]利用對相位長短延時的方法，短的時間延遲可以用于消除頻率模糊，長的延遲用于獲得較高的頻率分辨率的方法來提高測頻精度，但是這種方法受噪聲影響較大，在頻帶邊緣處測頻誤差較大。文獻(xiàn)[3]給出了最大似然頻率估計(jì)算法，該算法能夠達(dá)到頻率估計(jì)的下限：克拉美－羅限（CRB），但該算法計(jì)算量大，不利于工程實(shí)現(xiàn)。文中在分析Rife算法的基礎(chǔ)上，從FPGA實(shí)現(xiàn)的角度出發(fā)，提出了一種改進(jìn)的Rife算法，仿真和實(shí)測結(jié)果均表明該算法能在短數(shù)據(jù)條件下實(shí)現(xiàn)對信號頻率的快速、精確測量，具有較高的工程應(yīng)用價值。

1 Rife算法分析

文獻(xiàn)[4－5]采用最大譜線及其相鄰的次大譜線進(jìn)行插值，用插值的結(jié)果進(jìn)行頻率估計(jì)，稱這種算法為Rife算法。Rife算法的原理如下：

設(shè)接收的信號并被采樣后為：s（n）＝x（n）＋w（n），其中x（n）為接收到的純信號，w（n）為混入信號中的高斯白噪聲序列，其均值為0，方差為σ2。設(shè)x（n）＝aej（2πfc／fs＋φ0），其中a為信號幅度；fc為接收信號頻率；φ0為信號的初相；fs為信號采樣頻率。將采樣后的信號s（n）做N點(diǎn)FFT運(yùn)算，最大譜線值記為｜X（k0）｜，次大譜線值記為｜X（k0＋r）｜。則Rife算法所得到的頻率估計(jì)式為：

圖1 Rife算法示意圖

其中，當(dāng)｜X（k0＋1）｜≥｜X（k0－1）｜時，r＝1；相反，當(dāng)｜X（k0＋1）｜＜｜X（k0－1）｜時，r＝－1。圖1為r＝1時的情況。 由于0≤所以插值的偏移量在[k0－1／2，k0＋1／2]范圍內(nèi)。將偏移量與r的值代入式（1）即可估計(jì)出輸入信號的頻率值。

2 改進(jìn)Rife算法分析

Rife算法通過判斷FFT峰值兩側(cè)譜線幅度的大小來確定插值的方向，通過計(jì)算來選取偏移量大小。由此可見利用FPGA實(shí)現(xiàn)時，需要進(jìn)行除法運(yùn)算，而除法運(yùn)算需要占用較多的FPGA邏輯單元，運(yùn)算周期長，限制了測頻速度。

經(jīng)過仿真發(fā)現(xiàn)：將輸入信號加窗后，F(xiàn)FT 最大值兩側(cè)的兩個次大值的功率差G＝｜G（k0＋1）｜－｜G（k0－1）｜（其中G代表功率差）與在[k0－1／2，k0＋1／2]偏移量保持線性關(guān)系（其中k0為FFT幅度最大值所對應(yīng)的下標(biāo)索引值），如圖2所示。因此只要求出任意一段的線性斜率a，并且知道FFT的最大值下標(biāo)索引值k0，利用功率差值G就可以由式（2）確定出輸入信號的頻率。

圖2 FFT下標(biāo)索引值與功率差值關(guān)系圖

為了實(shí)現(xiàn)上述算法，首先需要求出斜率a的大小，在具體實(shí)現(xiàn)中可以利用Matlab計(jì)算出其結(jié)果，利用其結(jié)果便可以在硬件中實(shí)現(xiàn)改進(jìn)Rife算法測頻。具體步驟如下：

1）對采樣序列進(jìn)行加窗，加窗的長度與做FFT的點(diǎn)數(shù)相同；

2）加窗后的序列做FFT運(yùn)算，取幅度峰值｜X（k0）｜兩側(cè)的兩個值｜X（k0＋1）｜和｜X（k0－1）｜；

3）將｜X（k0＋1）｜和｜X（k0－1）｜轉(zhuǎn)換成對應(yīng)的功率值，并求出功率差：G＝｜G（k0＋1）｜－｜G（k0－1）｜，記錄G與k0；

4）在待測輸入信號頻率區(qū)間內(nèi)，輸入信號每次按照相同的步進(jìn)頻率重復(fù)1）～3），直至結(jié)束；

5）選取任意一段處于[k0－1／2，k0＋1／2]之間的頻率部分，根據(jù)3）記錄的G與k0的關(guān)系做出圖2，并求出該段斜率的倒數(shù)，即a值的大?。?/p>

6）利用5）求出的a和FFT的幅度峰值的下標(biāo)索引值k0，以及功率差G，代入式（2）便可以求出輸入信號的頻率值。

3 仿真分析

設(shè)置仿真參數(shù)為：輸入信號頻率范圍：fc＝227.5＋Δfc，Δfc取[－fs／2 N，fs／2 N]頻率間隔為0.1MHz，采樣頻率為240MHz，采樣點(diǎn)數(shù)為64點(diǎn)，信號初始相位φ0取[0，2π]隨機(jī)數(shù)，噪聲均值為0，方差σ2＝0.5。對每個頻率點(diǎn)做1000次Monte Caro仿真實(shí)驗(yàn)，當(dāng)輸入信號為實(shí)信號時，對應(yīng)的克拉美－羅限（CRB）由式（3）確定[7]。

圖3為SNR＝6dB時Rife算法和改進(jìn)Rife算法的測頻誤差仿真結(jié)果，可以看出在｜Δfc｜≤1／3時，Rife算法的測頻誤差較大，最大均方誤差為0.7MHz，而此時改進(jìn)的Rife算法的測頻誤差較小。相反當(dāng)｜Δfc｜＞1／3時，改進(jìn)Rife算法的測頻誤差變大，最大均方誤差為0.5MHz。

圖4為不同輸入信號信噪比條件下，改進(jìn)Rife算法和Rife算法的測頻性能仿真結(jié)果。其中SNR變化范圍為：0～20dB，每個信噪比下做1000次Monte－Caro試驗(yàn)，每次試驗(yàn)信號輸入頻率為fc＝277.5＋Δfc（Δfc取[－fs／2 N，fs／2 N]內(nèi)的隨機(jī)數(shù)）。由圖可以看出在信噪比變大時，Rife算法與改進(jìn)Rife算法的測頻精度都在提高，但改進(jìn)Rife算法的測頻誤差始終小于Rife算法。

圖3 SNR＝6dB時Rife算法和改進(jìn)Rife算法測頻誤差比較

圖4 不同輸入信噪比下Rife算法與改進(jìn)Rife算法測頻誤差比較

4 改進(jìn)Rife算法的FPGA實(shí)現(xiàn)

將改進(jìn)的Rife算法在Altera公司的EP3SL200F1152I4FPGA中進(jìn)行了實(shí)現(xiàn)。選用的窗函數(shù)為海明窗，長度為64點(diǎn)，量化后存入ROM中。利用計(jì)數(shù)器產(chǎn)生讀ROM地址，并在每個時鐘上升沿將ROM中所對應(yīng)的窗函數(shù)順次讀出。在Quartus II中調(diào)用乘法器來完成輸入信號與窗函數(shù)的乘積運(yùn)算。結(jié)果送給FFT IP核，F(xiàn)FT點(diǎn)數(shù)為64點(diǎn)，輸出的結(jié)果轉(zhuǎn)換成幅度并查表轉(zhuǎn)換成功率，取出功率峰值處的下標(biāo)索引k0，并將索引值為k0－1和k0＋1處的功率值調(diào)用減法器做差，得到功率的差值。由功率差值做地址在FPGA中進(jìn)行查表求出在一個區(qū)間內(nèi)的頻率偏移量，同時利用功率峰值下標(biāo)索引k0及斜率a就可以計(jì)算出輸入信號的頻率。

除了主要模塊部分外，還包括控制FFT模塊進(jìn)數(shù)時序以及延時控制等部分。本次選取窗函數(shù)量化位數(shù)為11位，F(xiàn)FT模塊輸入方式為streaming模式。

4.1 實(shí)測數(shù)據(jù)與仿真對比結(jié)果

信號源采用Agilent E4438C，利用噪聲源NC6110A對輸入信號加入相應(yīng)功率的噪聲，ADC采樣速度為240MHz。輸入信號頻率為277.5 MHz，輸入信噪比SNR分別為0dB、3dB、6dB、9dB，脈沖寬度為2μs，對脈沖前沿樣64點(diǎn)進(jìn)行頻率估計(jì)。利用signaltapII采樣4K輸出結(jié)果，利用Matlab軟件計(jì)算其均方誤差，實(shí)驗(yàn)結(jié)果如表1所示。由表1可以看出，實(shí)測均方誤差稍大于理論均方誤差，其原因是：對頻率進(jìn)行量化會引入量化誤差，同時ADC采樣不理想也會引入一定的誤差，對實(shí)驗(yàn)結(jié)果產(chǎn)生一定的影響。

表1 實(shí)測結(jié)果

4.2 FPGA資源消耗

在Quartus II中建立工程，選取FPGA芯片為EP3SL200F1152I4，經(jīng)過綜合并布線后消耗總資源為4%，各部分資源消耗情況如表2所示。

表2 改進(jìn)的Rife算法資源消耗表

經(jīng)過實(shí)測，單獨(dú)編譯一個64點(diǎn)FFT IP核所消耗的資源為2%，可見改進(jìn)的Rife算法占用資源很少。

5 結(jié)論

針對Rife算法測量誤差大、硬件實(shí)現(xiàn)速度慢的特點(diǎn)，提出了一種改進(jìn)的Rife測頻算法。從工程應(yīng)用的角度，對算法FPGA實(shí)現(xiàn)的關(guān)鍵步驟進(jìn)行了詳細(xì)闡述。該算法能夠?qū)崿F(xiàn)對低信噪比、短數(shù)據(jù)條件下的未知信號進(jìn)行快速、準(zhǔn)確測頻，在0dB時測頻誤差小于0.6MHz／σ，處理速度可達(dá)240MHz，適合應(yīng)用于電子對抗、偵察等領(lǐng)域。

[1]王曉君，陳禾，仲順安.改進(jìn)的脈沖雷達(dá)信號加窗DFT頻譜校正算法[J].北京理工大學(xué)學(xué)報(bào)，2008 ，28（2）：164－167.

[2]趙亮，程少云.寬帶數(shù)字接收機(jī)測頻方法的研究[J].微計(jì)算機(jī)信息，2008，24（19）：292－294.

[3]Rife D C，Boorstyn R R.Single－tone parameter estimation from discrete time observations[J].IEEE Trans.Inform.Theory，1974，20（5）：591－598.

[4]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：197－228.

[5]王旭東，劉渝，鄧振淼.基于修正Rife算法的正弦波頻率估計(jì)及FPGA實(shí)現(xiàn)[J].系統(tǒng)工程與電子技術(shù)，2008，30（4）：621－624.

[6]Altera.FFT MegaCore Function v8.0user Guide[EB／OL].http：／／www.altera.com.

[7]齊國清.離散實(shí)正弦信號參數(shù)估計(jì)的Cramer Rao方差下限[J].數(shù)據(jù)采集與處理，2003，18（2）：151－155.