劉志鵬, 田楨熔, 毛 重, 張國毅
(1.空軍航空大學 信息對抗系,吉林 長春 130022;2.空軍航空大學 飛行訓練基地第一飛行團機務大隊,吉林 長春 130022;3.空軍航空大學 訓練部信息管理中心,吉林 長春 130022)
采用傳統的基于5大參數對截獲輻射源信號進行分選和識別的方法適用于雷達數量少、脈沖流密度低、信號樣式簡單、信號參數固定的電磁環境,而隨著現代復雜體制雷達大量出現,其性能下降明顯。這主要因為新體制雷達的參數如RF、PRI等具有跳變乃至偽隨機捷變的能力,加之雷達數量的增多,使不同信號在傳統參數域交疊嚴重。解決思路之一就是針對復雜的雷達信號本身提取出新的分選和識別特征,為輻射源的分選和識別提供新的依據。近年來提出的瞬時頻率特征、復雜度特征及小波包特征等都是典型的基于脈內分析的分選特征,但都存在著信噪比要求高、計算復雜等局限。
維 格 納 分 布(Wigner-Ville Distribution,WVD)作為一種典型的Cohen類時頻分析方法,廣泛地應用于非平穩信號的分析與處理,是一種重要的雙線性時頻分布,具有分辨率高、能量集中以及優良的時頻邊緣特性等品質[1-2]。但WVD為二維分布,需要對其降維并進一步提取出更為有效的分選特征。
文中從信號的WVD特征著手,首先利用自適應短時傅里葉窗(Adaptive Short-Time Fourier Transform,ASTFT)抑制WVD的交叉項干擾[3],然后提取出WVD的主脊線,最后借鑒描述序列復雜度的符號熵方法,計算主脊線的符號熵特征作為新的分選特征。通過各種仿真及實際信號的測試可知,新特征在不同信噪比下性能穩定,且具有良好的可區分性。
信號s(t)的WVD定義為[4]:

式中:τ——時間差變量;
ω——頻率。
在時頻和頻域上對信號WVD變換進行積分,得到信號的能量E,即:

從上式可知,WVD的物理意義為信號能量在時域和頻域中的分布[5]。
注意式(1)中信號s(t)在右端出現了兩次,因而WVD屬于雙線性時頻分布。當利用WVD分析多分量信號時,結果會出現交叉項,從而使WVD的解釋變得困難,交叉項成了WVD的噪聲。ASTFT譜抑制的方法具有較好的時頻邊緣分辨率,且能有效抑制交叉項,文中利用ASTFT譜抑制來對WVD性能進行改善,并利用改進的WVD進行脈內特征的提取。
若交叉項與自項進行內積運算,其分量的正負部分將相互中和,其結果接近零。即

于是,利用rA(t,τ)與rs(t,τ)進行互相關,即可以消除rs(t,τ)中的交叉項rC(t,τ),再對τ進行傅里葉變換,從而得到消除交叉項WVDs,cross(t,ω)的WVD分布WVDs,Mod(t,ω)

由于信號的自項WVDs,auto(t,ω)不能預先知道,因此考慮利用自適應STFT的時頻譜ASTFTs(t,ω)來代替,即令

于是有

選取BPSK信號說明利用ASTFT抑制WVD交叉項的效果,信號采用13位巴克碼,歸一化后頻率為0.1,采樣點數為512,信噪比為20dB,對該信號進行時頻分析的結果如圖1所示。

圖1 13位巴克碼信號的時頻分析
圖1(a)為原始WVD,盡管有較好的能量聚集性和較高的分辨率,但產生的交叉項干擾十分明顯;圖1(b)為利用ASTFT譜抑制交叉項的WVD的效果,其分辨率增高,毛刺減少,不僅保留了WVD較好的能量聚集性和較高的分辨率,而且有效抑制了交叉項。
由于所求的WVD為三維分布幅度譜,即使一個64個采樣點的信號,也需要4 096個像素點進行表示,如果使用原始數據直接用作分選特征,顯然維數太大,因此必須對原始WVD分布特征進行降維處理。
因為WVD反映了信號在頻域和時域的能量分布,WVD圖像上的峰值點即為能量的主要分布點,連接每一采樣時刻的幅度峰值點即形成了一條表示能量分布的曲線,把這條曲線段稱為WVD能量主脊(WVD Main Ridge,WVDMR)。由WVD變換的特點可知,信號中混雜的噪聲經互相關計算后能量較低,因此從能量高點提取出的WVDMR具有良好的抗噪性,幾乎完全反應了信號本身的特點。分別作WVDMR在時間-頻率平面和幅度-時間平面的投影,得到WVDMR1和WVDMR2,則WVDMR1和WVDMR2分別描述了WVD主要能量分布在時間-頻率平面和幅度-時間平面的走向,因而可以很好地反映出WVD的特征,其本質上是WVD圖像的二值化。
6種典型雷達信號在信噪比10dB下,采用點數為512點時的WVDMR1和WVDMR2分別如圖2和圖3所示。

圖2 6種典型信號的WVDMR1


圖3 6種典型信號的WVDMR2
其中,CON信號的歸一化頻率為0.2;LFM信號的起始歸一化頻率為0.1,終止頻率為0.15;NFLM信號基于正切波形,歸一化帶寬為0.15,時間副瓣電平因子為2.5;Costas信號的起止歸一化頻率分別為0.2和0.4,調頻序列為(10 9 7 3 6 1 2 4 8 5);BPSK采用13位巴克碼,頻點為0.1;Frank碼信號每個碼組有64個相位,初始歸一化頻點為0.2。
可見,不同脈內調制的信號WVDMR1和WVDMR2具有較為明顯的差異,一般情況下,因為WVDMR1反映的是能量在時頻平面的分布走勢,頻率編碼信號和相位編碼信號的曲線相比復雜度更高,而線性調頻信號和常規信號的曲線較為規律;WVDMR2體現了能量峰值隨時間變化的趨勢,不同調制類型的信號曲線震蕩頻率和幅度不同,表現的復雜程度也不同。考慮到WVDMR1和WVDMR2基于不同類型的信號在不規則度上表現出較大相異性,文中將討論提取特征脊線的熵特征作為分選新特征補充。
由于不同類型的信號在脊線形狀的復雜度上表現出較大差異,為了度量不同信號脊線形狀的復雜程度,引入符號化時間序列的概念。符號化時間序列分析是一種新的信息分析方法,由符號動力學理論、信息理論和混沌時間分析理論綜合演變發展而來,具有運算簡捷、對噪聲不敏感以及抑制測量設備干擾等優點,現已在故障診斷、多相流測量以及股票分析等領域得到廣泛運用[6-7]。
符號化分析的基本思想是將原始數據狀態空間劃分為一系列區間,每個區間對應不同的符號,通過數據落入的不同區間,為其分配對應的符號,同一區間的數據用相同符號表示,即在幾個可能值上對時間序列進行重采樣,從而把原始數據序列空間變換到離散的符號空間上進行處理。該過程可看作是“粗粒化”的過程,可以較好地捕獲大尺度特征,從而降低噪聲對時間序列的影響[8]。時間序列符號化分析可歸結為兩部分:時間符號序列轉換和統計學定量分析。
2.1.1 時間符號序列轉換
設原始時間序列為X,長度為n,即X={xi|i=1,2,…,n},先對X按式(6)進行分段集成降維重采樣。

從而得到X′={xj|j=1,2,…,N},N?n,令τ=,τ即為采樣時延。
引入m=r+1維符號集合Q={0,1,…,r},依照式(7)將原始數據序列分配對應符號。

式中:z1,z2,…,zr——原始數據劃分區間閾值;
sj——每個區間所對應的符號。
于是得到符號序列S={s1,s2,…,sN}。但是m進制的符號序列S仍無法對原始數據的復雜程度做出有效描述,因此采用將符號序列編碼為十進制序列來捕捉序列的局部和總體特征。

式中:k=1,2,…,N-L+1。
將子序列s′(k)轉化為十進制數序列Sd(k),轉稱為字。

其中,s′(k)=s(k+i-1),轉換過程如圖4所示。

圖4 時間序列轉化為符號序列
此時L=3,將長度為26的二進制序列轉化為長度為24的由8個十進制字構成的符號序列。
最終,原始數據序列轉換成長度為N-L+1的十進制符號序列SD={Sd(k)|k=1,2,…,NL+1}。此時的符號序列用最簡單的方式描述了原始數據系統的動力學特征和其軌跡的復雜程度。
2.1.2 符號序列的統計學分析
由符號序列轉化的過程可知,符號序列不同的字代表不同的子序列結構[9],因此,我們可以通過對序列中字的統計學特征來描述出系統的結構特點。文中利用修正的香農熵來對符號序列字進行統計學處理,其定義為

式中:N′——符號序列中出現的總字符數量;
pi——第i個字出現的概率。
易見,0≤Hs≤1。
2.2.1 時間序列的符號熵特征
由以上分析可知,將時間序列符號化,利用符號熵特征可有效地描述系統的復雜程度,可以用來作為WVD的脊線特征,但需確定時延參數τ,符號集合維數m和字長L。時延參數τ決定了符號序列的采樣尺度,若τ取值過小,則失去了符號化方法大尺度反映總體特征的能力,同時放大了噪聲的干擾;τ值取的過大,則會造成符號化方法捕捉的特征“失真”。τ值的選取依據離散混沌系統相空間重構中確定時延所采用的自相關方法,其表示式為

式中:μ——序列均值;
σ2——序列的方差。
一般取自相關函數下降到1~1/e時對應的τ值。
符號集合維數m和子序列長度L決定了描述系統特征的精細程度和待分析結構尺度的大小。在工程實踐中,一般通過多次試驗不同參數來尋找m和L的最優或近優取值。
2.2.2 WVD脊線符號熵特征提取方法
綜上所述,WVD脊線熵特征提取方法總結如下:
1)選取τ參數進行重采樣。對信號進行n點采樣,按頻率f0進行歸一化,對預處理的信號進行改進WVD變換;然后提取WVD圖像峰值點,每個點對應的歸一化頻率構成WVDMR1時間序列,其對應的幅度序列構成WVDMR2時間序列,將得到的時間序列進行重采樣,得到N點的時間序列。
2)選取m和L,將重采樣的時間序列轉換為十進制字符號序列。劃分符號區間,確定區間閾值,將時間序列轉換為符號序列,提取子序列,并完成十進制字轉化,形成十進制字符號序列。
3)統計各個符號字出現頻度,計算信號的WVDMR1和WVDMR2特征符號熵。
仿真環境:Windows7,AMD CPU A6,4GB內存,編程工具為MATLAB7.10.1。
選擇1.2中采用的6種典型的雷達信號類型進行仿真實驗:常規信號CON、線性調頻信號LFM、非線性調頻信號NLFM、Costas編碼信號、BPSK信號及Frank編碼信號,采樣點數均為512,信號其他參數設置同1.2。首先產生分析樣本,對每一種類型的雷達信號在-5~30dB的范圍內每隔5dB產生20個不同初相的信號,提取每個信號的WVDMR1和WVDMR2特征,每種特征對應160個樣本。
3.1.1 WVDMR1特征的符號熵
選取不同符號化參數,對上述樣本經1 000次Monte Carlo仿真實驗可以得出,對于信號的WVDMR1特征,當取τ=4,m=6,L=8時,符號熵在不同信噪比下最為穩定,且不同信號間的區分度最大。
在τ=4,m=6,L=8情況下,6種典型信號在信噪比-5~30dB范圍內變化時的符號熵取值情況如圖5所示。

圖5 6種典型信號的WVDMR1符號熵
由圖中可直觀看出,常規信號(CON)的WVDMR1特征由于近似為一條水平直線,其符號熵最小,接近于0;而采用Costas編碼信號的WVDMR1特征規律性最差,其對應的符號熵值也最大;不同類型信號的符號熵值分布曲線保持一定間隔,互不交疊,且幾乎平行于橫軸,證明在參數取值為τ=4,m=6,L=8時,信號WVDMR1的符號熵可作為有效的分選特征。
3.1.2 WVDMR2特征的符號熵
對于WVDMR2特征,通過1 000次Monte Carlo實驗得出了當τ=8,m=3,L=5時,符號熵特征在不同信噪比下最為穩定,且不同信號間的區分度最大。6種典型信號在信噪比-5~30dB范圍內變化時的WVDMR2符號熵變化值情況如圖6所示。

圖6 6種典型信號的WVDMR2符號熵
可以看出,常規信號(CON)的WVDMR2特征較為規律,因而其符號熵最小,接近于0.2;而Costas編碼信號的WVDMR2特征規律性最差,其符號熵值最大;不同類型信號的符號熵值分布曲線互不交疊,與橫軸近似平行,具有較好的可分性。
3.1.3 新特征統計特性
通過MATLAB仿真得到不同參數的上述6種典型信號,WVDMR1的符號化參數選取τ=4,m=6,L=8,WVDMR2選取τ=8,m=3,L=5,分別進行1 000次Monte Carlo實驗,得到新特征的統計特性見表1。

表1 WVD符號熵的統計特性
綜合以上仿真實驗分析可知,適當選取符號化參數,可以從改進維格納分布脊線特征WVDMR1,WVDMR2得到穩定且可分性較強的符號熵特征,利用這些特征可以較好地區分不同脈內調制類型的信號,且新特征對噪聲不敏感。
實驗數據采用實際偵察數據,經人工分析,頻率歸一化后的信號類型及相應的參數見表2。

表2 實驗數據參數
對實際信號進行重采樣,提取WVD主脊線WVDMR1,WVDMR2并計算符號熵,參數設置同2.2節,得到WVDMR1,WVDMR2的符號熵二維分布如圖7所示。
從圖7可明顯看出新特征良好的可分性,幾乎不會產生交疊。繼續采用不同信號,利用最近鄰分類器進行分類,200次Monte Carlo實驗后得到的分類準確率見表3。

圖7 實際信號的符號熵二維分布

表3 對不同信號分類正確率
從表中可以看出:對實際信號利用新特征進行分類的準確率可達96.5%以上,CON的平均識別率可達99.2%;NLFM和BPSK的分類正確率稍低,這是二者因為噪聲干擾會在分類邊界有微量交疊,但平均正確率也分別到達了96.5%和97.0%。實驗結果證明了改進WVD脊線符號熵特征作為新特征在實際分選中的有效性。
從信號的維格納分布(WVD)入手,提出了一種有效描述時頻能量特征的方法。首先,對信號進行WVD變換,并利用ASTFT譜對信號的WVD進行加窗處理,在保證時頻分辨率的同時有效消除交叉項;但無法直接利用時頻圖像作為分選特征,須對其進行簡化。沿著時間軸方向提取每一時刻對應的能量最高點,形成WVD主脊曲線(WVDMR),對WVDMR分別做時間-頻率平面和時間-能量平面的投影,形成特征曲線WVDMR1和WVDMR2;計算WVDMR1和WVDMR2符號熵,把主脊線的符號熵作為描述時頻能量分布的新特征,從而有效地描述了信號的脈內調制特征。仿真和實際數據實驗表明,新特征具有良好的抗噪性、可分性及穩定性。
[1] 林炎,張友益.超級WVD對多分量LFM信號參數的估計[J].現代雷達,2014,36(1):47-51.
[2] 張靜,李梁,廖曦,等.基于WVD的雷達信號瞬時頻率提取方法研究[J].魯東大學學報:自然科學版,2013,29(1):38-41.
[3] 程發斌,湯寶平,鐘佑明.利用ASTFT譜有效抑制WVD交叉項的方法[J].電子與信息學報,2008,12(10):2299-2302.
[4] 朱曉華.雷達信號分析與處理[M].北京:國防工業出版社,2011:109-193.
[5] 宋宇,游海龍,翁新武,等.基于希爾伯特-黃變換的信號方法[J].長春工業大學學報,2015,36(5):374-378.
[6] 張雨,胡篤慶.基于符號樹信息熵的機械振動瞬態信號特征提取[J].國防科技大學學報,2003,25(4):79-81.
[7] Moon F C.Chaotic and fractal dynamics[M].New York:John Wiley and Sons,1992.
[8] Tang X Z,Tracye R,Brown R.Symbol statistics and spatio-temporal systems[J].Physical Review D,1997,10(2):253-261.
[9] Abarbanel H D I.Analysis of observed chaotic data[M].New York:Springer-Verlag,1996.