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

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

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

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

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

非彈道式；臨近空間；機動模式；跟蹤濾波

1　引言

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

2　臨近空間非彈道式機動模式

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

圖1　地心慣性坐標系示意圖

根據運動方程可以求解得出：

式中：

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

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

圖2　等高條件下飛行高度與飛行速度關系圖

根據目標維持等高飛行，得到三個方向目標真實加速度為：

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

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

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

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

圖6　等高條件下經度與緯度關系示意圖

此種機動模式主要存在于目標巡航段，這種機動模式可以使目標橫向機動能力達到最大，有效提高目標機動能力，對于規避敵方攔截起到至關重要的作用，符合高超聲速飛行器大范圍規避式機動特點。

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

圖7　等動壓條件下飛行高度示意圖

圖8　等動壓條件下飛行速度示意圖

圖9　等動壓條件下地面射程示意圖

式中，κ為玻爾茲曼常數，ρ為大氣密度，г為目標質心與地心距離。

參數C可通過跟蹤濾波進行確定。上式可作為辨別等動壓飛行的基本判據之一。圖10所示給出了實際飛行彈道與2lnV=κг+c擬合結果的對比。

圖10　等動壓條件下飛行速度與飛行高度示意圖

圖11　動壓隨飛行高度變化示意圖

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

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

等動壓模式還可得到彈道傾角隨飛行時間的對應關系：

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

圖12　彈道傾角隨飛行時間變化關系圖

因此，對于等動壓飛行段，可近似解析求出飛行高度、飛行速度、彈道傾角、地面射程變化規律，對進行非彈道式機動目標軌跡跟蹤具有重要意義。

該種機動模式由于考慮到飛行器結構防熱以及發動機工作條件要求，一般在飛行上升段使用該模式。

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

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

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

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

該模式主要用于飛行巡航段和滑翔段，以使高超聲速飛行器達到射程最遠，機動能力最強。

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

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

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

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

該模式主要用于飛行巡航段，一般在距離敵方目標較遠的位置飛行時可以采用該模式。

3　目標彈道方程及跟蹤濾波器設計

3.1目標彈道方程

由于目標機動模式為側向擺式機動，高度和速度基本保持不變，設置目標位置坐標為（X,Y,Z），速度為（X,Y,Z），加速度為（X,Y,Z）。

目標彈道方程為：

當高度在20km-30km時，目標的最大控制過載為3；當高度在30km-40km時，目標的最大控制過載為2；當高度在40km以上時，目標的最大控制過載小于1。

目標在地心慣性系的OIZI軸方向可以作機動突防，采用擺式機動，具體用不同頻率的正弦變化的加速度組合來描述，例如：

式中，wz1、wz2、wz3代表不同的機動頻率，Az1、Az2、Az3代表各個機動頻率相對應的最大機動幅值。

3.2目標跟蹤濾波器設計

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

式中，γx和γy分別代表OIXI軸方向和OIyI軸方向上目標機動時間常數的倒數，waI和waI分別代表零均值高斯白噪聲。

目標在慣性坐標系的Z軸方向上的多頻率組合機動用正弦形式的機動來描述。設：

式中，Atz代表目標加速度的最大值，wtz代表正弦目標加速度的頻率，ξtz代表正弦型目標加速度的初始相位。對（12）式求兩次導數得到：

由（12）式可得：

可以看出，目標機動模型只與機動頻率有關，與機動幅值和初始相位無關。

設目標位置在發射點慣性系中的投影為ra，ra，ra，目標速度在發射點慣性系中的投影為va，va，va，則可以建立起如下三組運動方程：

對于系統，定義狀態變量Xa-[tavaaaI]，則可以得到如下狀態方程：

其中：

對上式以一定的周期進行離散化后，得到 ：

則離散化后的狀態方程為：

目標位置由地面目標跟蹤裝置測量得到，則測量方程為：

該測量信息數據更新率較低，設為0.1s更新一次。

這樣根據上述模型，設計如下異步Kalman濾波器，其預報方程為：

其中，Qa為模型預報誤差協方差矩陣。濾波器的測量修正方程為：

同理，對目標在慣性坐標系的Y軸方向上的運動，即系統也采用這一方法設計濾波器，對其位置、速度和加速度進行估計。

然后，應用Kalman濾波對以上系統進行估計。其預報方程為：

當沒有測量信息時，只運行預報方程，有測量信息時，則同時運行預報方程和測量修正方程。

4　仿真分析

4.1仿真條件

為了簡化起見，假定地面有兩個無源雷達站對目標進行定位跟蹤。設定基站和目標的相對位置關系示意圖如圖19所示。試驗條件：兩個雷達站的位置坐標分別為，單位：km，雷達站采樣周期為T=0.01s。兩個雷達站測得的方位角分別為。設定角度，距離測量精度分別為：，俯仰角分別為，量測信息為，目標距離雷達站距離

圖19　雙站跟蹤定位目標飛行器示意圖

假定目標在三維空間中做側向擺式機動，目標運動時間從82s到 312s，總計運動230s。目標初始位置為目標初始速度為本次仿真利用卡爾曼濾波算法、擴展卡爾曼濾波算法、無跡卡爾曼濾波算法、粒子濾波算法對三軸方向加速度估計及誤差進行了仿真分析，總共進行了100次蒙特卡羅仿真試驗。

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

計算時設置如下參數：

設x，y，z三個方向的方差預測矩陣分別為：

4.2仿真結果

經過仿真得到臨近空間目標飛行器隨時間變化的運動軌跡如圖20所示：

圖20　目標運動軌跡示意圖

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

a）常值攻角飛行，傾斜角為常數

三個方向加速度估計及誤差如圖21-23所示：

圖21　x方向加速度估計及誤差

圖22　y方向加速度估計及誤差

圖23　z方向加速度估計及誤差

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

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

b）等高飛行，橫向機動最大

三個方向加速度估計及誤差如圖24-26所示：

圖24　x方向加速度估計及誤差

圖25　y方向加速度估計及誤差

圖26　　z方向加速度估計及誤差

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

從圖24-26和表2可以看出，三個方向粒子濾波算法加速度估計誤差均比其余三種濾波算法的誤差要小，說明在該種機動模式下使用粒子濾波算法能夠得到較高的跟蹤濾波精度，可以實現較高精度跟蹤臨近空間目標飛行器。

c）等動壓飛行，傾斜角為常數

三個方向加速度估計及誤差如圖27-29所示：

圖27　x方向加速度估計及誤差

圖28　y方向加速度估計及誤差

圖29　z方向加速度估計及誤差

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

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

三個方向加速度估計及誤差如圖30-32所示：

圖31　y方向加速度估計及誤差

圖32　z方向加速度估計及誤差

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

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

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

5　結論

本文針對臨近空間機動目標跟蹤問題，使用了基于KF、EKF、UKF、PF的跟蹤濾波算法。仿真結果表明在四種臨近空間非彈道式機動模式下，粒子濾波算法的加速度估計精度均明顯優于其余三種跟蹤濾波算法。

［1］孫鵬，楊建軍.臨近空間平臺防空反導作戰運用仿真研究［J］.彈箭與制導學報，2011，31（2）：14-16.

［2］孫鵬，楊建軍.臨近空間中繼制導反TBM關鍵技術分析［J］.裝備指揮技術學院學報，2011，22（2）：61-65.

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

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

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

［6］付夢印，鄧志紅，張繼偉，編著.Kalman濾波理論及其在導航系統中的應用［M］. 北京：科學出版社，2003，10.

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

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

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