孫照強(qiáng)，王志貴，孟 飛，李陸雨，于 中，陳 燕

（北京無線電測量研究所，北京 100854）

0 引 言

彈道導(dǎo)彈是現(xiàn)代高技術(shù)戰(zhàn)爭中的重要作戰(zhàn)武器，為了實(shí)現(xiàn)有效的攻擊，彈道導(dǎo)彈常采用以下突防方式：飽和式、多波次攻擊，多彈頭突防，釋放輕、重誘餌等。為了實(shí)現(xiàn)對彈道導(dǎo)彈的有效防御，地面跟蹤制導(dǎo)雷達(dá)必須具備跟蹤多批次目標(biāo)的能力以及單批次團(tuán)目標(biāo)跟蹤能力。另外，為了對來襲目標(biāo)進(jìn)行識別并對制導(dǎo)攔截彈進(jìn)行攔截，又要求雷達(dá)具備較高的跟蹤精度。然而，由于雷達(dá)資源有限，當(dāng)雷達(dá)跟蹤目標(biāo)數(shù)目較多時(shí)，需要降低跟蹤數(shù)據(jù)率以保證對當(dāng)前照射區(qū)域內(nèi)所有目標(biāo)的穩(wěn)定跟蹤，這就出現(xiàn)了跟蹤數(shù)據(jù)率和跟蹤精度之間的矛盾。例如，地面制導(dǎo)雷達(dá)當(dāng)前時(shí)刻對10批目標(biāo)進(jìn)行跟蹤，而每批目標(biāo)團(tuán)內(nèi)又有10個(gè)目標(biāo)，此時(shí)需要對100個(gè)目標(biāo)同時(shí)進(jìn)行濾波處理，若采用常規(guī)的濾波算法，計(jì)算量將急劇上升，不利于工程實(shí)現(xiàn)。因此，迫切需要研究一種在低數(shù)據(jù)率下具有高精度跟蹤能力、且計(jì)算量適中的彈道目標(biāo)濾波方法。

彈道導(dǎo)彈從發(fā)射點(diǎn)到落點(diǎn)的整個(gè)軌跡通常分為兩個(gè)階段：主動段和被動段［1］。在導(dǎo)彈整個(gè)飛行過程中，其會受重力、推力、氣動阻力、地球自轉(zhuǎn)偏向力的公共作用，其運(yùn)動特征表現(xiàn)為高度非線性，對其跟蹤的核心問題是設(shè)計(jì)精確的非線性濾波器，對其運(yùn)動狀態(tài)進(jìn)行估計(jì)。高精度彈道導(dǎo)彈跟蹤需要滿足兩個(gè)條件：精確的目標(biāo)運(yùn)動模型以及性能優(yōu)良的濾波方法［2］。

在彈道導(dǎo)彈運(yùn)動模型方面，動力學(xué)建模方法是目前公認(rèn)的較為準(zhǔn)確的一種描述目標(biāo)飛行過程中受力情況的方法［3-5］。在對被動段彈道導(dǎo)彈跟蹤中，使用較多的是球形模型［6-7］，然而當(dāng)雷達(dá)采樣率較低并且彈道導(dǎo)彈目標(biāo)射程較遠(yuǎn)時(shí)，常用的球形模型會無法滿足高精度跟蹤需求，此時(shí)需要采用基于橢球J2修正的重力加速度模型［8-10］。

在彈道導(dǎo)彈濾波方法方面，一般采用一些對于非線性系統(tǒng)具有良好適應(yīng)性的濾波算法［11-13］，其中具有代表性的算法有擴(kuò)展卡爾曼濾波（extended Kalman filter，EKF）［14-15］、無跡卡爾曼濾波（unscented Kalman filter，UKF）［16-18］和粒子濾波（particle filter，PF）［19-21］。EKF算法采用泰勒級數(shù)展開方法將非線性濾波過程轉(zhuǎn)化為一個(gè)近似線性問題，工程實(shí)現(xiàn)簡單，應(yīng)用范圍廣。UKF算法采用Sigma點(diǎn)進(jìn)行非線性逼近，可以取得三階泰勒展開精度［22］，但是其濾波穩(wěn)定性易受中心采樣點(diǎn)的權(quán)值影響，產(chǎn)生較大的波動［23］。PF方法一般用于非高斯系統(tǒng)，但是其計(jì)算量較大，無法進(jìn)行工程應(yīng)用。文獻(xiàn)［24］對彈道目標(biāo)跟蹤的4種濾波器的性能進(jìn)行了比較，然而，目標(biāo)運(yùn)動模型采用了分段常加速模型，對彈道導(dǎo)彈運(yùn)動過程描述不夠精確，導(dǎo)致濾波精度受限。文獻(xiàn)［25］針對彈道系數(shù)己知的再入彈道目標(biāo)進(jìn)行了跟蹤，對不同濾波器性能進(jìn)行了仿真比較，但是其僅對目標(biāo)運(yùn)動過程中受到的重力與空氣阻力進(jìn)行了建模，未考慮地球自轉(zhuǎn)帶來的影響，目標(biāo)運(yùn)動模型不夠精確。文獻(xiàn)［26］提出了一種迭代線性化彈道導(dǎo)彈被動段跟蹤方法，在混合坐標(biāo)系下完成了彈道導(dǎo)彈的跟蹤濾波。文獻(xiàn)［27］利用多模型對彈道導(dǎo)彈在不同飛行階段的運(yùn)動狀態(tài)進(jìn)行描述，并采用粒子濾波方法對其進(jìn)行跟蹤濾波，但是其具有較大的運(yùn)算量。文獻(xiàn)［28］針對彈道導(dǎo)彈目標(biāo)提出一種基于高斯粒子濾波的狀態(tài)依賴交互多模型方法，對彈道導(dǎo)彈3個(gè)飛行階段進(jìn)行跟蹤，取得了較好的效果。

綜上所述，考慮到彈道導(dǎo)彈跟蹤時(shí)數(shù)據(jù)率和精度之間的矛盾，以及目前所用模型不夠精確的問題，本文針對被動段飛行的彈道導(dǎo)彈跟蹤展開研究。首先建立了彈道導(dǎo)彈被動段的精確跟蹤模型，隨后詳細(xì)推導(dǎo)了基于彈道運(yùn)動方程的EKF濾波過程，使雷達(dá)在較低數(shù)據(jù)率下仍然保持較高的跟蹤精度，同時(shí)具有較小的計(jì)算量，最后通過仿真實(shí)驗(yàn)對所提的運(yùn)動模型及濾波方法進(jìn)行了驗(yàn)證，結(jié)果表明本文所提算法計(jì)算量小，且跟蹤精度滿足當(dāng)前制導(dǎo)雷達(dá)跟蹤需求。

1 導(dǎo)彈被動段質(zhì)心運(yùn)動方程

本文主要研究導(dǎo)彈在被動段的跟蹤過程，為此本節(jié)首先給出彈道目標(biāo)在被動段的運(yùn)動模型。

考慮到雷達(dá)測量值為雷達(dá)大地直角坐標(biāo)系下的目標(biāo)徑向距離、方位角、俯仰角，因此將導(dǎo)彈運(yùn)動方程建立在雷達(dá)大地直角坐標(biāo)系下，這樣可以減少坐標(biāo)變換過程，提高濾波效率。為了提高低數(shù)據(jù)率濾波精度，本文將考慮導(dǎo)彈在被動段所受到的諸多作用力，以盡可能建立與實(shí)際運(yùn)動情況相符的運(yùn)動模型。

如圖1所示，以地心OE為坐標(biāo)原點(diǎn)建立地心慣性坐標(biāo)系。雷達(dá)位于O點(diǎn)，P點(diǎn)為彈道上的一點(diǎn)，其在雷達(dá)大地直角坐標(biāo)系中的坐標(biāo)為（x，y，z），將點(diǎn)P的地心矢徑r表示為r＝R0＋ρ，其模r＝r，其中R0為地球半徑，ρ為雷達(dá)距P點(diǎn)的徑向距離。

圖1 彈道上點(diǎn)P和雷達(dá)點(diǎn)的地心矢徑關(guān)系圖Fig.1 Relation between the trajectory point P and geocentric arrow diameter of the radar

設(shè)導(dǎo)彈相對雷達(dá)大地直角坐標(biāo)系的速度V＝（V x，V y，V z）T，其模V＝V。在地面雷達(dá)看來，該階段導(dǎo)彈受到地心引力、離心慣性力、科氏慣性力、空氣阻力的共同影響，其在雷達(dá)大地直角坐標(biāo)系中的導(dǎo)彈質(zhì)心動力學(xué)方程［29-30］為

和gωe分別是引力加速度g在地心矢徑方向和地球自轉(zhuǎn)角速度方向的分量；?為點(diǎn)P的地心緯度值；J2為攝動常數(shù)；ωe為地球自轉(zhuǎn)角速度；μ為引力常數(shù)；ρ為目標(biāo)位置的大氣密度值，按照美國標(biāo)準(zhǔn)大氣76模型［31］計(jì)算；β為目標(biāo)的質(zhì)阻比。

2 基于彈道運(yùn)動方程的EKF濾波算法模型

雷達(dá)量測值為雷達(dá)大地球坐標(biāo)系下的R、Az、E。為了估計(jì)目標(biāo)質(zhì)阻比，將質(zhì)阻比表示為β＝β0eγ，β0為設(shè)置的初始質(zhì)阻比大小，γ為一待估分量。將導(dǎo)彈的運(yùn)動方程建立在雷達(dá)大地直角坐標(biāo)系下，狀態(tài)變量為

導(dǎo)彈在雷達(dá)大地直角坐標(biāo)系中的狀態(tài)微分方程為

為了應(yīng)用EKF的濾波算法，需要求出狀態(tài)函數(shù)和量測函數(shù)的Jacobi矩陣。

狀態(tài)函數(shù)Jacobi矩陣為

當(dāng)?shù)玫綘顟B(tài)函數(shù)和量測函數(shù)的Jacobi矩陣后，結(jié)合EKF［32］過程即可得到最終的基于彈道運(yùn)動方程的EKF算法，詳細(xì)濾波過程如下：

（1）狀態(tài)預(yù)測

式中：Q（k）為系統(tǒng)噪聲方差陣；Φ（k）為狀態(tài)轉(zhuǎn)移矩陣，在本文中可表示為

（2）狀態(tài)估計(jì)

式中：增益K（k＋1）可表示為

式中：R（k＋1）為量測噪聲方差陣。

本文研究的基于彈道運(yùn)動方程的EKF算法，采用了基于橢球地球模型的J2修正重力加速度模型，同時(shí)考慮了該被動段導(dǎo)彈受地心引力、離心慣性力、科氏慣性力、空氣阻力的共同影響，建立了精確的導(dǎo)彈質(zhì)心動力學(xué)方程。將該精確模型與非線性EKF濾波結(jié)合從而實(shí)現(xiàn)了彈道目標(biāo)的跟蹤精度。

3 BEKF濾波算法的仿真和驗(yàn)證

3.1 彈道目標(biāo)運(yùn)動場景

仿真用的彈道目標(biāo)其射程為1 000 km，地面雷達(dá)觀測的徑向距離、方位角、俯仰角以及全速度變化曲線如圖2～圖5所示，目標(biāo)再入速度最大值接近3 km／s。

圖2 目標(biāo)距離與時(shí)間的變化情況Fig.2 Variation of target range with time

圖3 目標(biāo)方位角與時(shí)間的變化情況Fig.3 Variation of target azimuth with time

圖4 目標(biāo)仰角與時(shí)間的變化情況Fig.4 Variation of target pitch angle with time

圖5 目標(biāo)全速度與時(shí)間的變化情況Fig.5 Variation of target velocity with time

3.2 實(shí)驗(yàn)驗(yàn)證

采用上述目標(biāo)場景，設(shè)定目標(biāo)質(zhì)阻比初值為6 000 kg／m2。設(shè)雷達(dá)距離、方位和俯仰量測噪聲服從相互獨(dú)立的零均值正態(tài)分布，其方差根據(jù)雷達(dá)威力、目標(biāo)散射界面等計(jì)算得到。為了驗(yàn)證本文提出算法的有效性和優(yōu)點(diǎn)，在雷達(dá)跟蹤數(shù)據(jù)率為5 Hz情況下分別對傳統(tǒng)的基于常加速模型的EKF算法、基于彈道運(yùn)動方程的UKF算法以及本文所提的基于彈道運(yùn)動方程的EKF算法進(jìn)行比較，其中基于彈道運(yùn)動方程的UFK方法采樣點(diǎn)數(shù)為13，尺度參數(shù)設(shè)為0.1。以上3種方法分別進(jìn)行30次蒙特卡羅仿真，結(jié)果如圖6～圖8所示。

圖6 距離的均方根誤差曲線Fig.6 Range root mean square error

圖7 方位角的均方根誤差曲線Fig.7 Azimuth root mean square error

圖8 俯仰角的均方根誤差曲線Fig.8 Pitch angle root mean square error

從實(shí)驗(yàn)結(jié)果中可以看出：

（1）本文提出的基于彈道運(yùn)動方程的EKF算法，濾波精度和收斂速度明顯優(yōu)于傳統(tǒng)的基于常加速模型的EKF算法，與基于彈道運(yùn)動方程的UKF算法相當(dāng)；

（2）本文所提的算法在濾波精度上與基于彈道運(yùn)動方程的UKF算法相當(dāng)，考慮到計(jì)算時(shí)間和計(jì)算穩(wěn)定性的雙重優(yōu)勢，基于彈道運(yùn)動方程的EKF更具有工程應(yīng)用價(jià)值。

4 結(jié) 論

針對地面跟蹤雷達(dá)在低數(shù)據(jù)率下對多目標(biāo)的高精度需求，本文提出了基于彈道運(yùn)動方程的EKF算法，基于精確的彈道導(dǎo)彈被動段質(zhì)心運(yùn)動方程，推導(dǎo)了EKF濾波過程。通過與傳統(tǒng)的基于常加速模型的EKF算法和基于彈道運(yùn)動方程的UKF算法比較，驗(yàn)證了基于彈道運(yùn)動方程的EKF具有低數(shù)據(jù)率下濾波精度高、計(jì)算量小等優(yōu)點(diǎn)，解決了地面跟蹤雷達(dá)實(shí)際中遇到的問題，具有工程應(yīng)用價(jià)值。