張 瑋,王 平

(1.海軍工程大學， 武漢 430000; 2.解放軍92038部隊， 山東 青島 266041)

1 引言

跳頻[1-3](frequency hopping,FH)通信屬于擴頻通信的一種，跳頻信號可以在較寬的帶寬內(nèi)隨時間按偽隨機規(guī)律跳變，具有突出的抗干擾、抗截獲、抗多徑能力，因而在軍事領(lǐng)域和民用領(lǐng)域應用廣泛。在通信對抗領(lǐng)域中，對跳頻通信實施高效通信干擾的前提獲取跳頻通信信號參數(shù)，而如今用頻設備日益增多，電磁環(huán)境日趨復雜，獲取跳頻通信信號參數(shù)也變得更加困難，因此對跳頻信號參數(shù)估計逐漸成為通信對抗領(lǐng)域研究的熱點問題[4]。

目前，用于跳頻信號參數(shù)估計的方法主要有3類，圖像處理法[5-6]，稀疏重構(gòu)法[7-8]，時頻分析法[9-11]，而目前時頻分析法是跳頻信號分析的主要方法。時頻分析法是通過時間和頻率2個變量來描述在不同時間條件下不同頻率的信號能量強度的方法[11]。文獻[12]提出一種基于平滑偽維格納分布的參數(shù)估計方法，具有較高的精度，但是存在交叉干擾。文獻[13]提出一種基于小波分解和希爾伯特黃的跳頻檢測算法，但是小波算法的低頻部分分辨率較低，各頻率分量強度并不相同，不便于進行閾值分割。文獻[14]提出一種同步壓縮的方法，通過對時頻頻譜進行壓縮得到新的時頻頻譜，但是存在噪聲魯棒性差的缺點。文獻[15]提出一種慢跳頻信號的檢測方法，根據(jù)各信道輻射計檢測結(jié)果判斷是否存在跳頻信號，其缺點是抗噪性能差。文獻[16]提出了一種根據(jù)信號和背景噪聲頻域自相關(guān)統(tǒng)計量檢測的算法，但是這種算法只能在背景噪聲中完成跳頻信號檢測，在存在其他干擾時算法失效。文獻[17]提出一種跳頻周期估計的方法，通過對時頻矩陣中信號強度最大值的提取得到周期性函數(shù)，對函數(shù)進傅里葉變換運算，得到的峰值頻率即為跳變速率，取倒數(shù)得到跳頻周期的估計值，但是當存在定頻信號干擾時算法性能大幅度下降。文獻[10]通過提取時頻脊線的方法進行參數(shù)估計，利用迭代法去除噪聲，通過k-means聚類的方法進行參數(shù)估計，抗噪聲性能較好，但無法解決定頻信號與跳頻信號發(fā)生頻率碰撞的問題。文獻[18]改進能量對消的方法，能夠在完成對跳頻和定頻信號檢測的同時減少計算量，但是在突發(fā)干擾和低信噪比的情況下性能較差。

為解決以上問題，本文提出一種基于時頻分析的算法。通過短時傅里葉變換得到時頻圖像，利用能量對消去除定頻信號干擾，采用全局閾值法對圖像進行二值化處理，通過形態(tài)學濾波消除各類干擾，然后，通過對時頻脊線中頻率跳變時刻處理獲得跳頻周期的估計值，通過k-means聚類算法實現(xiàn)對頻率集的估計。設計仿真實驗并通過實測數(shù)據(jù)測試算法可行性。

2 跳頻信號數(shù)學模型

跳頻信號的頻率可隨時間偽隨機跳變，其數(shù)學模型可進行如下定義[19-20]：

(1)

(2)

其中:T1為起跳時刻；Th為跳頻周期；fk為第k個跳頻時隙的跳頻頻率；N為觀察時間內(nèi)頻率跳變數(shù)量；rectTh為門函數(shù)，其定義如式(2)所示。

3 跳頻信號特征參數(shù)估計方法

算法具體步驟如下：

1) 采用短時傅里葉變換對信號進行時頻變換，得到時頻圖像；

2) 采用能量對消的方法消除定頻信號干擾；

3) 采用全局閾值的方法處理時頻圖像，得到二值化圖像，其中閾值的計算采用OTSU算法；

4) 采用形態(tài)學濾波的方法去除信號毛刺、掃頻干擾以及猝發(fā)信號干擾，彌合裂縫、填補空洞，獲得高清晰度的時頻圖像；

5) 畫出時頻脊線，通過對頻率跳變時刻求解一階差分方程，可以獲得跳頻周期的估計值；采用k-means聚類算法，將頻率值的聚類結(jié)果按照時間順序進行排列，可以獲得頻率集的估計值。

算法總體流程框圖如圖1所示。

圖1 算法流程框圖

3.1 時頻變換

通過時頻圖像對跳頻信號進行分析是參數(shù)估計的效手段之一[21]，時頻圖像可以清晰地表示跳頻信號的頻率跳變情況。常見的時頻變換方法有STFT、Wigner-Ville分布、平滑偽Wigner-Ville分布、小波變換、希爾伯特黃變換等。本文采用短時傅里葉變換(STFT)進行跳頻信號的時頻變換，即：

(3)

式(3)中:s(t)表示要處理的信號;h(t)表示窗函數(shù);h*(t)表示窗函數(shù)的共軛函數(shù)。

為方便計算機計算，可將信號離散化處理：

3.2 能量對消

將時頻矩陣與頻率分量相減，可以得到對消矩陣I(n,k)：

由于在接收信號時間段內(nèi)，定頻信號始終存在，因此定頻信號的時頻能量特征隨時間變化不大；跳頻信號的頻率隨時間變化較快，每一頻率在觀測時間內(nèi)存在時間較短，因此其對應頻率上的均值遠小于其能量值，因此采用對消的方法可以有效去除定頻干擾保留跳頻信號。

3.3 二值化處理

為獲取清晰地時頻圖像，可將對消矩陣進行二值化處理，如式(4)所示，遍歷對消圖像矩陣，大于閾值μ則賦值為1，小于等于閾值μ則賦值為0，即可得到二值化圖像矩陣。

(4)

對于閾值的選取有多種算法，OTSU算法[23]、Bernsen算法[24]、Niblack算法[25]等，本文采取最大類間方差法(OTSU)進行圖像二值化處理，OTSU算法屬于全局閾值法，是一類基于整個圖像通過計算選取合適的閾值的方法[26]。

通過OTSU算法計算閾值通常需要兩步：

1) 計算圖像灰度均值，選取閾值L∈[0,255]，小于等于L的像素點所占比例為w0，均值為u0；大于L的像素點所占比例為w1，均值為u1，則圖片的灰度均值u為

u=w0×u0+w1×u1

其中p(i)為灰度為i的像素所占比例。

2) 計算類間方差：

u(L)=w0×(u0-u)2+w1×(u1-u)2

(5)

式(5)中，分別取L∈[0,255]，當u取最大值時，分離效果最好，此時的L值為最佳閾值。

3.4 形態(tài)學濾波

二值化處理后的圖像仍然存在一些干擾信號，如掃頻信號、突發(fā)信號以及未清除干凈的噪聲點，在圖像上表現(xiàn)為孤立的或與有用圖像相連的點或線，可以通過形態(tài)學濾波進一步清除干擾，提取有用信號。形態(tài)學濾波的主要運算包括腐蝕、膨脹、開運算和閉運算，通過設定結(jié)構(gòu)元素消除圖像中的孤立點以及填補空洞[27-28]。開運算為先腐蝕運算、再膨脹運算，可定義為：

I1(n,k)=I(n,k)°b=(I(n,k)⊕b)⊙b)

(6)

式(6)中:I(n,k)表示待處理圖像的時頻矩陣; ° 表示開運算; ⊙表示膨脹運算; ⊕表示腐蝕運算;b表示結(jié)構(gòu)元素。

閉運算為先膨脹運算、再腐蝕運算，定義為

I1(n,k)=I(n,k)·b=(I(n,k)⊙b)⊕b)

(7)

式(7)中，·表示閉運算。

在運算過程中結(jié)構(gòu)元素起到重要作用，需要保留的圖形要大于結(jié)構(gòu)元素，需要消除的圖像要小于結(jié)構(gòu)元素。本文中形態(tài)學處理流程如下：

1) 構(gòu)造小尺寸矩形元素，對時頻矩陣進行閉運算，用以彌合裂縫、填充空洞。本文中矩形元素尺寸為3×80。

2) 構(gòu)造水平的線形結(jié)構(gòu)元素，對圖像進行開運算，用以消除掃頻和突發(fā)信號的干擾，從而完成形態(tài)學濾波處理。用于降噪處理的線性結(jié)構(gòu)元素尺寸可以適當放大，應大于矩形元素的水平尺寸，小于跳頻信號每一簇的水平尺寸，因此本文中構(gòu)造的線形元素尺寸為400，能夠有效保留跳頻信號、去除孤立的噪聲點和掃頻干擾。

3.5 參數(shù)估計

3.5.1跳頻周期估計

提取完成形態(tài)學濾波的時頻圖的頻率脊線，對其進行修正，即可獲得跳頻信號的跳頻圖案。提取頻率跳變時刻，組成新的數(shù)組：

T={T1,T2,T3,T4,…,Tn}

對該數(shù)組求一階差分方程，可以得到新的數(shù)組，為各頻率持續(xù)時間的估計值，對新的數(shù)組求尾切平均數(shù)，即除去最大值和最小值求均值，即可得到跳變周期的估計值。

3.5.2頻率集估計

通過k-means聚類算法[29-31]，可以完成對跳頻頻率的估計。

聚類算法首先要預估聚類數(shù)量，在上一步中，可以通過跳變次數(shù)完成對聚類數(shù)量的估計，將圖像點集劃分為k個集合Z={Y1,Y2,Y3,…,Yk}，構(gòu)造目標函數(shù):

(8)

式(8)中:Ci為聚類i中所有點集的均值；Xj為聚類i中的點；Ni為聚類i中點的個數(shù)。

當均值與聚類中心不一致時，均值代替成為新的聚類中心，繼續(xù)計算，直至均值與聚類中心相等，此時均值的頻率坐標即為跳頻頻率估計值。

4 仿真實驗

在高斯白噪聲的條件下，產(chǎn)生的跳頻信號，信噪比為-1 dB，采樣率為1 000 kHz，仿真頻率集[100、150、300、400、200、450、350、470、300、200]kHz，信號長度為100 000點。干擾信號包括兩個定頻信號，其頻率分別為250 kHz和350 kHz，以及一個掃頻信號，掃頻范圍為250～300 kHz，掃頻周期為20 ms。其中350 kHz的干擾信號與跳頻信號的其中一個頻點相重合。

圖2為通過對信號進行短時傅里葉變換產(chǎn)生的時頻圖像，窗函數(shù)為長度225的Hamming窗，在圖像中可以清晰得看到跳頻信號和干擾信號的時頻圖譜。

圖2 接收信號的時頻圖像

圖3為時頻能量對消后的時頻圖譜，在圖像中可以明顯看到，定頻干擾信號已經(jīng)消除，而與定頻干擾信號重合的跳頻信號依然保留了下來。

圖3 對消后的時頻圖像

圖4為二值化處理后的圖像，經(jīng)過全局閾值法能保留大部分信號能量，而大部分能量集中于跳頻信號，因此處理后的二值化圖像消除了大部分白噪聲，保存了大部分的有用信號。

通過開運算、閉運算對二值化圖像進行處理，可以有效去除掃頻信號，去除突發(fā)信號，彌合細小裂縫，得到圖5。提取脊線后得到的跳頻圖案如圖6所示。

圖4 二值化的時頻圖像

圖5 形態(tài)學處理后的時頻圖像

圖6 跳頻圖案

為分析本算法的抗噪聲性能，對測量誤差進行如下規(guī)定：

(16)

對本文的同一仿真信號分別在[-3,5]dB的信噪比條件下進行200次蒙特卡洛仿真實驗，計算參數(shù)估計的誤差。同時，分析在高斯白噪聲、定頻干擾和掃頻干擾背景下本文算法與文獻[5]和文獻[10]的性能對比，其中圖7為3種算法跳頻信號頻率集估計誤差，圖8為3種算法跳頻信號跳頻周期估計誤差。

圖7 頻率集估計誤差

圖8 跳頻周期估計誤差

從仿真結(jié)果來看，本文算法明顯優(yōu)于2個對比算法。由圖7可知：本文算法在跳頻頻率集估計與2個對比算法的估計精度在同一數(shù)量級上，但本算法有更高的精度；由圖8可知，本算法在跳頻周期估計誤差總體在10-7數(shù)量級上，在信噪比大于-1 dB時，跳頻周期估計性能有明顯提升，取得了更好的效果。本算法可以消除各類干擾影響并在較低信噪比條件下完成跳頻信號的參數(shù)估計，獲得更高的估計精度。

5 實測數(shù)據(jù)分析

為驗證算法實用性，采集大疆Mavic Air2型無人機遙控器信號進行算法性能驗證。無人機的遙控器信號頻段為2.4～2.486 GHz，信道帶寬為2 MHz，信號采集采用USRP+GNU Radio平臺的方式實現(xiàn)。由于硬件平臺性能限制，最大采樣率為50 Msps，因此采集信號前先通過頻譜儀觀察信號頻率范圍，再設置采樣中心頻率，以確保能夠覆蓋大部分信號。中心頻率設置為2.42 GHz，采樣率設置為50 Msps，截取567 005個數(shù)據(jù)點，時長11.3 ms的信號進行算法驗證。

無人機遙控器實測信號的時頻圖像如圖9所示，經(jīng)算法處理后的時頻圖像如圖10所示。

圖9 實測信號時頻圖像

圖10 經(jīng)算法處理后的時頻圖像

對實測信號的參數(shù)估計結(jié)果如表1所示。由于實測信號的跳頻周期可變，因此算法對跳頻周期的估計失效。USRP對射頻信號實施了下變頻，得到的頻率參數(shù)估計結(jié)果為基帶信號頻率，由于中心頻率為2.42 GHz，將基帶信號頻率與中心頻率相加，可以得到跳頻信號頻率估計值。

表1 實測信號參數(shù)估計結(jié)果Table 1 Estimation results of real signal parameters

各頻率間隔約為2 MHz，符合說明書中對信道帶寬的描述。通過實測數(shù)據(jù)測試驗證，本算法可以實現(xiàn)跳頻信號參數(shù)估計，具有實用性。

6 結(jié)論

本文提出一種基于時頻分析的跳頻信號參數(shù)估計算法，通過設計仿真實驗檢驗了本算法能夠在較低信噪比條件下具有較高的估計精度，通過實測數(shù)據(jù)的測試表明本算法具有實用性。