基于動態(tài)模態(tài)分解的彈道目標平動補償與微動特征提取方法

摘 要：針對彈道目標平動導致微動特征難以準確提取的問題，提出一種基于動態(tài)模態(tài)分解（dynamic mode decomposition，DMD）的彈道目標平動補償與微動特征提取方法。首先，在彈道目標微動回波建模的基礎上，對目標的慢時間-距離像序列進行微多普勒（micro Doppler，m D）特征曲線分離;其次，將分離后的數(shù)據(jù)向量移位堆疊構建為增廣數(shù)據(jù)矩陣，并對其進行DMD;然后，利用分解后的模態(tài)幅值對各模態(tài)進行排序，結合損失函數(shù)等信息選取主要模態(tài);同時，利用主要模態(tài)中的零頻率模態(tài)完成彈道目標的平動補償，從其他主要模態(tài)中提取出自旋頻率和錐旋頻率等微動特征信息;最后，對基于DMD的彈道目標平動補償與微動特征提取方法進行性能分析與對比實驗，驗證了所提方法的可行性和穩(wěn)健性。

關鍵詞： 動態(tài)模態(tài)分解; 彈道目標; 微多普勒; 平動補償; 特征提取

中圖分類號： TN 95

文獻標志碼： ADOI：10.12305/j.issn.1001 506X.2025.02.12

Translational compensation and micro motion feature extraction method of

ballistic targets based on dynamic mode decomposition

LI Kaiming^1，2，*， DAI Xiaonan^1，3， ZHANG Yuanpeng4， YAO Jiawen5， LUO Ying^1，2

（1. College of Information and Navigation， Air Force Engineering University， Xi’an 710077， China;

2. Collaborative Innovation Center of Information Sensing and Understanding，

Air Force Engineering University， Xi’an 710077， China; 3. Unit 95806 of the PLA， Beijing 100076， China;

4. The First Department， Air Force Early Warning Academy， Wuhan 430019， China;

5. Beijing Research Institute of Telemetry， Beijing 100094， China）

Abstract：To address the problem of difficulty in accurately extracting micro motion features caused by ballistic target translational motion， a ballistic target translational motion compensation and micro motion feature extraction method based on dynamic mode decomposition （DMD） is proposed. Firstly， based on the modeling of micro motion echoes of ballistic targets， the slow time range profile sequence of the target is subjected to micro Doppler （m D） feature curve separation. Secondly， the separated data vectors are shifted and stacked to construct an augmented data matrix， and DMD is performed. Then， the decomposed mode amplitudes are used to rank each mode， and the main mode is selected combined with information such as the loss function. The zero frequency mode in the main modes is utilized to achieve translational compensation of ballistic targets， and micro motion feature information such as spin frequency and cone rotation frequency from other main modes are extracted. Finally， performance analysis and comparative experiments are conducted on the DMD based ballistic target translational compensation and micro motion feature extraction methods， verifying the feasibility and robustness of the proposed methods.

Keywords：dynamic mode decomposition （DMD）; ballistic target; micro Doppler （m D）; translational compensation; feature extraction

0 引 言

彈道導彈目標飛行速度快且伴有自旋、進動等微動^［1-3］形式，彈頭與誘餌、燃料箱等的微動特征存在明顯差異^［4-7］，這種獨特的微動特征差異為彈道目標識別提供了有效途徑^［8-10］。但彈道導彈平動會導致目標微多普勒（micro Doppler， m D）曲線產(chǎn)生變化，給后續(xù)的特征提取帶來困難^［11-12］。文獻［12］利用經(jīng)驗模態(tài)分解（empirical mode decomposition， EMD）得到的趨勢項完成彈道目標的平動補償，并提取出自旋頻率、錐旋頻率等微動特征。文獻［13］采用延遲共軛乘法估計目標的平動加速度，通過搜索補償信號頻譜峰值估計出自旋頻率等微動參數(shù)。文獻［14］基于m D對稱抵消效應，對旋轉對稱目標進行平動補償，但對非對稱目標的平動補償不再適用。文獻［15］基于高階模糊函數(shù)實現(xiàn)進動目標的平動補償，但求解雷達回波的高階矩需要較長觀測時間。文獻［16］基于高階模糊函數(shù)、延遲共軛相乘及時頻分析相結合的方法完成多目標的平動參數(shù)估計，但該方法受噪聲影響較大。

近年來，Schmid^［17］提出一種利用動態(tài)系統(tǒng)的特征值對流場進行預測和穩(wěn)定性分析的方法，即動態(tài)模態(tài)分解（dynamic mode decomposition， DMD）算法。該算法的本質是尋求流場動力學的降階模型，通過對流場快照進行特征分解，提取出能夠表征流場結構的低階模態(tài)及其對應的特征值頻率，且頻率具有單一的特點，并能夠基于頻率對流場進行排序，從而觀察不同頻率的流動結構對流場的影響程度和貢獻度^［18］。目前，DMD算法被廣泛應用于流體力學^［19-20］、機器人^［21］、電磁信號分析^［22］、疾病控制^［23-24］等多個領域。Rowley等^［25］討論DMD算法與Koopman算子之間的聯(lián)系，指出DMD算法可以近似表示Koopman模態(tài)和特征值，具有描述非線性流動中的觀測量（如壓強、速度等）隨時間演化的能力，即DMD的過程就是一個線性估計過程，通過線性假設實現(xiàn)非線性估計^{［18，26］}。

考慮到彈道目標的微動回波通常包含自旋、錐旋等不同的微動頻率分量，而m D信號在時頻域或慢時間-距離域具有非平穩(wěn)、非線性的特點，并且DMD算法可以得到的模態(tài)頻率是唯一的，因此非常適合被應用于DMD算法、對彈道目標的m D信號進行特征分解，提取出其中的微動頻率分量。綜合以上分析，本文將DMD算法應用于彈道目標的平動補償與微動特征提取。首先，基于寬帶雷達進行彈道目標微動回波建模，分離出m D曲線，對分離的m D曲線組成的數(shù)據(jù)向量移位堆疊構建增廣數(shù)據(jù)矩陣，并進行DMD。然后，利用分解模態(tài)的幅值對各個模態(tài)進行排序，結合損失函數(shù)、模態(tài)能量等信息選取具有較大貢獻度的主要模態(tài)，重構出m D曲線，利用主要模態(tài)中的零頻率模態(tài)完成彈道目標的平動補償，從其他主要模態(tài)中提取出自旋頻率、錐旋頻率等微動特征。最后，仿真分析DMD算法在被應用于彈道目標微動特征提取時的性能，并與集合EMD（ensemble EMD， EEMD）^［27］、變分模態(tài)分解（variational mode decomposition， VMD）^［28］等微動特征提取方法進行對比分析，驗證了所提方法在彈道目標微動特征提取方面的性能優(yōu)勢。

1 彈道目標微動回波建模

通常，彈道導彈目標自旋以保持姿態(tài)穩(wěn)定，導彈頭體分離后受到擾動而發(fā)生進動或章動。圖1為有翼錐體彈頭與雷達位置關系的示意圖，雷達坐標系（或全局坐標系）為（U，V，W），雷達靜止于坐標原點O;目標等效散射中心主要由錐頂散射點A以及尾翼散射點B和C構成。（X，Y，Z）為參考坐標系，與雷達坐標系平行，坐標原點為目標質心O′。設初始時刻雷達到目標質心O′的徑向距離為R0，在雷達坐標系中的方位角和俯仰角分別為α和β;坐標系（x，y，z）為目標本地坐標系，坐標原點也為O′。假設彈道目標在參考坐標系中以角速度ωs繞目標對稱軸O′z做自旋運動，同時還以角速度ωc繞空間某定向軸O′N做錐旋運動?？紤]到彈道目標在中段飛行時速度相對比較平穩(wěn)，因此采用勻加速運動模型來近似目標中段的飛行過程，設目標的徑向速度和徑向加速度分別為v和a。

假設雷達發(fā)射寬帶線性調頻信號，則目標上任一散射點的回波可表示為

sr（tk，tm）=σ·recttk－τ（tm）TP·

expj2πfc（tk－τ（tm））+12μ（tk－τ（tm））2（1）

式中：σ為散射點的散射系數(shù);tk為快時間;tm為慢時間;TP為脈沖寬度;fc為載頻;μ為調頻斜率;τ（tm）=2R（tm）/c為tm時刻散射點到雷達的距離時延，c為光速。

考慮到彈道目標的運動過程可以分解為目標質心的平動和散射點繞目標質心的微動兩個部分，因此散射點到雷達的徑向距離可表示為

R（tm）=RTM（tm）+RM（tm）（2）

式中：RTM（tm）=R0+vtm+at2m/2為平動項;RM（tm）為微動項。

假設精確的平動補償已完成，則以目標質心O′點的回波信號作為參考信號，其表達式為

sref（tk，tm）=recttk－τref（tm）Tref·

expj2πfc（tk－τref（tm））+12μ（tk－τref（tm））2（3）

式中：Tref為參考信號脈寬;τref（tm）=2Rref（tm）/c為O′點的回波時延，Rref（tm）=R0+vtm+at2m/2為目標質心到雷達的徑向距離。

進一步，通過將該參考信號與目標回波信號進行Dechirp處理，并做關于快時間tk的傅里葉變換，然后去除剩余視頻相位（residual video phase， RVP）項和包絡“斜置”項，可得

Sd（f，tm）=σ·TP sincTPf+2μcΔR（tm）·

exp－j4πfccΔR（tm）（4）

式中：ΔR（tm）為散射點到參考點O′的徑向距離，ΔR（tm）=R（tm）－Rref（tm）=RM（tm）。而sinc函數(shù)的峰值位于f=－2μΔR（tm）/c處，則f通過換算可以轉化為散射點到O′的徑向距離ΔR（tm）。考慮到彈道目標在空間中的整體微動，ΔR（tm）將隨慢時間周期性地變化，可反映出散射點的微動特性。

在實際中，雖然可以通過雷達測距和測速分別得到R0和v的測量值R′0和v′，但是測量值和真值總會存在一定的誤差，此時O′點到雷達的距離可表示為R′ref（tm）=R′0+v′tm，經(jīng)過脈壓后的距離像即為

S′d（f，tm）=σ·TPsincTPf+2μcΔR′（tm）·

exp－j4πfccΔR′（tm）（5）

式中：ΔR′（tm）為散射點到O′點的徑向距離變化，ΔR′（tm）=R（tm）－R′ref（tm）=RM（tm）+R′TM（tm），且R′TM（tm）=（R0－R′0）+（v－v′）tm+at2m/2。

設參考坐標系中目標自旋角速度矢量為ωs=［ωsx，ωsy，ωsz］T，此自旋運動對應的旋轉矩陣^［4］為

Rspinning=I+ω^′ssin（Ωstm）+ω^′2s（1-cos（Ωstm））（6）

式中：I為單位矩陣;Ωs=ωs;ω^′s為自旋斜對稱矩陣^［4］。

同理，若目標錐旋角速度矢量為ωc=［ωcx，ωcy，ωcz］T時，錐旋運動對應的旋轉矩陣即為

Rconing=I+ω^′csin（Ωctm）+ω^′2c（1-cos（Ωctm））（7）

式中：Ωc=ωc;ω^′c為錐旋斜對稱矩陣。

對于錐頂散射點A，由于其位于自旋軸上，其微動形式僅表現(xiàn)為錐旋運動，則有

RM（tm）=［RconingRinitrA］Tn（8）

式中：Rinit為歐拉旋轉矩陣;rA為錐頂散射點A在目標本地坐標系中對應的矢量;n為雷達視線方向（line of sight，LOS）上的單位矢量。

結合式（5）和式（8）可知，ΔR′（tm）=－（1/2）×fc/μ，且由于目標微動的影響，距離像峰值位置將隨ΔR′（tm）而周期性變化。此時，對于進動目標，錐頂散射點A只發(fā)生錐旋。分析式（8）可知，ΔR′（tm）右邊中的微動項RM（tm）隨慢時間做正弦或余弦運動，運動周期即為錐旋運動的周期，R′TM（tm）則是目標散射點的平動項，需要進行平動補償。

對于尾翼散射點B或C，以B點為例，其微動形式表現(xiàn)為進動，即自旋與錐旋的合成運動，因此有

RM（tm）=［RconingRspinningRinitrB］Tn（9）

式中：rB為尾翼散射點B在目標本地坐標系中對應的矢量。

同理，由式（5）和式（9）可見，尾翼散射點B的距離像峰值位置同樣受到RM（tm）的調制，且隨慢時間tm周期性變化。將式（6）和式（7）代入式（9）并化簡可知，與散射點A的m D效應不同，尾翼散射點B由進動引起的m D效應在慢時間-距離平面上的曲線已經(jīng)不再是簡單的正弦形式，而是多個正弦（或余弦）頻率的復合形式，主要包括Ωs+Ωc、Ωs、Ωc和|Ωs-Ωc|這4種頻率分量^［4］，但該曲線的變化仍然具有周期性，周期為自旋周期和錐旋周期的最小公倍數(shù)^［12］。

由以上分析可知，彈道目標散射點的回波在圖像域（慢時間和距離所構成的二維平面）均表現(xiàn)為周期性，即目標散射點的慢時間-距離像序列，無論是錐頂散射點A或是尾翼散射點B或C，散射點在各時刻的距離像峰值位置都隨慢時間呈現(xiàn)周期變化。進一步，采用文獻［4］所提方法對慢時間-距離像序列進行預處理，首先對|S′d（f，tm）|進行高斯平滑濾波和二值化，然后進行骨架提取，可搜索分離出每個散射點對應的m D曲線，此時每條m D曲線的數(shù)據(jù)所組成的向量仍然包含各散射點的微動特征，即m D曲線隨慢時間表現(xiàn)出非平穩(wěn)、非線性的周期性變化，而DMD算法能夠描述這種非線性變化隨時間的演化歷程。基于此，可采用DMD方法進行平動補償和微動特征的提取。

2 平動補償與微動特征提取

2.1 利用增廣數(shù)據(jù)矩陣實現(xiàn)DMD

通過在圖像域對目標的慢時間-距離像序列進行曲線分離，可得到目標散射點在N個慢時間的數(shù)據(jù)，進一步將N個慢時間的數(shù)據(jù)按照回波時間先后順序排列后形成一個數(shù)據(jù)向量［x1，x2，…，xN］，記第i個慢時間的數(shù)據(jù)為xi（i=1，2，…，N），任意兩個數(shù)據(jù)之間的時間間隔為Δt，即脈沖重復間隔（pulse repetition interval， PRI）。則定義兩個快照序列為

X^N－11=［x1，x2，…，xN－1］（10）

XN2=［x2，x3，…，xN］（11）

將DMD算法應用于雷達目標m D信號處理時，每條m D曲線的數(shù)據(jù)所組成的向量是一個行向量。如果對數(shù)據(jù)向量直接進行DMD，只能得到單個實特征值，而不是共軛成對的復特征值，因此無法完全捕獲到周期振蕩的m D曲線^{［23，29］}。為此，本文將數(shù)據(jù)向量移位堆疊為增廣數(shù)據(jù)矩陣，通過增加矩陣行數(shù)，周期振蕩的m D曲線中的頻率分量就能被完全捕獲^［29］。具體而言，可對式（10）構建如下增廣數(shù)據(jù)矩陣：

Xaug=

x1x2…xN－h(huán)

x2x3…xN－h(huán)+1

xhxh+1…xN－1

（12）

式中：h為移位堆疊的次數(shù)。

同理，可對式（11）移位堆疊形成增廣數(shù)據(jù)矩陣：

X′aug=

x2x3…xN－h(huán)+1

x3x4…xN－h(huán)+2

xh+1xh+2…xN

（13）

假設在較短時間內，目標散射點的微動回波數(shù)據(jù)滿足線性相關，即存在線性算子A，則有

X′aug=AXaug（14）

可以看出，這個過程就是一個線性估計過程，即通過線性假設實現(xiàn)非線性估計。進一步，對增廣數(shù)據(jù)矩陣進行DMD，由于DMD算法得到的模態(tài)、頻率、增長率/衰減率等信息與矩陣A的特征分解密切相關，但是矩陣A通常是一個h×h的高維矩陣，難以進行求解。因此，在實際處理中需要對其進行降維處理，尋求一個秩為r（rh）的低維矩陣A～∈C^r×r來近似表示高維矩陣A，則DMD算法的具體步驟如下。

步驟 1 對增廣數(shù)據(jù)矩陣Xaug進行奇異值分解：

Xaug=UΣV（15）

式中：U和V為酉矩陣，U∈C^h×r，V∈C^{（N－h(huán)）×r}，且滿足UU=I，VV=I;*表示復共軛轉置;Σ為對角陣，Σ∈C^r×r，對角線上有r個奇異值（σ1，σ2，…，σr）。

步驟 2 求解相似矩陣A～：

將式（15）代入式（14）可得

A=X′augVΣ^－1U（16）

則相似矩陣A～可表示為

A～=UAU=UX′augVΣ^－1（17）

步驟 3 對A～進行特征分解：

A～W=WΛ（18）

式中：W的列是特征向量;Λ是包含對應特征值λi的對角陣。如果特征值落在單位圓內，表示該模態(tài)穩(wěn)定，反之則不穩(wěn)定，可以計算該特征值的對數(shù)形式：

ωi=ln λiΔt（19）

ωi的實部代表相應DMD模態(tài)的增長率/衰減率，虛部決定了模態(tài)的頻率：

fi=im（ωi）2π（20）

式中：im（·）表示取ωi的虛部。

步驟 4 通過W和對角陣Λ重構矩陣A的特征分解，Λ中的λi為A的特征值，Φ的列Φi為A的特征向量（即DMD模態(tài)）^［30］：

Φ=X′augVΣ^－1W（21）

可通過前k階模態(tài)進行重構，則重構的數(shù)據(jù)矩陣可表示為

XDMD（t）≈∑ki=1Φiexp（ωit）bi=Φexp（Ωt）b（22）

式中：t表示時間;Ω=diag（ω）為一個對角陣，對角陣中的元素為ωi;bi是每個模態(tài)的初始幅值;b為由bi組成的向量。

如果將初始時刻的數(shù)據(jù)，即Xaug的第1列Xaug1代入式（22），可以得到Xaug1=Φb，則模態(tài)幅值b為

b=Φ+Xaug1（23）

式中：Φ+是矩陣Φ的偽逆矩陣。按照模態(tài)幅值b從大到小對各個模態(tài)重新進行排序，從而得到排序后的DMD模態(tài)^［31］，采用模態(tài)幅值排序的方法可以按照各個模態(tài)對m D曲線的貢獻度和影響程度進行排列，通過前幾個主要模態(tài)可獲得m D曲線的主要頻率分量，且可進一步通過前幾個主要模態(tài)實現(xiàn)m D曲線重構。

為確定移位堆疊次數(shù)h的最優(yōu)值，在避免移位堆疊次數(shù)過多浪費計算資源的同時，確保DMD完全分解并準確捕獲到周期振蕩的m D曲線，提取出其包含的自旋、錐旋頻率等微動特征。本文采用文獻［31］定義的損失函數(shù)來選擇最佳的堆疊次數(shù)：

Floss=100×Xaug－XDMD（t）FXaugF（24）

式中：·F為Frobenius范數(shù);XDMD（t）為重構的數(shù)據(jù)矩陣。

式（24）表明，增廣數(shù)據(jù)矩陣與Xaug數(shù)據(jù)重構越接近，損失函數(shù)的值越小，說明數(shù)據(jù)重構效果越好，散射點的微動回波分解效果越好。若隨著堆疊次數(shù)的增加，損失函數(shù)趨于一個恒定值，說明已經(jīng)能夠完全捕獲到周期振蕩，原則上繼續(xù)堆疊不會再增加DMD算法的準確性。

2.2 平動補償與m D曲線頻率提取

由式（5）可知，m D曲線受到ΔR′（tm）的調制，ΔR′（tm）包含微動項RM（tm）和平動項R′TM（tm），且由式（20）可以看出，DMD得到的每個模態(tài)都具有頻率單一的特點，能夠反映出微動信號所包含的頻率分量。尾翼散射點B對應的m D曲線，主要包括微動項RM（tm）調制產(chǎn)生的Ωs+Ωc、Ωs、Ωc和|Ωs-Ωc|這4種頻率分量和平動項R′TM（tm）的頻率分量，平動項的頻率分量為零。對m D曲線進行DMD并按照式（23）定義的幅值b排序后，前5個模態(tài)的幅值相比其他模態(tài)的幅值是最大的，這5個模態(tài)的頻率包含了0 Hz和Ωs+Ωc、Ωs、Ωc、|Ωs-Ωc|這5種主要頻率分量，其中0 Hz頻率分量即對應平動分量。對m D曲線進行DMD可提取出自旋頻率、錐旋頻率等微動特征以及需要進行平動補償?shù)牧泐l率模態(tài)，其余幅度較小的模態(tài)則是在距離維采樣中引入的量化誤差，導致分離出的m D曲線并非是理想和光滑的。

無論是對錐頂散射點A還是對尾翼散射點B或C而言，其平動項都是相同的，對m D曲線DMD得到的零頻率模態(tài)也都是相同的，因此對其中任一散射點的m D曲線，利用DMD得到的零頻率模態(tài)就能完成平動補償。具體而言，慢時間-距離像的每一列為目標的一次距離像，利用零頻率模態(tài)，采用循環(huán)移位方法，對距離像進行移位補償處理即可完成平動補償。

綜上分析，基于DMD的彈道目標平動補償與微動特征提取方法可以總結為如圖2所示的流程圖。首先，對分離出的m D曲線移位堆疊構建增廣數(shù)據(jù)矩陣，并進行DMD，對分解得到的模態(tài)按照模態(tài)幅值的大小進行排序;進一步判斷損失函數(shù)是否趨于恒定值，如果不是，繼續(xù)增加堆疊的次數(shù)，直到完成DMD;然后，再次使用損失函數(shù)選取模態(tài)數(shù)量，如果損失函數(shù)不趨于恒定值，則繼續(xù)增加模態(tài)數(shù)量，直至損失函數(shù)趨于恒定值;最后，選取模態(tài)中的零頻率模態(tài)進行平動補償，通過其他主要模態(tài)提取自旋頻率、錐旋頻率等微動特征，同時利用選取的主要模態(tài)可重構出目標散射點的m D曲線。

3 仿真分析

3.1 算法可行性分析與仿真驗證

雷達及彈道導彈目標仿真參數(shù)設置如表1所示。首先，在不考慮噪聲和雜波的條件下，進動彈頭目標慢時間-距離像序列的預處理結果如圖3所示，其中m D曲線分離結果如圖3（c）所示。可見，目標散射點在慢時間-距離平面（圖像域）的m D曲線呈現(xiàn)出正弦或近似正弦的變化規(guī)律，且由于平動項的存在，所有m D曲線整體呈現(xiàn)出向下的趨勢，給后續(xù)的微動特征提取帶來困難。

以散射點B為例，對分離出的m D曲線組成的數(shù)據(jù)向量移位堆疊后形成增廣數(shù)據(jù)矩陣。為了確定最優(yōu)堆疊次數(shù)，結合式（24）分析損失函數(shù)與堆疊次數(shù)的關系。從圖4中可以看出，針對分離后的目標散射點微動數(shù)據(jù)，堆疊次數(shù)為550時損失函數(shù)為0.88%，當堆疊次數(shù)為600時損失函數(shù)迅速下降到0.16%，隨著堆疊次數(shù)增加到約700，損失函數(shù)已經(jīng)趨于一個恒定值，約為0.06%左右，說明此時m D曲線的周期振蕩已經(jīng)被完全捕獲，繼續(xù)增加堆疊次數(shù)不會增加DMD算法的準確性和精確度，反而會帶來計算效率的下降。

進一步可以看到，根據(jù)式（23）定義的模態(tài)幅值對分解后的各個模態(tài)進行排序后，模態(tài)的特征值分布情況如圖5（a）所示，特征值均成對分布，因為其特征值為共軛的復數(shù)對，所以圖5中的橫軸均為特征值的實部（即re（λ））;同時，每一對模態(tài)（即每兩階模態(tài)）可以看作是單個模態(tài)，圖中模態(tài)幅值越大，特征值對應的圓圈顏色越深。由圖5（b） 模態(tài)特征值分布局部放大圖可以看出，模態(tài)排序后前10階模態(tài)的幅值是最大的（其中最中間的一個深色圓圈為第1～2階模態(tài)的重合），且特征值基本位于單位圓上，說明各個模態(tài)接近于穩(wěn)定狀態(tài)。再次利用式（24）選取模態(tài)數(shù)量，確保得到需要進行平動補償?shù)牧泐l率模態(tài)以及含有微動頻率的其他主要模態(tài)。圖6為DMD模態(tài)數(shù)量與損失函數(shù)的關系?？梢钥闯?，前10階模態(tài)的損失函數(shù)恒定于0.03%左右，已經(jīng)能夠反映m D曲線的主要頻率構成，因此通過前10階模態(tài)也能夠較好地實現(xiàn)m D曲線的重構。

下面對各個模態(tài)的頻率做進一步分析。圖7為DMD模態(tài)頻率與幅值的關系圖，雖然圖中第1～2階模態(tài)幅值最大，但是頻率為0。結合圖8（a）可以看到，第1～2階模態(tài)即對應目標的平動項，與圖3中m D曲線的基線走勢一致;第3～4階模態(tài)的頻率為3.5 Hz，對應目標的Ωs+Ωc頻率分量;第5～6階模態(tài)的頻率為3 Hz，對應目標的自旋頻率分量，即Ωs分量;第7～8階模態(tài)的頻率為0.500 2 Hz，對應目標的錐旋頻率分量，即Ωc分量;第9～10階模態(tài)的幅值最小，頻率為2.504 Hz，對應目標的|Ωs-Ωc|分量。在圖8（b）中，第3～10階模態(tài)基本呈周期性的簡諧運動，對應圖7中相應模態(tài)的頻率?？梢钥吹?，每個模態(tài)的頻率是單一的，對應目標回波中一個確定的頻率分量，這是DMD算法的一個優(yōu)點?？梢姡珼MD可以準確提取出彈道目標的自旋頻率、錐旋頻率等微動特征。

模態(tài)能量也是反映模態(tài)重要性的指標之一，表2為前i（i=2，4，6，8，10）階模態(tài)能量和在模態(tài)總能量中的占比?？梢钥闯觯?階模態(tài)能量和接近于模態(tài)總能量的100%，前2階模態(tài)對應的是彈道目標的平動項，反映了目標主體的多普勒回波。第3～10階的模態(tài)能量相對較小，反映了微動散射點的m D回波信號。表3為去除前2階模態(tài)后，前i（i=2，4，6，8）階模態(tài)（即原來的第3～10階模態(tài)）能量和在剩余模態(tài)總能量中的占比?？梢?，這8階模態(tài)能量和接近于100%，占據(jù)了目標微動散射點回波能量的絕大部分，這8階模態(tài)表示的是進動時彈道目標的Ωs+Ωc、Ωs、Ωc和|Ωs-Ωc|這4種頻率分量。

進一步，對目標散射點的微動回波在慢時間-距離域進行重構。如圖9所示，利用前10階模態(tài)可以較好地重構出目標散射點的m D曲線，與DMD前的m D曲線十分接近，可見前10階主要模態(tài)基本包含了目標微動散射點的回波信息，能夠體現(xiàn)出m D曲線的頻率構成分量。

下面對A、B和C這3個散射點的3條m D曲線各自進行DMD，可得到每條曲線的第1階和第2階模態(tài)的疊加結果，對每條m D曲線的平動項取平均后進行平動補償，補償結果如圖10所示?？梢姡絼訋淼膍 D曲線向下傾斜并得到抑制，取得了較好的補償效果，便于后續(xù)微動特征的準確提取。

3.2 噪聲條件下算法穩(wěn)健性分析與對比實驗

為驗證本文方法在低信噪比條件下的穩(wěn)健性，在目標回波中加入信噪比為-15 dB的白噪聲。限于文章篇幅，堆疊次數(shù)及模態(tài)數(shù)量選取等不再贅述，選取規(guī)則與第3.1節(jié)一致。

下面重點分析低信噪比條件下基于DMD算法的平動補償及微動特征提取。從圖11（a）的特征值分布情況可以看出，加噪后分解得到的特征值仍然成對分布; 從圖11（b）特征值分布情況的局部放大圖可以看出，雖然有些特征值受噪聲影響，偏離了單位圓（即出現(xiàn)了增長或衰減的情況），但前10階模態(tài)的特征值仍然位于單位圓上，這表明前10階模態(tài)基本包含了目標m D曲線的主要頻率特征。

從圖12（a）的模態(tài)頻率與幅值關系圖可以看出，第1～2階模態(tài)的幅值仍然最大，頻率為0 Hz，同樣對應彈道目標的平動項;從圖12（b）的模態(tài)頻率與幅值關系的局部放大圖可以看出，第3～10階模態(tài)仍然對應目標回波的Ωs+Ωc、Ωs、Ωc和|Ωs-Ωc|這4種頻率分量，其頻率位置仍然較為準確。從圖13中DMD后的前10階模態(tài)可以看出，當信噪比為-15 dB時，DMD后的第1～2階模態(tài)仍然能夠較好地反映出目標的平動項。同樣，第3～10階模態(tài)仍呈周期性運動，包含了目標散射點的自旋、錐旋頻率分量?？梢姡诘托旁氡葪l件下，DMD仍然可以準確地獲取目標m D曲線所包含的微動特征。

圖14為低信噪比條件下的平動補償結果?？梢姡诘托旁氡葧r，平動補償仍取得了較好的效果，表明DMD算法具有較好的抗噪性能。

為進一步驗證應用DMD算法進行平動補償及微動特征提取的準確性，在信噪比為-15 dB的條件下分別對目標在慢時間-距離域的回波采用EEMD、VMD等微動特征提取算法進行對比分析與仿真驗證，仿真參數(shù)設置與第3.1節(jié)相同。

圖15給出采用EEMD算法的分解及平動補償結果，從圖15（a）可以看到對散射點A的m D曲線EEMD的結果。由于在慢時間-距離像序列的預處理步驟中，距離維采樣引入了量化誤差，導致分離出的m D曲線不再是連續(xù)光滑的，不可避免地引入了高頻分量。同時，由于EEMD算法的分解結果一般是從高頻到低頻排列，所以首先得到的分量IMF1和IMF2是高頻分量，分量IMF3的頻率為0.5 Hz，對應目標的錐旋頻率分量Ωc，分量IMF4對應平動項。圖15（b）為對散射點B的m D曲線的EEMD結果，分量IMF3頻率為3.5 Hz，對應目標的Ωs+Ωc分量，而3 Hz的自旋頻率Ωs分量和2.5 Hz的|Ωs-Ωc|分量與3.5 Hz相比較為接近，難以準確進行分解。同時，從分量IMF4中可以看到，0.5 Hz的錐旋頻率Ωc分量也難以準確獲取?？梢?，雖然EEMD算法是EMD的改進算法，在一定程度上克服了模態(tài)混疊現(xiàn)象，但仍不能完全消除模態(tài)混疊現(xiàn)象，分量IMF5即為平動項。圖15（c）為使用EEMD得到的趨勢項進行平動補償?shù)慕Y果，由于端點效應的存在，從平動補償結果中可以看到，兩側的微動曲線由于受到端點效應的影響而產(chǎn)生了形變。

圖16為采用VMD算法的分解及平動補償結果。圖16（a）為對散射點A的m D曲線VMD的結果，其中分量IMF1為平動項，分量IMF2的頻率為0.5 Hz，對應目標的錐旋頻率分量Ωc，其他分量為高頻分量。圖16（b）為對散射點B的m D曲線進行VMD的結果，分量IMF1仍為平動項，分量IMF2的頻率為3.5 Hz，對應目標的Ωs+Ωc分量，分量IMF3的頻率為0.5 Hz，對應目標的錐旋頻率分量Ωc，分量IMF4的頻率為3 Hz，對應自旋頻率分量Ωs，其他分量為高頻分量。圖16（c）為對使用VMD得到的趨勢項進行平動補償?shù)慕Y果。由于VMD算法同樣也存在端點效應，兩側的微動曲線受其影響產(chǎn)生了形變，但受影響程度相比于EEMD算法而言較小。

下面在同一計算機配置條件下對EEMD、VMD和DMD這3種算法的運算時間進行比較，計算機硬件配置條件為Intel i5 8265U處理器，內存為8G，Windows 10操作系統(tǒng)。以對散射點B回波的模態(tài)分解為例，3種算法的運算時間如表4所示?？梢?，EEMD算法的運算時間最短，VMD算法運算時間最長，而DMD算法運算時間介于兩者之間，約為6.57 s。

同樣以散射點B為例，綜合比較3種算法的性能?？梢?，雖然EEMD算法在運算時間上具有優(yōu)勢，耗時相對最短，但是該算法分解得到的頻率分量是最少的，僅包含Ωs+Ωc分量，平動補償結果受端點效應的影響最大。相比于EEMD算法，VMD算法獲取的頻率分量相對較多，平動補償結果受端點效應的影響相對較小，但運算時間相對是最長的。而本文所提DMD算法能夠準確提取4種頻率分量，且在低信噪比條件下依然運行相對穩(wěn)健，運算時間介于EEMD、VMD兩種算法之間，平動補償效果較前兩種算法更佳。

4 結束語

本文提出基于DMD的彈道目標平動補償與微動特征提取方法，建立相應的數(shù)學模型，并通過仿真分析所提方法的可行性、穩(wěn)健性與運算效率，可以得出以下結論。

（1） 相比于EEMD和VMD等傳統(tǒng)算法，DMD算法能夠更好地實現(xiàn)彈道目標的平動補償。

（2） 在描述周期性的觀測量（如m D曲線）隨時間演化的過程中，對m D曲線的DMD不僅可以準確提取到目標回波中的主要頻率分量，且每個模態(tài)對應的頻率具有單一的特點，相比于EMD類方法（包括EMD、EEMD等），未產(chǎn)生模態(tài)混疊現(xiàn)象，有利于準確地提取出目標的微動頻率分量。

（3） 相比于EEMD算法和VMD算法，基于DMD的彈道目標微動特征提取方法的耗時介于兩種算法之間。隨著硬件水平的提升和快速算法的改進，預期可以進一步降低運算量，提升算法的實時性。

基于以上分析與仿真實驗可知，基于DMD的彈道目標平動補償與微動特征提取方法是可行的，且在低信噪比條件下仍然能夠較為準確地提取出目標的微動頻率特征，體現(xiàn)出較好的穩(wěn)健性，后續(xù)針對DMD的快速算法及實測數(shù)據(jù)條件下的平動補償與微動特征提取值得進一步深入研究。

參考文獻

［1］CHEN V C， LI F Y， HO S S， et al. Analysis of micro Doppler signatures［J］. IEE Proceedings Radar Sonar and Navigation， 2003， 150（4）： 271-276.

［2］CHEN V C， LI F Y， HO S S， et al. Micro Doppler effect in radar phenomenon， model and simulation study［J］. IEEE Trans.on Aerospace and Electronic Systems， 2006， 42（1）： 2-21.

［3］CHEN V C. Doppler signatures of radar backscattering from objects with micro motions［J］. IET Signal Processing， 2008， 2（3）： 291-300.

［4］ZHANG Q， LUO Y， CHEN Y A. Micro Doppler characteristics of radar targets［M］. London： Elsevier Press， 2016.

［5］ZHANG Y P， ZHANG Q， KANG L， et al. End to end recognition of similar space cone cylinder targets based on complex va lued coordinate attention networks［J］. IEEE Trans.on Geoscience and Remote Sensing， 2022， 60： 5106214.

［6］ZHANG Y P， ZHANG L， KANG L， et al. Space target classification with corrupted HRRP sequences based on temporal spatial feature aggregation network［J］. IEEE Trans.on Geoscience and Remote Sensing， 2023， 61： 5100618.

［7］ZHANG Y P， XIE Y， KANG L， et al. Feature level fusion recognition of space targets with composite micromotion［J］. IEEE Trans.on Aerospace and Electronic Systems， 2024， 60（1）： 934-951.

［8］LEE J I， KIM N， MIN S， et al. Space target classification improvement by generating micro Doppler signatures considering incident angle［J］. Sensors， 2022， 22（4）： 1653.

［9］HANIF A， MUAZ M， HASAN A， et al. Micro Doppler based target recognition with radars： a review［J］. IEEE Sensors Journal， 2022， 22（4）： 2948-2961.

［10］WANG Z H， LUO Y， LI K M， et al. Micro Doppler parameters extraction of precession cone shaped targets based on rotating antenna［J］. Remote Sensing， 2022， 14（11）： 2549.

［11］TIAN X D， BAI X R， XUE R H， et al. Fusion recognition of space targets with micro motion［J］. IEEE Trans.on Aerospace and Electronic Systems， 2022， 58（4）： 3116-3125.

［12］羅迎， 柏又青， 張群， 等. 彈道目標平動補償與微多普勒特征提取方法［J］. 電子與信息學報， 2012， 34（3）： 602-608.

LUO Y， BAI Y Q， ZHANG Q， et al. Translational motion compensation and micro Doppler feature extraction of ballistic targets［J］. Journal of Electronics amp; Information Technology， 2012， 34（3）： 602-608.

［13］GU F F， FU M H， LIANG B S， et al. Translational motion compensation and micro Doppler feature extraction of space spinning targets［J］. IEEE Geoscience and Remote Sensing Letters， 2018， 15（10）： 1550-1554.

［14］LI J Q， HE S S， FENG C Q， et al. Method for compensating translational motion of rotationally symmetric target based on local symmetry cancellation［J］. Journal of Systems Engineering and Electronics， 2017， 28（1）： 36-39.

［15］ZHUO Z Y， DU L， LU X F， et al. Ambiguity function based high order translational motion compensation［J］. IEEE Trans.on Aerospace and Electronic Systems， 2023， 59（2）： 2013-2019.

［16］馮存前， 李江， 黃大榮， 等. 彈道中段不同平動多目標的平動參數(shù)估計方法［J］. 電子與信息學報， 2021， 43（3）： 564-571.

FENG C Q， LI J， HUANG D R， et al. Estimation method of translational parameters for different translational of ballistic targets in midcourse［J］. Journal of Electronics amp; Information Technology， 2021， 43（3）： 564-571.

［17］SCHMID P J. Dynamic mode decomposition of numerical and experimental data［J］. Journal of Fluid Mechanics， 2010， 656： 5-28.

［18］寇家慶， 張偉偉， 高傳強. 基于POD和DMD方法的跨聲速抖振模態(tài)分析［J］. 航空學報， 2016， 37（9）： 2679-2689.

KOU J Q， ZHANG W W， GAO C Q. Modal analysis of transonic buffet based on POD and DMD method［J］. Acta Aeronautica et Astronautica Sinica， 2016， 37（9）： 2679-2689.

［19］SCHMID P J. Dynamic mode decomposition and its variants［J］. Annual Review of Fluid Mechanics， 2022， 54： 225-254.

［20］ROWLEY C W， DAWSON S T M. Model reduction for flow analysis and control［J］. Annual Review of Fluid Mechanics， 2017， 49： 387-417.

［21］BERGER E， SASTUBA M， VOGT D， et al. Estimation of perturbations in robotic behavior using dynamic mode decomposition［J］. Advanced Robotics： the International Journal of the Robotics Society of Japan， 2015， 29（5/6）： 331-343.

［22］鄭建擁， 魏光輝. 基于多分辨率動態(tài)模態(tài)分解的電磁信號時頻-能量分析［J］. 系統(tǒng)工程與電子技術， 2022， 44（5）： 1468-1474.

ZHENG J Y， WEI G H. Time frequency energy analysis of electromagnetic signals based on multi resolution dynamic modal decomposition［J］. Systems Engineering and Electronics， 2022， 44（5）： 1468-1474.

［23］BRUNTON B W， JOHNSON L A， OJEMANN J G， et al. Extracting spatial temporal coherent patterns in large scale neural recordings using dynamic mode decomposition［J］. Journal of Neuroscience Methods， 2016， 258： 1-15.

［24］PROCTOR J L， ECKHOFF P A. Discovering dynamic patterns from infectious disease data using dynamic mode decomposition［J］. International Health， 2015， 7（2）： 139-145.

［25］ROWLEY C W， MEZIC I， BAGHERI S， et al. Spectral analysis of nonlinear flows［J］. Journal of Fluid Mechanics， 2009， 641： 115-127.

［26］寇家慶， 張偉偉. 動力學模態(tài)分解及其在流體力學中的應用［J］. 空氣動力學報， 2018， 36（2）： 163-179.

KOU J Q， ZHANG W W. Dynamic mode decomposition and its applications in fluid dynamics［J］. Acta Aerodynamica Sinica， 2018， 36（2）： 163-179.

［27］DRAGOMIRETSKIY K， ZOSSO D. Variational mode decomposition［J］. IEEE Trans.on Signal Processing， 2013， 62（3）： 531-544.

［28］WU Z H， HUANG N E. Ensemble empirical mode decomposition： a noise assisted data analysis method［J］. Advances in Adaptive Data Analysis， 2009， 1（1）： 1-41.

［29］KUTZ J N， BRUNTON S L， BRUNTON B W， et al. Dynamic mode decomposition： data driven modeling of complex systems［M］. Philadelphia： Society for Industrial and Applied Mathematics，2016.

［30］TU J H， ROWLEY C W， LUCHTENBERG D M， et al. On dynamic mode decomposition： theory and applications［J］. Journal of Computational Dynamics， 2014， 1（2）： 391-421.

［31］JOVANOVIC M R， SCHMID P J， NICHOLS J W. Sparsity promoting dynamic mode decomposition［J］. Physics of Fluids， 2014， 26（2）： 024103.

作者簡介

李開明（1982—），男，副教授，碩士研究生導師，博士，主要研究方向為雷達成像與目標識別。

代肖楠（1991—），男，助理工程師，碩士，主要研究方向為雷達成像與目標識別。

張袁鵬（1984—），男，講師，博士，主要研究方向為雷達信號處理、雷達目標識別。

姚佳文（1995—），女，助理工程師，碩士，主要研究方向為雷達信號處理。

羅 迎（1984—），男，教授，博士研究生導師，博士，主要研究方向為雷達成像與目標識別。