孫 雷，李 冬，易東云

（1.國防科學(xué)技術(shù)大學(xué)理學(xué)院，湖南 長沙410073；2.南京炮兵學(xué)院，江蘇 南京211132）

1 引言

衛(wèi)星主要采用星載紅外傳感器，通過被動(dòng)探測(cè)，獲得彈道導(dǎo)彈的到達(dá)角信息，實(shí)現(xiàn)對(duì)目標(biāo)發(fā)射的監(jiān)視、跟蹤、彈道估計(jì)等功能［1，2］。彈道導(dǎo)彈的飛行過程分為主動(dòng)段和被動(dòng)段，主動(dòng)段包括垂直上升段、轉(zhuǎn)彎段、過渡段和穩(wěn)定段。紅外傳感器只能獲取導(dǎo)彈目標(biāo)的到達(dá)角和輻射強(qiáng)度信息，觀測(cè)情形如圖1所示。

衛(wèi)星在只有角度測(cè)量信息下對(duì)于彈道導(dǎo)彈主動(dòng)段的跟蹤屬于無源跟蹤，導(dǎo)彈目標(biāo)的狀態(tài)方程和測(cè)量方程都是非線性的，并且衛(wèi)星與目標(biāo)距離較遠(yuǎn)，可觀測(cè)性較弱。這些都對(duì)主動(dòng)段彈道跟蹤問題帶來了較大的困難，因此尋找一種快速收斂且具有較高精度的跟蹤算法十分關(guān)鍵。

Figure 1 Obervation of the trajectory form the satellite圖1 衛(wèi)星對(duì)導(dǎo)彈主動(dòng)段的觀測(cè)

擴(kuò)展卡爾曼濾波 EKF（Extended Kalman Fil－ter）是常用的非線性濾波算法，但在實(shí)際應(yīng)用中它存在明顯的缺陷。無跡卡爾曼濾波 UKF（Unscented Kalman Filter）是目前受到廣泛使用的非線性濾波算法，通過選取一些樣點(diǎn)更精確地給出狀態(tài)變量經(jīng)過非線性變化后的均值和方差，減小了線性化帶來的誤差，避免了EKF算法在求解Jacobi矩陣帶來的Fisher矩陣病態(tài)問題。但是，標(biāo)準(zhǔn)的UKF算法存在收斂速度慢、精度不高的問題。一種易于實(shí)現(xiàn)的提高非線性近似精度的方式是使用迭代測(cè)量更新的方式，這在EKF上已經(jīng)得到了成功的應(yīng)用［3］。與迭代EKF的思想類似，可以使用狀態(tài)的估計(jì)值代替預(yù)測(cè)值來重新進(jìn)行無跡變換UT（Unscented Transformation），得到新的濾波參數(shù)，并再次利用觀測(cè)值來改善估計(jì)狀態(tài)，從而提高狀態(tài)估計(jì)的精度，已有文獻(xiàn)提出迭代UKF的思想和易于實(shí)現(xiàn)的迭代測(cè)量更新的UKF算法［4，5］。

本文通過分析彈道導(dǎo)彈在主動(dòng)段的受力情況，引入一種基于重力轉(zhuǎn)彎模型的狀態(tài)方程。根據(jù)基于迭代測(cè)量更新的UKF思想，研究其在衛(wèi)星觀測(cè)信息下的迭代更新過程，從而建立基于迭代UKF的主動(dòng)段跟蹤算法。仿真結(jié)果表明，迭代UKF算法與標(biāo)準(zhǔn)UKF算法相比，收斂速度更快、跟蹤精度更好，是一種適應(yīng)性更強(qiáng)的導(dǎo)彈目標(biāo)主動(dòng)段跟蹤算法。

2 主動(dòng)段彈道估計(jì)模型

主動(dòng)段為彈道導(dǎo)彈飛行的第一個(gè)階段，其強(qiáng)機(jī)動(dòng)特性和變質(zhì)量過程使得主動(dòng)段彈道動(dòng)力學(xué)特性十分復(fù)雜。主動(dòng)段時(shí)間并不長，一般在幾十秒到幾百秒范圍內(nèi)，但是主動(dòng)段的跟蹤對(duì)于被動(dòng)段的狀態(tài)估計(jì)以及導(dǎo)彈的射程和落點(diǎn)卻起著決定性作用。

對(duì)主動(dòng)段彈道準(zhǔn)確地建模是實(shí)現(xiàn)彈道導(dǎo)彈跟蹤的前提，主動(dòng)段的跟蹤模型分為導(dǎo)彈的主動(dòng)段狀態(tài)模型和預(yù)警衛(wèi)星觀測(cè)模型。

2.1 主動(dòng)段彈道模型

導(dǎo)彈在主動(dòng)段的受力情況十分復(fù)雜，主要包括推力、重力、氣動(dòng)阻力等。對(duì)用于主動(dòng)段彈道估計(jì)的目標(biāo)狀態(tài)方程已有相關(guān)的研究［6，7］，其中一種為基于重力轉(zhuǎn)彎假設(shè)的狀態(tài)方程。重力轉(zhuǎn)彎假設(shè)是指假設(shè)導(dǎo)彈在主動(dòng)段飛行過程中，攻角很小，幾乎為零，推力與阻力均沿速度方向，重力是使彈道彎曲的唯一因素，因此稱為重力轉(zhuǎn)彎。本文使用的目標(biāo)狀態(tài)方程就是重力轉(zhuǎn)彎模型。

設(shè)導(dǎo)彈在地心慣性坐標(biāo)系ECI（Earth Cen－tered Intertial）下的位置速度矢量分別為r和v，則r＝ ［x y z］T，v＝˙r＝ ［˙x ˙y ˙z］T。導(dǎo)彈受到的力可以分為沿重力方向的力和沿速度方向的力，在重力轉(zhuǎn)彎假設(shè)條件下，導(dǎo)彈彈道參數(shù)的主動(dòng)段動(dòng)力學(xué)方程可記為：

其中，μ為地球引力常數(shù)，α為導(dǎo)彈推力加速度大小。考慮實(shí)際中大氣透過率的影響，紅外傳感器在導(dǎo)彈飛出稠密大氣層后才能發(fā)現(xiàn)目標(biāo)，因此這里忽略氣動(dòng)加速度的影響。

下面分析推力加速度大小α隨時(shí)間變化的情況。設(shè)導(dǎo)彈推力大小為常值T，則推力加速度可以寫成：

其中，m（t）為t時(shí)刻導(dǎo)彈的瞬時(shí)質(zhì)量，σ為發(fā)動(dòng)機(jī)秒耗量，t0為初始時(shí)刻。

上式對(duì)時(shí)間求導(dǎo)后可以得到：

結(jié)合公式（1）和公式（4），就可以得到基于重力轉(zhuǎn)彎模型的主動(dòng)段狀態(tài)方程：

目 標(biāo) 的 狀 態(tài) 記 為：X ＝ ［x，y，z，˙x，˙y，˙z，α，∑］T，則主動(dòng)段動(dòng)力學(xué)模型（5）可以表示為：

其中，X為目標(biāo)的狀態(tài)，w為零均值不相關(guān)的狀態(tài)噪聲。根據(jù)給定的初始條件，采用數(shù)值積分的方法對(duì)公式（6）求解就可以得到導(dǎo)彈在主動(dòng)段某一時(shí)刻的位置和速度。

2.2 衛(wèi)星角度觀測(cè)方程

在衛(wèi)星測(cè)量坐標(biāo)系 UEN（Up East North）中，角測(cè)量方程定義如下：

其中，ψ和γ分別是目標(biāo)在UEN坐標(biāo)系下的方位角和俯仰角，如圖2所示。

Figure 2 Satellite observation model圖2 衛(wèi)星觀測(cè)模型

觀測(cè)模型基于目標(biāo)的角度測(cè)量信息，設(shè)在Tk時(shí)刻目標(biāo)在ECI坐標(biāo)系下的狀態(tài)變量為Xk，其中位置變量為rk，衛(wèi)星的位置變量記為sk，根據(jù)ECI坐標(biāo)系與UEN坐標(biāo)系的關(guān)系，Tk時(shí)刻衛(wèi)星的視線測(cè)量矢量zk可以表示為：

3 迭代UKF濾波算法

UKF算法是基于UT變換和Kalman濾波技術(shù)的一種濾波算法［8，9］，其通過特別選取一些樣點(diǎn)，更好地近似隨機(jī)變量經(jīng)過非線性變化后的均值和方差。將系統(tǒng)圍繞狀態(tài)估計(jì)重新進(jìn)行UT變換計(jì)算濾波參數(shù)，然后再一次利用觀測(cè)值改善對(duì)狀態(tài)的估計(jì)，這就是迭代 UKF——IUKF（Iterated UKF）的思想。

3.1 UT變換

UT變換是一種計(jì)算隨機(jī)變量在經(jīng)歷非線性變換后的統(tǒng)計(jì)特性的方法?？紤]一個(gè)隨機(jī)變量x（維數(shù)為L），通過一個(gè)非線性函數(shù)y＝f（x）得到隨機(jī)變量y，假設(shè)x具有均值ˉx與協(xié)方差Px，為計(jì)算y的統(tǒng)計(jì)特性，首先根據(jù)下式：得到2L＋1個(gè)Sigma采樣點(diǎn)，其中λ＝φ2（L＋κ）－L為尺度參數(shù)，常數(shù)φ決定這些Sigma采樣點(diǎn)在xˉ附近的擴(kuò)展范圍；κ是另一個(gè)尺度因子，通常取為0或3－L；（是矩陣方根的第i列，矩陣方根可以通過Cholesky分解得到。將Sigma采樣點(diǎn)代入非線性函數(shù)：

通過下式計(jì)算y的均值和協(xié)方差矩陣：

其中，權(quán)Wi由下式給出：其中，β為狀態(tài)分布參數(shù)，對(duì)于高斯分布，β＝2是最優(yōu)的。如果狀態(tài)變量是單變量，則最佳的選擇是β＝0。

3.2 迭代UKF（IUKF）濾波算法

將UT變換應(yīng)用于Kalman濾波，就可以得到UKF濾波算法。將UKF濾波算法狀態(tài)估計(jì)結(jié)果Xk重新進(jìn)行UT變換，計(jì)算濾波的參數(shù)，再次利用觀測(cè)值改善對(duì)狀態(tài)的估計(jì)結(jié)果，這樣就得到了迭代UKF算法。

對(duì)于非線性的狀態(tài)方程和測(cè)量方程：

對(duì)于本文研究的導(dǎo)彈主動(dòng)段的跟蹤問題，上式中的狀態(tài)方程就是2.1節(jié)中的方程（6），測(cè)量方程就是2.2節(jié)中的方程（8）。對(duì)于狀態(tài)變量Xk的估計(jì)，首先由k時(shí)刻的UKF算法計(jì)算得到狀態(tài)估計(jì)和協(xié)方差分別為＾Xk和P－k，由＾Xk和Pk生成的新采樣點(diǎn)可重新求得濾波器的增益Kk，再利用觀測(cè)值對(duì)狀態(tài)估計(jì)進(jìn)行改善，具體步驟如下：

（1）初始化。

（2）計(jì)算采樣點(diǎn)。

（3）時(shí)間更新。

（4）測(cè)量更新。

其中，Rv為系統(tǒng)狀態(tài)噪聲方差，Rn為系統(tǒng)測(cè)量噪聲方差，xα＝ ［xTvTnT］。

（5）重新計(jì)算采樣點(diǎn)，更新估計(jì)結(jié)果。

Step 2 定義：

Step 4 j＝j(luò)＋1，g＝η·g，重新計(jì)算采樣點(diǎn)：

將公式（18）代入公式（16）和公式（17）式重新計(jì)算濾波參數(shù)，得到新的估計(jì)結(jié)果：

轉(zhuǎn)向Step 2。

最后將公式（19）作為k時(shí)刻的濾波狀態(tài)和協(xié)方差估計(jì)結(jié)果。

在上述IUKF算法流程中，為滿足迭代要求，一般取0＜η≤1。

與UKF算法相比，IUKF算法增加了迭代測(cè)量更新運(yùn)算，每次迭代增加的運(yùn)算量主要是樣點(diǎn)的非線性變換?？傮w來說，IUKF算法比UKF算法增加了一些計(jì)算量，但其計(jì)算復(fù)雜度仍然為同一量級(jí)。由于不同的非線性變換形式帶來的計(jì)算量有較大不同，這將影響IUKF算法與UKF算法之間計(jì)算量的差別。

4 仿真實(shí)驗(yàn)

在仿真計(jì)算機(jī)CPU為P4雙核、2.7GHz、內(nèi)存為1GB的情況下，使用Matlab工具進(jìn)行仿真計(jì)算。仿真實(shí)驗(yàn)中，測(cè)量數(shù)據(jù)間隔Ts取0.5s，視線測(cè)量誤差σLOS分別取30μrad和60μrad。目標(biāo)為某2 000km彈道導(dǎo)彈，主動(dòng)段飛行時(shí)間為100s，首次觀測(cè)時(shí)間為第20s。選取兩顆地球同步軌道衛(wèi)星觀測(cè)。Monte－Carlo次數(shù)為50次。圖3是測(cè)量數(shù)據(jù)間隔為0.5s、σLOS為30μrad情況下UKF算法和IUKF算法對(duì)目標(biāo)速度和位置的跟蹤誤差仿真結(jié)果。圖4是測(cè)量數(shù)據(jù)間隔為0.5s、σLOS為60 μrad情況下UKF算法和IUKF算法仿真結(jié)果。

Figure 3 Simulation results of UKF and IUKF whenδLOS ＝30μrad圖3 σLOS取30μrad時(shí)UKF與IUKF仿真結(jié)果

Figure 4 Simulation results of UKF and IUKF whenδLOS ＝60μrad圖4 σLOS取60μrad時(shí)UKF與IUKF仿真結(jié)果

從圖3和圖4的仿真結(jié)果可以看出，IUKF算法對(duì)新的觀測(cè)值反應(yīng)更迅速，在相同的仿真條件下，IUKF算法比UKF算法收斂速度更快、收斂精度更高。

5 結(jié)束語

本文圍繞衛(wèi)星在僅有角度測(cè)量信息下彈道導(dǎo)彈的主動(dòng)段跟蹤問題，對(duì)常用的非線性濾波UKF算法進(jìn)行改進(jìn)，提出了IUKF算法。采用重力轉(zhuǎn)彎模型建立主動(dòng)段狀態(tài)方程，結(jié)合UKF算法與IUKF算法進(jìn)行仿真實(shí)驗(yàn)。實(shí)驗(yàn)表明，本文提出的方法收斂速度更快、收斂精度更高，適用于基于衛(wèi)星角度測(cè)量信息下導(dǎo)彈目標(biāo)主動(dòng)段跟蹤的應(yīng)用。

［1］ Danis J.Space－based tactical missile launch parameter estimation［J］.IEEE Transactions on AES，1993，29（2）：412－424.

［2］ Li Yi－cong，Kirubarajan T，Bar－Shalom Y，et al.Trajectory and launch point estimation for ballistic missiles from boost phase LOS measurements［C］∥Proc of IEEE Aerospace Conference，1999：425－442.

［3］ Bellaire R L，Kamen E W，Zabin S M.New nonlinear iterated filter with applications to target tracking［C］∥Proc of the SPIE’95，1995：240－251.

［4］ Van Zandt J R.Boost phase tracking with an unscented filter［C］∥Proc of the SPIE’02，2002：263－274.

［5］ Zhan Rong－h(huán)ui，Wang Ling，Wan Jian－wei.Research on robust algorithm for single observer passive target tracking［J］.Signal Processing，2007，23（3）：464－468.（in Chinese）

［6］ Li X R，Jilkov V P.Survey of maneuvering target tracking part II：Motion models of ballistic and space target［J］.IEEE Transactions on Aerospace and Electronic Systems，2010，46（1）：96－119.

［7］ Li Dun.Target location and prediction with space early warning system［D］.Changsha：National University of Defense Technology，2001.（in Chinese）

［8］ Julier S J，Uhlmann J K，Durrant－Whyte H F.A new method for the nonlinear transformation of means and covariances in filters and estimatimators［J］.IEEE Transactions on Automatic Control，2000，45（3）：477－482.

［9］ Wan E A，Merwe R，Nelson A T.Dual estimation and the unscented transformation［M］∥Advances in Neural Information Processing Systems 12，Cambridge：MIT Press，2000：666－672.

附中文參考文獻(xiàn)：

［5］ 占榮輝，王玲，萬建偉.穩(wěn)健的單站無源目標(biāo)跟蹤算法研究［J］.信號(hào)處理，2007，23（3）：464－468.

［7］ 李盾.空間預(yù)警系統(tǒng)對(duì)目標(biāo)的定位與預(yù)報(bào)［D］.長沙：國防科學(xué)技術(shù)大學(xué)，2001.