張 瑋,王 平,解西坤,3
(1.海軍工程大學(xué) 電子工程學(xué)院,武漢 430033;2.中國(guó)人民解放軍92038部隊(duì),山東 青島 266041;3.中國(guó)人民解放軍91202部隊(duì),遼寧 葫蘆島 125000)
跳頻[1-2](Frequency Hopping,FH)通信是最常用的擴(kuò)頻通信的方式之一,其信號(hào)的載波頻率可隨時(shí)間按偽隨機(jī)規(guī)律跳變,具有較好的抗干擾、抗截獲、抗多徑能力,在軍事和民用領(lǐng)域都有廣泛的應(yīng)用。獲得準(zhǔn)確的跳頻信號(hào)參數(shù)是通信對(duì)抗領(lǐng)域的重要內(nèi)容[3],因此,跳頻參數(shù)估計(jì)逐漸成為信息對(duì)抗領(lǐng)域的熱點(diǎn)問(wèn)題。
目前跳頻信號(hào)參數(shù)估計(jì)算法可主要分為時(shí)頻分析法[4]和非時(shí)頻分析法[5],而時(shí)頻分析法是當(dāng)前跳頻信號(hào)分析的最常用的方法。時(shí)頻分析法是通過(guò)時(shí)間和頻率兩個(gè)維度來(lái)描述能量密度變化的方法。文獻(xiàn) [6]提出了一種結(jié)合圖像處理的跳頻信號(hào)參數(shù)估計(jì)算法,采用一種迭代匹配算法抑制復(fù)雜的毛刺干擾信號(hào)。文獻(xiàn)[7]提出了一種基于最大值的慢跳頻信號(hào)檢測(cè)算法,根據(jù)統(tǒng)計(jì)各信道輻射強(qiáng)度判斷是否存在跳頻信號(hào)。文獻(xiàn)[8]以傳統(tǒng)能量對(duì)消方法為基礎(chǔ)進(jìn)行改進(jìn),能夠在完成信號(hào)檢測(cè)的同時(shí)具有較好的魯棒性和更小的計(jì)算量。但是,上述三種算法在低信噪比時(shí)性能較差。文獻(xiàn)[9]提出了一種改進(jìn)時(shí)頻脊線的跳頻參數(shù)估計(jì)算法,抗噪聲性能較好。文獻(xiàn)[10]提出了一種提出基于有噪模板的互相關(guān)監(jiān)測(cè)方案,能夠解決快跳頻系統(tǒng)檢測(cè)難算法性能低的問(wèn)題。但是,上述兩種算法在含有其他突發(fā)干擾時(shí)性能較差。文獻(xiàn)[11]采用修正時(shí)頻圖和二次估計(jì)的方法實(shí)現(xiàn)跳頻頻率和跳變時(shí)刻的估計(jì),具有較好的估計(jì)精度和抗噪聲性能,但是僅能在高斯白噪聲的條件下完成,存在其他干擾時(shí)算法性能下降較快。
為解決以上問(wèn)題,本文提出一種基于自相關(guān)和時(shí)頻分析的跳頻參數(shù)估計(jì)算法。首先,通過(guò)基于能量檢測(cè)的分段自相關(guān)算法對(duì)信號(hào)進(jìn)行預(yù)處理預(yù)處理,利用短時(shí)傅里葉變換(Short-Time Fourier Transform,STFT)得到時(shí)頻圖像;然后,采用二值化和形態(tài)學(xué)濾波進(jìn)行降噪提取;最后,通過(guò)聚類(lèi)算法完成跳頻信號(hào)參數(shù)估計(jì)。
在觀測(cè)時(shí)間T內(nèi)接收到N個(gè)跳頻簇,則跳頻信號(hào)的數(shù)學(xué)模型可進(jìn)行如下定義:
(1)
(2)
式中:rectTh為門(mén)函數(shù);τ0為起跳時(shí)刻;Th為跳頻周期;fk為第k個(gè)跳頻簇的跳頻頻率;N是觀察時(shí)間內(nèi)頻率跳變數(shù)量。
接收端信號(hào)r(t)中除了包含跳頻信號(hào)s(t)外,還包括干擾信號(hào)J(t),可表示為
r(t)=s(t)+J(t)。
(3)
自相關(guān)是信號(hào)在不同時(shí)刻的數(shù)值相關(guān)性的一種描述,是信號(hào)時(shí)域處理的重要方法[12]。自相關(guān)函數(shù)定義如下:
(4)
式中:x(t)為待處理信號(hào);τ為自相關(guān)函數(shù)的延時(shí)值;T為信號(hào)周期。
根據(jù)自相關(guān)函數(shù)的性質(zhì),隨著延時(shí)τ的增加,高斯噪聲的自相關(guān)值隨之快速衰減,并逐漸趨近于零,而具有周期特性的信號(hào)自相關(guān)數(shù)值不會(huì)隨著延時(shí)τ的增加而衰減,其周期特性會(huì)有效保留,因此可以采用自相關(guān)的方法對(duì)信號(hào)實(shí)施預(yù)處理,實(shí)現(xiàn)信號(hào)降噪。
跳頻信號(hào)不是周期信號(hào),但是其每一簇信號(hào)都是周期性的,因此需要對(duì)信號(hào)進(jìn)行分段,算法的具體流程如下:
1)基于能量檢測(cè)的時(shí)域分段
通過(guò)能量檢測(cè)判斷跳頻信號(hào)數(shù)量,若已知跳頻周期,時(shí)域分段的長(zhǎng)度應(yīng)與跳頻周期長(zhǎng)度相等,此時(shí)自相關(guān)運(yùn)算效果最佳。在盲估計(jì)的情況下,設(shè)計(jì)一種時(shí)域分段算法,具體的執(zhí)行流程如下:
輸入:步長(zhǎng)T0
輸出:段長(zhǎng)L
1 初始化i=1,從信號(hào)起點(diǎn)截取段長(zhǎng)L=i×T0的信號(hào)x1(n),此時(shí)段長(zhǎng)L=T0;
2 求取截取信號(hào)x1(n)的功率譜X1(m),求取X1(m)的最大值;
3 設(shè)置i=i+1;
4 重復(fù)步驟1~2,截取段長(zhǎng)為L(zhǎng)=i×T0的信號(hào)xi(n),計(jì)算截取信號(hào)功率譜Xi(m)的最大值;
5 判斷Xi(m)的最大值是否繼續(xù)升高,如果剛好不再升高,則輸出此時(shí)的段長(zhǎng)L=i×T0,如果最大值繼續(xù)升高,則重復(fù)步驟4,直至最大值不再升高,然后輸出此時(shí)的段長(zhǎng)L=i×T0。
信號(hào)截取時(shí)刻的功率譜如圖1和圖2所示。觀察圖像可知,當(dāng)截取信號(hào)長(zhǎng)度為跳頻周期時(shí),功率譜圖的最大值為第一簇跳頻信號(hào)的頻點(diǎn),截取長(zhǎng)度大于跳頻周期時(shí),第一簇信號(hào)頻點(diǎn)的信號(hào)強(qiáng)度不再變化,而新截取的部分第二簇信號(hào)所在頻點(diǎn)的信號(hào)強(qiáng)度隨著截取長(zhǎng)度變長(zhǎng)而變大,因此當(dāng)功率譜最大值不再發(fā)生變化時(shí)的信號(hào)長(zhǎng)度即為時(shí)域分段的段長(zhǎng)。
圖1 第一簇信號(hào)功率譜
圖2 過(guò)量截取信號(hào)功率譜
2)自相關(guān)運(yùn)算
使用公式(4)對(duì)每一段信號(hào)進(jìn)行自相關(guān)運(yùn)算,得到每一段削弱噪聲的時(shí)域信號(hào)。
3)信號(hào)拼接
將每一段信號(hào)拼接起來(lái),可以得到完整的預(yù)處理信號(hào)。
因?yàn)槎虝r(shí)傅里葉變換具有較好的時(shí)頻聚集性,因此通過(guò)STFT的方式進(jìn)行時(shí)頻變換,可得到信號(hào)的時(shí)頻圖像。
圖3為接收端信號(hào)未經(jīng)預(yù)處理直接時(shí)頻變換得到的時(shí)頻圖像,設(shè)置信噪比為-7 dB,圖4為經(jīng)過(guò)預(yù)處理降噪后的時(shí)頻圖像,可以看出經(jīng)處理后消除了大部分高斯白噪聲影響。在圖4中可以看到豎直的細(xì)線形狀的噪聲,是因?yàn)楫?dāng)τ=0時(shí)噪聲的自相關(guān)值最大,每?jī)蓷l細(xì)線的間隔為選取的步長(zhǎng),掃頻信號(hào)具有周期特性,但是由于其周期大于步長(zhǎng),因此在一個(gè)步長(zhǎng)內(nèi)的掃頻信號(hào)并無(wú)周期特征,自相關(guān)運(yùn)算也消除了掃頻信號(hào)帶來(lái)的干擾。
圖3 接收端信號(hào)時(shí)頻圖像
圖4 預(yù)處理后信號(hào)時(shí)頻圖像
經(jīng)過(guò)時(shí)頻變換后,對(duì)時(shí)頻圖像進(jìn)行降噪提取,主要有以下兩個(gè)步驟:
1)二值化濾波
為獲得更精確的參數(shù)估計(jì)值,需要獲取清晰的時(shí)頻圖像,因此需要對(duì)時(shí)頻圖像作進(jìn)一步降噪處理。通常采用閾值法對(duì)圖像進(jìn)行二值化,大于閾值μ的點(diǎn)賦值為1,小于等于閾值 的點(diǎn)賦值為0。計(jì)算過(guò)程定義如下:
(5)
閾值μ可采用迭代法[9]來(lái)計(jì)算。首先選取初始閾值μ1,由時(shí)頻矩陣STFT的最大值STFTMax和最小值STFTMin決定:
μ1=STFTMax+STFTMin。
(6)
然后以μ1為閾值將時(shí)頻矩陣分為STFT1和STFT2兩個(gè)部分,其中,STFT1表示大于閾值μ1的部分,STFT2表示小于閾值μ1的部分,分別計(jì)算各部分的均值,表示為s1和s2:
(7)
式中:N1為STFT1中點(diǎn)的個(gè)數(shù);N2為STFT2中點(diǎn)的個(gè)數(shù)。
根據(jù)下式計(jì)算得到新的閾值μ2:
μ2=(s1+s2)/2。
(8)
進(jìn)行迭代運(yùn)算不斷獲得新的閾值,直至μk+1=μk時(shí),結(jié)束運(yùn)算并獲得最終閾值,并將閾值μk+1代入公式(5)進(jìn)行降噪,得到二值化圖像。圖5為經(jīng)過(guò)二值化處理后的時(shí)頻圖像。
圖5 二值化圖像
2)形態(tài)學(xué)濾波
在干擾較強(qiáng)的情況下,二值化處理不能將所有干擾處理干凈,同時(shí)也會(huì)有部分信號(hào)被當(dāng)作干擾清除而導(dǎo)致跳頻簇的時(shí)頻圖像不完整,形態(tài)學(xué)濾波可以較好地解決這個(gè)問(wèn)題。形態(tài)學(xué)濾波是通過(guò)設(shè)定不同的結(jié)構(gòu)元素對(duì)圖像進(jìn)行運(yùn)算進(jìn)而獲取圖像的局部特征[13]。膨脹運(yùn)算、腐蝕運(yùn)算、開(kāi)運(yùn)算和閉運(yùn)算是形態(tài)學(xué)濾波的基本運(yùn)算,其中,開(kāi)運(yùn)算可以消除干擾、平滑圖像,閉運(yùn)算可以填補(bǔ)空洞、彌補(bǔ)圖像。開(kāi)運(yùn)算為先進(jìn)行腐蝕運(yùn)算再進(jìn)行膨脹運(yùn)算,可做如下定義:
Ft1(n,k)=Ft(n,k)°b=(Ft(n,k)⊕SE)⊙SE)。
(9)
式中:Ft(n,k)表示待處理圖像的時(shí)頻矩陣;°表示開(kāi)運(yùn)算;⊙表示膨脹運(yùn)算;⊕表示腐蝕運(yùn)算;SE表示結(jié)構(gòu)元素。
閉運(yùn)算為先進(jìn)行膨脹運(yùn)算,再進(jìn)行腐蝕運(yùn)算,定義為
Ft1(n,k)=Ft(n,k)·b=(Ft(n,k)⊙SE)⊕SE)。
(10)
式中:·表示閉運(yùn)算。
構(gòu)造結(jié)構(gòu)元素是形態(tài)學(xué)濾波的關(guān)鍵環(huán)節(jié),結(jié)構(gòu)元素應(yīng)大于需要消除的部分、小于需要保留的部分,這樣才能獲得較好地濾波效果。因此,本文中形態(tài)學(xué)濾波的處理流程和結(jié)構(gòu)元素構(gòu)造策略如下:
①對(duì)時(shí)頻矩陣進(jìn)行閉運(yùn)算,用于彌合裂縫、填補(bǔ)空洞,因此構(gòu)造小尺寸的矩形元素,本文中的矩形元素尺寸為3×100。
②對(duì)時(shí)頻矩陣進(jìn)行開(kāi)運(yùn)算,用于平滑圖像,消除未處理干凈的噪聲點(diǎn),完成形態(tài)學(xué)濾波處理。跳頻簇在時(shí)頻矩陣中表現(xiàn)為水平狀態(tài),因此需要構(gòu)造水平線型元素,尺寸小于跳頻簇的水平尺寸,大于噪聲點(diǎn)的水平尺寸,本文中的線性元素尺寸為300。
圖6為經(jīng)過(guò)形態(tài)學(xué)濾波處理后的時(shí)頻圖像,可以看到基本完成跳頻信號(hào)的提取,為接下來(lái)的參數(shù)估計(jì)做好了準(zhǔn)備。
圖6 形態(tài)學(xué)濾波后的時(shí)頻圖像
采用聚類(lèi)算法思想完成信號(hào)參數(shù)估計(jì),重點(diǎn)在于對(duì)各個(gè)跳頻簇聚類(lèi)中心的計(jì)算。計(jì)算流程如下:
1)采用8連通標(biāo)記算法完成對(duì)時(shí)頻矩陣的區(qū)域劃分,總計(jì)k簇跳頻信號(hào)。圖7為區(qū)域聯(lián)通標(biāo)記圖像。
圖7 區(qū)域聯(lián)通標(biāo)記時(shí)頻圖像
2)在跳頻信號(hào)完成聯(lián)通標(biāo)記后,所有數(shù)據(jù)點(diǎn)被標(biāo)記賦值,假設(shè)時(shí)頻矩陣內(nèi)數(shù)值為1的點(diǎn)的數(shù)量為N,提取時(shí)頻矩陣內(nèi)所有數(shù)值為1的點(diǎn)構(gòu)成集合P={P1(t1,f1),P2(t2,f2),…,PN(tN,fN)},假設(shè)該跳頻信號(hào)共計(jì)k個(gè)跳頻簇,按照如下公式計(jì)算第i跳跳頻簇聚類(lèi)中心:
(11)
3)所有聚類(lèi)中心構(gòu)成集合C={C1(t1,f1),C2(t2,f2),…,Ck(tk,fk)},此時(shí)輸出各聚類(lèi)中心的頻率坐標(biāo)作為跳頻頻率參數(shù)估計(jì)值:
(12)
T(i)=(ti+1-ti)|i∈[1,k-1],
(13)
(14)
綜上所述,本文算法具體步驟如下:
Step1 通過(guò)自相關(guān)算法對(duì)接收端信號(hào)進(jìn)行降噪預(yù)處理,可去除白噪聲和不含有周期特性的干擾信號(hào)。
Step2 采用STFT算法對(duì)信號(hào)進(jìn)行時(shí)頻變換,得到時(shí)頻矩陣。
Step3 通過(guò)二值化和形態(tài)學(xué)濾波,完成對(duì)時(shí)頻圖像的降噪提取。
Step4 采用聚類(lèi)算法的思想而完成跳頻信號(hào)參數(shù)估計(jì)。
算法總體流程框圖如圖8所示。
圖8 算法流程
為分析在不同信噪比條件下的算法性能,對(duì)測(cè)量誤差做如下規(guī)定:
(15)
以下實(shí)驗(yàn)結(jié)果均是200次蒙特卡洛實(shí)驗(yàn)結(jié)果的平均值。
設(shè)置采樣率為1 000 kHz,頻率集為[410,350,330,400,360,440,340,430,370,350]kHz,跳頻周期為10 ms;設(shè)置掃頻干擾,掃頻范圍250~300 kHz,掃頻周期為20 ms;信號(hào)長(zhǎng)度為100 000點(diǎn),仿真時(shí)長(zhǎng)100 ms,比較在不同強(qiáng)度高斯白噪聲條件下本文算法與文獻(xiàn)[6]和文獻(xiàn)[9]算法的性能。
3種算法對(duì)跳頻頻率的估計(jì)誤差如圖9所示,對(duì)跳頻周期的估計(jì)誤差如圖10所示。由圖9和圖10可見(jiàn),本文算法在跳頻頻率估計(jì)和跳頻周期估計(jì)方面都優(yōu)于對(duì)比算法,而本文算法在信噪比低于-10 dB時(shí)性能大幅度下降,而文獻(xiàn)[6]算法在低于-4 dB時(shí)性能開(kāi)始衰減,文獻(xiàn)[9]算法在低于-6 dB時(shí)性能開(kāi)始衰減,因此本文算法在抗噪聲性能方面更加穩(wěn)定。
圖9 跳頻頻率的估計(jì)誤差
圖10 跳頻周期的估計(jì)誤差
由于查閱的資料中,很少有采用自相關(guān)運(yùn)算對(duì)信號(hào)進(jìn)行預(yù)處理,因此本實(shí)驗(yàn)將本算法有無(wú)自相關(guān)預(yù)處理作為變量之一,檢驗(yàn)預(yù)處理過(guò)程對(duì)算法性能的影響。跳頻信號(hào)的參數(shù)同3.1節(jié),有無(wú)預(yù)處理的參數(shù)估計(jì)誤差如圖11所示。
圖11 有無(wú)預(yù)處理的參數(shù)估計(jì)誤差
由圖11可見(jiàn),添加自相關(guān)預(yù)處理后,參數(shù)估計(jì)精度有不同程度提高;無(wú)預(yù)處理時(shí),跳頻頻率估計(jì)性能在信噪比低于-7 dB后衰減較快,跳頻周期估計(jì)無(wú)明顯變化,由此可以判斷自相關(guān)預(yù)處理可以使跳頻參數(shù)估計(jì)在抗噪聲性能方面有明顯提升。
本文提出了一種基于自相關(guān)和時(shí)頻分析的跳頻參數(shù)估計(jì)算法,并通過(guò)仿真實(shí)驗(yàn)對(duì)算法性能進(jìn)行了驗(yàn)證。本文算法通過(guò)對(duì)接收端信號(hào)進(jìn)行自相關(guān)運(yùn)算,有效消除高斯白噪聲等無(wú)周期特性的噪聲及干擾,通過(guò)STFT對(duì)信號(hào)進(jìn)行時(shí)頻變換,采用二值化和形態(tài)學(xué)濾波的方法完成降噪提取,最后通過(guò)聚類(lèi)算法完成跳頻參數(shù)估計(jì)。仿真實(shí)驗(yàn)表明,本文算法在低信噪比時(shí)具有較高的估計(jì)精度,同時(shí)驗(yàn)證了自相關(guān)運(yùn)算可較大幅度提升算法性能。