任季中,趙 倩
(1.中國電子科技集團公司第五十四研究所,河北石家莊050081;2.國家廣電總局951臺,河北石家莊050701)
脈沖信號是現代雷達主要采用的信號形式,脈沖信號頻率測量是雷達偵察中不可或缺的環節,對雷達對抗起著重要的作用。數字化處理是雷達對抗系統發展的趨勢之一,常用的數字測頻方法包括過零點檢測法、相位差分法、快速傅里葉變換(FFT)法和現代譜估計法[1]。其中FFT法工程可實現性強,實時性好,且適用于寬帶偵收,因此在工程中得到廣泛應用。
本文以時寬較短(0.2~1μs)的正弦波脈沖信號為研究對象,分析了傳統FFT測頻法的不足之處,從工程應用角度分析了提高測頻精度的改進方法,并提出了基于FPGA的全數字實現流程。
信號x(t)經過數字化采樣后為x(n),n=0,1,2,…,N-1,為對其進行頻譜分析,進行離散傅里葉變換(DFT),將信號從時域轉換到頻域,如式(2)所示:

DFT實現時采用的快速算法即是FFT,經FFT處理后,信號的頻率分辨率為:

式中,fs為采樣率,設信號的時寬為T,則信號的點數為T×fs,信號的頻率分辨率可表示為:

可見,FFT測頻的頻率分辨率只與信號時寬有關,根據譜線的最大值來換算信號的頻率,如果信號的頻率正好落在一根譜線上,得到的頻率測量結果是準確的,而在多數情況下,信號頻率落在兩根譜線之間,由最大值譜線位置反映的頻率不再準確,最大測頻誤差為Δf/2。
脈沖是雷達最常采用的信號形式,根據需要,雷達有時會采用脈內帶調制的信號類型,例如相位編碼、線性調頻等,對于此類復雜信號可采用各種信號處理方法將其轉化為普通正弦波信號,因此正弦波脈沖的測頻方法具有通用性。根據上文分析結果,對于時寬較長的脈沖,采用FFT測頻法易于實現較高測頻精度,滿足設備指標要求。但是對于短脈沖,例如一個0.2μs寬的脈沖,根據式(3),理論能達到的測頻精度只有2.5MHz,難以滿足偵察要求。
補零[1]是指在進行FFT運算之前在時域數據的尾部添加一些零,并使總的時域數據點數保持為2的冪次方。由于補零不增加任何新的信息,所以并不改變頻譜形狀和頻率分辨率,補零只是在原始點數的FFT結果中內插了一些頻率分量。對于點數較少的FFT結果,在大多數情況下,從中找到峰值比較困難,也很難觀察到頻譜的細微結構。而補零之后,功率譜的峰值位置可以較清晰的顯露出來,有助于提高對主瓣峰值頻率分量進行精確定位的能力,由此提高測頻精度。
補零技術的缺點是額外增加了處理量,補零越多,處理時間也就越長。此外,對于存在噪聲的情況,補零也不能改善信噪比,存在頻譜峰值點定位錯誤的可能,造成測頻誤差增大。
插值FFT估計頻率方法利用真正的頻譜峰值兩側的2根FFT譜線,求其幅度比值,建立一個以修正頻率為變量的方程,解方程得到修正頻率值,對FFT最大譜線位置進行校正,以實現對信號頻率更高精度的估計,如圖1所示。相比上節補零的方法,不必增加FFT的長度以及由此帶來的運算處理量,只需從FFT結果中找出兩個點就足夠。

圖1 矩形窗頻譜函數
在圖1中插值頻率校正即求出矩形窗譜主瓣中心與相鄰譜線的橫坐標差,對于譜線位置x、x+1,其矩形窗譜函數為sinc函數[2],表示為f(x),頻譜值為yx、yx+1,矩形窗譜函數和頻譜值已知,可構成一方程如下:

在圖1中,sinc函數以峰值橫坐標為零點,頻率修正值δ=-x,只要根據式(4)求解出x,即可得到頻率修正值。
對矩形窗譜函數歸一化,求模可得:

帶入式(4),得到:

式中,α=yx/yx+1。實際應用中,已知FFT譜峰最大值位置k1,相鄰次大值位置k2,頻率分辨率Δf,利用修正頻率值校正頻率可得[3]:

當k2=k1+1時,取加號;k2=k1-1時,取減號。
以上對插值FFT頻率估計法進行了理論分析,實際應用中,不可避免的會有背景噪聲,本小節將在加性高斯白噪聲背景下,通過仿真分析插值FFT頻率估計法的性能。
設定仿真參數,信號采樣率fs為1280MHz,脈沖寬度0.2μs,頻率分別設f1為102.4MHz,f2為100.4MHz,按照10dB信噪比加入高斯白噪聲[4]。
以信號頻率f1進行仿真,連續測頻1000次,仿真結果如圖2所示。由圖可知,最大測頻誤差不超過300kHz。

圖2 測頻誤差變化圖
以信號頻率f2進行仿真,連續測頻1000次,仿真結果如圖3所示。由圖3可知,最大測頻誤差超過1MHz。

圖3 測頻誤差變化圖
由以上結果易知,噪聲背景下的插值法測頻誤差與頻率位置的選取有關,準確的說,是與實際頻率位置偏離FFT譜線的距離,即與頻率修正值δ大小有關。一般情況下,FFT幅度最大值k1和相鄰次大值k2都位于矩形窗函數的主瓣內,當實際頻率位置位于k1、k2中間附近時,信號向兩邊泄漏的能量都較多,在一定信噪比下,使得k1、k2電平均大于噪聲電平,確保了k2位置不會找錯,這對應了圖2的情況。而當δ值接近0時,較多信號能量集中在k1處,k2處幅度較小,而最大譜線相鄰另一側的幅值k3由于受噪聲影響,與k2幅度接近,因此會造成最大譜線相鄰的次大譜線位置找錯,導致式(7)中加或減符號錯誤,使得測頻結果出現較大誤差,對應了圖3的情況??梢?,在噪聲背景下,插值FFT測頻法有局限性,即只有在δ值大于某一閾值時,才能達到較理想的測頻精度[5]。
為抑制頻譜泄漏,進行FFT之前常對采樣數據進行加窗處理。抑制泄漏的同時,加窗會使得頻譜主瓣加寬[6]。對于插值FFT法求頻率,無論頻譜最大值偏離實際FFT譜線距離遠近,最大值及其相鄰兩側譜線都被包含在主瓣之內,在一定信噪比條件下,次大值不會趨近于噪聲電平,使得抗噪聲性能增強。
加窗后頻率校正值仍隨k1、k2幅度大小變化,但變化規律不再依據sinc函數,文獻[7]給出了幾種窗函數對應的頻率校正計算公式,當選用漢寧(Hanning)窗時,計算式較易于實現。對采樣數據加Hanning窗,利用k1和k2的比值α帶入窗函數,經推導可得:

利用α估計頻率修正值δ的解析式如下:

校正頻率的方法如式(10)所示。
設定仿真參數,信號采樣率、脈沖寬度不變,仍按照10dB信噪比加入高斯白噪聲。連續測頻1 000次,頻率f1仿真結果如圖4所示,頻率f2仿真結果如圖5所示。

圖4 測頻誤差變化圖

圖5 測頻誤差變化圖
由仿真結果可知,最大測頻誤差不超過500kHz。加窗處理后,在常規信噪比條件下,次大值方向錯誤的概率大大降低,由此造成的頻率估計誤差已可以忽略。
加漢寧窗插值FFT測頻的實現框圖如圖6所示。整個算法可在一片FPGA中實現,采樣數據進入FPGA后,與漢寧窗數值相乘,漢寧窗值可預先存儲在FPGA內ROM中,以查表方式讀出[8]。加窗后的數據進入FFT模塊進行流水處理,得到信號的頻譜結果,對頻譜結果進行峰值搜索,并與檢測門限比較,判斷是否存在信號,當頻譜峰值大于檢測門限時,找出峰值位置相鄰幅度較大的譜線位置,按照式(8)經過插值換算,得到頻率估計值[9]。

圖6 加窗插值FFT測頻實現框圖
式(10)中存在除法計算,實現時可將除法轉化為先對除數求倒數,再與被除數相乘的過程,利用FPGA中豐富的RAM資源,求倒計算利用查表完成[10]。除此之外,運算只由常規加、乘組成,便于FPGA實現。
某寬帶偵察接收機,指標要求適應脈沖寬度0.2~1000μs,測頻誤差不大于500kHz。實現時信號檢測與頻率測量由FPGA硬件完成,算法采用定點實現,頻率的分辨率設為15.625kHz。測頻結果送出至軟件顯示,誤差單位為kHz,取整。根據要求設置信號幅度在接收機實測靈敏度以上3dB,頻率選擇在1001~1003mHz和200kHz步進,脈沖寬度分別設為1μs、0.5μs和0.2μs。測試結果如表1所示。

表1 雷達信號測頻精度測試結果
可見在不同頻率、不同脈寬時測頻最大誤差均小于500kHz,滿足指標要求。
論述了一種易于工程實現的脈沖信號實時測頻算法,與傳統方法相比可以達到更高的測頻精度。經過試驗證明,可以滿足目前常規雷達偵察接收機的指標要求,可應用于目標為脈沖信號的電子對抗系統,具有較高的應用價值。
[1]TSUI J.寬帶數字接收機[M].北京:電子工業出版社,2002.
[2]程佩青.數字信號處理教程[M].北京:清華大學出版社,2001.
[3]劉建林,陳兵.一種高精度單頻信號頻率估計算法[J].無線電工程,2011,41(4):26-28.
[4]萬永革.數字信號處理的MATLAB實現[M].北京:科學出版社,2012.
[5]齊國清,賈欣樂.插值FFT估計正弦信號頻率的精度分析[J].電子學報,2004,32(4):625-629.
[6]劉麗格,魏利輝,馬媛.一種高速跳頻信號實時載頻測量方法[J].無線電工程,2011,41(7):30-32.
[7]丁康,謝明,王志杰.離散頻譜的幅值、相位和頻率的校正方法及誤差分析[J].動態分析與測試技術,1996,14(1):10-29.
[8]吳爭,李輝.高速并行FFT的FPGA實現[J].無線電工程,2013,43(10):16-18.
[9]趙雷鳴.一種跳頻MSK信號檢測算法及FPGA實現[J].無線電通信技術,2011,37(1):56-58,64.
[10]Uwe Meyer-Baese.數字信號處理的FPGA實現[M].劉凌,譯.北京:清華大學出版社,2006.