FFT幅相聯(lián)合的快速高精度頻率估計方法

馬羅文，劉寧，胡心語，倪麗花

（1.中國人民解放軍95865部隊，北京 102200；2.中國人民解放軍32802部隊，北京 100082；3.電子科技大學(xué)，四川 成都 611731）

0 引言

在雷達［1-3］、通信［4-5］、聲吶［6］等領(lǐng)域?qū)π盘栴l率進行快速、準(zhǔn)確的測量具有重要意義。目前，常用的頻率估計算法按照工作原理可以分為：計數(shù)法、相位法、快速傅里葉變換算法、子空間分解類算法等。其中，單頻信號［7］的估計問題一般采用計數(shù)法和相位法；而多頻信號估計問題一般采用的方法為子空間分解類算法和快速傅里葉變換算法。子空間分解類算法［8-9］中典型代表為多重信號分類（Multiple Signal Classification，MUSIC）算法［10］和旋轉(zhuǎn)不變技術(shù)估計信號參數(shù)（Estimating Signal Parameter Via Rotational Invariance Techniques，ESPRIT）算法［11］，此類子空間分解算法具有較高的估計精度，但由于需要搜索幅度譜峰值位置或特征值分解導(dǎo)致計算量相對較大。快速傅里葉變換算法［12］只能對有限點序列進行處理，會造成能量泄露，同時存在柵欄效應(yīng)，當(dāng)信號頻譜峰值譜線不能與主瓣的中心重合時，信號真實頻率落于主瓣內(nèi)的兩根譜線之間，便會產(chǎn)生較大的頻率測量誤差［13-16］。

由于快速傅里葉變換的幅度譜頻率搜索間隔與信號采樣數(shù)或快速傅里葉變換的長度成反比，為了降低快速傅里葉變換法幅度譜頻率搜索間隔的量化誤差對頻率測量精度的影響，常采用填零的方法增加信號采樣數(shù)或快速傅里葉變換的長度。但是，快速傅里葉變換的計算量與信號采樣數(shù)或快速傅里葉變換長度成正比，這種提高測頻精度的方法也會使頻率測量的計算量成倍地增加，導(dǎo)致其難以應(yīng)用于實時性要求較高的場景。

高精度和實時性要求的矛盾一直都是頻率估計在高速全數(shù)字解調(diào)器中應(yīng)用的瓶頸。因此，研究一種既保證高精度，又滿足低計算量的頻率估計算法，具有重要的應(yīng)用前景［17］。本文提出同時利用信號采樣的快速傅里葉變換的幅度譜峰值位置和頻率搜索間隔量化誤差校正值確定信號頻率的方法，能實現(xiàn)高精度快速測頻的目的。

1 信號模型

暫不考慮噪聲時，單頻復(fù)諧波信號表示為

式中：A為復(fù)正信號的幅度；f0為頻率，Hz；θ0為初相，rad。在不做特殊說明的情況下，本文基于該信號展開算法理論推導(dǎo)與實驗分析。

經(jīng)過N點采樣后的離散信號表示為

對x(n)進行N點快速傅里葉變換，得到

幅度最大值處的離散頻率索引值為m，根據(jù)X[k]幅度最大值處的位置m可以得到信號頻率的估計值為由于柵欄效應(yīng)的存在，當(dāng)信號頻譜峰值譜線不能與主瓣的中心重合時，信號真實頻率落于主瓣內(nèi)的兩根譜線之間，便會產(chǎn)生較大的頻率測量誤差。

因此，本文利用信號相位信息，消除了傳統(tǒng)方法僅使用幅度信息受到的搜索間隔的限制，提出了FFT幅相聯(lián)合的快速高精度頻率估計方法。

2 FFT幅相聯(lián)合的高精度頻率估計方法

首先對前N-1個信號采樣x(n)，n=0，1，2，…，N-2，采用快速傅里葉變換，得到頻譜序列y1(k)，為

然后經(jīng)過一位滑窗平移得到后N-1個信號采樣x(n)，n=1，2，…，N-1，并對其采用快速傅里葉變換，得到頻譜序列y2(k)，為

由于頻譜序列y1(k)和y2(k)是對具有相同頻率的采樣信號進行FFT運算得到的，那么y1(k)和y2(k)幅度譜的峰值位置相同。根據(jù)頻譜序列y1(k)和y2(k)，得到幅度譜z(k)為

根據(jù)式（6），信號頻率估計值所對應(yīng)的序列號，即幅度譜的譜峰位置m為

根據(jù)式（3），式（4）和式（5）在幅度譜的譜峰位置m處的比值為

FFT幅相聯(lián)合的高精度頻率估計f0為

式中：angle表示復(fù)數(shù)的相位角。

由式（9）可知，利用了頻譜序列y1(m)和y2(m)的相位差對頻率進行估計，頻率的估計結(jié)果不再受限于頻率搜索間隔的量化誤差，從而實現(xiàn)高精度頻率估計。

3 FFT幅相聯(lián)合的快速高精度頻率估計方法

根據(jù)式（4）和式（5）可知，為了求解y1(k)與y2(k)，需要進行兩次N-1點FFT變換，對此，提出FFT幅相聯(lián)合的快速高精度頻率估計方法，利用迭代，僅需要一次FFT變換求解y1(k)，然后通過對y1(k)進行簡單的加法和乘法迭代運算就可以得到y(tǒng)2(k)，避免了FFT變換重復(fù)計算的問題可實現(xiàn)快速高精度頻率估計。

由式（4）和式（5）得

式中：v為信號采樣的尾首樣本差。

根據(jù)頻譜序列y1(k)和信號采樣的尾首樣本差v，即可快速得到幅度譜z(k)為

由式（7）知信號頻率估計值所對應(yīng)的序列號m。

進一步利用頻譜序列y1(k)和y1(k)+v在幅度譜峰位置m處的相位差，消除FFT算法受限于量化間隔的影響。根據(jù)式（8）和式（10），在幅度譜的譜峰位置m處，可以得到

頻率搜索間隔的量化誤差校正值u為

得到FFT幅相聯(lián)合的快速高精度頻率估計f0為

上述算法流程圖如圖1所示。

圖1 FFT幅相聯(lián)合快速高精度算法流程圖Fig.1 Flow chart of fast high precision algorithm based on FFT amplitude and phase combination

綜上，F(xiàn)FT幅相聯(lián)合的快速高精度頻率估計算法，一方面通過使用滑窗迭代算法，只需要進行一次FFT變換就可以得到頻譜序列y1(k)和y2(k)，減少了冗余計算；另一方面，利用幅度譜峰y1(m)和y1(m)+v的相位差對頻率搜索間隔量化誤差進行校正，達到了高精度快速測頻的目的。

4 計算復(fù)雜度分析

為了求解y1(k)與y2(k)，需要進行2次N-1點FFT變換，計算復(fù)雜度為

由式（10）可知，通過一次FFT變換求解y1(k)，然后對y1(k)進行簡單的加法和乘法迭代運算就可以得到y(tǒng)2(k)，從而避免了FFT變換重復(fù)計算的問題，因此，與傳統(tǒng)FFT算法相比，本文所提的FFT幅相聯(lián)合的快速高精度頻率估計算法的計算復(fù)雜度節(jié)省了O((N-1)log(N-1)-2)，從而實現(xiàn)快速高精度頻率估計。

5 仿真結(jié)果與分析

為了驗證FFT幅相聯(lián)合的快速高精度頻率估計算法的性能，進行多次蒙特卡洛仿真實驗，并與經(jīng)典的FFT算法進行比較。測量信號頻率時，F(xiàn)FT算法只利用幅度譜的譜峰位置測量信號頻率，F(xiàn)FT幅相聯(lián)合的快速高精度頻率估計算法不僅利用了幅度譜的譜峰位置，還利用了頻率搜索間隔的量化誤差校正值。

仿真過程中，信號采樣數(shù)分四種情況，即N=32，64，128，256，信號采樣中包括歸一化頻率（相對于采樣頻率）等于0.2140，信噪比為18 dB的信號，以及歸一化頻率等于0.1070、干信比為11 dB的干擾。

當(dāng)N=32，64，128，256時，頻譜序列y1(k)和y2(k)的幅度譜及其相位如圖2至圖5所示?？梢姡l譜序列y1(k)和y2(k)的幅度譜y1，y2的譜峰位置相同，驗證了通過式（6）估計譜峰位置是合理的；且隨著采樣點數(shù)增加，這一個特征越明顯。

圖2 N=32時，y1(k)和y2(k)的幅度譜及相位Fig.2 Amplitude spectrums and phases of y1(k)and y2(k) for N=32

圖3 N=64時，y1(k)和y2(k)的幅度譜及相位Fig.3 Amplitude spectrums and phases of y1(k)r and y2(k) for N=64

圖5 N=256時，y1(k)和y2(k)的幅度譜及相位Fig.5 Amplitude spectrums and phases of y1(k)and y2(k) for N=256

當(dāng)N= 32，64，128，256時，F(xiàn)FT算法與FFT幅相聯(lián)合的快速高精度頻率估計算法的信號歸一化測頻均方根誤差的累計概率分布函數(shù)（CDF）如圖6所示?？梢?，F(xiàn)FT幅相聯(lián)合的快速高精度頻率估計算法在信號采樣數(shù)N= 32時的CDF接近FFT算法在信號采樣數(shù)N= 256時的CDF，達到了高精度測頻的目的。

圖4 N=128時，y1(k)和y2(k)的幅度譜及相位Fig.4 Amplitude spectrums and phases of y1(k)and y2(k) for N=128

圖6 不同采樣數(shù)，F(xiàn)FT算法與FFT幅相聯(lián)合快速算法的歸一化測頻誤差的CDF及局部圖Fig.6 CDFs of normalized frequency measurement error of FFT algorithm and fast FFT amplitude and phase combination algorithm with different sampling numbers

由仿真結(jié)果可見，不同采樣數(shù)下，F(xiàn)FT算法的信號歸一化測頻均方根誤差分別為0.0091，0.0044，0.0023和0.0011；FFT幅相聯(lián)合的快速高精度頻率估計算法的信號歸一化測頻均方根誤差分別為0.0014，0.0007，0.0003和0.0002?？梢缘贸?，與FFT算法相比，利用同樣的采樣信號，當(dāng)N=32，64，128，256時，F(xiàn)FT幅相聯(lián)合的快速高精度頻率估計算法的歸一化均方根誤差分別降低了84.15%，84.57%，85.24%和84.73%，F(xiàn)FT幅相聯(lián)合的快速高精度頻率估計算法在信號采樣數(shù)N=32時的測頻精度接近FFT算法在信號采樣數(shù)N=256時的測頻精度，達到了高精度測頻的目的。

6 結(jié)論

本文提出了一種FFT幅相聯(lián)合的快速高精度頻率估計算法，一方面，只需要進行一次FFT變換就可以得到兩個頻譜序列，從而減少了冗余計算；另一方面，利用兩個頻譜序列的相位差對頻率搜索間隔量化誤差進行校正，達到了高精度快速測頻的目的。實驗結(jié)果和分析表明，由于同時利用了幅度譜的譜峰位置與頻率搜索間隔的量化誤差校正值，該算法在信號采樣數(shù)N=32時的測頻精度接近FFT算法在信號采樣數(shù)N=256時的測頻精度，適用于需要快速、高精度測頻的應(yīng)用場合。未來的研究還可以利用更多的頻譜序列，進一步提高頻率估計的精度。