低通濾波與卡爾曼濾波相結合的制導律識別

杜潤樂，劉佳琪，李志峰，牛振紅，張 力

(試驗物理與計算數學國家級重點實驗室(北京航天長征飛行器研究所), 北京 100076)

低通濾波與卡爾曼濾波相結合的制導律識別

杜潤樂，劉佳琪，李志峰，牛振紅，張 力

(試驗物理與計算數學國家級重點實驗室(北京航天長征飛行器研究所), 北京 100076)

為利用觀測數據對執行追蹤任務的非合作目標飛行器進行制導律識別，采用一種通用的隱式制導函數對制導行為進行統一描述，并且給出了一種采用時變制導系數的模型來逼進隱式制導函數.在此模型基礎上，利用低通濾波與卡爾曼濾波相結合的頻域分離自適應卡爾曼濾波方法對觀測數據進行實時處理分析，最終得到目標飛行器的制導律識別結果.仿真結果表明，在目標飛行器采用比例導引律或Bang-Bang變結構制導律的情況下，識別算法能夠從觀測數據出發，對目標飛行器的制導律完成有效和快速識別.

制導律識別；隱式制導函數；時變制導系數；低通濾波；卡爾曼濾波

在一些特定的背景下，為了提高我方飛行器的生存能力，需要設計回避目標飛行器的制導律.這種回避制導律需要知道關于目標飛行器的制導律、參數和狀態的完整信息.因此，需要一種針對目標飛行器制導律、參數和狀態的辨識估計方法，在我方飛行器與目標飛行器相遇過程中進行辨識.從目標飛行器外部測量數據出發，對其制導律進行辨識，對于飛行器制導規律的優化設計、效能評估、有效觀測、軌跡預測等方面具有重要的理論和應用價值.

目標飛行器可以從已有的制導律中選取一種制導律，例如，比例導引律、滑模制導律、最優制導律等.制導律識別問題則是從目標飛行器外部觀測數據出發，對未知的目標飛行器制導律進行辨識和參數估計.

由于制導律的多樣性，其識別問題的難點在于構造統一的通用數學模型，對制導行為進行描述.Grove等[1]在認定目標飛行器采用比例導引律的前提下，以比例導引系數作為未知量建立觀測數據方程，通過求解制導控制量方程得到制導律的參數.Lin等[2]指出，對于未知的制導律可以假設為比例導引律，通過調整比例導引系數來擬合任意未知的制導律.Horton[3]則使用卡爾曼濾波器和神經網絡相結合的方法對飛行器自身的空氣動力學參數的進行實時識別.Li等[4]歸納總結了針對運動目標軌跡的跟蹤方法，討論了利用動力學模型對機動目標的建模和軌跡跟蹤方法.在其他類似的研究中，Raviadra等[5]使用卡爾曼濾波器對目標飛行器的軌跡進行識別.Kunwar等[6]研究了移動目標攔截過程中的軌跡規劃問題.

綜合上述研究結果可以看出，對目標飛行器制導行為的識別以及對識別結果的利用目前需要進一步深入研究.本文利用了通用的隱式制導函數方法對制導行為模式進行統一描述，并且給出了采用時變制導系數的近似模型來逼進隱式制導函數.在該數學模型的基礎上，利用低通濾波與卡爾曼濾波相結合的頻域分離自適應卡爾曼濾波方法(frequency divided adaptive Kalman filter, FD-AKF)對觀測數據進行實時處理分析，并得到目標飛行器的制導律識別結果.

由于采用了低通濾波與卡爾曼濾波相結合的方法，測量信息中的噪聲得到了更好地抑制.為了在動態噪聲或測量噪聲統計特性未知情況下設計卡爾曼濾波器，Yang等[7]提出了需要預知殘差向量和預測向量的理論協方差矩陣自適應參數調整方法.Huang等[8]提出的實時切換預定模型以適應噪聲參數變化的方法.Kim等[9]采用自適應衰減濾波器使用遺忘因子對統計信息進行在線分析增強濾波方法的穩定性.Karavasilis等[10]通過動態調整狀態切換矩陣的方法，完成濾波器對輸入的主動適應和調整.Chai等[11]整合了神經網絡的學習適應能力和濾波器的估計識別能力.Kim等[12]利用自適應卡爾曼濾波器校正小波神經網絡的學習速率，以保證神經網絡的穩定收斂.Wang等[13]將Sage-Husa自適應濾波器與卡爾曼濾波器整合，利用自適應濾波器設定濾波器的噪聲水平參數.Yang[14]則指出自適應卡爾曼濾波器的核心是對狀態噪聲和測量噪聲方差的實時估計.

已有研究表明，對預測噪聲和測量噪聲的方差水平識別需要同時進行，才能夠有效地保證算法的收斂性[14]，但是這種做法會帶來大量的計算，而基于窗口方法和遺忘因子設計的濾波器，其中的大量參數需要根據實際情況和經驗值才能確定，不利于算法和方法的推廣.本文提出的低通濾波與卡爾曼濾波相結合的FD-AKF方法利用了系統中有效信號與測量噪聲在頻域上的不同分布特性，使用低通濾波器快速確定測量噪聲方差水平，同時保持了有效信號不受影響.

1 制導律識別問題

三維空間中的制導問題如圖1所示.其中：慣性參考坐標系為末制導初始時刻的視線坐標系相對于慣性空間固化不動；M為我方飛行器，其在慣性參照系中的坐標為x,y,z；T為非合作目標飛行器，其在慣性參照系里面的坐標為xt,yt,zt；符號v、vt分別為我方飛行器和目標飛行器的速度向量.目標飛行器根據兩者的相對運動關系和自身的制導律，調整自身的加速度，從而達到命中我方飛行器的目的.

根據外在的測量參數到機動加速度之間的映射關系，未知的制導律就可以簡化并抽象成一個隱式的待識別的制導函數[3]：

(1)

圖1 制導問題

分析現有制導律設計結果可得出，制導過程中目標飛行器的機動加速度可以近似為與接近速度和視線轉率的乘積成比例關系，所以識別過程中引入參數Ny和Nz.在目前的研究階段，采取簡化的數學模型將函數關系假設為

(2)

制導律的識別具體化為未知的等效參數Ny和Nz的識別，它們甚至可以是變化的.

2 制導律識別算法

2.1 基于擴展卡爾曼濾波的識別算法

2.1.1 狀態方程的建立

選取狀態向量為

(3)

則不難寫出系統的狀態方程：

式中：uxt、uyt、uzt分別為目標飛行器制導律控制下的加速度指令；τ為目標飛行器過載回路的時間常數；ωNy、ωNz分別為狀態遷移引入的動態噪聲.該方程可以寫作狀態方程：

2.1.2 觀測方程

濾波器的觀測向量由飛行器的位置、速度、目標飛行器的位置、速度、機動加速度構成.上述觀測結果向量，實際工程上可以通過組網定位等技術對目標飛行器進行測量，并且利用卡爾曼濾波的方法對其機動加速度進行實時處理.本文中將忽略具體的實現過程，假設測量結果能夠直接獲得，作為識別算法的原始測量內容參與識別.制導過程中需要的相對運動參數，如相對位置和視線轉率等參數可以借助原始測量值間接獲得，對制導行為的識別轉化成對復雜等效制導參數的識別.

式中：εx、εy、εz分別為位置向量的觀測噪聲；εdx、εdy、εdz分別為速度向量的觀測噪聲；εxt、εyt、εzt分別為目標飛行器位置向量的觀測噪聲；εdxt、εdyt、εdzt分別為目標飛行器速度向量的觀測噪聲；εax、εay、εaz分別為目標飛行器機動加速度的觀測噪聲.

2.1.3 擴展卡爾曼濾波算法(EKF)

針對系統(2)，推廣卡爾曼濾波算法(EKF)[5,7,9]為

在系統(1)中

其中:

式中，rx=x-xt，ry=y-yt，rz=z-zt，vrx=vx-vxt，vry=vy-vyt，vrz=vz-vzt.[15]

觀測方程(3)以觀測矩陣形式描述，則

2.2 過程噪聲的在線識別

由于制導律辨識模型具有不確定性，通過對狀態噪聲和觀測噪聲協方差矩陣的實時處理來補償模型的不確定性.因為目標飛行器有主動的機動行為，特別的，當目標飛行器制導行為中的等效參數變化復雜時，系統狀態的預測上會產生較大的誤差.

已有研究中Sage-Husa算法[14]和其他改進的算法[16-17]中，都是以此為出發點，對狀態噪聲和測量噪聲進行方差水平的識別.根據觀測值與處理值之間的差值，就可以獲得當前一次迭代的狀態噪聲和測量噪聲估計值，然后使用窗口法則或者漸消法則，完成方差水平的計算.

2.2.1 低通濾波的理論依據

根據采樣定理，卡爾曼濾波器的迭代計算頻率一定要大于有效信號頻率的2倍.由于狀態噪聲產生于每一步預報計算，觀測噪聲產生于每一步觀測測量，他們的特征頻率必然接近于卡爾曼濾波的計算頻率和觀測頻率.系統中的狀態噪聲和測量噪聲的能量密度集中在頻譜的高端部分.由于有效信號和噪聲的能量分別分布在低頻段和高頻段，可以通過低通數字濾波的方法降低噪聲的影響.

經過快速傅里葉變換，由圖2單邊帶頻域能量分布中可以看出，有效的系統信號與噪聲信號在頻域出現了明顯的分離.低頻部分代表有效數據的變化，而高頻部分則來自于狀態噪聲或者觀測噪聲.通過低通濾波，在預測輸出的數據序列中提出有效數據，同時將預測輸出與低通濾波以后的差值作為噪聲的估計，進行后續的計算.

圖2 模型預測向量的頻域分析結果

同理，Zk和Vk也可以通過低通濾波的方法進行分離.

2.2.2 噪聲向量的低通濾波分離法

濾波過程中，由于Φk,k-1Xk-1的頻率低于低通濾波截止頻率，所以有LPF(Φk,k-1Xk-1)≈Φk,k-1Xk-1.與此同時，Wk的頻率高于低通濾波截止頻率，所以有LPF(Wk)≈0.那么，對迭代過程中得到的預測數據進行低通濾波處理，將會得到：

LPF(Xk)=LPF(Φk,k-1Xk-1+Wk)≈ LPF(Φk,k-1Xk-1)≈Φk,k-1Xk-1.

于是Wk≈Xk-LPF(Xk).

同理，Vk≈Zk-LPF(Zk)，低通濾波器的迭代計算公式LPF：yt=yt-1*(1-α)+xt*α，其中α為低通濾波參數，如果對應截止頻率為f的時候，α的計算公式為α=f*2πdt，dt為采樣周期.

2.2.3 方差水平的估計

Wk和Vk的近似序列獲得以后，可以采用窗口方法或者漸消濾波方法計算方差水平.本文提出了再次低通濾波的方法對Wk和Vk的方差水平進行識別.Wk和Vk雖然本身是高頻的，但是其方差水平的變化頻率卻遠低于其自身的頻率，所以文中采用了低通濾波的方法對其方差水平進行估計.

2.2.4 Q和R的條件數保護

在低通濾波方法中，Q、R的正定性能夠得到保證.但是隨著濾波器的迭代計算，系統的狀態噪聲會逐漸縮減.導致Q矩陣雖然是主對角陣，但是個別分量的數值特別的小，以至于矩陣的條件數非常大，影響后續求逆陣過程中的數值計算精度.本文借鑒已有的保護策略，為Q、R設定了條件數保護措施，即

Q=Q+0.5*I，

R=R+0.5*I.

保護措施的特點是疊加一個基礎的保護矩陣，而不是使用條件判斷.當Q的值大于保護矩陣數值的時候，保護矩陣會被真實數據所掩蓋；當Q的值小于保護矩陣數值的時候，保護矩陣保證Q的分量不小于設定的保護范圍.

3 結果和分析

仿真的過程中假設：1)忽略任何外在場的影響，目標飛行器具有理想的軌控發動機，能夠產生連續推力.2)我方飛行器和目標飛行器均具有理想的姿態控制系統，保持飛行器的姿態主軸與速度向量相重合.3)執行識別任務的飛行器被目標飛行器追蹤過程中，能夠獲得目標飛行器動力學測量信息.

本文將以Sage-Husa估計器為基礎的自適應Kalman濾波識別算法，與文中提出的FD-AKF進行對比以檢驗該方法的有效性.Sage-Husa估計器算法中窗口深度設定為5，使用RAE方式對預測噪聲和測量噪聲的方差進行估計.FD-AKF只有一個低通濾波參數需要設定，由于待識別目標的運動特征頻率不會超過10 Hz，所以算法的低通濾波截止頻率設定為小于10 Hz的一個經驗值，2.5 Hz.

仿真過程中在多種不同的制導律配置情況下，對制導律識別的效果進行驗證.

情況1 當目標飛行器使用固定的制導律參數N=4.83的時候，等效制導系數的識別結果如圖3、4所示.

圖3、4中，上部的藍色曲線為俯仰通道制導參數Ny的理論值，紅色曲線為識別結果.下部的紫色曲線為偏航通道制導參數Nz的理論值，青色曲線為識別結果.可以看出，對于固定的制導系數，兩種算法都能夠完成很好的識別，結果的誤差水平都在1%以內.但是，在仿真后期，隨著目標飛行器機動加速度逐漸衰減，Sage-Husa算法出現了發散.

圖3 FD-AKF對固定制導參數的識別

圖4 Sage-Husa算法對固定制導參數的識別

情況2 目標飛行器等效制導行為參數隨著時間推移以正弦規律變化.實際中可能不會出現，但是這種虛擬的制導律可以驗證識別算法對處于變化中的制導規律的識別能力.

圖5、6中可以看出，當制導參數發生動態變化的時候，FD-AKF算法和Sage-Husa算法都能夠保持跟蹤能力.但是，在仿真后期，由于目標飛行器機動加速度逐漸衰減，Sage-Husa算法出現了發散.

圖5 FD-AKF對時變制導參數的識別結果

情況3 目標飛行器等效制導行為參數發生突變.

為進一步說明濾波器對于變化的制導律的識別能力，設計了制導參數會發生突變的制導律，用以檢驗識別算法跟蹤制導律變化的能力.

圖6 Sage-Husa對時變制導參數的識別結果

圖7、8中可以看出，當制導參數發生動態切換的時候，FD-AKF算法和Saga-Husa算法都能夠保持跟蹤能力.但是，Sage-Husa在制導參數切換的時刻，識別結果出現了明顯的跳變.而且在仿真后期，隨著目標飛行器的機動加速度逐漸衰減，Sage-Husa算法出現了發散，而FD-AKF則能夠保持穩定.

圖7 FD-AKF對跳變制導參數的識別結果

圖8 Sage-Husa對跳變制導參數的識別結果

從仿真結果中可以看出，當目標飛行器的制導規律不發生變化的時候，FD-AKF和Sage-Husa算法都能夠很好的完成目標行為規律的識別，制導律參數的收斂速度很快，而且最終的相對誤差穩定在1%以內.當目標飛行器的制導行為規律發生變化的時候，文中提出的FD-AKF在收斂速度和跟隨變化的性能上明顯優于Sage-Husa算法，原因在于FD-AKF算法可以直接從預測數據和測量數據出發，分離出高頻預測噪聲和測量噪聲的同時，保留了低頻率的有效數據變化規律.所以，目標的制導行為規律發生變化的時候，FD-AKF算法能夠及時響應，正確跟蹤其變化過程.

無論是FD-AKF還是Sage-Husa算法，當目標飛行器的過載低于1 g以后，制導行為的識別效果明顯變差.

情況4 目標飛行器采取BANG-BANG變結構制導律.

為進一步檢驗FD-AKF對目標飛行器制導行為的辨識效果，假設目標飛行器采用Bang-Bang變結構制導律.由于飛行器的姿態控制系統令飛行器體坐標系跟蹤視線坐標系，因此，設計制導律為：

式中：ay1為沿彈體坐標系oy1軸方向上的加速度；az1為沿彈體坐標系oz1軸方向上的加速度；M=const.>0為理想的軌控發動機推力作用在導彈上產生的加速度.Bang-Bang變結構制導律由于采用了非線性的控制策略，在切換面附件遷移的過程中，等效時變制導參數會發生高頻率的跳變.由于這種跳變的存在，FD-AKF對Bang-Bang制導律等效制導參數識別過程中，低通濾波的截止頻率調高到50 Hz.

俯仰、偏航通道辨識的等效制導參數分別如圖9、10所示.

圖9 FD-AKF對Bang-Bang變結構制導律識別

圖9、10中，上部的藍色和紫色曲線為等效制導參數Ny、Nz的理論值，理論值是根據控制加速度和接近速度視線轉率換算得來的.紅色和青色曲線為對應的識別結果.可以看出，FD-AKF算法對等效時變制導參數的識別明顯優于Sage-Husa算法的識別結果.

圖10 Sage-Husa對于Bang-Bang變結構制導律的識別

圖11、12分別是圖9、10俯仰通道的局部放大.可以看出，對于藍色曲線表示出的等效制導系數理論值在0～12之間的變化，FD-AKF識別結果在0～12之間跟隨，而Sage-Husa的識別結果的最大、最小值分別為7和2.偏航通道的局部放大圖可以得到類似的結果.

圖11 FD-AKF對于Bang-Bang變結構制導律的識別結果(局部放大)

圖12 Sage-Husa對于Bang-Bang變結構制導律的識別結果(局部放大)

4 結 論

1)通用制導模型中引入時變的等效制導系數解決了數學模型不統一的問題.在此基礎上，低通濾波和卡爾曼濾波相結合的頻域分離自適應卡爾曼濾波算法(FD-AKF)可以完成對制導律的在線識別.仿真實驗結果表明，在飛行器機動的情況下，本文提出的制導律識別方法能夠從比例導引和Bang-Bang變結構制導律的觀測數據出發，對目標飛行器的制導律完成有效的快速識別.

2)在本文工作的基礎上，可以進一步通過實驗驗證制導律FD-AKF識別算法對于更復雜制導律的識別效果.影響制導律識別效果的各種因素也需要進一步確認.

[1] GROVES G W, BLAIR W D, HOFFMAN S. Some concepts for trajectory prediction for ship self defense[J]. Proceedings of the American Control Conference, 1995, 6: 4121-4126. DOI: 10.1109/ACC.1995.532707.

[2] LIN C S, RAETH P G. Prediction of missile trajectory[C]//Proceedings of the IEEE International Conference on Systems, Man and Cybernetics Intelligent Systems for the 21st Century. [S.l.]: IEEE, 1995: 2558-2563. DOI: 10.1109/ICSMC.1995.538167.

[3] HORTON M P. Real-time identification of missile aerodynamics using a linearised kalman filter aided by an artificial neural network[J]. IEEE Transactions on Control Theory and Applications, 1997, 144(4):299-308. DOI: 10.1049/ip-cta:19971125.

[4] LI X R, JILKOV V P. Survey of maneuvering target tracking. Part II: motion models of ballistic and space targets[J]. IEEE Transactions on Aerospace & Electronic Systems, 2010, 46(1): 96-119. DOI: 10.1109/TAES.2010.5417150.

[5] RAVINDRA V C, BAR-SHALOM Y, WILLETT P. Projectile identification and impact point prediction[J]. IEEE Transactions on Aerospace & Electronic Systems, 2010, 46(4):2004-2021. DOI: 10.1109/TAES.2010.5595610.

[6] KUNWAR F, SHERIDAN P K, BENHABIB B. Predictive guidance-based navigation for mobile robots: a novel strategy for target interception on realistic terrains[J]. Journal of Intelligent & Robotic Systems, 2010, 59(3/4): 367-398.DOI:10.1007/s10846-010-9401-3.

[7] YANG Yuanxi, GAO Weiguang. An optimal adaptive kalman filter[J]. Journal of Geodesy, 2006, 80(4):177-183. DOI: 10.1007/s00190-006-0041-0.

[8] HUANG Guanwen, ZHANG Qin. Realtime estimation of satellite clock offset using adaptively robust Kalman Filter with classified adaptive factors[J]. GPS Solutions, 2012, 16(4): 531-539. DOI 10.10007/s10291-012-0254-z.

[9] KIM K H, JEE G I, PARK C G, et al. The stability analysis of the adaptive fading extended Kalman filter using the innovation covariance[J]. International Journal of Control Automation & Systems, 2009, 7(1): 49-56. DOI: 10.1007/s12555-009-0107-x.

[10]KARAVASILIS V, NIKOU C, LIKAS A. Visual tracking by adaptive kalman filtering and mean shift[M]. Berlin, Heidelberg: Springer-Verlag,2010:153-162. DOI: 10.1007/978-3-642-12842-4_19.

[11]CHAI Lin, YUAN Jianping, FANG Qun, et al. Neural network aided adaptive kalman filter for multi-sensors integrated navigation[M]. Berlin, Heidelberg: Springer-Verlag, 2004:381-386. DOI: 10.1007/978-3-540-28648-6_60.

[12]KIM K J, JIN B P, CHOI Y H. The adaptive learning rates of extended kalman filter based training algorithm for wavelet neural networks[M]. Berlin, Heidelberg: Springer-Verlag, 2006:327-337. DOI:10.1007/11925231_31.

[13]WANG Sisi, WANG Lijun. Maneuvering target tracking based on adaptive square root cubature kalman filter algorithm[M]. Berlin, Heidelberg: Springer-Verlag, 2012:901-908. DOI: 10.1007/978-3-642-34528-9_95.

[14]YANG Yuanxi. Adaptively robust Kalman Filters with applications in navigation [M]. Heidelberg: Springer-Verlag, 2010. DOI 10.1007/978-3-642-11741-1_2.

[15]周荻，尋的導彈新型導引規律[M]. 北京：國防工業出版社，2012. ZHOU Di. New guidance laws for homing missile[M]. Beijing: National Defense Industry Press, 2012.

[16]LIN S G. Assisted adaptive extended Kalman filter for low-cost single-frequency GPS/SBAS kinematic positioning[J]. GPS Solutions, 2015, 19(2):215-223. DOI: 10.1007/s10291-014-0381-9.

[17]ZHAO Yan, GAO Shesheng, ZHANG Jing, et al. Robust predictive augmented unscented Kalman filter[J]. International Journal of Control Automation & Systems, 2014, 12(5):996-1004. DOI: 10.1007/s12555-013-0048-2.

(編輯 張 紅)

A LPF enhanced adaptive Kalman filter for guidance law recognition

DU Runle, LIU Jiaqi, LI Zhifeng, NIU Zhenhong, ZHANG Li

(State Key Laboratory of Science and Technology on Test Physics and Numerical Mathematics (Beijing Institute of Space Long March Vehicle), Beijing 100076, China)

The recognition of the guidance law of a non-cooperative target vehicle is very useful for the counteraction design of the pursuing vehicle. A generic implicit guidance function is presented to describe all known guidance laws, and a practical model with time-varying guidance parameters is used to approximate the generic implicit guidance function. Based on this model, a LPF enhanced Adaptive Kalman Filter is proposed to recognize the guidance law in real-time. Simulation results show that the proposed algorithm can recognize the guidance law fast and effectively when the non-cooperative vehicle employs either the Proportional Navigation guidance law or the Bang-Bang guidance law.

guidance law recognition; implicit guidance function; time-varying guidance parameter; low pass filter; Kalman filter

10.11918/j.issn.0367-6234.201506116

2015-06-29

國家自然科學基金(61471023)

杜潤樂(1975—)，女，高級工程師

劉佳琪，liujiaqi_business@163.com

TP29

A

0367-6234(2017)04-0066-07