基于多種濾波算法跟蹤臨近空間非彈道式目標(biāo)

北京電子工程總體研究所　秦　雷　謝曉瑛　李君龍

基于多種濾波算法跟蹤臨近空間非彈道式目標(biāo)

北京電子工程總體研究所秦雷謝曉瑛李君龍

近些年來(lái)各國(guó)在臨近空間高超聲速飛行器領(lǐng)域有了迅猛發(fā)展，其中以X-51A和HTV-2高超聲速飛行器飛行試驗(yàn)相繼成功為代表，它們多采用非彈道式機(jī)動(dòng)飛行方式，飛行速度快，具有較高的升阻比，且在大氣層內(nèi)長(zhǎng)時(shí)間飛行，其運(yùn)動(dòng)軌跡往往呈現(xiàn)出“跳躍”特征，因此存在非彈道式目標(biāo)加速度估計(jì)與軌跡跟蹤預(yù)報(bào)難的問(wèn)題。本文介紹了臨近空間目標(biāo)四種典型的非彈道式機(jī)動(dòng)模式，然后對(duì)目標(biāo)彈道方程及跟蹤濾波器設(shè)計(jì)進(jìn)行了介紹，最后使用四種跟蹤濾波方法對(duì)目標(biāo)加速度進(jìn)行估計(jì)。仿真結(jié)果得出針對(duì)以上四種非彈道式機(jī)動(dòng)模式粒子濾波算法效果最好，該算法保證了目標(biāo)跟蹤精度在允許范圍之內(nèi)，Matlab仿真結(jié)果驗(yàn)證了該算法的有效性。

非彈道式；臨近空間；機(jī)動(dòng)模式；跟蹤濾波

1　引言

目標(biāo)［1-2］的彈道估計(jì)是目標(biāo)跟蹤的一個(gè)重要方面。由于臨近空間高超聲速目標(biāo)具有飛行速度快，大范圍機(jī)動(dòng)等特點(diǎn)，而且臨近空間環(huán)境極其惡劣，飛行彈道受到各種擾動(dòng)和氣動(dòng)阻力的影響，使得狀態(tài)和量測(cè)方程是高度非線性的。所以臨近空間機(jī)動(dòng)式彈道濾波問(wèn)題［3-5］是一個(gè)復(fù)雜的非線性濾波問(wèn)題。跟蹤系統(tǒng)的核心是濾波算法，目前大多數(shù)用于跟蹤定位臨近空間機(jī)動(dòng)目標(biāo)是基于卡爾曼濾波算法［6-7］，但是使用該算法會(huì)造成跟蹤濾波精度變差，不能獲得目標(biāo)狀態(tài)的良好估計(jì)，存在非彈道式目標(biāo)加速度估計(jì)與軌跡跟蹤預(yù)報(bào)難的問(wèn)題。因此為了提高目標(biāo)跟蹤濾波器的跟蹤效果與穩(wěn)定性，本文通過(guò)使用卡爾曼濾波算法、擴(kuò)展卡爾曼濾波算法、無(wú)跡卡爾曼濾波算法、粒子濾波算法分別對(duì)非彈道式機(jī)動(dòng)目標(biāo)加速度進(jìn)行了仿真估計(jì)，通過(guò)蒙特卡洛仿真試驗(yàn)得出粒子濾波算法跟蹤精度較高，跟蹤誤差較小，且在允許范圍之內(nèi)，算法比較穩(wěn)定和有效，產(chǎn)生了很好的跟蹤濾波效果，在工程實(shí)踐中有一定的實(shí)用價(jià)值。

2　臨近空間非彈道式機(jī)動(dòng)模式

在地心慣性坐標(biāo)系中，原點(diǎn)OI位于地心上，OIyI指向攔截導(dǎo)彈發(fā)射點(diǎn)，OIxI軸位于攔截導(dǎo)彈射擊平面內(nèi)，與OIyI軸垂直，指向目標(biāo)為正，OIzI按右手定則確定，如圖1所示。根據(jù)地心慣性坐標(biāo)系中彈道傾角及彈道偏角計(jì)算公式得到OIyI方向速度為vyI=0，q代表彈道傾角。

圖1　地心慣性坐標(biāo)系示意圖

根據(jù)運(yùn)動(dòng)方程可以求解得出：

式中：

因此，如果，總可以找到參數(shù)密切相關(guān)。圖2所示給出了不同升力系數(shù)取值情況下維持等高飛行的最小速度隨高度的變化情況。從圖中可以看出，在20km以下，維持等高飛行的最小速度隨高度的升高變化不大，且受升力系數(shù)不確定性的影響相對(duì)較小。，使得目標(biāo)維持等高飛行。速度與目標(biāo)飛行高度為維持等高飛行，需增大飛行攻角，若在最大許用攻角前提下仍然無(wú)法滿足V>Vmin，則等高飛行難以維系，飛行高度將迅速降低?？紤]到由于側(cè)向力主要由升力系數(shù)和傾斜角產(chǎn)生，且近似正比于CLsinσ，為確保側(cè)向機(jī)動(dòng)最大，σ應(yīng)盡可能大。

不失一般性，σ的選取應(yīng)使得在最大情況下，下式取最小值：

圖2　等高條件下飛行高度與飛行速度關(guān)系圖

根據(jù)目標(biāo)維持等高飛行，得到三個(gè)方向目標(biāo)真實(shí)加速度為：

假定初始速度為6km/s，初始時(shí)刻飛行高度為50km，分析得出的典型飛行彈道如圖3-6所示：

圖3　等高條件下飛行高度示意圖

圖4　等高條件下飛行速度示意圖

圖5　等高條件下地面射程示意圖

圖6　等高條件下經(jīng)度與緯度關(guān)系示意圖

此種機(jī)動(dòng)模式主要存在于目標(biāo)巡航段，這種機(jī)動(dòng)模式可以使目標(biāo)橫向機(jī)動(dòng)能力達(dá)到最大，有效提高目標(biāo)機(jī)動(dòng)能力，對(duì)于規(guī)避敵方攔截起到至關(guān)重要的作用，符合高超聲速飛行器大范圍規(guī)避式機(jī)動(dòng)特點(diǎn)。

假定初始速度為6km/s，初始時(shí)刻飛行高度為50km，分析得出的典型飛行彈道如圖7-9所示：

圖7　等動(dòng)壓條件下飛行高度示意圖

圖8　等動(dòng)壓條件下飛行速度示意圖

圖9　等動(dòng)壓條件下地面射程示意圖

式中，κ為玻爾茲曼常數(shù)，ρ為大氣密度，г為目標(biāo)質(zhì)心與地心距離。

參數(shù)C可通過(guò)跟蹤濾波進(jìn)行確定。上式可作為辨別等動(dòng)壓飛行的基本判據(jù)之一。圖10所示給出了實(shí)際飛行彈道與2lnV=κг+c擬合結(jié)果的對(duì)比。

圖10　等動(dòng)壓條件下飛行速度與飛行高度示意圖

圖11　動(dòng)壓隨飛行高度變化示意圖

從圖11中可以看出，在動(dòng)壓變化不大的飛行段，上述規(guī)律吻合較好。

另外，通過(guò)相關(guān)近似分析，在阻力系數(shù)變化不大的情況下，還可近似認(rèn)為彈道傾角近似與速度平方成反比，這也可以作為判斷等動(dòng)壓模式的依據(jù)。從圖12所示可以看出，這一結(jié)論在等動(dòng)壓飛行段也較為吻合。

等動(dòng)壓模式還可得到彈道傾角隨飛行時(shí)間的對(duì)應(yīng)關(guān)系：

式中，κ為玻爾茲曼常數(shù)，?為升阻比。

圖12　彈道傾角隨飛行時(shí)間變化關(guān)系圖

因此，對(duì)于等動(dòng)壓飛行段，可近似解析求出飛行高度、飛行速度、彈道傾角、地面射程變化規(guī)律，對(duì)進(jìn)行非彈道式機(jī)動(dòng)目標(biāo)軌跡跟蹤具有重要意義。

該種機(jī)動(dòng)模式由于考慮到飛行器結(jié)構(gòu)防熱以及發(fā)動(dòng)機(jī)工作條件要求，一般在飛行上升段使用該模式。

設(shè)定升阻比?=L/D，假定初始速度為6km/s，初始時(shí)刻飛行高度為90km，通過(guò)計(jì)算得到最大升阻比時(shí)的攻角，并保持該攻角、常值傾斜角飛行。分析得出的典型飛行彈道如圖13-15所示：

圖13　最大升阻比條件下飛行高度變化示意圖

圖14　最大升阻比條件下飛行速度變化示意圖

圖15　最大升阻比條件下地面射程變化示意圖

該模式主要用于飛行巡航段和滑翔段，以使高超聲速飛行器達(dá)到射程最遠(yuǎn)，機(jī)動(dòng)能力最強(qiáng)。

針對(duì)AHW滑翔彈頭，升阻比?=L/D，假定初始速度為6km/s，初始時(shí)刻飛行高度為90km，分析得出的典型飛行彈道如下圖16-18所示：

圖16　常值攻角和傾斜角條件下飛行高度變化示意圖

圖17　常值攻角和傾斜角條件下飛行速度變化示意圖

圖18　常值攻角和傾斜角條件下地面射程變化示意圖

該模式主要用于飛行巡航段，一般在距離敵方目標(biāo)較遠(yuǎn)的位置飛行時(shí)可以采用該模式。

3　目標(biāo)彈道方程及跟蹤濾波器設(shè)計(jì)

3.1目標(biāo)彈道方程

由于目標(biāo)機(jī)動(dòng)模式為側(cè)向擺式機(jī)動(dòng)，高度和速度基本保持不變，設(shè)置目標(biāo)位置坐標(biāo)為（X,Y,Z），速度為（X,Y,Z），加速度為（X,Y,Z）。

目標(biāo)彈道方程為：

當(dāng)高度在20km-30km時(shí)，目標(biāo)的最大控制過(guò)載為3；當(dāng)高度在30km-40km時(shí)，目標(biāo)的最大控制過(guò)載為2；當(dāng)高度在40km以上時(shí)，目標(biāo)的最大控制過(guò)載小于1。

目標(biāo)在地心慣性系的OIZI軸方向可以作機(jī)動(dòng)突防，采用擺式機(jī)動(dòng)，具體用不同頻率的正弦變化的加速度組合來(lái)描述，例如：

式中，wz1、wz2、wz3代表不同的機(jī)動(dòng)頻率，Az1、Az2、Az3代表各個(gè)機(jī)動(dòng)頻率相對(duì)應(yīng)的最大機(jī)動(dòng)幅值。

3.2目標(biāo)跟蹤濾波器設(shè)計(jì)

由于目標(biāo)在OIXI軸方向上只受重力加速度分量-μX/r3以及擾動(dòng)加速度aaI的作用，而目標(biāo)在地心慣性系的OIyI軸方向上的控制加速度aaI的作用是保持其飛行高度不變，在這兩個(gè)方向上目標(biāo)不會(huì)做較強(qiáng)烈的機(jī)動(dòng)突防。因此用Singer模型來(lái)分別描述目標(biāo)加速度這兩個(gè)分量的變化，即：

式中，γx和γy分別代表OIXI軸方向和OIyI軸方向上目標(biāo)機(jī)動(dòng)時(shí)間常數(shù)的倒數(shù)，waI和waI分別代表零均值高斯白噪聲。

目標(biāo)在慣性坐標(biāo)系的Z軸方向上的多頻率組合機(jī)動(dòng)用正弦形式的機(jī)動(dòng)來(lái)描述。設(shè)：

式中，Atz代表目標(biāo)加速度的最大值，wtz代表正弦目標(biāo)加速度的頻率，ξtz代表正弦型目標(biāo)加速度的初始相位。對(duì)（12）式求兩次導(dǎo)數(shù)得到：

由（12）式可得：

可以看出，目標(biāo)機(jī)動(dòng)模型只與機(jī)動(dòng)頻率有關(guān)，與機(jī)動(dòng)幅值和初始相位無(wú)關(guān)。

設(shè)目標(biāo)位置在發(fā)射點(diǎn)慣性系中的投影為ra，ra，ra，目標(biāo)速度在發(fā)射點(diǎn)慣性系中的投影為va，va，va，則可以建立起如下三組運(yùn)動(dòng)方程：

對(duì)于系統(tǒng)，定義狀態(tài)變量Xa-[tavaaaI]，則可以得到如下?tīng)顟B(tài)方程：

其中：

對(duì)上式以一定的周期進(jìn)行離散化后，得到 ：

則離散化后的狀態(tài)方程為：

目標(biāo)位置由地面目標(biāo)跟蹤裝置測(cè)量得到，則測(cè)量方程為：

該測(cè)量信息數(shù)據(jù)更新率較低，設(shè)為0.1s更新一次。

這樣根據(jù)上述模型，設(shè)計(jì)如下異步Kalman濾波器，其預(yù)報(bào)方程為：

其中，Qa為模型預(yù)報(bào)誤差協(xié)方差矩陣。濾波器的測(cè)量修正方程為：

同理，對(duì)目標(biāo)在慣性坐標(biāo)系的Y軸方向上的運(yùn)動(dòng)，即系統(tǒng)也采用這一方法設(shè)計(jì)濾波器，對(duì)其位置、速度和加速度進(jìn)行估計(jì)。

然后，應(yīng)用Kalman濾波對(duì)以上系統(tǒng)進(jìn)行估計(jì)。其預(yù)報(bào)方程為：

當(dāng)沒(méi)有測(cè)量信息時(shí)，只運(yùn)行預(yù)報(bào)方程，有測(cè)量信息時(shí)，則同時(shí)運(yùn)行預(yù)報(bào)方程和測(cè)量修正方程。

4　仿真分析

4.1仿真條件

為了簡(jiǎn)化起見(jiàn)，假定地面有兩個(gè)無(wú)源雷達(dá)站對(duì)目標(biāo)進(jìn)行定位跟蹤。設(shè)定基站和目標(biāo)的相對(duì)位置關(guān)系示意圖如圖19所示。試驗(yàn)條件：兩個(gè)雷達(dá)站的位置坐標(biāo)分別為，單位：km，雷達(dá)站采樣周期為T(mén)=0.01s。兩個(gè)雷達(dá)站測(cè)得的方位角分別為。設(shè)定角度，距離測(cè)量精度分別為：，俯仰角分別為，量測(cè)信息為，目標(biāo)距離雷達(dá)站距離

圖19　雙站跟蹤定位目標(biāo)飛行器示意圖

假定目標(biāo)在三維空間中做側(cè)向擺式機(jī)動(dòng)，目標(biāo)運(yùn)動(dòng)時(shí)間從82s到 312s，總計(jì)運(yùn)動(dòng)230s。目標(biāo)初始位置為目標(biāo)初始速度為本次仿真利用卡爾曼濾波算法、擴(kuò)展卡爾曼濾波算法、無(wú)跡卡爾曼濾波算法、粒子濾波算法對(duì)三軸方向加速度估計(jì)及誤差進(jìn)行了仿真分析，總共進(jìn)行了100次蒙特卡羅仿真試驗(yàn)。

利用公式（34）計(jì)算出斜距估計(jì)誤差：

計(jì)算時(shí)設(shè)置如下參數(shù)：

設(shè)x，y，z三個(gè)方向的方差預(yù)測(cè)矩陣分別為：

4.2仿真結(jié)果

經(jīng)過(guò)仿真得到臨近空間目標(biāo)飛行器隨時(shí)間變化的運(yùn)動(dòng)軌跡如圖20所示：

圖20　目標(biāo)運(yùn)動(dòng)軌跡示意圖

使用四種濾波算法分別對(duì)臨近空間飛行器可能存在的四種非彈道式機(jī)動(dòng)模式進(jìn)行跟蹤濾波，分別得到三軸方向加速度估計(jì)及誤差如下圖所示：

a）常值攻角飛行，傾斜角為常數(shù)

三個(gè)方向加速度估計(jì)及誤差如圖21-23所示：

圖21　x方向加速度估計(jì)及誤差

圖22　y方向加速度估計(jì)及誤差

圖23　z方向加速度估計(jì)及誤差

表1　三軸方向加速度均方根誤差比較

從圖21-23和表1可以看出，三個(gè)方向粒子濾波算法加速度估計(jì)誤差均比其余三種濾波算法的誤差要小，說(shuō)明在該種機(jī)動(dòng)模式下粒子濾波算法跟蹤濾波精度更高，跟蹤效果更好。

b）等高飛行，橫向機(jī)動(dòng)最大

三個(gè)方向加速度估計(jì)及誤差如圖24-26所示：

圖24　x方向加速度估計(jì)及誤差

圖25　y方向加速度估計(jì)及誤差

圖26　　z方向加速度估計(jì)及誤差

表2　三軸方向加速度均方根誤差比較

從圖24-26和表2可以看出，三個(gè)方向粒子濾波算法加速度估計(jì)誤差均比其余三種濾波算法的誤差要小，說(shuō)明在該種機(jī)動(dòng)模式下使用粒子濾波算法能夠得到較高的跟蹤濾波精度，可以實(shí)現(xiàn)較高精度跟蹤臨近空間目標(biāo)飛行器。

c）等動(dòng)壓飛行，傾斜角為常數(shù)

三個(gè)方向加速度估計(jì)及誤差如圖27-29所示：

圖27　x方向加速度估計(jì)及誤差

圖28　y方向加速度估計(jì)及誤差

圖29　z方向加速度估計(jì)及誤差

從圖27-29和表3可以看出，三個(gè)方向粒子濾波算法加速度估計(jì)誤差均比其余三種濾波算法的誤差要小，說(shuō)明在該種機(jī)動(dòng)模式下使用粒子濾波算法能夠得到最高的跟蹤濾波精度，使用該算法狀態(tài)估計(jì)跟蹤效果最好。

d）最大升阻比飛行，傾斜角為常數(shù)

三個(gè)方向加速度估計(jì)及誤差如圖30-32所示：

圖31　y方向加速度估計(jì)及誤差

圖32　z方向加速度估計(jì)及誤差

表4　三軸方向加速度均方根誤差比較

從圖30-32和表4可以看出，三個(gè)方向粒子濾波算法加速度估計(jì)誤差均比其余三種濾波算法的誤差要小，說(shuō)明在該種機(jī)動(dòng)模式下使用粒子濾波算法能夠得到最高的跟蹤濾波精度，具有一定的參考作用。

綜上所述，從四種機(jī)動(dòng)模式的仿真結(jié)果來(lái)看，由于“等高飛行，橫向機(jī)動(dòng)最大”機(jī)動(dòng)模式要求攻角取最大值，彈道傾角為0，目標(biāo)的機(jī)動(dòng)性相對(duì)較高，而卡爾曼濾波算法主要是在線性高斯情況下利用最小均方誤差準(zhǔn)則獲得目標(biāo)的動(dòng)態(tài)估計(jì)，該算法濾波得到的加速度估計(jì)誤差最大，跟蹤濾波效果最差，甚至出現(xiàn)了加速度反向的情況。粒子濾波算法相比擴(kuò)展卡爾曼濾波算法和無(wú)跡卡爾曼濾波算法來(lái)說(shuō)，該算法的適應(yīng)環(huán)境最為廣泛，更適合于復(fù)雜的非高斯環(huán)境。所以經(jīng)過(guò)PF算法跟蹤濾波后的加速度估計(jì)精度最好，跟蹤濾波誤差最??；當(dāng)以“常值攻角和常值傾斜角”機(jī)動(dòng)模式飛行時(shí)，由于該模式運(yùn)動(dòng)較為平緩，所以使用四種濾波算法得到的均方根誤差相差不多，并且數(shù)值較小，沒(méi)有出現(xiàn)加速度反向的情況，但是相比來(lái)說(shuō)，PF算法的均方根誤差最?。划?dāng)以“等動(dòng)壓飛行，常值傾斜角”的機(jī)動(dòng)模式飛行時(shí)，由于該模式相比其他模式機(jī)動(dòng)性最強(qiáng)，因此使用卡爾曼濾波算法時(shí)加速度估計(jì)誤差最大，雖然使用EKF算法、UKF算法、PF算法的均方根誤差相比其他機(jī)動(dòng)模式要大，但是PF算法的加速度估計(jì)誤差還是最小，說(shuō)明該算法在強(qiáng)機(jī)動(dòng)性情況下仍然可以保持相對(duì)較好的估計(jì)精度；當(dāng)以“最大升阻比飛行，常值傾斜角”的機(jī)動(dòng)模式飛行時(shí)，該機(jī)動(dòng)模式仍然保持較強(qiáng)的機(jī)動(dòng)性，從y向加速度估計(jì)誤差較大可以看出，PF算法的加速度估計(jì)精度最好，說(shuō)明該算法具有一定的工程實(shí)用價(jià)值。

5　結(jié)論

本文針對(duì)臨近空間機(jī)動(dòng)目標(biāo)跟蹤問(wèn)題，使用了基于KF、EKF、UKF、PF的跟蹤濾波算法。仿真結(jié)果表明在四種臨近空間非彈道式機(jī)動(dòng)模式下，粒子濾波算法的加速度估計(jì)精度均明顯優(yōu)于其余三種跟蹤濾波算法。

［1］孫鵬，楊建軍.臨近空間平臺(tái)防空反導(dǎo)作戰(zhàn)運(yùn)用仿真研究［J］.彈箭與制導(dǎo)學(xué)報(bào)，2011，31（2）：14-16.

［2］孫鵬，楊建軍.臨近空間中繼制導(dǎo)反TBM關(guān)鍵技術(shù)分析［J］.裝備指揮技術(shù)學(xué)院學(xué)報(bào)，2011，22（2）：61-65.

［3］常思江，王中原，牛春峰.基于卡爾曼濾波的彈箭飛行狀態(tài)估計(jì)方法［J］.彈道學(xué)報(bào)，2010，22（3）：94-98.

［4］徐明友.彈道濾波引論［J］.南京理工大學(xué)學(xué)報(bào)，1988（3）：25-31.

［5］修觀，王良明，楊榮軍.線性彈道模型建立與仿真［J］.海軍工程大學(xué)學(xué)報(bào)，2010，22（2）：84-96.

［6］付夢(mèng)印，鄧志紅，張繼偉，編著.Kalman濾波理論及其在導(dǎo)航系統(tǒng)中的應(yīng)用［M］. 北京：科學(xué)出版社，2003，10.

［7］宋文堯，張牙編著.卡爾曼濾波［M］.北京：科學(xué)出版社，1991，1.

謝曉瑛（1962—），北京人，研究員，現(xiàn)就職于北京電子工程總體研究所。

秦雷【通訊作者】（1987—），吉林長(zhǎng)春人，博士，工程師，主要研究方向：導(dǎo)航、制導(dǎo)與控制。