周曉青，李合生，陶榮輝，蔡英武

(中國工程物理研究院 電子工程研究所，四川 綿陽 621900)

在現(xiàn)代戰(zhàn)爭中，數(shù)字接收機是電子偵察系統(tǒng)的重要組成部分，它在廣大的作戰(zhàn)區(qū)域中準(zhǔn)確、快速、全面地獲得各種目標(biāo)信息的作用巨大，其地位是不可取代的。它可以捕獲、定位、分析和識別作戰(zhàn)區(qū)域內(nèi)的雷達輻射源，從而獲取有用信息。其中電磁波的頻率無疑是非常重要的一個參數(shù)，它反映了雷達的功能和用途、雷達的頻率捷變頻范圍和譜寬，是信號分選和威脅識別的重要參數(shù)之一[1]。為了實時、準(zhǔn)確地獲得該參數(shù)，需要對載波頻率進行瞬時測量。在20世紀(jì)五六十年代，瞬時測頻技術(shù)IFM(Instantaneous Frequency Measurement)應(yīng)運而生，瞬時測頻主要分為模擬瞬時測頻和數(shù)字瞬時測頻兩種。其中，模擬瞬時測頻設(shè)備簡單易于實現(xiàn)，但測頻精確度不是很高，并且穩(wěn)定性差，容易受環(huán)境的影響，該技術(shù)較早在工程中得到應(yīng)用。數(shù)字瞬時測頻具有穩(wěn)定、精確度高、實時性好等優(yōu)點。90年代中期以后，隨著高速集成電路和數(shù)字算法的發(fā)展，數(shù)字瞬時測頻DIFM(Digital Instantaneous Frequency Measurement)技術(shù)得到飛速的發(fā)展，并且有逐漸取代模擬瞬時測頻的趨勢。

數(shù)字瞬時測頻即用數(shù)字信號處理的方法對信號的載頻進行測量[2]?，F(xiàn)代電子戰(zhàn)爭的主要特點為：(1)寬頻帶、大視場；(2)瞬時信號檢測、測量和高速信號處理[1]。所以要求數(shù)字測頻必須可以對大帶寬的信號進行高速、高精度的測量。隨著數(shù)字接收機的發(fā)展，對數(shù)字測頻算法的研究也越來越多，這些方法各有優(yōu)缺點，需要根據(jù)實際的應(yīng)用場合選擇使用。

本文主要按照算法所利用的信號的不同信息分量，對常用的一些算法進行了分類介紹，并給出了算法的基本思想以及目前比較常用的改進算法，對其中的直接計數(shù)法、相位推算法、頻率推算法、傅里葉變換法及譜估計等方法進行了MATLAB仿真?？偨Y(jié)了各算法的主要優(yōu)缺點及適用范圍，為設(shè)計者在工程實用中選擇合適的算法提供幫助。

數(shù)字測頻算法種類繁多，根據(jù)算法所利用的信號不同信息分量，數(shù)字測頻算法分為：利用幅度信息的測頻算法、利用相位信息的測頻算法、利用頻率信息的測頻算法、利用功率譜信息的測頻算法4類。下面將分別對這4類算法進行詳細的分析研究。

1 利用幅度信息的測頻算法

利用信號的幅度信息進行測頻的典型算法為直接計數(shù)法[3-6](也稱為過零計數(shù)法、過零檢測法)，是通過計算一定數(shù)量的信號周期總共占有多長時間來推算信號頻率。其原理框圖如圖1所示。

圖1 直接計數(shù)法原理框圖

設(shè)采樣量化的時鐘頻率為fs，計算信號的n個半周期，即計算出一段信號極性變化n次的時間，總共的信號時間為m個時鐘周期，則信號的周期和頻率分別為：

由于計數(shù)結(jié)果總是整數(shù)，可能發(fā)生的最大誤差為1個時鐘周期，則可能的最大周期誤差為2/nfs?；蛘哂孟鄬φ`差表述，測量信號周期和頻率的相對誤差均為1/m，當(dāng)m較小時，計數(shù)法的測頻誤差較大，要提高m的數(shù)值，只能是提高采樣頻率和信號存在時間。

圖2給出了直接計數(shù)法在有噪聲和無噪聲情況下的比較。圖 2(a)、圖2(c)分別為當(dāng)信噪比等于 10 dB時的待測信號和直接計數(shù)法的計數(shù)點；圖 2(b)和圖2(d)分別為沒有噪聲時的待測信號和直接計數(shù)法的計數(shù)點?？梢钥闯觯?dāng)沒有噪聲時，計數(shù)點清晰可見，可完全滿足頻率計算的要求。但在信噪比為10 dB時，已經(jīng)很難利用計數(shù)法對待測信號進行頻率測量了。

為了進一步提高直接計數(shù)法的測頻精度，參考文獻[4]提出可以采用雙門限整形的方式來降低在較低信噪比情況下的計數(shù)誤差；參考文獻[3]提出引入等效時鐘來減小量化誤差。由于直接計數(shù)法本身的抗噪聲性能較差，所以這些方法雖然對算法有所改進，但性能仍有待進一步提高。

由以上分析可以看出，直接計數(shù)法的主要優(yōu)點為：測頻方法簡單、運算量小、速度快、工程實現(xiàn)的體積小、重量輕。而其主要缺點：抗噪聲性能差、測頻精度較低，并且只能用于單載頻測量。

2 利用相位信息的測頻算法

在信號處理中，信息調(diào)制一般在信號的幅度、相位和頻率上，在這三者中，將有用信息調(diào)制在幅度上的信號抗噪聲性能一般是最差的，所以希望其能夠從相位和頻率中提取所需要的有用信息。

相位推算法[7-9]先把信號幅度轉(zhuǎn)換為信號相位，再采用相位信息來計算頻率值，如圖3所示。實踐證明，相位推算法可以得到較高的測量精度。

圖3 相位推算法原理框圖

假設(shè)被測量的信號具有如下形式：

其中，A為信號幅度因子，ω為角頻率，φ為相位，記錄m個連續(xù)的采樣點為：

式中，T為采樣周期，經(jīng)過一系列推算[7]，可得信號頻率為：

其中 φi為：

式中，n為所使用的相位點數(shù)，fs為采樣率。

用MATLAB對相位推算法進行了仿真，當(dāng)采樣頻率fs=500 MHz、脈沖寬度 τ=2 μs、信噪比 SNR=10 dB、量化位數(shù)為 8 bit時，對 100 MHz～150 MHz的頻率范圍，每隔1 MHz進行1次測頻，測得的頻率誤差如圖4所示。

從仿真結(jié)果可以看出，測頻誤差最大不會超過3 kHz，達到了較高的測頻精度。同時，相位推算法的測頻實時性和抗噪聲性能均較好。但缺點在于測頻算法較為復(fù)雜，需要的硬件資源較多。

因此，相位推算法適用于對測頻精度及實時性要求較高，且信號信噪比較高，有較多硬件資源的情況。

3 利用頻率信息的測頻算法

利用頻率信息進行測頻，包括在時域上和頻域上進行測頻兩類算法。在時域上對頻率進行推算的典型算法為頻率推算法，而在頻域上對頻率進行推算的典型算法為傅里葉變換法。其中，對傅里葉變換法進行的研究較多，發(fā)展較快，目前有很多改進算法，同時也帶動了對短時傅里葉變換以及小波變換在瞬時測頻方面的研究。

3.1 頻率推算法

頻率推算法[2，8-10]又稱為曲線擬合法，是在假設(shè)信號為一個單一載頻的情況下導(dǎo)出的一種算法，其原理框圖如圖5所示。

假設(shè)被測量的信號具有如下形式：

如果采樣周期為T，則信號變成如下序列：

則可推出：

采樣率fs=1/T，則信號頻率與采樣率之間的關(guān)系：

在已知采樣率情況下，只需要3個采樣點就可以得出信號的頻率。但必須強調(diào)的是，這一算法的前提是信號必須是單一載頻的。

如果要通過式(8)求取信號頻率，必須滿足2個條件：(1)si≠0，而且可以求取反余弦函數(shù)；(2)比較容易滿足求取反余弦函數(shù)(可以利用CORDIC算法來求取)。第1個條件比較困難，在量化精度較低的情況下，較小的信號往往也會被量化為0，因此在工程上往往采用其改進算法，即用3點進行測頻變?yōu)橛枚帱c進行測量[8]。

令 yi=si+1+si-1，xi=2si，K=cos(ωT)，有：

對于某一個被測頻率，K=cos(ωT)是一個常數(shù)，則(9)式實際上是一個直線方程，直線的斜率為K。采用線性擬合，求一條直線，使這些點到直線的距離的平方和最小，這樣就得到一個優(yōu)化的K值，再用這個值去計算信號的頻率。

本文對頻率推算法進行了MATLAB仿真，當(dāng)采樣頻率 fs=500 MHz、脈沖寬度 τ=2 μs、信噪比 SNR=10 dB、量化位數(shù)為 8 bit時，對 100 MHz～150 MHz的頻率范圍，每隔1 MHz進行1次測頻，測得的頻率誤差結(jié)果如圖6(a)所示； 當(dāng)待測頻率 f=200 MHz、脈沖寬度 τ=2 μs、信噪比SNR=10 dB、量化位數(shù)為 8 bit時，對 400 MHz～800 MHz的采樣頻率范圍，每隔1 MHz進行1次測頻，測得的頻率誤差結(jié)果如圖6(b)所示。

從仿真結(jié)果可以看出，頻率推算法的測頻精度受采樣頻率的影響較大。在相同條件下，用頻率推算法進行測頻的最大誤差不超過300 kHz，測頻精度明顯差于相位推算法，但仍具有工程應(yīng)用價值，同時瞬時性較好。另外，由于頻率推算法運算量大，所以在硬件實現(xiàn)時會耗費較多的硬件資源。

因此，頻率推算法適用于對瞬時性及測頻精度要求稍高，同時有較高的信噪比和充裕的硬件資源的場合。

3.2 傅里葉變換法

用傅里葉變換法來求解信號的頻率是當(dāng)前研究最多的一種算法。 傅里葉變換法[2，9，11-13]是將一段信號進行傅里葉變換，將時域信息轉(zhuǎn)化為頻域信息，在頻域上的峰值點即對應(yīng)該信號的頻率。其測頻精度依賴于采樣點數(shù)的多少。

傅里葉變換法尤其適用于信噪比較低的情況，它的測頻精度隨著信噪比的降低并沒有明顯的降低。在時域上，當(dāng)信噪比為0 dB時，傅里葉變換法的測頻基本上還是正確的。此外，傅里葉變換法還可以對復(fù)雜信號進行頻譜分析。本文對傅里葉變換法在不同信噪比的測評效果進行了仿真，當(dāng)采樣頻率為10 MHz、FFT點數(shù)為2 048點、量化位數(shù)為8 bit時，仿真結(jié)果如圖7所示。

圖7(a)、(b)中信號都是頻率為 1 MHz和2 MHz的多載頻信號，圖7(c)、(d)分別為其所對應(yīng)的頻譜。其中，圖7(a)的信號疊加了噪聲，信噪比為 0 dB，而圖 7(b)的信號沒有疊加任何噪聲，信噪比為無窮大。從仿真結(jié)果可以看出，盡管信號疊加的噪聲較大，但是在信號頻譜上仍然可以較好地分辨出信號譜線，這說明傅里葉變換法具有很好的抗噪聲性能。同時，傅里葉變換法還具有測量多載頻信號頻率的功能，圖中可以清楚看到，1 MHz和2 MHz分別對應(yīng)譜線。正是這些優(yōu)點使得傅里葉變換法成為當(dāng)今數(shù)字測頻應(yīng)用最為廣泛的一種算法。

但是，傅里葉變換法也有非常明顯的缺點：測頻精度及瞬時性較差，其測頻精度完全是依靠采樣點數(shù)的多少，采樣點數(shù)越多精度越高。這就要求要提高信號的采樣率，而這一點對于高頻信號來說非常困難，所以要求通過其他的方法來獲得測頻精度的提高。

但總體上來說，傅里葉變換法比較適用于信噪比低，對瞬時性和測頻精度要求不太高，具有足夠的硬件資源的場合。

4 利用功率譜信息的測頻算法

譜估計法測頻是通過求解信號的功率譜來確定信號的頻率，它對多載頻信號的分辨能力較強，對信噪比的要求較低，主要包括經(jīng)典功率譜估計法和現(xiàn)代譜估計法。經(jīng)典功率譜估計主要包括周期圖法(又稱直接法)和BT(Blackman-Tukey)法(又稱為間接法或相關(guān)圖法)。現(xiàn)代譜估計主要包括參數(shù)模型法(包括常用AR模型、最大熵法、最小熵法和線性預(yù)測法等)、最大似然法、自適應(yīng)法、高階譜估計等。下面將主要對經(jīng)典譜估計中的周期圖法和現(xiàn)代譜估計中的最大熵法進行分析。

4.1 周期圖法

周期圖法[2，14-17]是將平穩(wěn)隨機序列 x(n)的 N點觀察數(shù)據(jù)xN(n)視為能量有限的信號，直接將離散信號xN(n)進行傅里葉變換來求取功率譜估計的。即首先由采樣序列 x(0)，x(1)，…，x(N-1)，求得：

再由下式求得功率譜：

由于快速傅里葉變換(FFT)的提出，使DFT算法的運算量顯著降低。因此，周期圖法得到廣泛應(yīng)用。該方法的測頻范圍由采樣率決定，測頻精度和頻率分辨率由采樣率和采樣序列長度決定。其頻率分辨率為：fs/(2N)，fs為采樣率，N為序列長度。在有限的信號樣本長度下，為提高測頻精度可使用內(nèi)插法。

周期圖法最大的優(yōu)點在于可以對多載頻信號進行測頻，并且抗噪聲性能較好。對信號x(n)=sin(2πf1n)+cos(2πf2n)+w(n)進 行 了 仿 真 ， 其 中 ，f1為 0.2 Hz，f2為0.3 Hz，w(n)為白噪聲，信噪比為 0 dB，測頻結(jié)果如圖 8所示。

從圖8中可以看出，在頻率0.2 Hz和0.3 Hz處功率譜有2個峰值，說明信號中含有0.2 Hz和0.3 Hz的周期成分。并且，采樣點數(shù)的增加并沒有使周期圖的方差減小，而是采樣點數(shù)越多，周期圖上的隨機起伏反而越密集和劇烈。用有限長樣本序列的FFT來表示隨機序列的功率譜是一種近似算法，必然存在著誤差。為了減小誤差，使功率譜估計更加平滑，可以采用分段求平均的方法，即 Welch法[18]和 Bartlett法[19]，當(dāng)前對 Welch法的研究更為廣泛。參考文獻[17]采用Welch算法對在4種窗函數(shù)情況下進行仿真，結(jié)合得到的結(jié)果，給出了選取合適的窗函數(shù)應(yīng)注意的一般原則。

從以上分析可以看出，周期圖法適用于需要測量多載頻信號頻率、且信噪比較低、對測頻精度和測頻實時性要求不高的情況。

4.2 最大熵法

由于周期圖法功率譜估計需要對信號序列“截斷”進行處理，所以必然存在誤差。 最大熵功率譜估計[2，14-16，19-20]的目的是最大限度地保留截斷后丟失的其他信號的信息，使估計譜的熵最大。

如果數(shù)據(jù)為 x(0)，x(1)，…，x(N-1)，其算法步驟為：

(1)計算預(yù)測誤差功率的初始值P0和前、后向預(yù)測誤差的初始值f0(n)。

并令 m=1。

(2)求反射系數(shù)。

式中，*表示取共軛。

(3)計算前向預(yù)測濾波器系數(shù)。

(4)計算預(yù)測誤差功率。

(5)計算濾波器輸出

(6)令 m←m+1，并重復(fù)步驟(2)～步驟(5)，直到預(yù)測誤差功率Pm不再明顯變小。

最大熵法最大的優(yōu)點在于可以較清晰地分辨相隔較近的兩個信號，而這一點正是周期圖法欠缺的，并且測頻精度較好。對信號 x(n)=sin(2πf1t)+sin(2πf2t)+sin(2πf3t)+w(n)進行了最大熵法測頻仿真，其中，f1為50 Hz、f2為 60 Hz、f3為 70 Hz、w(n)為白噪聲、信噪比為 10 dB，仿真結(jié)果如圖9所示。

在圖9中，可以非常清楚地分辨出3個頻率，這是最大熵法相較于周期圖法的優(yōu)勢所在，但是最大熵法仍然存在分辨率較低、頻偏等問題。參考文獻[20]提出了一種改進算法，使得其分辨率和運算速度有了一定的提高。參考文獻[21]提出了一種修正算法，該方法在階數(shù)估計時引入收斂因子，從而更為有效地估計了階數(shù)，這種算法減小了頻偏、提高了分辨率。

所以，最大熵法適用于需要測量頻率相隔較近的多載頻信號頻率，且信噪比較低、對測頻精度和測頻實時性要求不高的情況。

通過以上的分析研究，對各種算法性能的比較如表1所示。

表1 各種算法性能比較

由表1可以得出：在信噪比非常差的情況下，可以優(yōu)先考慮采用傅里葉變換法；在多載頻信號測頻時，可以考慮譜估計法；在對測頻精度要求較高的場合可以優(yōu)先考慮相位推算法；在對實時性要求較高，同時信號具有很高信噪比的情況下，可以采用直接計數(shù)法；在對實時性要求高，但信號的信噪比不太高時，可以采用相位推算法。

各種數(shù)字測頻算法都有其不同的優(yōu)缺點，需要在實際應(yīng)用時根據(jù)不同情況加以選擇使用。本文對這些算法加以總結(jié)，給出了各自的特點，其目的是為設(shè)計者在今后的實際工程應(yīng)用中選擇合適的算法提供參考。

[1]趙國慶.雷達對抗原理[M].西安：西安電子科技大學(xué)出版社，2005.

[2]王新春.高速寬帶數(shù)字測頻技術(shù)研究[D].西安：西安電子科技大學(xué)，2005.

[3]李相平.一種數(shù)字計數(shù)瞬時測頻的方法研究[J].國外電子測量技術(shù)，2006，25(5)：7-8.

[4]胡來招.數(shù)字瞬時測頻-直接計數(shù)法測頻[J].電子對抗，2004(5)：44-49.

[5]魏鐘記，劉峰，汪銘東.數(shù)字測頻方法及實現(xiàn)[J].電子測量技術(shù)，2006，29(4)：28-29.

[6]張金剛，孟小芬，唐志凱，等.計數(shù)式瞬時測頻精度分析[J].國外電子測量技術(shù)，2005(5)：10-11.

[7]喻旭偉.相位推算法瞬時測頻研究[J].電子信息對抗技術(shù)，2007，22(6)：6-9.

[8]胡來招.數(shù)字瞬時測頻-相位推算法測頻[J].電子對抗，2005(3)：36-41.

[9]王偉旭.用相位推算法實現(xiàn)瞬時測頻[D].成都：電子科技大學(xué)，2007.

[10]何偉.新型寬帶數(shù)字接收機[D].成都：電子科技大學(xué)，2004.

[11]胡來招.數(shù)字瞬時測頻-傅里葉變換法測頻[J].電子對抗，2005(1)：42-47.

[12]陳文靜.一種基于DFT快速頻率檢測的算法與實現(xiàn)[J].電子信息對抗技術(shù)，2006，21(3)：3-6.

[13]王峰.基于離散傅里葉變換的數(shù)字高速測頻[J].數(shù)據(jù)采集與處理，2007，22：370-373.

[14]皇甫堪.現(xiàn)代數(shù)字信號處理[M].北京：電子工業(yè)出版社，2004.

[15]萬建偉.信號處理仿真技術(shù)[M].長沙：國防科技大學(xué)出版社，2008.

[16]張賢達.現(xiàn)代信號處理[M].北京：清華大學(xué)出版社，2003.

[17]余訓(xùn)鋒.改進周期圖法功率譜估計中的窗函數(shù)仿真分析[J].計算機仿真，2008(3):111-114.

[18]WELCH P D.The use of fast fourier transform for the estimation of power spectra:a method based on time averaging over short，modified periodograms[J].IEEE Trans on Audio and Electroacoustics， 1967，15(2):70-73.

[19]BARTLETT M S.Properties of sufficiency and statistical tests[J].Proceedings of the Royal Society of London.Series A， Mathematicaland PhysicalSciences，1937，160(901):268-282.

[20]高勤.一種改進的 Burg算法[J].指揮技術(shù)學(xué)院學(xué)報，2000(4):81-84.

[21]羅豐.基于 Burg算法的短序列譜估計研究[J].電子科技大學(xué)學(xué)報，2005(10):724-728.