混合噪聲條件下的目標(biāo)被動(dòng)定位算法

陳林秀， 郝明瑞， 趙佳佳

(復(fù)雜系統(tǒng)控制與智能協(xié)同技術(shù)重點(diǎn)實(shí)驗(yàn)室， 北京 100074)

0 引言

基于只測(cè)角體制的飛行器以其導(dǎo)引頭作用距離遠(yuǎn)、不主動(dòng)發(fā)射電磁波、隱蔽性好等優(yōu)點(diǎn)，已成為打擊雷達(dá)等輻射源目標(biāo)的主戰(zhàn)武器。為應(yīng)對(duì)目標(biāo)雷達(dá)突然關(guān)機(jī)，飛行器制導(dǎo)系統(tǒng)設(shè)計(jì)的一個(gè)關(guān)鍵問(wèn)題就是如何通過(guò)導(dǎo)引頭探測(cè)到的信息對(duì)目標(biāo)進(jìn)行被動(dòng)定位，即基于角度觀測(cè)信號(hào)的被動(dòng)定位濾波算法設(shè)計(jì)。

現(xiàn)有的目標(biāo)定位濾波算法大多采用卡爾曼濾波(KF)算法，基本KF算法是一種線性最小方差估計(jì)。但是在被動(dòng)定位系統(tǒng)中，不論在直角坐標(biāo)系中還是極坐標(biāo)系中，目標(biāo)定位系統(tǒng)都是一個(gè)強(qiáng)非線性系統(tǒng)。近年來(lái)出現(xiàn)了眾多關(guān)于非線性濾波的理論和方法，其中Arasaratnam等[1]于2009年提出的容積卡爾曼濾波(CKF)算法在濾波精度、算法穩(wěn)定性和計(jì)算量等方面表現(xiàn)較為突出。研究成果表明該算法在處理高維濾波問(wèn)題時(shí)，可獲得比擴(kuò)展卡爾曼濾波(EKF)、中心差分卡爾曼濾波(CDKF)、無(wú)跡卡爾曼濾波(UKF)等非線性估計(jì)算法更高的濾波精度、穩(wěn)定性和相對(duì)低的計(jì)算復(fù)雜度[2-5]。然而，CKF算法和基本線性濾波方法一樣，均假定過(guò)程噪聲及觀測(cè)噪聲為高斯白噪聲，實(shí)際上在被動(dòng)定位應(yīng)用中，過(guò)程噪聲來(lái)源于運(yùn)動(dòng)學(xué)建模，其建模準(zhǔn)確度較高；而觀測(cè)噪聲受電磁波傳播空間、目標(biāo)所處環(huán)境、電磁干擾以及天線罩等因素的影響，往往是含有時(shí)變有色噪聲、跳變?cè)肼暤榷喾N特性的混合噪聲，在此情況下，若仍將其當(dāng)作白噪聲對(duì)待則可能會(huì)導(dǎo)致濾波精度下降甚至濾波器發(fā)散[6-10]。文獻(xiàn)[11]推導(dǎo)了有色噪聲條件下的EKF算法，該算法結(jié)構(gòu)簡(jiǎn)單，但存在復(fù)雜的雅克比矩陣計(jì)算問(wèn)題，且濾波精度不高，受參數(shù)影響較大。文獻(xiàn)[12]利用歷元噪聲的相關(guān)性特征構(gòu)建多步相關(guān)的噪聲協(xié)方差陣，通過(guò)線性變換得到了改進(jìn)的狀態(tài)協(xié)方差和增益矩陣。該算法能減弱有色噪聲對(duì)系統(tǒng)的擾動(dòng)，但是過(guò)程噪聲與觀測(cè)噪聲之間的協(xié)方差計(jì)算相當(dāng)復(fù)雜，當(dāng)考慮相關(guān)步數(shù)增多時(shí)算法計(jì)算量急劇增大。目前，被動(dòng)定位背景下的濾波估計(jì)算法大多基于EKF或UKF濾波算法，且對(duì)傳感器觀測(cè)噪聲的適應(yīng)性還有待提升。

本文以基本CKF算法為基礎(chǔ)，在建立的被動(dòng)定位濾波模型基礎(chǔ)上，探索能夠適應(yīng)混合噪聲的CKF算法，以期進(jìn)一步提升算法對(duì)混合觀測(cè)噪聲的適應(yīng)性和魯棒性。

1 被動(dòng)定位問(wèn)題的濾波模型

1.1 狀態(tài)方程的建立

選取制導(dǎo)坐標(biāo)系Oxyz下飛行器- 目標(biāo)相對(duì)位置矢量、相對(duì)速度矢量、目標(biāo)加速度矢量在其3個(gè)軸上的投影分量(Δx、Δy、Δz)、(Δvx、Δvy、Δvz)、(atx、aty、atz)作為狀態(tài)向量，表示為

X=[ΔxΔyΔzΔvxΔvyΔvzatxatyatz]T，

式中：X表示系統(tǒng)的狀態(tài)向量；下標(biāo)t表示目標(biāo)。則系統(tǒng)狀態(tài)方程表示為

(1)

式中：A為狀態(tài)的轉(zhuǎn)移矩陣；G為飛行器加速度分配矩陣；U表示飛行器加速度矢量在制導(dǎo)坐標(biāo)系3個(gè)軸上的分量，U=[axayaz]T；Γ為系統(tǒng)過(guò)程噪聲分配矩陣；W表示過(guò)程噪聲，W=[wxwywz]T為零均值高斯分布白噪聲向量，且各分量間彼此相互獨(dú)立，wx、wy、wz分別為制導(dǎo)坐標(biāo)系3個(gè)軸上的過(guò)程噪聲。

考慮艦船目標(biāo)的弱機(jī)動(dòng)性，本文中目標(biāo)模型采用Singer加速度模型[13]，用α表示目標(biāo)機(jī)動(dòng)時(shí)間常數(shù)的倒數(shù)，即目標(biāo)機(jī)動(dòng)頻率。則相關(guān)矩陣表示為

(2)

G=[03×3-I3×303×3]T，

(3)

Γ=[03×303×3I3×3]T.

(4)

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

xk=Fk|k-1xk-1+Bk-1uk-1+Γk|k-1Wk-1，

(5)

式中：xk表示系統(tǒng)在離散時(shí)間k的狀態(tài)向量，k表示離散的時(shí)間量；Fk|k-1表示離散時(shí)間k-1到k的狀態(tài)轉(zhuǎn)移矩陣，

(6)

T為濾波周期即仿真步長(zhǎng)；Bk-1表示離散時(shí)間k-1的飛行器加速度分配矩陣，

(7)

uk-1表示離散時(shí)間k-1的飛行器加速度矢量在制導(dǎo)坐標(biāo)系3個(gè)軸上的分量；Γk|k-1表示離散時(shí)間k-1到k的系統(tǒng)過(guò)程噪聲分配矩陣，

(8)

1.2 觀測(cè)方程的建立

被動(dòng)傳感器可獲得目標(biāo)視線高低角qf和視線方位角qh，如圖1所示。

圖1 被動(dòng)傳感器所測(cè)信息Fig.1 Information measured by passive sensor

qf、qh的理論計(jì)算公式為

(9)

式中：hf(x)為計(jì)算視線高低角理論值的非線性函數(shù)；hh(x)為計(jì)算視線方位角理論值的非線性函數(shù)。

考慮傳感器的觀測(cè)誤差，可建立如下觀測(cè)方程：

(10)

2 非線性濾波被動(dòng)定位算法

2.1 基本CKF算法

基本CKF算法是依據(jù)貝葉斯濾波器的遞推過(guò)程，通過(guò)假設(shè)各概率密度函數(shù)均服從高斯分布獲得高斯域的貝葉斯濾波框架，在高斯濾波框架上應(yīng)用三度容積規(guī)則近似得到的，文獻(xiàn)[1]給出了其具體推導(dǎo)過(guò)程。此處，結(jié)合第1節(jié)建立的被動(dòng)濾波定位模型，給出被動(dòng)跟蹤問(wèn)題中的基本CKF過(guò)程。

考慮如(5)式和(10)式建立的非線性隨機(jī)系統(tǒng)，寫(xiě)成一般形式如(11)式所示：

(11)

在此條件下，CKF算法的計(jì)算流程如下：

2) 時(shí)間更新：

k|k-1=Fk|k-1k-1+Bk-1uk-1，

(12)

(13)

3) 觀測(cè)更新：

①計(jì)算Cubature點(diǎn)：

(14)

(15)

②傳播Cubature點(diǎn)：

Zi,k|k-1=h(Xi,k|k-1)，

(16)

(17)

互協(xié)方差陣Pxz,k|k-1和觀測(cè)預(yù)測(cè)協(xié)方差陣Pzz,k|k-1為

(18)

(19)

于是，可以得到觀測(cè)矩陣為

Hk=((Pk|k-1)-1Pxz,k|k-1)T，

(20)

濾波增益為

Kk=Pxz,k|k-1(Pzz,k|k-1)-1，

(21)

以及后驗(yàn)狀態(tài)估計(jì)值及其協(xié)方差陣為

k=k|k-1+Kk(zk-k|k-1)，

(22)

由以上濾波過(guò)程可知，獲得Pxz,k|k-1之后，可結(jié)合Pk|k-1計(jì)算得到觀測(cè)矩陣Hk，于是非線性系統(tǒng)中的濾波問(wèn)題可轉(zhuǎn)化至線性系統(tǒng)中進(jìn)行分析和處理。

2.2 有色噪聲條件下的CKF算法

2.2.1 建立有色噪聲模型

傳感器觀測(cè)噪聲中，有色隨機(jī)噪聲為主要分布噪聲，因此本節(jié)假設(shè)觀測(cè)噪聲vk為有色隨機(jī)噪聲特征，在基本CKF算法的基礎(chǔ)上推導(dǎo)得到有色噪聲條件下的容積卡爾曼濾波(CKF-CMN)算法。傳感器觀測(cè)噪聲模型中對(duì)有色噪聲的描述為：正態(tài)分布的白噪聲通過(guò)時(shí)間常數(shù)為τ的1階濾波器后的輸出(標(biāo)準(zhǔn)差為σ)表示為[σ,τ]。根據(jù)歐拉法求解微分方程可以獲得有色噪聲序列，同時(shí)可以將有色噪聲寫(xiě)成如下形式：

vk=(1-T/τ)·vk-1+(T/τ)·rk，

(23)

式中：rk表示1階濾波器輸入的正態(tài)分布白噪聲，rk的方差為Rr.將(23)式寫(xiě)成(24)式所示的1階馬爾可夫形式：

vk=φvk-1+ξk，

(24)

式中：φ=1-T/τ；ξk=(T/τ)·rk，ξk的方差Rk=(T/τ)2·Rr.

2.2.2 有色噪聲的白化

考慮與(11)式等價(jià)的線性系統(tǒng)：

(25)

式中：wk-1與(11)式中的Γk|k-1Wk-1等價(jià)，vk=φvk-1+ξk.為了使CKF算法在有色觀測(cè)噪聲條件下仍然保持優(yōu)良的濾波性能，首先對(duì)有色噪聲進(jìn)行白化，構(gòu)造新的觀測(cè)量zk-φzk-1，經(jīng)等價(jià)替換及整理得

zk-φzk-1=(HkFk|k-1-φHk-1)xk-1+
HkBk-1uk-1+Hkwk-1+ξk.

(26)

令

(27)

(28)

(29)

則觀測(cè)方程可以改寫(xiě)為

(30)

(31)

(32)

(33)

通過(guò)以上分析可知，將有色噪聲白化后，系統(tǒng)的觀測(cè)噪聲模型便符合基本KF器的觀測(cè)噪聲為白噪聲這一假設(shè)，然而與基本KF假設(shè)不同的是，此時(shí)觀測(cè)噪聲和系統(tǒng)噪聲之間存在一定的相關(guān)性。下面考慮采用待定系數(shù)去相關(guān)法對(duì)上述相關(guān)性進(jìn)行消除。

2.2.3 基于待定系數(shù)法的噪聲去相關(guān)方法

完成有色觀測(cè)噪聲白色化后，由(30)式可將(25)式等效為(34)式所示的一般形式：

(34)

(34)式中的狀態(tài)方程改寫(xiě)為

(35)

(36)

(37)

選擇

(38)

(39)

(40)

綜上所述，可得到CKF-CMN算法計(jì)算過(guò)程如下：

3)根據(jù)(23)式～(40)式，構(gòu)造新的觀測(cè)方程和狀態(tài)方程，按以下步驟對(duì)濾波過(guò)程進(jìn)行更新：

濾波增益：

(41)

狀態(tài)估計(jì)及其協(xié)方差：

(42)

2.3 ACKF算法

針對(duì)混合觀測(cè)噪聲中的跳變?cè)肼曇约坝猩肼暱赡軙r(shí)變給CKF-CMN算法帶來(lái)的影響，在CKF-CMN算法基礎(chǔ)上設(shè)計(jì)自適應(yīng)濾波算法，實(shí)時(shí)對(duì)混合噪聲的噪聲協(xié)方差進(jìn)行估計(jì)，并對(duì)觀測(cè)噪聲中的有害信息進(jìn)行剔除。

2.3.1 對(duì)噪聲協(xié)方差的估計(jì)

(43)

(44)

(45)

2.3.2 有害觀測(cè)信息的剔除

H0：觀測(cè)信息有效；

H1：觀測(cè)信息無(wú)效，剔除。

采用檢驗(yàn)統(tǒng)計(jì)量Dk

(46)

3 數(shù)學(xué)仿真

圖2 飛行器與目標(biāo)運(yùn)動(dòng)軌跡Fig.2 Moving trajectories of aircraft and target

本文仿真時(shí)間為100 s，仿真步長(zhǎng)為0.01 s，蒙特卡洛仿真次數(shù)N為50次。假設(shè)被動(dòng)傳感器測(cè)角噪聲為時(shí)變有色噪聲和跳變?cè)肼暤幕旌显肼?。?duì)于有色噪聲[σ,τ](參考2.2節(jié)，σ、τ分別表示其標(biāo)準(zhǔn)差和時(shí)間常數(shù))。某次仿真中其變化特性如圖3所示，圖中Tc表示噪聲方差發(fā)生變化的時(shí)刻，在Tc=50 s處，有色噪聲特性由[0.5°，0.02 s]變?yōu)閇1.0°，0.02 s]，跳變?cè)肼暦禐?.5°，周期為5 s，某次仿真中其變化特性如圖4所示，二者的混合噪聲如圖5所示。

圖3 時(shí)變有色噪聲Fig.3 Time-varying colored noise

圖4 跳變?cè)肼旻ig.4 Jumping noise

圖5 被動(dòng)傳感器觀測(cè)噪聲隨時(shí)間變化Fig.5 Passive sensor observation noise vs. time

在觀測(cè)噪聲為混合噪聲條件下采用基本CKF、CKF-CMN以及ACKF算法在同一條件下進(jìn)行對(duì)比仿真，結(jié)果如圖6～圖10所示，圖中縱軸為目標(biāo)位置分量的均方根誤差(RMSE)，即x軸、y軸和z軸方向目標(biāo)位置真實(shí)值和估計(jì)值的均方根誤差(RMSEx、RMSEy、RMSEz)，RMSEx、RMSEy、RMSEz根據(jù)(47)式計(jì)算得到：

圖6 x軸方向被動(dòng)定位估計(jì)誤差Fig.6 Estimated errors of passive positioning in x direction

圖10 z軸方向被動(dòng)定位估計(jì)誤差(放大圖)Fig.10 Estimated errors of passive positioning in z direction (enlarged view)

(47)

圖7 x軸方向被動(dòng)定位估計(jì)誤差(放大圖)Fig.7 Estimated errors of passive positioning in x direction (enlarged view)

圖8 y軸方向被動(dòng)定位估計(jì)誤差Fig.8 Estimated errors of passive positioning in y direction

圖9 z軸方向被動(dòng)定位估計(jì)誤差Fig.9 Estimated errors of passive positioning in z direction

取圖中第1～100 s的數(shù)據(jù)進(jìn)行計(jì)算，得到x軸、y軸、z軸3個(gè)方向RMSE的平均值如表1所示。

由圖6～圖10可知：當(dāng)被動(dòng)定位系統(tǒng)的觀測(cè)數(shù)據(jù)中存在時(shí)變有色噪聲及跳變?cè)肼晻r(shí)，基本KF算法的精度和適用性受到明顯制約，難以實(shí)現(xiàn)高精度的目標(biāo)位置估計(jì)。在觀測(cè)量包含混合噪聲時(shí)，3種被動(dòng)定位算法的精度依次為CKF算法

表1 3種濾波定位算法精度分析

4 結(jié)論

本文從混合觀測(cè)噪聲條件下的被動(dòng)定位精度和濾波穩(wěn)定性提升需求出發(fā)，在分析CKF算法的基礎(chǔ)上，通過(guò)重新構(gòu)建觀測(cè)方程和狀態(tài)方程提高所建立的濾波器對(duì)混合噪聲的適應(yīng)性。得出主要結(jié)論如下：

1)針對(duì)混合噪聲中的有色噪聲使得基本CKF算法定位精度下降的問(wèn)題，設(shè)計(jì)了CKF-CMN算法。通過(guò)觀測(cè)重構(gòu)、待定系數(shù)去相關(guān)的思想，構(gòu)建新的濾波狀態(tài)方程和觀測(cè)方程，將有色噪聲轉(zhuǎn)變?yōu)榘自肼?，使得濾波估計(jì)精度得到提升。

2)針對(duì)有色噪聲的時(shí)變性以及跳變?cè)肼暿沟脼V波性能受損的問(wèn)題，在CKF-CMN算法的基礎(chǔ)上設(shè)計(jì)了ACKF算法，采用實(shí)時(shí)噪聲方差估計(jì)和有害信息剔除法使得被動(dòng)定位精度得到進(jìn)一步提升。

3)50次蒙特卡洛仿真實(shí)驗(yàn)結(jié)果表明，ACKF算法的被動(dòng)濾波定位精度得到提升且魯棒性更好，基本CKF、CKF-CMN、ACKF 3種算法的定位精度依次為CKF算法