一種基于航向估計(jì)的多平臺純方位目標(biāo)轉(zhuǎn)向機(jī)動檢測算法

李巖昊,袁富宇

(江蘇自動化研究所,江蘇 連云港 222061)

純方位目標(biāo)跟蹤(Bearing-Only Tracking, BOT)是通過觀測到的目標(biāo)方位序列估計(jì)目標(biāo)狀態(tài)[1],在軍事領(lǐng)域和民用領(lǐng)域發(fā)揮著重要作用,例如水下監(jiān)視[2]、空中目標(biāo)跟蹤以及無人機(jī)(UAV)路徑規(guī)劃[3-4]。

針對BOT領(lǐng)域,國內(nèi)外學(xué)者對跟蹤模型以及非線性濾波兩個方面進(jìn)行了廣泛研究。在建立跟蹤模型研究方面,Bar-Shalom等人基于廣義偽貝葉斯(GPB)算法[5],提出了一種通過馬爾可夫狀態(tài)轉(zhuǎn)移矩陣來進(jìn)行模型交互的交互多模型算法(Interactive Multiple Model, IMM)[6],相較于單一模型,該算法能夠有效跟蹤機(jī)動形式多變的目標(biāo),但存在自適應(yīng)能力較差的問題。X R Li等人在此基礎(chǔ)上做出改進(jìn),提出了基于可變模型集的變結(jié)構(gòu)多模型算法(Variable Structure Multiple Model, VSMM)[7],該算法能夠在跟蹤過程中自適應(yīng)調(diào)整模型集,提高了跟蹤效率。在處理非線性濾波問題方面,有基于貝葉斯濾波的Kalman濾波及其改進(jìn)算法[8-10],但這類濾波器對初值選取敏感,容易產(chǎn)生濾波發(fā)散的問題。粒子濾波算法[11-12]是通過蒙特卡洛方法來求解貝葉斯估計(jì)中的最優(yōu)估計(jì),在非線性、非高斯系統(tǒng)中被廣泛應(yīng)用,但系統(tǒng)的后驗(yàn)概率密度需要通過選取大量樣本來近似,導(dǎo)致算法復(fù)雜度高,不宜應(yīng)用于復(fù)雜多變的跟蹤場景。此外,其他學(xué)者所提出的基于方位預(yù)測誤差序列均值變化的重復(fù)利用序列概率比檢測[13-14]、批處理方法如Taylor級數(shù)法[15-16]、最小二乘法[17-18]和極大似然估計(jì)[19]等算法,給BOT問題的研究提供了更加多元化的思路。Taylor級數(shù)要素解算法迭代收斂速度快,跟蹤遠(yuǎn)距離勻速直航目標(biāo)時(shí)截?cái)嗾`差小[15],被應(yīng)用到各種被動定位系統(tǒng)中。

本文以Taylor級數(shù)要素解算為基礎(chǔ),針對目標(biāo)轉(zhuǎn)向機(jī)動,提出了基于兩段運(yùn)動要素聯(lián)合解算、兩段運(yùn)動要素獨(dú)立解算兩種解算模型以及基于估計(jì)相鄰段航向差序列的機(jī)動檢測算法,并對機(jī)動時(shí)刻進(jìn)行估計(jì)。仿真結(jié)果表明,兩種解算模型下的機(jī)動檢測算法能夠有效地對轉(zhuǎn)向機(jī)動目標(biāo)進(jìn)行機(jī)動檢測。

1 基于雙平臺Taylor級數(shù)要素解算

水下雙平臺觀測形式如圖1所示,平臺1、平臺2坐標(biāo)為(xo,1,yo,1),(xo,2,yo,2),目標(biāo)在機(jī)動前tk時(shí)刻的運(yùn)動狀態(tài)為

Xk=(xk,yk,VTx,VTy)′+Wkk=1,…,n

(1)

其中,(xk,yk)′為目標(biāo)的位置坐標(biāo),(VTx,VTy)′為目標(biāo)的速度分量(假定目標(biāo)轉(zhuǎn)向機(jī)動前勻速直航),Wk為過程噪聲:Wk=(wx,k,wy,k,wVTx,k,wVTy,k)′k=1,…,n。

圖1 雙平臺觀測Fig.1 Dual-platform observation

目標(biāo)位置與速度的關(guān)系為:

xk=x0+VTx·Δtkk=1,…,n

(2)

yk=y0+VTy·Δtkk=1,…,n

(3)

(x0,y0)′為目標(biāo)初始位置坐標(biāo),Δtk為初始時(shí)刻到tk時(shí)刻的時(shí)間間隔,tk=tk-t0,k=1,…,n,(wx,k,wy,k)′,(wVTx,k,wVTy,k)′為過程噪聲的距離分量和速度分量,取其為高斯白噪聲,過程噪聲距離分量與過程噪聲速度分量的關(guān)系為wx,k=(wVTx,1+…+wVTx,k)·Δt,wy,k=(wVTy,1+…+wVTy,k)·Δt,k=1,…,n。Δt為采樣間隔,由于位置坐標(biāo)可由(2)、(3)式推算出來,解算的目標(biāo)運(yùn)動要素為

X0=(x0,y0,VTx,VTy)′

(4)

平臺1、平臺2的方位序列量測方程為

(5)

(xk-1+wx,k-1,yk-1+wy,k-1)為目標(biāo)在tk-1時(shí)刻的位置,經(jīng)化簡后,(5)式等價(jià)于

i=1,2k=1,…n

(6)

其中,σBi,k為量測噪聲,φBi,k為體現(xiàn)在量測方位上的過程噪聲,給定初值X0,對兩平臺觀測方位在X0處進(jìn)行一階Taylor展開:

Bi,k≈Bi,k(X0)+Bi,k(X0)·dX+φBi,k+σBi,k

i=1,2k=1,…n

(7)

其中,

dX=X-X0

(8)

i=1,2k=1,…n

(9)

將(6)式代入(9)式中,得

(10)

令

(11)

Hi,k=Bi,k(X0)i=1,2k=1,…,n

(12)

則(7)式可改寫為

(13)

(13)式就轉(zhuǎn)化為最小二乘濾波問題,采用線性加權(quán)最小二乘法可證:

i=1,2k=1,…,n

(14)

更新迭代:

X1=X0+dX

(15)

2 基于相鄰段航向差序列的目標(biāo)機(jī)動檢測

目標(biāo)機(jī)動形式為轉(zhuǎn)向機(jī)動。目標(biāo)在機(jī)動時(shí)刻前后勻速直航,量測方位噪聲是高斯白噪聲。

2.1 兩段運(yùn)動要素聯(lián)合解算

假設(shè)目標(biāo)發(fā)生轉(zhuǎn)向機(jī)動,機(jī)動時(shí)刻未知,目標(biāo)機(jī)動前后運(yùn)動要素聯(lián)合解算模型如下:

令機(jī)動前后目標(biāo)速度為:V1=(VTx1,VTy1)′,V2=(VTx2,VTy2)′,機(jī)動時(shí)刻為tm,需要解算的目標(biāo)運(yùn)動要素由(4)式中的四維矢量擴(kuò)維成六維矢量:

X=(x0,y0,VTx1,VTy1,VTx2,VTy2)′

(16)

疊加過程噪聲后,目標(biāo)相對于兩平臺的方位序列量測方程

i=1,2k=1,…,n;tk

(17)

(18)

其中,Δtk=tk-t0,Δt1k=tm-t0,Δt2k=tk-tm,機(jī)動時(shí)刻前,(10)式改為

(19)

機(jī)動時(shí)刻后,(9)式改為

(20)

2.2 兩段運(yùn)動要素獨(dú)立解算

對于目標(biāo)機(jī)動前后的要素解算,另一種方法是把機(jī)動前段和機(jī)動后段獨(dú)立出來,首先估計(jì)出目標(biāo)機(jī)動前的狀態(tài),推算出目標(biāo)在發(fā)生機(jī)動時(shí)的狀態(tài)(位置和速度),作為機(jī)動后段的初始狀態(tài),再對機(jī)動后段進(jìn)行獨(dú)立解算。目標(biāo)機(jī)動前的初始狀態(tài)為

X1=(x0,y0,VTx1,VTy1)′

(21)

目標(biāo)機(jī)動后的初始狀態(tài)為:

X2=(xm,ym,VTx2,VTy2)′

(22)

xm=x0+VTx1·(tm-t0)

(23)

ym=y0+VTy1·(tm-t0)

(24)

其中,tm為機(jī)動時(shí)刻,(xm,ym)′為待估計(jì)的目標(biāo)在第2段的初始位置。

2.3 純方位目標(biāo)機(jī)動檢測

本論文探討的機(jī)動檢測法是基于航向差序列變化進(jìn)行分析的,假定目標(biāo)勻速直航或進(jìn)行轉(zhuǎn)向機(jī)動,但機(jī)動時(shí)刻tm未知,如圖2所示。

圖2 機(jī)動時(shí)刻未知的目標(biāo)航路Fig.2 Unknown-maneuvering target situation

(25)

其中,β表示航向差門限,g0表示置信水平門限。N(|ΔCTl|>β)表示航向差絕對值大于門限β的對應(yīng)的假定機(jī)動點(diǎn)個數(shù)。

2.4 機(jī)動時(shí)刻估計(jì)

(26)

(27)

(28)

3 仿真分析

3.1 評價(jià)指標(biāo)

考慮如下評價(jià)指標(biāo):

2)虛警率:把目標(biāo)的勻速直航檢測成機(jī)動的概率,計(jì)算模式同1)。

3)機(jī)動時(shí)刻的估計(jì)精度:在正確檢測的條件下,統(tǒng)計(jì)估計(jì)機(jī)動時(shí)刻的精度,用估計(jì)機(jī)動時(shí)刻均方誤差σΔtm表示。

4)機(jī)動檢測延遲時(shí)間:是指檢測到目標(biāo)機(jī)動的時(shí)間與實(shí)際目標(biāo)發(fā)生機(jī)動的時(shí)間差,用均方誤差σΔ(tk-tm)表示,其中,tk表示檢測到機(jī)動的當(dāng)前時(shí)刻。

3.2 算例航路參數(shù)

以大地直角坐標(biāo)系為例,觀測平臺1、平臺2坐標(biāo)設(shè)為P1(0,0),P2(15,0)(單位:km),目標(biāo)初始位置(0,D0)(下文用目標(biāo)初距D0表示);目標(biāo)機(jī)動前的速度為VT1,機(jī)動后的速度為VT2;目標(biāo)機(jī)動前的航向?yàn)镃T1,機(jī)動后的航向?yàn)镃T2,機(jī)動轉(zhuǎn)向角α(順時(shí)針為正,逆時(shí)針為負(fù));過程噪聲φB,量測噪聲σB;兩平臺量測到的方位角B1、B2如圖3所示(圖3假設(shè)目標(biāo)初始航向?yàn)?0°);取β=10°,g0=75%,采樣間隔為1 s,量測噪聲σB=1°,過程噪聲φB=1 kn,噪聲類型均為高斯白噪聲,目標(biāo)機(jī)動方式為轉(zhuǎn)向機(jī)動。

仿真計(jì)算航路參考數(shù)值如表1所示,共有3×3×3×5×10=1 350種組合:

圖3 目標(biāo)航路Fig.3 Target situation

表1 仿真參考數(shù)值

迭代初值選取:目標(biāo)初距D0在(6,30)km區(qū)間內(nèi)隨機(jī)選取,速度VT1在(8,30)節(jié)內(nèi)隨機(jī)選取,初始舷角QT在(-180°,180°)內(nèi)隨機(jī)選取。

3.3 仿真結(jié)果及分析

1)勻速直航虛警率

以目標(biāo)初始航向90°,目標(biāo)速度18 kn,目標(biāo)初始位置分別為30、40、50(單位:km),量測噪聲σB=1°,過程噪聲φB=1 kn,目標(biāo)總航行時(shí)間12 min,算法分別在T=5 min、T=6 min、T=7 min時(shí)刻啟動,對目標(biāo)進(jìn)行轉(zhuǎn)向機(jī)動檢測。每組航路仿真10 000次,統(tǒng)計(jì)目標(biāo)勻速直航時(shí)基于雙平臺Taylor級數(shù)要素解算算法的迭代收斂次數(shù)、基于兩種解算模型的目標(biāo)機(jī)動檢測算法的虛警率。迭代收斂依照(13)式中‖dX‖2<10-7來判斷,表2為統(tǒng)計(jì)結(jié)果。

目標(biāo)勻速直航時(shí),對于目標(biāo)初始距離30 km、速度18 kn、量測噪聲1°、過程噪聲1 kn的航路下,解算10 000 次航路,結(jié)果全部迭代收斂,算法的虛警率為0%,隨著目標(biāo)初始距離增大,兩種解算模式下的機(jī)動檢測算法均出現(xiàn)虛警,且虛警率與算法啟動時(shí)刻相關(guān),以目標(biāo)初距50 km的航路為例,算法啟動時(shí)刻T=5 min、T=6 min、T=7 min時(shí),聯(lián)合解算對應(yīng)的機(jī)動檢測虛警率分別為15.11%、9.25%、6.25%;獨(dú)立解算對應(yīng)的機(jī)動檢測虛警率分別為7.69%、6.69%、4.69%。原因是算法啟動時(shí)間越長,收集的量測數(shù)據(jù)越多,解算精度越高,虛警率越低。

2)機(jī)動檢測率

為驗(yàn)證算法對于機(jī)動角度的適應(yīng)范圍,統(tǒng)計(jì)在仿真航路參數(shù)范圍內(nèi)目標(biāo)機(jī)動檢測延遲時(shí)間,包括在目標(biāo)機(jī)動后兩分鐘內(nèi)以及機(jī)動后三分鐘內(nèi)的機(jī)動檢測率,其中目標(biāo)機(jī)動時(shí)刻tm=6 min,目標(biāo)總觀測時(shí)間tn=12 min。算法啟動時(shí)刻T=6 min。機(jī)動轉(zhuǎn)角α在15°～60°范圍內(nèi)選取。表3為統(tǒng)計(jì)結(jié)果。

表2 勻速直航解算虛警率Tab.2 Convergence and false alarm rate of constant speed direct flight

表3 不同機(jī)動轉(zhuǎn)角下目標(biāo)機(jī)動檢測率Tab.3 Target maneuver detection rate at different maneuvering angles

在機(jī)動轉(zhuǎn)角大于等于30°的航路下,兩種解算模型下算法在機(jī)動后兩分鐘內(nèi)的機(jī)動檢測率大于90°。仿真結(jié)果表明算法適用于機(jī)動轉(zhuǎn)角大于等于30°的航路。

3)航路參數(shù)對機(jī)動檢測結(jié)果影響

統(tǒng)計(jì)目標(biāo)轉(zhuǎn)向機(jī)動時(shí),目標(biāo)初始位置、目標(biāo)速度、機(jī)動轉(zhuǎn)角大小對目標(biāo)機(jī)動檢測延遲時(shí)間和目標(biāo)機(jī)動時(shí)刻估計(jì)精度的影響。

圖4、圖5表示目標(biāo)初始距離對兩種解算模型下機(jī)動檢測算法的機(jī)動檢測延遲時(shí)間和機(jī)動時(shí)刻估計(jì)精度的影響,對應(yīng)仿真的目標(biāo)航路參數(shù):目標(biāo)初距D0在(30,50)km范圍變化,目標(biāo)速度18 kn、機(jī)動轉(zhuǎn)角45°、初始弦角90°、量測噪聲σB=1°;過程噪聲φB=1 kn。

圖4 目標(biāo)初距對機(jī)動檢測延遲時(shí)間的影響Fig. 4 Influence of target initial distance on the delay time of maneuver detection

圖5 目標(biāo)初距對機(jī)動時(shí)刻估計(jì)精度的影響Fig. 5 Influence of target initial distance on the accuracy of maneuvering time estimation

圖6、圖7表示目標(biāo)速度對兩種解算模型下機(jī)動檢測算法的機(jī)動檢測延遲時(shí)間和機(jī)動時(shí)刻估計(jì)精度的影響。對應(yīng)仿真的目標(biāo)航路參數(shù):目標(biāo)速度VT1、VT2在(10,24)kn范圍變化,目標(biāo)初距為30 km,機(jī)動轉(zhuǎn)角45°,初始弦角90°,量測噪聲1°,過程噪聲1 kn。

圖6 目標(biāo)速度對機(jī)動檢測延遲時(shí)間的影響Fig. 6 Influence of target speed on maneuver detection delay time

圖8、圖9表示機(jī)動轉(zhuǎn)角對兩種解算模型下機(jī)動檢測算法的機(jī)動檢測延遲時(shí)間和機(jī)動時(shí)刻估計(jì)精度的影響。對應(yīng)仿真的目標(biāo)航路參數(shù):機(jī)動轉(zhuǎn)角α在(-60°,60°)內(nèi)范圍變化,目標(biāo)初距30 km,目標(biāo)速度18 kn,初始弦角90°,量測噪聲1°,過程噪聲1 kn。

圖8 機(jī)動轉(zhuǎn)角對機(jī)動檢測延遲時(shí)間的影響Fig. 8 Influence of maneuver angle on maneuver detection delay time

目標(biāo)航路對機(jī)動檢測算法的影響主要表現(xiàn)在:隨著目標(biāo)與觀測平臺初始距離增大,機(jī)動檢測延遲時(shí)間增大,且出現(xiàn)虛警,同時(shí),機(jī)動時(shí)刻估計(jì)精度降低;隨著目標(biāo)速度降低,觀測到的方位序列變化不明顯,這導(dǎo)致目標(biāo)機(jī)動檢測時(shí)間增大,機(jī)動時(shí)刻估計(jì)精度降低。目標(biāo)轉(zhuǎn)向方向同樣對機(jī)動檢測結(jié)果存在影響,圖9中,左舷側(cè)機(jī)動和右舷側(cè)機(jī)動時(shí),機(jī)動時(shí)刻估計(jì)結(jié)果相差較大,尤其是正負(fù)30°轉(zhuǎn)角處相差最大。這是由于在向右舷側(cè)轉(zhuǎn)向機(jī)動時(shí),目標(biāo)靠近觀測平臺,相對于觀測平臺的距離更近,方位變化明顯;目標(biāo)向左舷側(cè)轉(zhuǎn)向機(jī)動則遠(yuǎn)離觀測平臺,方位變化并沒有右舷機(jī)動明顯,為驗(yàn)證這一觀點(diǎn),選定對稱航路進(jìn)行仿真實(shí)驗(yàn),令目標(biāo)初始位置位于觀測平臺中央,即目標(biāo)初始位置(x0,y0)=(10,30)km,目標(biāo)機(jī)動前后速度VT1=VT2=18 kn,初始航向CT1=180°,過程噪聲фB=1 kn,量測噪聲σB=1°,機(jī)動時(shí)刻tm=6 min,統(tǒng)計(jì)機(jī)動轉(zhuǎn)角|α|=30°～60°時(shí),機(jī)動檢測延遲時(shí)間和機(jī)動時(shí)刻估計(jì)精度。每組航路仿真10 000 次,機(jī)動檢測結(jié)果圖10、圖11所示。

圖10 機(jī)動轉(zhuǎn)角對機(jī)動檢測延遲時(shí)間的影響Fig.10 Influence of measurement error on the delay time of maneuver detection

圖11 機(jī)動轉(zhuǎn)角對機(jī)動時(shí)刻估計(jì)精度的影響Fig.11 Influence of measurement error on the accuracy of maneuvering time estimation

仿真結(jié)果表明,對稱航路下,左舷機(jī)動和右舷機(jī)動檢測結(jié)果大致相同,二者的機(jī)動檢測延遲時(shí)間均低于90 s,機(jī)動時(shí)刻估計(jì)精度均低于15 s。

4 結(jié)束語

本文針對純方位條件下的變向機(jī)動目標(biāo),根據(jù)目標(biāo)航向變化規(guī)律進(jìn)行機(jī)動檢測,基于反推方位均方誤差進(jìn)行機(jī)動時(shí)刻估計(jì)。該算法針對目標(biāo)初距30～50 km、速度10～25 kn,過程噪聲速度分量1 kn,量測噪聲1°,轉(zhuǎn)向角不小于30°的航路,算法全部迭代收斂,機(jī)動檢測延遲時(shí)間在2 min以內(nèi),機(jī)動時(shí)刻估計(jì)誤差小于60 s,機(jī)動檢測算法的虛警率低于10%,仿真結(jié)果表明該算法可以有效地對水面及水下轉(zhuǎn)向機(jī)動目標(biāo)進(jìn)行機(jī)動檢測及跟蹤。