基于自適應(yīng)滑動(dòng)時(shí)間窗的四維航跡預(yù)測(cè)*

朱代武 劉 豪 黃邦菊 史繼龍 陳澤暉

（1.中國(guó)民用航空飛行學(xué)院圖書館 廣漢 618307）（2.中國(guó)民用航空飛行學(xué)院空中交通管理學(xué)院 廣漢 618307）

1 引言

隨著航空器數(shù)量的劇增，使得在其運(yùn)行過(guò)程中出現(xiàn)空域資源緊張的問(wèn)題，進(jìn)而造成終端區(qū)和航路擁堵，為此美國(guó)和歐洲分別提出下一代空中交通管理系統(tǒng)-NextGen（Next Generation）和 SESAR（Single European Sky ATM Research）［1～3］。上述系統(tǒng)的核心皆是基于航跡的空域運(yùn)行，即航行方式實(shí)現(xiàn)從基于傳感器導(dǎo)航到基于性能導(dǎo)航的轉(zhuǎn)變。同時(shí)航跡預(yù)測(cè)是輔助管制員決策的關(guān)鍵，因此航跡預(yù)測(cè)成為未來(lái)空管自動(dòng)化決策的中堅(jiān)力量。

目前對(duì)航跡預(yù)測(cè)主要采用以下兩種方式：一是基于卡爾曼波的無(wú)參數(shù)估計(jì)。常用的是Lymperopoulos［4］等采用蒙特卡洛法簡(jiǎn)介對(duì)航跡進(jìn)行評(píng)估及Wu［5］等采用的基于歷史數(shù)據(jù)的多元線性回歸法。二是基于航空器性能或氣動(dòng)學(xué)模型的方法。常見(jiàn)的是Warren等采用航空器運(yùn)動(dòng)方程，利用混雜系統(tǒng)遞推求解航跡模型。但上述兩種方法經(jīng)常受到航空器機(jī)動(dòng)性能多樣性、大氣環(huán)境及單一模型濾波的使用限制，因此本文提出基于自適應(yīng)滑動(dòng)時(shí)間窗的四維航跡預(yù)測(cè)方法，通過(guò)滑動(dòng)時(shí)間窗估計(jì)下一時(shí)刻的狀態(tài)矢量，進(jìn)而提高航跡預(yù)測(cè)的精度。

2 預(yù)測(cè)模型

在航空器實(shí)際運(yùn)行過(guò)程中，需通過(guò)ADS-B等雷達(dá)設(shè)備及時(shí)獲取并進(jìn)行狀態(tài)評(píng)估：航班號(hào)、位置信息、運(yùn)動(dòng)狀態(tài)及二次應(yīng)答機(jī)編碼?？柭鼮V波是一個(gè)對(duì)動(dòng)態(tài)序列進(jìn)行線性最小方差估計(jì)的算法，適合用來(lái)解決目標(biāo)狀態(tài)跟蹤及位置識(shí)別等問(wèn)題［6～7］。其四維航跡預(yù)測(cè)模型的結(jié)構(gòu)示意圖如圖1所示（Visio2015版）。

圖1 四維航跡預(yù)測(cè)模型的結(jié)構(gòu)示意圖

圖2 航空器地速分解示意圖

假設(shè)航空器的地速為V，航向?yàn)棣?，?duì)其進(jìn)行速度分解為X軸的速度分量Vx，Y軸的速度分量為Vy即［8］：

卡爾曼濾波根據(jù)系統(tǒng)的狀態(tài)方程和觀測(cè)方程，估計(jì)動(dòng)態(tài)系統(tǒng)下一時(shí)刻的狀態(tài)。即：

其中X（m）為系統(tǒng)狀態(tài)矢量，表征在m時(shí)刻點(diǎn)的所有的運(yùn)動(dòng)向量；A（m|m-1）為描述目標(biāo)體運(yùn)動(dòng)的狀態(tài)轉(zhuǎn)移方程；L（m|m-1）為干擾轉(zhuǎn)移矩陣；w（m）表示運(yùn)動(dòng)模型的系統(tǒng)噪聲，是p維零均值的白噪聲限；Z（m）為所觀測(cè)向量，描述m時(shí)刻的系統(tǒng)觀測(cè)值；H（m）為所觀測(cè)矩陣；v（m）為運(yùn)動(dòng)估計(jì)過(guò)程中所產(chǎn)生的噪聲。

系統(tǒng)噪聲w（m）與觀測(cè)噪聲v（m）相互獨(dú)立存在的高斯噪聲，其統(tǒng)計(jì)特性為

其中Q（m）是系統(tǒng)噪聲w（m）的對(duì)稱非負(fù)矩陣，R（m）為所觀測(cè)到的總噪聲v（m）的對(duì)稱方差矩陣。

假設(shè)時(shí)間更新?tīng)顟B(tài)初始分布的均值和協(xié)方差矩陣分別為X0|0=X0和∑0|0=∑0，先驗(yàn)分布p（Xm|H0：m-1）在m時(shí)刻的均值和協(xié)方差定義為Xm|m-1和∑m|m-1，后驗(yàn)分布p（Xm|H0：m）在m時(shí)刻的均值和協(xié)方差定義為Xm|m和∑m|m［9］，其中先驗(yàn)分布和后驗(yàn)分布又被稱為誤差方差矩陣，可表示為P（m|m-1）。K（m）為濾波增益矩陣，對(duì)于m=1，…M，卡爾曼濾波可表示為時(shí)間更新和觀測(cè)更新的動(dòng)態(tài)方程［10］：

時(shí)間更新方程負(fù)責(zé)推進(jìn)時(shí)間軸進(jìn)而預(yù)測(cè)下一狀態(tài)的矢量值和誤差協(xié)方差估計(jì)值，為下一刻狀態(tài)構(gòu)建預(yù)評(píng)估模型；觀測(cè)方程對(duì)其矢量及協(xié)方差估計(jì)偏差結(jié)合構(gòu)造校驗(yàn)后的估計(jì)。根據(jù)時(shí)間和狀態(tài)更新方程預(yù)測(cè)四維航跡從m-1時(shí)刻到m時(shí)刻的狀態(tài)［10］。

3 滑動(dòng)時(shí)間窗口構(gòu)建

由于卡爾曼濾波本質(zhì)上是將非線性問(wèn)題線性處理，將會(huì)導(dǎo)致濾波單位精度的下降甚至分散。同時(shí)超前時(shí)間預(yù)測(cè)窗體過(guò)大會(huì)使預(yù)測(cè)精度降低，預(yù)測(cè)窗口過(guò)小航空器無(wú)法及時(shí)進(jìn)行沖突識(shí)別及后續(xù)采取避讓措施，因此根據(jù)飛機(jī)性能、航路環(huán)境建立自適應(yīng)滑動(dòng)窗體模型，對(duì)四維航跡進(jìn)行高精度預(yù)測(cè)，其滑動(dòng)模型如圖3所示。

圖3 自適應(yīng)時(shí)間窗體滑動(dòng)模型

滑動(dòng)時(shí)間窗體可定義為當(dāng)t為超前預(yù)測(cè)步數(shù)時(shí)，從當(dāng)前時(shí)刻m到超前時(shí)刻m+t之間的滑動(dòng)時(shí)間窗Tw=［Tm，Tm+t］［11］。通過(guò)找尋合適的預(yù)測(cè)時(shí)間窗體Tw即可求解超前預(yù)測(cè)步數(shù)t的過(guò)程，就移動(dòng)目標(biāo)與飛行沖突給出預(yù)測(cè)時(shí)間窗體的預(yù)測(cè)方法。

設(shè)ADS-B傳感器的預(yù)測(cè)上限為tmax，超前預(yù)測(cè)步數(shù)應(yīng)當(dāng)滿足t≤tmax。由卡爾曼可知加速度隨機(jī)誤差服從正態(tài)分布并且相互獨(dú)立，即εk～N（0，δ2），速度與加速度的估計(jì)值為［vxm，vym，axm，aym］，以y方向的矢量為例，根據(jù)對(duì)應(yīng)狀態(tài)運(yùn)動(dòng)方程可得出：

整理可得：

由正態(tài)分布的可累加性得：

設(shè)uy，m+t為vy，m+t的估計(jì)值，利用［uxm，uym，axm，aym］進(jìn)行航跡矢量信息預(yù)測(cè)，可知當(dāng)t＞1時(shí)（對(duì)X方向同樣適用），并可估算出航空器移動(dòng)的平均速率u：

令t=tmax，并利用概率分布產(chǎn)生的隨機(jī)矩陣即可得出滑動(dòng)時(shí)間窗體內(nèi)航空器的平均速率。同時(shí)采用自適應(yīng)調(diào)整噪聲協(xié)方差P（m|m-1）以減少外界干擾及信號(hào)丟失，引入自適應(yīng)矩陣Gm，表述如下：

式中Sm=diag｛s1，s2，s3…sm｝，m為對(duì)應(yīng)的航空器矢量維數(shù)，根據(jù)上述系統(tǒng)方程的定義及參數(shù)物理關(guān)系，可得狀態(tài)轉(zhuǎn)移的矩陣A（m|m-1）和觀測(cè)矩陣H（m）。

假定航空器的飛行高度滿足Hmin≤h≤hmax，根據(jù)航空器的機(jī)動(dòng)性能及環(huán)境限制，可構(gòu)造約束條件：

式中Δθ為航向改變角，由于航空器機(jī)動(dòng)過(guò)載限制設(shè)其最大轉(zhuǎn)角為Δθmax；tf為航空器相鄰航跡點(diǎn)的飛行時(shí)間；cmin為最短的航跡段長(zhǎng)度，即航空器在一個(gè)滑動(dòng)時(shí)間窗體內(nèi)的最小直線距離；cmax為最大航程。

4 仿真分析

以ZSNB-ZUUU航路為監(jiān)測(cè)對(duì)象，以2018年7月份某航班CES9937（A320機(jī)型）的ADS-B實(shí)測(cè)數(shù)據(jù)為例，選取寧波櫟社國(guó)際機(jī)場(chǎng)的ARP（機(jī)場(chǎng)基準(zhǔn)點(diǎn)）作為空間坐標(biāo)系的原點(diǎn)。部分A320數(shù)據(jù)如表1、表2所示。

表1 CES9937航行時(shí)間位置圖

表2 CES9937航行速度高度圖

在仿真實(shí)驗(yàn)中令采樣周期控制在35±5s區(qū)間內(nèi)，觀測(cè)噪聲R（n）=0.01×I4×4，Q（n）=I1×1，P（0|0）=8×I6×6，系統(tǒng)噪聲參數(shù)為B（n|n-1）=0.01×I6×1，ADS-B掃描的等時(shí)間間隔為2s，因此航空器的航跡由離散的航跡點(diǎn)組成。航班CES9937的離場(chǎng)航跡及巡航航程實(shí)測(cè)航跡如圖4所示。

圖4 CES9937航空器ZSNB離場(chǎng)圖

從圖4可知由于信號(hào)干擾、阻擋等因素會(huì)出現(xiàn)數(shù)據(jù)丟失問(wèn)題，圓圈部分表示數(shù)據(jù)丟失，按時(shí)間間隔2s進(jìn)行m元線性擬合。從圖5可知CES9937航空器的真實(shí)值和觀測(cè)值通過(guò)卡爾曼波曲線進(jìn)行偏差分析，可知傳統(tǒng)卡爾曼濾波預(yù)測(cè)方法對(duì)于四維航跡的預(yù)測(cè)誤差較大。

圖5 CES9937航空器ZSNB-ZUUU航路圖

采樣周期T=2s，利用自適應(yīng)滑動(dòng)窗體模型的卡爾曼波并與傳統(tǒng)的卡爾曼波的仿真結(jié)果進(jìn)行仿真，通過(guò)利用機(jī)載傳感器1，2對(duì)航空器運(yùn)動(dòng)測(cè)量的航跡及濾波預(yù)測(cè)結(jié)果進(jìn)行偏差分析，其偏差如圖6所示，同時(shí)對(duì)自適應(yīng)滑動(dòng)窗體模型與傳統(tǒng)卡爾曼波X軸和Y軸的預(yù)測(cè)誤差進(jìn)行對(duì)比［12］。

圖6 基于傳統(tǒng)算法的卡爾曼濾波

由圖6、圖7及圖8可知當(dāng)T=2s時(shí)，利用傳統(tǒng)卡爾曼算法進(jìn)行四維航跡預(yù)測(cè)時(shí)，X軸和Y軸的預(yù)測(cè)誤差P（m|m-1）控制在±27m范圍區(qū)間，而利用基于自適應(yīng)滑動(dòng)時(shí)間窗的卡爾曼濾波控制在±22m范圍區(qū)間，對(duì)應(yīng)誤差如表3所示。

表3 傳統(tǒng)卡爾曼與基于時(shí)間窗卡爾曼對(duì)比表

圖7 基于自適應(yīng)滑動(dòng)時(shí)間窗的卡爾曼濾波

圖8 傳統(tǒng)卡爾曼濾波與自適應(yīng)滑動(dòng)時(shí)間窗的卡爾曼濾波對(duì)比

按照國(guó)際民航組織ICAO對(duì)航行通信要求，這里選取±50m作為導(dǎo)航基準(zhǔn)源，區(qū)域概率和航跡的平滑性成正比，區(qū)域概率越大則誤差越大，反之相反。通過(guò)對(duì)上述數(shù)據(jù)進(jìn)行分析，可知傳統(tǒng)卡爾曼濾波的區(qū)域概率為0.2916，基于滑動(dòng)時(shí)間窗的卡爾曼濾波區(qū)域概率為0.1000，相比可知改進(jìn)后的四維航跡區(qū)域誤差減少65.8%，因此改進(jìn)后的濾波發(fā)散性減弱，航跡的平滑性增強(qiáng)［10～13］，對(duì)應(yīng)的濾波增益矩陣K（m）系數(shù)差異減小，按照實(shí)驗(yàn)數(shù)據(jù)對(duì)比為表4所示。

表4 傳統(tǒng)卡爾曼與基于時(shí)間窗卡爾曼區(qū)域概率對(duì)比

以傳統(tǒng)卡爾曼波區(qū)域概率0.2916為全誤差基準(zhǔn)，則改進(jìn)后的誤差概率縮減為原來(lái)的（0.2916-0.1000）/0.2916=0.658，對(duì)應(yīng)誤差縮減率具體如表5所示。

表5 傳統(tǒng)卡爾曼與基于時(shí)間窗卡爾曼區(qū)域誤差對(duì)比

由于航空器在實(shí)際飛行過(guò)程中，對(duì)地航向及速度收到航空器機(jī)動(dòng)性能的影響不斷變化，預(yù)測(cè)誤差往往會(huì)上下波動(dòng)偏置［14］，跟蹤誤差也必然存在波動(dòng)偏置。通過(guò)觀測(cè)上述圖中的濾波，改進(jìn)后的滑動(dòng)時(shí)間窗能有效減少εk，在受到外界干擾的情況下其濾波性能夠很快收斂，魯棒性較好［15］。

5 結(jié)語(yǔ)

本文提出基于自適應(yīng)滑動(dòng)時(shí)間窗的四維航跡預(yù)測(cè)模型，解決了航空器四維航跡實(shí)時(shí)規(guī)劃問(wèn)題。通過(guò)對(duì)時(shí)間特性賦予的航空器機(jī)動(dòng)進(jìn)行評(píng)估［16］，引入自適應(yīng)矩陣Gm保持卡爾曼濾波中的混合高斯模型適應(yīng)外界干擾的多模態(tài)特點(diǎn)［17］，并在此基礎(chǔ)上利用矢量因素及設(shè)定的誤差矩陣P（m|m-1）求解實(shí)時(shí)四維航跡［18～20］。仿真結(jié)果證明基于自適應(yīng)滑動(dòng)時(shí)間窗的四維航跡預(yù)測(cè)模型是有效的，利用機(jī)載傳感器及ADS-B數(shù)據(jù)能夠?qū)嵤┚珳?zhǔn)四維航跡預(yù)測(cè)。