基于HCKS-EM的戰(zhàn)斗機蛇形機動目標跟蹤算法

盧春光， 周中良， 劉宏強， 寇添， 楊遠志

(空軍工程大學航空工程學院, 西安 710038)

蛇形機動是戰(zhàn)斗機飛行員在隱蔽接敵、機動規(guī)避及協(xié)同探測等戰(zhàn)術(shù)動作中經(jīng)常采用的戰(zhàn)術(shù)機動類型。在一定的空戰(zhàn)態(tài)勢下，飛行員進行蛇形機動反映了其戰(zhàn)術(shù)意圖，因此在空戰(zhàn)中如何快速而又準確地識別出目標的蛇形機動模式對于明確其戰(zhàn)術(shù)意圖以及評估當前的戰(zhàn)場態(tài)勢具有十分重要的意義。目前，機動模式的識別主要從機動模式的幾何特征(宏觀)以及運動參量特征(微觀)2個方面著手：文獻[1]通過將戰(zhàn)術(shù)機動進行分割并按時間序列進行編碼，采用隱馬爾可夫模型進行訓練、識別；文獻[2]通過提取目標航向角變化率、高度變化率等運動參量特征，設(shè)計了一個2級識別方案，采用模糊推理和時間自動機方法實現(xiàn)了機動模式的識別。對于蛇形機動目標而言，轉(zhuǎn)彎角速度是其最為關(guān)鍵的運動參量，因此通過載機雷達的量測數(shù)據(jù)精確辨識出轉(zhuǎn)彎角速度對于提高戰(zhàn)斗機蛇形機動模式的識別率具有重要意義。

目前，機動目標轉(zhuǎn)彎角速度辨識的問題已經(jīng)得到了許多學者的關(guān)注。文獻[3-5]將轉(zhuǎn)彎角速度作為目標狀態(tài)變量，通過重新構(gòu)造狀態(tài)方程，分別使用擴展卡爾曼濾波(EKF)/無跡卡爾曼濾波(UKF)、基于5階球徑容積規(guī)則的高階容積卡爾曼濾波(HCKF)等非線性濾波算法估計角速度和其他狀態(tài)變量，其缺點是角速度估計的精度依賴于目標初始狀態(tài)及協(xié)方差、過程噪聲等因素。文獻[6]根據(jù)角速度與目標速度方向角之間的物理關(guān)系，通過采用濾波算法對目標速度方向角進行估計，間接求出了角速度的估計值，該方法的缺點是采用方向角等中間變量間接地估計角速度，將引入新的估計誤差。文獻[7]將角速度當作系統(tǒng)的未知參數(shù)，采用期望最大化(EM)算法將角速度與目標狀態(tài)進行聯(lián)合估計與辨識，該方法的缺點是使用的是全量測數(shù)據(jù)，算法運行時間長，不能滿足空戰(zhàn)中實時識別的要求。

在實際空戰(zhàn)中，目標的角速度是未知且時變的，并且與目標狀態(tài)相互耦合，因此在進行聯(lián)合狀態(tài)估計與角速度辨識時，為了獲得角速度辨識的解析解，需要解除這種耦合關(guān)系。對在傳統(tǒng)意義上的狀態(tài)估計與參數(shù)辨識中，通常認為量測噪聲之間是相互獨立的，但是在實際跟蹤過程中，雷達的采樣頻率過高[8]，使得相鄰量測噪聲之間的相關(guān)性往往不能夠被忽略，因此，在對角速度辨識和目標狀態(tài)估計的過程中必須考慮這種相關(guān)性。針對上述問題，本文基于EM算法框架提出一種帶有色量測噪聲的聯(lián)合估計與辨識算法,并引入滑窗思想，以提高辨識的實時性：E-step利用當前估計的角速度以及經(jīng)過帶有色量測噪聲的高階容積卡爾曼平滑(HCKS)獲得的后驗平滑概率密度和聯(lián)合概率密度，計算完備數(shù)據(jù)似然函數(shù)的條件期望；M-step通過使期望最大化進而更新獲得下一次迭代的角速度估計量。

1 問題描述

蛇形機動可以分解成若干個具有不同轉(zhuǎn)彎角速度的圓弧轉(zhuǎn)彎機動，因此蛇形機動目標轉(zhuǎn)彎角速度的辨識問題轉(zhuǎn)化成圓弧轉(zhuǎn)彎機動的角速度辨識問題。在有色量測噪聲背景下，假設(shè)目標在二維平面中運動，轉(zhuǎn)彎機動模型狀態(tài)方程及量測方程為

(1)

f(xk,Ωk)=

vk=ψk,k-1vk-1+ζk-1

(2)

由式(1)可知，蛇形機動目標的轉(zhuǎn)彎角速度Ωk這一未知參數(shù)耦合在狀態(tài)轉(zhuǎn)移矩陣之中，為了獲取轉(zhuǎn)彎角速度的解析解，需要解除轉(zhuǎn)彎角速度與狀態(tài)轉(zhuǎn)移矩陣之間的非線性耦合關(guān)系，從而將狀態(tài)方程等價轉(zhuǎn)換成如式(3)形式：

xk+1=F(xk)θ+xk+wk

(3)

式中：

θ1=sin(ΩkT)/Ωk

θ2=(1-cos(ΩkT))/Ωk

θ3=1-cos(ΩkT)

θ4=sin(ΩkT)

由于在有色噪聲條件下，傳統(tǒng)的高斯近似濾波器和平滑器估計性能不佳，進而影響EM算法的辨識效果，所以為了實現(xiàn)目標狀態(tài)的精確估計和轉(zhuǎn)彎角速度的準確辨識，需要解除相鄰量測之間相互耦合的關(guān)系，實現(xiàn)有色量測噪聲的白化[10]。本文采取量測差分的方法實現(xiàn)有色量測噪聲的白化，重新構(gòu)造量測方程如下：

(4)

依據(jù)式(1)和式(2)可以得到具有一步狀態(tài)延遲的量測方程：

(5)

2 基于HCKS-EM的聯(lián)合估計與辨識算法

基于極大似然估計準則，本文提出了一種帶有色量測噪聲的聯(lián)合估計與辨識算法對未知角速度進行辨識，具體的算法框架如圖1所示。在第t次迭代時，通過將帶有色噪聲的量測序列經(jīng)過HCKS，獲得目標狀態(tài)的平滑估計量，并采用EM算法進行參數(shù)θ的辨識，辨識后的參數(shù)θ用于第t+1次的狀態(tài)估計，不斷迭代直到滿足設(shè)定的要求為止[12-13]。

圖1 聯(lián)合估計與辨識算法Fig.1 Joint estimation and identification algorithm

(6)

(7)

式中：

(8)

dxk-idxk-i-1

(9)

(10)

其中：l為滑動窗口的長度，其長度大小可依據(jù)實際需求自行設(shè)定。當l=k-1時，I1是一個常數(shù)，與待辨識量θ無關(guān)；當0≤l

pθ(xk-i|xk-i-1)～N(xk-i;F(xk-i-1)θ+xk-i-1,Q)

(11)

ψh(xk-i-1),R)

(12)

2.1 E-step

xk-i-1)TQ-1F(xk-i-1)·

(13)

式中：后驗平滑概率密度和基于滑窗內(nèi)數(shù)據(jù)的xk-i和xk-i-1的聯(lián)合概率密度滿足如下分布：

為了求解I2，首先定義如下2種積分函數(shù)：

(14)

(15)

然后分別將F(xi)和xi按xi的分量進行分解，可得

Γ(xk-i-1)和Δ(xk-i,xk-i-1)可分別表示為

(16)

(17)

(18)

(19)

則

(20)

(21)

2.2 M-step

(22)

當I2取得極大值時滿足

(23)

則可求得參數(shù)θ的迭代表達式為

(24)

2.3 帶有色量測噪聲的HCKS算法的設(shè)計

眾所周知，平滑是在濾波的基礎(chǔ)之上進行的，因此首先設(shè)計一種帶有色量測噪聲的HCKF算法，算法如下：

算法1帶有色量測噪聲的高斯濾波算法框架。

步驟1狀態(tài)預測。

(25)

步驟2狀態(tài)修正。

(26)

(27)

帶有色量測噪聲的HCKS算法是一類前向后向平滑器，其主要包括2個部分：前向一步平滑和后向固定區(qū)間平滑，具體算法如下：

算法2帶有色量測噪聲的高斯平滑算法框架。

步驟1前向一步平滑。

(28)

步驟2后向固定區(qū)間平滑。

(29)

2.4 帶有色量測噪聲的HCKS算法的實現(xiàn)

本文中采用了Jia等[5]提出的5階球徑容積規(guī)則計算式(25)、式(27)和式(28)中的非線性高斯積分：

(30)

ξi=

(31)

式中：ei為n維單位向量，且其第i個元素為1。

(32)

容積點所對應的權(quán)重ωi為

(33)

可將本文算法總結(jié)如下：

算法3基于HCKS-EM的聯(lián)合估計與辨識算法

步驟1E-step。已知第t次迭代后參數(shù)θ的估計值為θt，由式(28)和式(29)可以求得后驗平滑概率密度和基于滑窗內(nèi)數(shù)據(jù)的xk-i和xk-i-1的聯(lián)合概率密度，進而可以求解得到I2的值。

遞歸：t=0，k←k+1，并將參數(shù)θ的初值設(shè)為上次迭代結(jié)束時的辨識值，進入下一輪迭代。

3 仿真分析

本文利用水平方向上的轉(zhuǎn)彎機動非線性動態(tài)模型，仿真出一條蛇形機動軌跡。假設(shè)機動目標在1～80 s以Ω1=-2(°)/s作轉(zhuǎn)彎運動，在k=81 s時轉(zhuǎn)彎角速度突變?yōu)棣?=2(°)/s并持續(xù)到300 s。設(shè)置轉(zhuǎn)彎角速度初始值為Ω0=-1(°)/s，采樣周期T=1 s，q1=0.1 m2/s3，過程噪聲wk的協(xié)方差為Q。

Q=diag(q1M,q1M)

通過載機雷達可以獲得目標與載機之間的相對距離r、方向角φ的信息，則可獲得系統(tǒng)的非線性量測方程為

基于HCKS-EM聯(lián)合估計與辨識算法和擴維法的初始狀態(tài)以及協(xié)方差分別為

擴維法初始狀態(tài)和協(xié)方差的分別設(shè)置為

100 mrad/s]

為了評估分析本文聯(lián)合估計和辨識算法的性能，設(shè)置了如下3組仿真：

仿真1窗口長度l分別設(shè)置為5、10、15、20，最大迭代次數(shù)均為5次，各執(zhí)行100次蒙特卡羅仿真。圖2為不同窗口長度下角速度辨識結(jié)果，圖3為不同窗口長度下角速度辨識均方根誤差(RMSE)。從圖2和圖3可以看出，隨著窗口長度的增大，該算法在第一階段收斂于真實值的時刻就越早，精度越高，并且越穩(wěn)定，但是當角速度發(fā)生突變時，對于突變的角速度反應的就越慢。并且從圖4和表1可以看出，窗口長度越大，該算法估計的目標狀態(tài)整體精度就越高，但是消耗的時間越長，這顯然是時間與精度之間的“博弈”問題，當窗口大于10時，由窗口長度即量測數(shù)據(jù)帶來的精度收益優(yōu)勢不太明顯，相反運算時間的增加確實比較可觀的，也就是說算法的執(zhí)行效率降低了。并且在角速度發(fā)生突變時刻及其附近，窗口長度為10時的估計效果明顯好于其他3個窗口。

圖2 不同窗口長度下轉(zhuǎn)彎角速度辨識結(jié)果Fig.2 Turn rate identification results with different window lengths

圖3 不同窗口長度下轉(zhuǎn)彎角速度辨識的均方根誤差Fig.3 RMSE of turn rate identification with different window lengths

圖4 不同窗口長度下位置和速度估計均方根誤差Fig.4 RMSE of position and velocity estimation with different window lengths

滑動窗口長度時間/s50.0465100.1048150.1222200.1649

仿真2窗口長度設(shè)置為10，最大迭代次數(shù)r分別設(shè)置為3、5、10、15，各執(zhí)行100次蒙特卡羅仿真。從圖5和圖6可以看出，隨著迭代次數(shù)的增加，該算法在初始時刻收斂于真實值的時刻就越早，并且對于角速度突變反應的也比較靈敏。從圖7目標狀態(tài)4個分量的RMSE可以看出，迭代次數(shù)越大，該算法估計的目標狀態(tài)整體精度就越高，但是根據(jù)表2中的運算時間可知花費的時間也會相應地增大。尤其是當最大迭代次數(shù)大于5時，最大迭代次數(shù)增加所帶來的精度收益越發(fā)的不明顯。

圖5 不同最大迭代次數(shù)下轉(zhuǎn)彎角速度辨識結(jié)果Fig.5 Turn rate identification result with different numbers of maximum iterations

圖6 不同迭代次數(shù)下轉(zhuǎn)彎角速度辨識的均方根誤差Fig.6 RMSE of turn rate identification with different numbers of iterations

仿真3將本文聯(lián)合估計與辨識算法窗口長度設(shè)置為10，最大迭代次數(shù)設(shè)置為5，將該算法與帶有色的UKF算法和帶有色的HCKF算法各執(zhí)行100次蒙特卡羅仿真。從圖8和圖9可以看出，本文EM算法收斂于角速度真實值的時刻比帶有色的UKF算法和帶有色的HCKF算法都要早，并且穩(wěn)定、精度高，但是當角速度發(fā)生突變時，對于突變的角速度反應的速度比較慢。從圖10目標狀態(tài)4個分量的RMSE可以看出，除了在角速度突變時刻及其鄰近時刻之外，本文EM算法估計的精度都要明顯高于帶有色的UKF算法和帶有色的HCKF算法估計的精度，帶有色的HCKF算法比帶有色的UKF算法估計的精度要高。

圖7 不同迭代次數(shù)下位置和速度估計均方根誤差Fig.7 RMSE of position and velocity estimation with different numbers of iterations

最大迭代次數(shù)時間/s30.049550.0827100.1599150.2577

圖8 EM算法與擴維法之間的對比Fig.8 Comparison between EM algorithm and augmentation method

圖9 轉(zhuǎn)彎角速度辨識均方根誤差Fig.9 RMSE of turn rate identification

圖10 位置和速度估計均方根誤差Fig.10 RMSE of position and velocity estimation

4 結(jié) 論

針對蛇形機動目標角速度辨識與目標狀態(tài)聯(lián)合估計的問題，基于EM算法框架提出一種帶有色量測噪聲的聯(lián)合估計與辨識算法。主要結(jié)論有：

1) EM算法設(shè)計：在E-Step，通過將完備似然函數(shù)進行分解，從而將其轉(zhuǎn)換成帶參數(shù)θ的解析表達式，提高了辨識與估計的精度，降低了由狀態(tài)擴維法擴維所帶來的復雜性；在M-Step，通過極大化完備似然函數(shù)求得高精度的解析解。

2) 帶有色量測噪聲的HCKS算法的設(shè)計：考慮到在實際空戰(zhàn)中，相鄰噪聲序列之間的相關(guān)性對于目標狀態(tài)估計與轉(zhuǎn)彎角速度辨識的影響，設(shè)計了帶有色量測噪聲的HCKS算法用于消除這種影響，提高了聯(lián)合辨識與估計的精度。