湯新民, 鄭鵬程

(南京航空航天大學(xué)民航學(xué)院, 江蘇 南京 211106)

0 引 言

隨著空域內(nèi)航空器密度的增加,管制員的負荷不斷增加,當(dāng)管制員負荷高于一定水平之后,空域內(nèi)產(chǎn)生運行不安全事故的可能性會大大增加。對于在航路運行的航空器,飛行員的工作內(nèi)容少、工作負荷低,因此,如果在航路飛行階段將部分間隔保持責(zé)任由管制員移交給飛行員,可以在降低管制負荷的同時保證運行安全。

將間隔保持責(zé)任轉(zhuǎn)移到飛行員的前提是該空域范圍內(nèi)的航空器具有感知周圍環(huán)境態(tài)勢的能力,在當(dāng)前及可預(yù)見的未來一段時間內(nèi),航空器能夠通過廣播式自動相關(guān)監(jiān)視(automatic dependent surveillance-broadcast, ADS-B)的接收功能(簡稱為ADS-B-IN)獲取周圍航空器的狀態(tài)信息,但目前沒有關(guān)于飛機廣播意圖信息的技術(shù)規(guī)范,所以難以獲取周圍航空器的意圖信息,而周圍航空器的意圖信息使得本機能夠預(yù)測他機在未來一段時間內(nèi)的軌跡,若在計算間隔保證策略時能夠獲得他機的意圖信息,將會極大地提升間隔保持策略的安全與效率。因此,需要根據(jù)航空器的歷史狀態(tài)信息推測航空器的意圖信息,并根據(jù)航空器的意圖信息進行外推,得到航空器在未來一段時間的外推航跡。

單運動模型跟蹤算法在跟蹤目標(biāo)發(fā)生機動時會產(chǎn)生跟蹤模型的失配問題,交互式多模型 (interacting multiple model, IMM)算法是在廣義偽貝葉斯基礎(chǔ)上提出的一種具有馬爾可夫轉(zhuǎn)移概率的算法,IMM算法考慮了目標(biāo)可能處于的多種運動狀態(tài),分別建立多個運動子模型,并考慮子模型之間的交互作用,以達到對機動目標(biāo)的自適應(yīng)跟蹤。

IMM算法的優(yōu)劣取決于模型集的選取是否恰當(dāng),目前常用的模型集包括勻速運動模型、勻加速運動模型、協(xié)同轉(zhuǎn)彎模型等,對機動目標(biāo)加入Singer模型可以在目標(biāo)發(fā)生機動時無需進行機動檢測,就能進行無時間滯后的機動目標(biāo)跟蹤。Singer模型是一個零均值的一階相關(guān)模型,該模型認為目標(biāo)在加速度的均值為0,然而目標(biāo)在機動過程中,下一時刻加速度的均值應(yīng)該是當(dāng)前時刻加速度的均值,因此,根據(jù)Singer模型進行改進得到的“當(dāng)前”統(tǒng)計模型被廣泛用于IMM算法運動模型集中。為了提高目標(biāo)跟蹤的精度,會盡量選取更多的模型加入到模型集中,從而使該算法能夠覆蓋盡可能多的運動情況,但是隨著模型數(shù)量的增加,在每個時間段內(nèi)的計算量會大量增加。此外,模型數(shù)量過多會導(dǎo)致模型間出現(xiàn)競爭問題,導(dǎo)致算法的性能下降。為了在保證跟蹤準(zhǔn)確率的情況下減小計算量,提出了一系列的變結(jié)構(gòu)模型集的IMM算法,例如自適應(yīng)網(wǎng)格算法、機動判別算法、有向圖切換算法等。針對多普勒雷達目標(biāo)跟蹤中非線性的特點,以IMM算法為基本框架,根據(jù)最佳線性無偏估計的卡爾曼濾波和序貫濾波器更新模型概率,最后的估計是序列濾波器輸出和模型概率的加權(quán)和,仿真結(jié)果表明該算法能實現(xiàn)較好的機動目標(biāo)跟蹤精度。針對復(fù)雜機動情況下跟蹤精度低、易發(fā)散的問題,文獻[21-22]對IMM算法中的濾波算法進行了改進,文獻[21]提出了一種基于增益矩陣和轉(zhuǎn)移概率矩陣實時動態(tài)調(diào)整思想的IMM強跟蹤平方空間卡爾曼濾波;文獻[22]提出了一種改進的IMM自適應(yīng)強跟蹤隨機加權(quán)容積卡爾曼濾波;文獻[23]提出了使用極限梯度提升(extreme gradient boosting, XGBoost)方法來替代最大似然估計方法計算目標(biāo)最終的狀態(tài)估計,從而最大化利用系統(tǒng)的先驗信息。針對傳統(tǒng)恒速度子模型和恒加速度子模型在降噪方面的不足,文獻[24]改變了狀態(tài)空間的結(jié)構(gòu),使用了動態(tài)狀態(tài)轉(zhuǎn)移矩陣,從而提高模型的收斂速度并降低噪聲影響。

當(dāng)目標(biāo)的觀測數(shù)據(jù)丟失,IMM算法采用的軌跡外推算法是將記錄最后時刻各模型的概率取值,并根據(jù)模型的馬爾可夫狀態(tài)轉(zhuǎn)移矩陣進行概率矩陣的外推,得到各個模型的概率取值后,各個模型交互得到該時刻目標(biāo)的狀態(tài)信息。

本文采用IMM算法,首先跟蹤航空器的飛行航跡,并辨識航空器子模型的概率分布,然后預(yù)測航空器關(guān)鍵運動參數(shù)的變化趨勢,最后根據(jù)末時航空器狀態(tài)和子模型概率分布對航空器進行合理的短期航跡外推。

1 問題描述和目標(biāo)運動模型

1.1 問題描述

機載ADS-B的發(fā)射功能(簡稱為ADS-B-OUT)可以將航空器的經(jīng)度、緯度、高度、速度、航向、爬升率等信息編碼后進行廣播,本機的ADS-B-IN設(shè)備接收到ADS-B報文后進行解碼,能夠獲取周圍他機的狀態(tài)信息,用于對周圍空域內(nèi)的他機進行追蹤監(jiān)視。對于機載自主間隔保持系統(tǒng)來說,周圍航空器狀態(tài)信息的準(zhǔn)確性至關(guān)重要,因此本文采用IMM算法結(jié)合卡爾曼濾波對ADS-B報文中的航空器狀態(tài)信息進行處理,辨識航空器的關(guān)鍵運動參數(shù),最后根據(jù)航空器的運動趨勢進行短期航跡的外推。

1.2 坐標(biāo)系的選擇

圖1 空間直角坐標(biāo)系中的航空器位置坐標(biāo)Fig.1 Aircraft position coordinates in rectangular coordinate system

圖2 空間直角坐標(biāo)系中的線性外推航跡與理想外推航跡的對比Fig.2 Comparison between linear and ideal extrapolation track in rectangular coordinate system

圖3 航空器在大地坐標(biāo)系中的位置示意圖Fig.3 Position of aircraft in geodetic coordinate system

1.3 目標(biāo)運動模型

航空器的運動狀態(tài)可以用如下的離散系統(tǒng)方程表示:

狀態(tài)方程:

(+1)=()()+()()

(1)

觀測方程:

()=()()+()

(2)

式中:()為狀態(tài)轉(zhuǎn)移矩陣;()為噪聲驅(qū)動矩陣;()∈為狀態(tài)方程白噪聲,其協(xié)方差矩陣為();()為觀測矩陣;()∈為觀測噪聲,其協(xié)方差矩陣為();()∈為系統(tǒng)在時刻的狀態(tài)向量;()為時刻航空器狀態(tài)的觀測向量。其中,()與()為互不相關(guān)的零均值高斯白噪聲。

在大地坐標(biāo)系中,選取的目標(biāo)航空器的運動狀態(tài)變量為

(3)

利用IMM算法跟蹤目標(biāo)航空器的關(guān)鍵在于選取適當(dāng)?shù)倪\動模型作為子模型,模型集應(yīng)當(dāng)能夠覆蓋航空器所有可能的運動狀態(tài)。本文的模型集中共選取3個方向上的6個運動子模型,如表1所示。

首先對巡航階段航空器運動狀態(tài)進行分析,判斷模型集能否滿足覆蓋性要求。

巡航階段飛機可能的運動狀態(tài)及模型集組合如表2所示。

表1 IMM算法模型集

表2 飛機可能的運動狀態(tài)與子模型組合的關(guān)系

根據(jù)分析,可以看出表1中的模型集能夠覆蓋巡航階段航空器的各種運行狀態(tài)。

1.3.1 勻角速度-零速運動模型

(1) 經(jīng)度方向,目標(biāo)做勻角速度運動時,其連續(xù)時間狀態(tài)方程為

(4)

(+1)=()()+()()

(5)

(2) 緯度方向,目標(biāo)的勻角速度運動模型與經(jīng)度方向的勻角速度運動模型一致,所以可得緯度方向的離散時間狀態(tài)方程為

(+1)=()()+()()

(6)

其中,

=[2,, 1]

(3) 高度方向,由于零速運動模型中,航空器高度方向的速度為0,所以離散時間狀態(tài)方程為

(+1)=()()+()()

(7)

其中,

132 “當(dāng)前”運動模型

“當(dāng)前”統(tǒng)計模型是對Singer模型的改進,Singer模型是零均值的時間相關(guān)模型,而實際上航空器下一時刻的加速度的均值是當(dāng)前時刻的加速度,隨機加速度在時間軸上仍然符合一階時間相關(guān)過程,即

(8)

(9)

加速度的方差為

(10)

(11)

(12)

(1) 經(jīng)度方向,目標(biāo)做變加速運動時,根據(jù)式(11)與式(12)可得其連續(xù)時間狀態(tài)方程為

(13)

相應(yīng)的離散形式的運動狀態(tài)方程為

(14)

其中,

=[2,, 1]

(2) 緯度方向,目標(biāo)的變加速運動模型與經(jīng)度方向的變加速運動模型一致,所以可得緯度方向的離散時間狀態(tài)方程為

(15)

其中,

=[2,, 1]

(3) 高度方向,目標(biāo)做變加速升降運動時,根據(jù)式(11)與式(12)可得其連續(xù)時間狀態(tài)方程為

(16)

相應(yīng)的離散形式的運動狀態(tài)方程為

(17)

其中,

=[2,, 1]

133 目標(biāo)運動解耦

本文采用大地坐標(biāo)系進行建模,在使用IMM算法跟蹤目標(biāo)航空器時,如果將經(jīng)度、緯度和高度3個方向上運動模型結(jié)合在一起,可能會導(dǎo)致航空器在3個方向上的運動產(chǎn)生耦合現(xiàn)象。子模型的耦合作用會導(dǎo)致3個方向上的“勻角速度-零速運動模型”的概率取值相同,3個方向上的“當(dāng)前”運動模型的概率取值也相同,這種耦合現(xiàn)象會導(dǎo)致跟蹤效果較差,因此本文將3個方向上的運動進行解耦。

具體的解耦方法如下:

(1) 經(jīng)度方向,對軌跡數(shù)據(jù)中的經(jīng)度數(shù)據(jù)使用經(jīng)度方向的“勻角速度運動模型”和“當(dāng)前”運動模型進行跟蹤;

(2) 緯度方向,對軌跡數(shù)據(jù)中的緯度數(shù)據(jù)使用緯度方向的“勻角速度運動模型”和“當(dāng)前”運動模型進行跟蹤;

(3) 高度方向,對軌跡數(shù)據(jù)中的高度數(shù)據(jù)使用高度方向的“零速運動模型”和“當(dāng)前”運動模型進行跟蹤。

將3個方向的跟蹤數(shù)據(jù)進行合并,即可得到IMM算法的跟蹤結(jié)果。

2 IMM算法

IMM算法是在廣義偽貝葉斯基礎(chǔ)上提出的一種具有馬爾可夫轉(zhuǎn)移概率的算法,其本質(zhì)是將前一時刻多個模型的輸出進行加權(quán)組合作為各個模型的當(dāng)前輸入,多個模型并行估計,最后得到組合狀態(tài)估計。

完整的IMM 循環(huán)由以下4 部分組成: 輸入交互、濾波、模型概率更新、輸出交互。

假設(shè)IMM算法模型集包含有個模型,在時刻,從模型轉(zhuǎn)移到模型(1≤,≤)服從一個給定狀態(tài)轉(zhuǎn)移概率的馬爾可夫矩陣:

={()|(-1)}

(18)

2.1 輸入交互

-1時刻,模型轉(zhuǎn)移到模型的概率為

(-1)=·(-1)

(19)

因此,-1時刻,由其他模型轉(zhuǎn)移到模型的總概率(歸一化常數(shù))為

(20)

模型到模型混合概率為

(21)

模型的混合狀態(tài)估計為

(22)

模型的混合協(xié)方差估計為

{(-1|-1)+

(23)

2.2 狀態(tài)濾波

本文采用卡爾曼濾波算法作為各個子模型的濾波器對輸入交互的結(jié)果進行濾波,首先根據(jù)狀態(tài)轉(zhuǎn)移方程計算模型的估計狀態(tài),然后將目標(biāo)的觀測數(shù)據(jù)作為先驗信息來對狀態(tài)估計進行修正,得到濾波后的狀態(tài)估計?？柭鼮V波的執(zhí)行步驟如下。

根據(jù)狀態(tài)轉(zhuǎn)移方程進行狀態(tài)一步預(yù)測:

(24)

一步預(yù)測的協(xié)方差矩陣為

(25)

卡爾曼濾波的增益矩陣為

()=(|-1)··[(|-1)+]

(26)

使用觀測數(shù)據(jù)作為先驗信息對目標(biāo)狀態(tài)進行更新:

(27)

協(xié)方差矩陣更新:

(|)=[-()()](|-1)

(28)

2.3 模型概率更新

對于第個模型,其似然函數(shù)為

(29)

模型的概率更新方程為

(30)

2.4 輸出交互

濾波器的總輸出為對各個模型進行濾波后的估計結(jié)果的加權(quán)值,加權(quán)的狀態(tài)估計為

(31)

加權(quán)后的協(xié)方差估計為

(32)

3 航空器航跡外推算法

IMM算法可以提供周圍他機的當(dāng)前狀態(tài)信息,但是他機在未來一段時間內(nèi)的航跡對于航空器自主間隔保持也是同等重要的,本節(jié)的主要內(nèi)容是根據(jù)航空器的歷史軌跡信息外推得到他機未來一段時間的航跡,為航空器自主間隔保持算法提供數(shù)據(jù)支持。

對于單一運動子模型的系統(tǒng),可以根據(jù)最近一段時間的狀態(tài)轉(zhuǎn)移矩陣使用式(24)進行狀態(tài)的一步預(yù)測,但對于含有多個子模型的系統(tǒng),狀態(tài)轉(zhuǎn)移矩陣難以確定。

(33)

(34)

式中:()是時刻子模型的概率。由于子模型之間的轉(zhuǎn)移概率服從馬爾可夫過程,因此可以通過-1時刻各個子模型的概率(-1)和得出()：

(35)

4 仿真

4.1 航跡跟蹤仿真

本文使用2021年3月3日,航班號為CCA1516的航空器的ADS-B軌跡數(shù)據(jù)進行濾波跟蹤。首先,將采集過的數(shù)據(jù)進行可視化,判斷軌跡數(shù)據(jù)是否連續(xù),然后將采集到的ADS-B數(shù)據(jù)進行數(shù)據(jù)清洗,然后使用分段3次 Hermite 插值將離散的軌跡點插值成時間間隔=1 s的采樣點,得到1 113個軌跡點。仿真采用的參數(shù)如下。

仿真得到的IMM跟蹤軌跡和觀測軌跡曲線圖如圖4所示。

圖4 航空器觀測航跡和IMM跟蹤軌跡曲線圖Fig.4 Observation track and IMM tracking curve of aircraft

IMM算法跟蹤誤差曲線圖如圖5所示,可以看出,大部分時間段內(nèi),IMM算法跟蹤的誤差小于50 m,僅在極少數(shù)位置會出現(xiàn)大于250 m的誤差,總體的跟蹤效果較好。

圖5 IMM軌跡跟蹤誤差曲線圖Fig.5 IMM trajectory tracking error curve

IMM跟蹤過程中,3個方向上,各個子模型的概率變化曲線圖如圖6所示。

圖6 子模型概率變化曲線圖Fig.6 Probability variation curve of sub model

4.2 航跡外推仿真

選取某一時刻航空器的狀態(tài)作為初始狀態(tài),外推時間為250 s,外推航跡與觀測航跡的對比如圖7(a)所示,外推航跡的誤差隨時間變化的關(guān)系如圖7(b)所示。

圖7 基于大地坐標(biāo)系的航跡外推結(jié)果Fig.7 Track extrapolation result based on geodetic coordinate system

第1.2節(jié)僅從理論上說明了空間直角坐標(biāo)系在航跡外推方面的不足,本節(jié)使用兩種參考系,針對同一初始狀態(tài)、使用相同的外推方法計算航空器的外推航跡?；诳臻g直角坐標(biāo)系的外推航跡與觀測航跡的對比如圖8(a)所示,外推航跡的誤差隨時間變化的關(guān)系如圖8(b)所示。

圖8 基于空間直角坐標(biāo)系的航跡外推結(jié)果Fig.8 Track extrapolation result based on spatial Cartesian coordinate system

基于大地坐標(biāo)系和空間直角坐標(biāo)系的航跡跟蹤誤差對比如圖9所示,在航跡跟蹤過程中,二者相對于觀測航跡的誤差均保持在較低的水平。

圖9 大地坐標(biāo)系和空間直角坐標(biāo)系下的跟蹤誤差對比Fig.9 Tracking error comparison between geodetic coordinate system and spatial cartesian coordinate system

基于大地坐標(biāo)系和空間直角坐標(biāo)系的航跡外推的誤差對比如圖10所示?？梢钥闯?隨著外推時間的增加,基于空間直角坐標(biāo)系的航跡外推的誤差在70 s的時候已經(jīng)增長到大約1 000 m,而此時基于大地坐標(biāo)系的航跡外推的誤差僅為50 m左右;在第250 s時,基于空間直角坐標(biāo)系的航跡外推的誤差超過了7 000 m,而基于大地坐標(biāo)系的航跡外推的誤差為700 m左右。仿真結(jié)果表明,基于大地坐標(biāo)系在航跡外推方面比基于空間直角坐標(biāo)系具有更好的性能。

圖10 大地坐標(biāo)系和空間直角坐標(biāo)系下的外推誤差對比Fig.10 Extrapolation error comparison between geodetic coordinate system and spatial cartesian coordinate system

5 結(jié) 論

本文研究了航路飛行階段的航跡跟蹤和短期航跡外推的問題:

(1) 針對航路運行狀態(tài)下的航空器,使用IMM算法結(jié)合卡爾曼濾波跟蹤航空器的運動軌跡,通過對大地坐標(biāo)系與空間直角坐標(biāo)系進行比較后,采用大地坐標(biāo)系來定位航空器的位置;

(2) 基于大地坐標(biāo)系,將航空器的運動狀態(tài)解耦成3個方向上的獨立運動,分別推導(dǎo)各個方向上可能的運動模型,建立航空器運動方程;

(3) 基于線性外推的思想,根據(jù)航空器末時狀態(tài)和子模型的概率進行航跡外推,同時對比大地坐標(biāo)系和空間直角坐標(biāo)系下航跡外推的性能。結(jié)果表明,在跟蹤性能相近的情況下,采用大地坐標(biāo)系進行航跡外推的性能優(yōu)于空間直角坐標(biāo)系。

本文提出的基于大地坐標(biāo)系的IMM算法的跟蹤效果在與基于空間直角坐標(biāo)系的IMM算法性能相近的情況下,其外推性能有了很大的提升,為航空器短期航跡外推提供了可靠的方法。未來的工作會基于本文提出的短期航跡外推算法進行航空器之間潛在沖突的探測。