李立元 栗蘋 李國(guó)林 張廣為
(北京理工大學(xué) 機(jī)電動(dòng)態(tài)控制重點(diǎn)實(shí)驗(yàn)室, 北京 100081)
近年來,在各代干擾機(jī)、微波武器、大功率信號(hào)等一系列的干擾下,無線電引信極有可能因此出現(xiàn)早炸或失效的情況。 太赫茲波的工作頻率避開了高功率微波武器、雷達(dá)、干擾機(jī)等常用的工作頻段,再加上太赫茲波在大氣中的衰減較大[1-2],主要適用于近程探測(cè),因此太赫茲波應(yīng)用于引信之上具有天然的優(yōu)勢(shì)[3]。
對(duì)地的無線電引信獲得的回波能量取決于地面散射系數(shù),而地面散射系數(shù)又取決于地面類型、落角等一系列因素[4-5]。 面對(duì)不同類型的地面,回波能量具有極大的差異。 而在很長(zhǎng)一部分中國(guó)邊境線的高原上,由于空氣的稀薄,同樣對(duì)回波能量的大小造成了影響。 綜上所述,獲取太赫茲信號(hào)在高原不同地貌下的回波特性,并實(shí)現(xiàn)對(duì)回波數(shù)據(jù)快速分類從而確定炸高,對(duì)于保證引信的毀傷效能顯得尤為重要。
因此,實(shí)地采集了高原各種地貌下的回波數(shù)據(jù),本文僅僅對(duì)灌木地形下的回波雙譜用4 種積分雙譜進(jìn)行降維與分類,從而在占用很少資源的情況下實(shí)現(xiàn)對(duì)引信炸高的識(shí)別。
高階譜是通過累積量函數(shù)(cum)來定義的,因此也被稱作累積量譜。 假如是一個(gè)零均值的k階平穩(wěn)過程,則將該過程的k階累積量定義為
式中:ω為傅里葉變換到頻域得到的頻率分量;功率譜、雙譜(三階譜)、三譜(四階譜)分別為k=2,3,4 時(shí)的特殊情況。
本文所用雙譜(三階累積量譜)為
雙譜是二維函數(shù),可以直接通過卷積神經(jīng)網(wǎng)絡(luò)(CNN)等方式直接進(jìn)行分類,但由于數(shù)據(jù)計(jì)算量過于龐大,嚴(yán)重限制了雙譜在引信上的應(yīng)用。因此,可以按照一定的規(guī)則對(duì)雙譜值進(jìn)行積分,將二維函數(shù)轉(zhuǎn)化為一維函數(shù)[6],從而大大減少分類時(shí)的難度及計(jì)算量,通過雙譜的對(duì)稱性,計(jì)算量則會(huì)進(jìn)一步減少。 根據(jù)選定的積分路徑,將積分雙譜分為以下4 種[7]。
1.2.1 徑向積分雙譜
如圖1 所示,在雙譜平面上,積分路徑取過原點(diǎn)的直線,得到的積分值即為徑向積分雙譜(radially integrated bispectrum, RIB)。
圖1 徑向積分雙譜Fig.1 Radially integrated bispectrum
RIB 具有一個(gè)明顯的缺陷:如果斜率的步長(zhǎng)選擇過大,那即使積分路徑上的所有ω2均取整,仍然會(huì)有點(diǎn)被遺漏,而這些點(diǎn)可能包含重要的信息。 如果斜率的步長(zhǎng)選擇過小,路徑上取整后的ω2必然還會(huì)有所重復(fù),從而造成一些雙譜值被重復(fù)積分,從而導(dǎo)致雙譜積分值特征不明顯。
1.2.2 軸向積分雙譜
如圖2 所示,在雙譜平面上,積分路徑取平行于ω1軸或平行于ω2軸的直線,得到的積分值即為軸向積分雙譜(axially integrated bispectrum,AIB)。AIB 使用了每一個(gè)雙譜值并確保雙譜值不會(huì)被重復(fù)使用,但其丟失了絕大部分相位信息。
圖2 軸向積分雙譜Fig.2 Axially integrated bispectrum
1.2.3 圓周積分雙譜
如圖3 所示,在雙譜平面上,積分路徑取以原點(diǎn)為圓心的圓周,得到的積分值即為圓周積分雙譜(circularly integrated bispectrum, CIB)。
圖3 圓周積分雙譜Fig.3 Circularly integrated bispectrum
同樣的,由于步長(zhǎng)的原因,CIB 在積分時(shí)必定會(huì)遺漏一部分雙譜值。
1.2.4 圍線積分雙譜
如圖4 所示,在雙譜平面上,積分路徑取以原點(diǎn)為中心的正方形,得到的積分值即為圍線積分雙譜(surrounding-line integral bispectrum, SIB)。SIB 同樣使用了每一個(gè)雙譜值并確保雙譜值不會(huì)被重復(fù)使用,同時(shí)保留了信號(hào)的尺度信息及部分相位信息[8-9]。
圖4 圍線積分雙譜Fig.4 Surrounding-line integral bispectrum
2.1.1 目標(biāo)亮點(diǎn)模型
如圖5(a)所示,以地面為目標(biāo)時(shí),高原上地面粗糙,檢測(cè)出的回波信號(hào)必然是多個(gè)目標(biāo)回波信號(hào)的疊加,因此可能導(dǎo)致無法識(shí)別出正確的炸高。 圖5(b)為實(shí)采數(shù)據(jù)頻譜,通過圖像可以明顯看出,頻率值具有多個(gè)峰值,說明檢測(cè)到了不同距離目標(biāo)的回波信號(hào)。
圖5 多目標(biāo)回波示意圖Fig.5 Schematic diagram of multi-target echo
大部分無線電引信以地面分布的對(duì)象作為目標(biāo),為了準(zhǔn)確地探測(cè)地面目標(biāo),必須建立一個(gè)有效的回波信號(hào)模型。
目標(biāo)回波是目標(biāo)對(duì)于信號(hào)源發(fā)射的信號(hào)做出的一種響應(yīng),目標(biāo)的固有散射特性和發(fā)射電磁波都會(huì)對(duì)回波信號(hào)造成影響[10]。 回波特性與入射的初始時(shí)間無關(guān),因此可以將目標(biāo)看做一個(gè)線性時(shí)不變系統(tǒng),回波就是對(duì)入射波的響應(yīng)。 可以將單個(gè)目標(biāo)的傳遞函數(shù)Hn(ω)定義為
Hn(ω)=An(ω)exp(iωτ+ i?) (4)
式中:An(ω)為幅度反射系數(shù),是關(guān)于頻率的函數(shù);τ為時(shí)延,是關(guān)于目標(biāo)反射距離的函數(shù),τ=2R/c,R為參考點(diǎn)到信號(hào)源的波程,c為光速;?為回波產(chǎn)生的相位差。
因此,函數(shù)可以視為一個(gè)多亮點(diǎn)疊加的模型,總的傳遞函數(shù)表示為
2.1.2 信號(hào)仿真
如圖6 所示,通過在Simulink 模型中添加2 個(gè)不同的距離,從而得到雙亮點(diǎn)接收信號(hào)的數(shù)字特征。 模型中默認(rèn)幅度反射系數(shù)為相同值。
圖6 雙亮點(diǎn)回波Simulink 模型Fig.6 Simulink model of double bright spot echo
如圖7 所示,能很明顯看出目標(biāo)亮點(diǎn)模型的仿真結(jié)果與圖5(b)中的實(shí)采信號(hào)包絡(luò)具有極高的相似度,一定程度上模擬出了地面目標(biāo)的信號(hào)特征,之后還可以通過添加地雜波進(jìn)行進(jìn)一步研究。
圖7 雙目標(biāo)亮點(diǎn)模型仿真結(jié)果Fig.7 Simulation results of double-target highlight model
如圖8 所示,將前端通過云臺(tái)綁定在標(biāo)志桿上,通過改變標(biāo)志桿的高度h及發(fā)射信號(hào)與地面的夾角,從而獲取不同的回波信號(hào)。 在本次實(shí)驗(yàn)中,高度h分別取2 m、3 m、4 m、5 m,發(fā)射信號(hào)與地面的夾角θ分別取20°、45°、65°、85°。
圖8 數(shù)據(jù)采集平臺(tái)Fig.8 Data collection platform
2.3.1 經(jīng)驗(yàn)?zāi)B(tài)分解
經(jīng)驗(yàn)?zāi)B(tài)分解(empirical mode decomposition,EMD)是Hilbert-Huang 變換(HHT)最重要的步驟,HHT 利用EMD,得到信號(hào)各階的本征模態(tài)函數(shù)(intrinsic mode function,IMF)。EMD 過程:首先,對(duì)原始信號(hào)進(jìn)行處理,求出各個(gè)局部極大值點(diǎn)。 然后,用插值法連接這些點(diǎn)組成上包絡(luò)線,同樣找到極小值點(diǎn)連接這些點(diǎn)組成下包絡(luò)線,用原始信號(hào)減去包絡(luò)線均值就得到了第一個(gè)residue,如果其不滿足IMF 定義(極值點(diǎn)和過零點(diǎn)數(shù)目至多相差一個(gè),且包絡(luò)線均值為0),就以該residue 作為原始信號(hào)繼續(xù)重復(fù)上面的操作,直到得到第一階IMF。 最后,用原始信號(hào)減去該IMF,得到新的信號(hào),作為初始信號(hào)再重復(fù)上面的篩選過程,直到所得到的剩余部分為單一信號(hào)或其值小于設(shè)定值時(shí),分解完成。
2.3.2 數(shù)據(jù)分解
將回波數(shù)據(jù)通過EMD 得到數(shù)據(jù)的各階IMF及殘余分量,對(duì)各階本征模態(tài)函數(shù)進(jìn)一步提取各種雙譜特征用于驗(yàn)證分類的準(zhǔn)確性。
圖9 和圖10 為初始信號(hào)及各階IMF 示意圖,其中采樣頻率為1 MHz,即采樣時(shí)間間隔為1 μs。 圖10 中縱坐標(biāo)無實(shí)際意義,僅為將各階本征模態(tài)函數(shù)直觀地表示在同一張圖中的相對(duì)參考位置。
圖9 初始信號(hào)Fig.9 Initial signal
圖10 各階IMF 示意圖Fig.10 Schematic diagram of IMF in each stage
為了更直觀地看出IMF 信號(hào)瞬時(shí)頻率的變化,每?jī)呻A提取一個(gè)IMF 信號(hào)的瞬時(shí)頻率。 從圖11可以明顯看出,各階IMF 分量的瞬時(shí)頻率隨著階數(shù)的增加而降低,因此提取不同距離回波信號(hào)的同階IMF 分量雙譜特征,再次進(jìn)行分類。 半數(shù)以上的當(dāng)多次分類結(jié)果為同一結(jié)果時(shí),結(jié)束對(duì)下一階IMF 分量進(jìn)行分類,否則持續(xù)此步驟,此時(shí)得到的信號(hào)為最終分類結(jié)果。 回波信號(hào)分類流程如圖12所示。
圖12 回波信號(hào)分類流程Fig.12 Flow chart of echo signal classification
2.4.1 灌木數(shù)據(jù)對(duì)比
在不同的高度上,回波信號(hào)具有不同的雙譜特征,如圖13 所示,分別為在θ=45°時(shí)距地面2,3,4,5 m 處的雙譜分布。
如圖13 所示,以f1、f2這2 個(gè)頻率作為坐標(biāo)軸來表示信號(hào)的偏斜程度在頻域上的分解,分別得到在θ=45°時(shí)距地面,2,3,4,5 m 處的雙譜分布。 可以明顯地看出,2,3,4,5 m 處的雙譜分布在幅值上具有明顯的區(qū)別,但對(duì)于工作平臺(tái)很小的引信來說,實(shí)時(shí)性是一個(gè)非常重要的指標(biāo),將二維數(shù)據(jù)轉(zhuǎn)換為一維數(shù)據(jù)能夠?qū)?shù)據(jù)量直接從x2(x為fft 次數(shù))降為x,因此,科學(xué)有效地提取出一維雙譜切片特征對(duì)于降低運(yùn)算速度具有極大的幫助。
圖13 θ =45°時(shí)不同高度雙譜分布Fig.13 Bispectrum distribution at different heights at θ =45°
2.4.2 四種積分雙譜特征
以θ=45°時(shí)4 個(gè)不同高度的4 種積分雙譜為例,得到的結(jié)果如圖14 所示。
圖14 θ =45°時(shí)不同積分雙譜分布Fig.14 Bispectrum distribution by different integral bispectrum at θ =45°
從圖14 中可以看出,將二維的雙譜特征按照4 種方式進(jìn)行積分并歸一化之后,仍然具有顯著的差異,保證了一維數(shù)據(jù)的分類處理仍具有很高的可行性。
目前,主要采用的分類方法有概率神經(jīng)網(wǎng)絡(luò)(PNN)[11]、支持向量機(jī)(SVM)[12-13]、卷積神經(jīng)網(wǎng)絡(luò)(CNN)[14-15],由于最鄰近節(jié)點(diǎn)算法(KNN)計(jì)算簡(jiǎn)便、計(jì)算速度快[16],符合引信實(shí)時(shí)性的要求,因此本文采用KNN 算法。
最鄰近節(jié)點(diǎn)算法也稱為k最近鄰算法,該方法簡(jiǎn)單且直觀,如圖15 所示。 每一個(gè)樣本x通過計(jì)算到其他樣本的距離,用與其最鄰近的k個(gè)樣本來進(jìn)行分類,將樣本x分為k個(gè)樣本中數(shù)量最多的類別。 由于同樣情況下回波信號(hào)的樣本集重疊很多,用在本文非常合適。 該方法無需訓(xùn)練,適合應(yīng)用于引信上。
圖15 最鄰近節(jié)點(diǎn)算法Fig.15 KNN algorithm
將4 種積分雙譜得到的4 個(gè)高度、4 個(gè)角度的數(shù)據(jù)各取30 組,利用KNN 對(duì)其進(jìn)行分類,得到結(jié)果如圖16 ~圖19 所示,黑色方框中的數(shù)據(jù)分別表示在對(duì)應(yīng)數(shù)據(jù)下的分類正確的數(shù)量及所占分類數(shù)據(jù)總量的比例。
圖16 RIB 分類結(jié)果Fig.16 RIB classification results
圖17 AIB 分類結(jié)果Fig.17 AIB classification results
圖18 CIB 分類結(jié)果Fig.18 CIB classification results
3.2.1 RIB 分類結(jié)果分析
對(duì)20°、45°、65°、85°這4 個(gè)角度120 組分類數(shù)據(jù)RIB 積分?jǐn)?shù)據(jù)進(jìn)行分類,最終可以得到4 次分類的準(zhǔn)確率分別為61.7%、61.7%、69.2%、66.7%,利用RIB 對(duì)480 組數(shù)據(jù)分類的總準(zhǔn)確率為64.8%。
3.2.2 AIB 分類結(jié)果分析
對(duì)20°、45°、65°、85°這4 個(gè)角度120 組分類數(shù)據(jù)AIB 積分?jǐn)?shù)據(jù)進(jìn)行分類,最終可以得到4 次分類的準(zhǔn)確率分別為89.2%、91.7%、80.0%、82.5%,利用AIB 對(duì)480 組數(shù)據(jù)分類的總準(zhǔn)確率為85.8%。
3.2.3 CIB 分類結(jié)果分析
對(duì)20°、45°、65°、85°這4 個(gè)角度120 組分類數(shù)據(jù)CIB 積分?jǐn)?shù)據(jù)進(jìn)行分類,最終可以得到4 次分類的準(zhǔn)確率分別為84.2%、80.8%、80.0%、80.0%,利用CIB 對(duì)480 組數(shù)據(jù)分類的總準(zhǔn)確率為81.3%。
3.2.4 SIB 分類結(jié)果分析
對(duì)20°、45°、65°、85°這4 個(gè)角度120 組分類數(shù)據(jù)SIB 積分?jǐn)?shù)據(jù)進(jìn)行分類,最終可以得到4 次分類的準(zhǔn)確率分別為81. 7%、64. 2%、60. 0%、53.3%,利用SIB 對(duì)480 組數(shù)據(jù)分類的總準(zhǔn)確率為64.8%。
顯然,利用KNN 對(duì)AIB 和CIB 得到的一維數(shù)據(jù)分類結(jié)果更為準(zhǔn)確,準(zhǔn)確率可以達(dá)到80% 以上,而對(duì)RIB 和SIB 得到的一維數(shù)據(jù)分類結(jié)果并不理想,分類的準(zhǔn)確率只有65%左右。
通過表1、表2 可以看出,利用經(jīng)EMD 對(duì)不同頻率的數(shù)據(jù)進(jìn)行特征提取后,分類結(jié)果有了明顯的改進(jìn),成功率最多可提升約8.3%,分類總成功率提升3.7%。
表1 分類結(jié)果匯總Table 1 Summary of classification results
表2 優(yōu)化雙譜分類結(jié)果Table 2 Optimization of bispectral classification results
本文主要圍繞雙譜分析在信號(hào)處理中的應(yīng)用,對(duì)實(shí)測(cè)高原地區(qū)灌木的回波進(jìn)行了分析,通過4 種不同的積分雙譜方式對(duì)雙譜特征進(jìn)行了降維,進(jìn)而對(duì)一維特征進(jìn)行了分類,分類結(jié)果表明:
1) 利用KNN 算法,得到了對(duì)2 m、3 m、4 m、5 m這4 個(gè)高度的分類成功率,最高可以達(dá)到74.1%。
2) 雖然在理論上SIB 算法利用了所有的點(diǎn)數(shù)和相位信息,但在實(shí)際數(shù)據(jù)分類中效果并不好。而AIB 與CIB 兩種算法準(zhǔn)確率較高,分別可以達(dá)到85.4% ~85.8%和81.5% ~81.3%,其中對(duì)于45°的數(shù)據(jù),AIB 分類效果更為卓越,可以達(dá)到91.7%。
3) 利用EMD 對(duì)數(shù)據(jù)進(jìn)行分解后,優(yōu)化分類結(jié)果有了明顯的提高,其中AIB 成功率最高可達(dá)92.5%。
4) 說明太赫茲波段的引信完全可以通過AIB 和CIB 兩種算法實(shí)現(xiàn)快速準(zhǔn)確的分類,利用HHT 優(yōu)化的雙譜數(shù)據(jù)對(duì)于提升分類識(shí)別結(jié)果的準(zhǔn)確度具有一定的作用。