基于改進(jìn)希爾伯特-黃變換的進(jìn)動目標(biāo)周期估計

金家偉，阮懷林，孫 兵

（國防科技大學(xué)電子對抗學(xué)院，合肥 230037）

0 引言

物體或物體結(jié)構(gòu)的微動可能會對返回的信號產(chǎn)生頻率調(diào)制，從而產(chǎn)生多普勒頻移的邊帶，稱為微多普勒效應(yīng)。彈頭目標(biāo)為了保持飛行的穩(wěn)定性，不僅會繞自身對稱軸作自旋運(yùn)動，還會受到?jīng)_擊力矩的作用，力矩消失后對稱軸將在平衡位置作圓錐運(yùn)動，即進(jìn)動。進(jìn)動是目標(biāo)保持穩(wěn)定的一種特殊微動，而誘餌一般不具備姿態(tài)控制裝置，因此，進(jìn)動是中段彈頭特有的微動方式。通過提取具有進(jìn)動錐體目標(biāo)的幾何特征，可以區(qū)分它是彈頭還是誘餌。

至今為止，由于時頻分析可將一維信號映射到二維時頻平面，從而展現(xiàn)信號的時頻信息，因此時頻分析技術(shù)被廣泛應(yīng)用于分析時變的微多普勒頻率特征，因?yàn)闀r頻分析可以通過將一維信號映射到二維時頻平面來同時產(chǎn)生信號的時頻信息。文獻(xiàn)［6］結(jié)合Gabor 時頻分布和變分模態(tài)分解，來估計目標(biāo)的自旋頻率和錐旋頻率；文獻(xiàn)［7］利用短時分?jǐn)?shù)階傅立葉變換分離肢體和軀干的微多普勒信號；文獻(xiàn)［8］采用基于加權(quán)迭代自適應(yīng)的時頻分析方法結(jié)合逆Radon 變換，分離重構(gòu)不同散射點(diǎn)的微多普勒分量。然而，在實(shí)際應(yīng)用中，大多數(shù)現(xiàn)有的時頻方法在提取微多普勒信息方面仍然存在一些不可避免的不足，包括邊界失真和干擾項(xiàng)。另外，當(dāng)目標(biāo)含有多個強(qiáng)散射中心時，微動目標(biāo)返回的雷達(dá)信號是多分量的。因此，所需的時頻分析方法應(yīng)該同時抑制交叉項(xiàng)和自動項(xiàng)。

近年來，一種新的分析非線性、非平穩(wěn)數(shù)據(jù)的方法——希爾伯特- 黃變換（hilbert-huang transform，HHT）越來越受到人們的重視。該方法的關(guān)鍵部分是經(jīng)驗(yàn)?zāi)B(tài)分解（empirical mode decomposition，EMD），其主要思想是將需要進(jìn)行處理的非線性非平穩(wěn)信號分解，得到一組頻率分量各有不同的本征模態(tài)函數(shù)（intrinsic mode functions，IMF）。通過希爾伯特變換（hilbert transform，HT），IMF 能夠產(chǎn)生隨時間變化的瞬時頻率，其最終表示形式是能量-頻率- 時間分布，稱為希爾伯特- 黃譜（hilberthuang spectrum，HHS）。希爾伯特-黃變換的優(yōu)點(diǎn)之一是希爾伯特- 黃譜不涉及頻率分辨率和時間分辨率的概念，而是瞬時頻率。因其良好的適應(yīng)性，已被廣泛應(yīng)用于微多普勒效應(yīng)的研究，文獻(xiàn)［12］使用了希爾伯特-黃變換對振動目標(biāo)進(jìn)行雷達(dá)微多普勒特征分析；文獻(xiàn)［13］提出了一種基于經(jīng)驗(yàn)?zāi)B(tài)分解的圓錐目標(biāo)雷達(dá)微動特征提取方法。然而，在實(shí)際應(yīng)用中，EMD 算法并不完美，同樣存在一些問題，如模態(tài)混疊、端點(diǎn)效應(yīng)等。

針對基于希爾伯特- 黃變換的進(jìn)動目標(biāo)微多普勒特征提取問題，由于雷達(dá)信號不可避免地含有較高水平的隨機(jī)噪聲，其瞬時頻率和瞬時幅度非常分散，常規(guī)希爾伯特-黃變換難以將目標(biāo)信號與噪聲信號分離。本文采用EMD 的衍生算法CEEMDAN先對雷達(dá)回波進(jìn)行分解，然后再進(jìn)行希爾伯特譜分析，以此提取進(jìn)動目標(biāo)的微多普勒周期。

1 目標(biāo)模型

空間進(jìn)動錐體目標(biāo)模型如圖1 所示，坐標(biāo)系為雷達(dá)坐標(biāo)系（U，V，W），雷達(dá)靜止于原點(diǎn)Q。O 為目標(biāo)質(zhì)心，以O(shè) 為原點(diǎn)、目標(biāo)對稱軸Oz 為z 軸建立目標(biāo)本體坐標(biāo)系O-xyz。以O(shè) 為原點(diǎn)建立參考坐標(biāo)系O-XYZ，以初始時刻與目標(biāo)對稱軸Oz、進(jìn)動軸OZ共面且垂直于OZ 的方向?yàn)閅 軸，X 軸根據(jù)右手準(zhǔn)則確定。目標(biāo)在平動的同時，以角速度ω繞對稱軸z軸作自旋運(yùn)動，同時以角速度ω繞OZ 軸錐旋運(yùn)動（ω和ω均采用參考坐標(biāo)系中的表達(dá)式），自旋軸和錐旋軸之間的夾角為進(jìn)動角θ。

圖1 目標(biāo)模型

在光學(xué)區(qū)，雷達(dá)目標(biāo)的整體散射特性通常可以等效為若干個散射中心的疊加。不失一般性，假設(shè)目標(biāo)的微動是周期性進(jìn)動，目標(biāo)等效為K 個散射中心，各散射中心各向同性，雷達(dá)發(fā)射的電磁波為連續(xù)單頻波，載頻為f。雷達(dá)向由i 個散射點(diǎn)組成的目標(biāo)發(fā)射電磁波，不考慮目標(biāo)平動帶來的影響，返回的基帶信號可表示為

2 改進(jìn)的希爾伯特-黃變換

2.1 希爾伯特-黃變換

2.2 完備總體經(jīng)驗(yàn)?zāi)B(tài)分解

Torres 等提出了一種對噪聲有自適應(yīng)分解特性的完備總體經(jīng)驗(yàn)?zāi)B(tài)分解（CEEMDAN），在每一階段的分解中都加入一個特定的白噪聲，再通過計算殘差得到該階段的模態(tài)。CEEMDAN 是一種噪聲輔助數(shù)據(jù)分析方法，所加入的白噪聲在整個時頻空間分布相對均勻，以提供時頻空間中濾波的參考尺度，而不同尺度的信號自動投影到由白噪聲建立的合適的參考尺度上。雖然每次試驗(yàn)產(chǎn)生的結(jié)果都會包含噪聲，但是當(dāng)試驗(yàn)的數(shù)量足夠時，其結(jié)果的總體平均當(dāng)中的噪聲將消失，即總體平均當(dāng)中所添加的白噪聲往往會相互抵消，從而只包含所分解得到的信號。

3 仿真與分析

為驗(yàn)證希爾伯特- 黃變換在進(jìn)動目標(biāo)微多普勒特征提取方面的效果，本文分別以EMD 和CEEMDAN 來進(jìn)行信號分解，然后進(jìn)行希爾伯特譜分析，并對比兩種方法的提取效果。

仿真條件：假設(shè)雷達(dá)工作在10 GHz，且雷達(dá)坐標(biāo)系中本體坐標(biāo)系原點(diǎn)O 的坐標(biāo)為（400，500，100）km，本地坐標(biāo)系和參考坐標(biāo)系之間的初始?xì)W拉角（x-y-z 序列）為（30°，60°，45°）。假設(shè)目標(biāo)繞z 軸旋轉(zhuǎn)，目標(biāo)上有兩個散射點(diǎn)：第1 散射點(diǎn)A 位于錐頂，在本體坐標(biāo)系中的坐標(biāo)為（0，0，1）m；第2 個散射點(diǎn)B 位于錐底的尾翼，在本體坐標(biāo)系中的坐標(biāo)為（0.5，0，-0.5）m。雷達(dá)照射時間為8 s。旋轉(zhuǎn)頻率為f=1 Hz，圓錐運(yùn)動頻率為f=0.5 Hz。假設(shè)環(huán)境噪聲為高斯白噪聲。

當(dāng)不存在噪聲時，基于EMD 和CEEMDAN 的進(jìn)動周期提取結(jié)果分別如圖2 和圖3 所示，圖2（a）和圖3（a）分別是兩種算法獲取的回波信號希爾伯特譜，圖2（b）和圖3（b）分別是兩種算法獲取的IMF前3 個分量的希爾伯特譜。

圖2 基于EMD 的估計結(jié)果（無噪聲）

希爾伯特譜是時間、頻率和能量的三維譜，當(dāng)沒有噪聲時，圖2（a）和圖3（a）都展現(xiàn)了回波信號時頻分布的細(xì)節(jié)，其中的頻率都是取模后的頻率，也就是正頻率。微多普勒信號經(jīng)過EMD 后都被分解在第1 個IMF 分量中，如圖2（b）所示；而微多普勒信號經(jīng)過CEEMDAN 后被分解在各個IMF 分量中，如圖3（b）所示。由圖2 和圖3 可看出，兩種算法在不存在噪聲時，都能夠準(zhǔn)確估計出進(jìn)動目標(biāo)的微多普勒周期為2 s。

圖3 基于CEEMDAN 的估計結(jié)果（無噪聲）

當(dāng)雷達(dá)回波信號中混雜著噪聲時，同樣的仿真條件下，假設(shè)信噪比SNR=10 dB，圖4（a）為基于EMD 的回波信號希爾伯特譜，圖4（b）是基于EMD得到的IMF 前3 個分量的希爾伯特譜。

圖4 基于EMD 的估計結(jié)果（SNR=10 dB）

從圖4 可以看出，隨著信噪比的降低，信號的希爾伯特譜中的微多普勒頻率逐漸淹沒在噪聲中，無法用以提取進(jìn)動目標(biāo)微多普勒周期，如圖4（a）所示。由于EMD 未能將微多普勒信號和噪聲分離開，而是都混在第1 個IMF 分量中，如圖4（b）所示，隨著信噪比SNR 的降低，將同樣無法通過第1 個IMF分量的希爾伯特譜來提取目標(biāo)信息。但在其他條件不變的情況下，改用CEEMDAN 對雷達(dá)回波進(jìn)行處理，然后進(jìn)行希爾伯特譜分析。圖5（a）為回波信號希爾伯特譜，圖5（b）為采用CEEMDAN 所獲取的IMF 前3 個分量的希爾伯特譜。

圖5 基于CEEMDAN 的估計結(jié)果（SNR=10 dB）

從圖5（a）中可以看出，當(dāng)信噪比持續(xù)降低時，與圖4（a）相似，其中的微多普勒信號伴隨著大量的噪聲，難以直接通過回波信號地希爾伯特譜來提取微多普勒信息。但通過CEEMAN 分解后所得到的一組IMF，噪聲頻率主要集中于第1 個IMF 分量中，雖然第1 個IMF 分量中也包含了一部分的微多普勒頻率信號，但更多的微多普勒信號被分解到了第2 個和第3 個分量中，從而達(dá)到將噪聲和目標(biāo)信號分離的目的，如圖5（b）所示，則從該圖中的第2 個和第3 個IMF 分量的希爾伯特譜中，可以準(zhǔn)確提取出目標(biāo)的微多普勒周期。

進(jìn)一步降低信噪比進(jìn)行對比分析，圖6（a）和圖6（b）分別給出了信噪比分別為5 dB 和1 dB 時，基于EMD 得到的IMF 前3 個分量的希爾伯特譜。

圖6 EMD 獲取的IMF 前3 個分量的希爾伯特譜

由于EMD 算法難以從噪聲分離出微多普勒信號，在其IMF 分量的希爾伯特譜中，即使信噪比為5 dB，仍然無法估計出進(jìn)動周期，如圖6（a）所示。但采用CEEMAN 時，圖7（a）和圖7（b）分別給出了信噪比分別為5 dB 和1 dB 時，基于CEEMDAN 得到的IMF 前3 個分量的希爾伯特譜。

圖7 CEEMDAN 獲取的IMF 前3 個分量的希爾伯特譜

將圖7 與圖6 對比可知，隨著信噪比的進(jìn)一步降低，即使信噪比低至1 dB，仍能通過其IMF 分量的希爾伯特譜提取出進(jìn)動周期，如圖7（b）所示。這也進(jìn)一步驗(yàn)證了在低信噪比環(huán)境中，以CEEMDAN改進(jìn)的希爾伯特-黃變換具備更加良好的性能。仿真中假設(shè)環(huán)境噪聲為高斯白噪聲，旨在驗(yàn)證CEEMDAN 具備比EMD 更好的抗噪性能，但在實(shí)際工程中的噪聲不一定是高斯白噪聲，算法具體的抗噪性能需要以實(shí)際雷達(dá)回波數(shù)據(jù)進(jìn)行驗(yàn)證，這是本文算法存在的局限性。

4 結(jié)論

本文將CEEMDAN 算法替換EMD 算法，以改進(jìn)希爾伯特-黃變換，用以提取進(jìn)動目標(biāo)的微多普勒周期。仿真實(shí)驗(yàn)表明，改進(jìn)的希爾伯特-黃變換能夠獲取進(jìn)動目標(biāo)微多普勒信號的希爾伯特譜，從而提取出目標(biāo)的進(jìn)動周期；相對于EMD 算法，CEEMDAN 算法能夠?qū)⑽⒍嗥绽招盘柡驮肼曅盘栠M(jìn)行分離，從而有效克服了常規(guī)希爾伯特-黃變換在噪聲環(huán)境中性能較差的不足。