章　濤，高　陽(yáng)，張成偉，吳仁彪

（1.中國(guó)民航大學(xué)天津市智能信號(hào)和圖像處理重點(diǎn)實(shí)驗(yàn)室，天津　300300；2.民航深圳空中交通管理站，深圳　518128）

基于KF聯(lián)合EKF參數(shù)辨識(shí)的短時(shí)4D航跡預(yù)測(cè)

章濤1，高陽(yáng)1，張成偉2，吳仁彪1

（1.中國(guó)民航大學(xué)天津市智能信號(hào)和圖像處理重點(diǎn)實(shí)驗(yàn)室，天津300300；2.民航深圳空中交通管理站，深圳518128）

提出了一種KF（Kalman filter）和EKF（efxtended Kalman filter）聯(lián)合算法辨識(shí)運(yùn)動(dòng)模型參數(shù)的4D航跡預(yù)測(cè)方法。該方法在等角航跡飛行模型的基礎(chǔ)上，運(yùn)用KF和EKF聯(lián)合算法辨識(shí)航空器的地速，以此計(jì)算航空器未來(lái)特征位置的過(guò)點(diǎn)時(shí)間。仿真實(shí)驗(yàn)結(jié)果表明，該方法采用航空器的經(jīng)緯度作為觀測(cè)變量分別更新的方法，可降低運(yùn)算復(fù)雜度，并能夠較精確地預(yù)測(cè)航空器等速巡航階段的短期飛行航跡。

卡爾曼濾波；擴(kuò)展卡爾曼濾波；參數(shù)辨識(shí)；等角航跡；4D航跡預(yù)測(cè)

隨著航空運(yùn)輸業(yè)的發(fā)展，傳統(tǒng)的飛行計(jì)劃結(jié)合空管動(dòng)態(tài)調(diào)配的空中交通管理方式已不能滿足空中交通的發(fā)展，空域資源有限和空域交通需求的矛盾日益突出[1-2]。因此，美國(guó)提出構(gòu)建下一代空中交通運(yùn)輸管理系統(tǒng)。高精度的4D航跡預(yù)測(cè)是下一代空中交通運(yùn)輸管理系統(tǒng)的技術(shù)基礎(chǔ)[3]。4D航跡通過(guò)控制航空器到達(dá)航路點(diǎn)的“時(shí)間窗”對(duì)航空器的位置進(jìn)行精確描述，從而有效適應(yīng)未來(lái)高密度、大流量、小間隔的空域情況。

基于4D航跡預(yù)測(cè)在下一代空中交通運(yùn)輸管理系統(tǒng)的重要性[4]，國(guó)內(nèi)外已對(duì)此展開(kāi)大量研究。4D航跡預(yù)測(cè)研究算法可以分為2種：①基于KF或數(shù)據(jù)挖掘等無(wú)參數(shù)的預(yù)測(cè)方法[5-7]。無(wú)參數(shù)的4D航跡預(yù)測(cè)模型簡(jiǎn)單、適用范圍廣，但可擴(kuò)展性差，不能很好地適應(yīng)管制、氣象環(huán)境變化等突發(fā)情況的干擾。②基于航空器動(dòng)力學(xué)或運(yùn)動(dòng)學(xué)模型的預(yù)測(cè)方法[8-10]。通過(guò)對(duì)航空器每個(gè)飛行階段進(jìn)行空氣運(yùn)動(dòng)學(xué)模型建模，推測(cè)航空器飛行環(huán)境（如風(fēng)速、風(fēng)向和溫度變化等氣象因素）對(duì)航跡的影響。但由于機(jī)型、載重等差異帶來(lái)的航空器飛行性能的差異，導(dǎo)致預(yù)測(cè)精度不理想?；跀?shù)據(jù)挖掘的無(wú)參數(shù)預(yù)測(cè)方法和基于航空器動(dòng)力學(xué)模型的預(yù)測(cè)方法多側(cè)重全程航跡的預(yù)測(cè)。但飛行器的飛行過(guò)程受到管制、氣象因素、空域流量和航空公司經(jīng)濟(jì)性等影響，往往會(huì)偏離飛行剖面，因此有必要對(duì)航跡進(jìn)行實(shí)時(shí)預(yù)測(cè)和修正。徐琴等[11]提出了基于參數(shù)辨識(shí)的短期4D航跡預(yù)測(cè)，利用實(shí)時(shí)航跡信息，采用極大似然準(zhǔn)則和牛頓-拉夫遜迭代算法對(duì)運(yùn)動(dòng)模型中的校正空速進(jìn)行辨識(shí)，進(jìn)一步推測(cè)航空器的過(guò)點(diǎn)時(shí)間[11]。然而，牛頓-拉夫遜算法迭代次數(shù)多、耗時(shí)長(zhǎng)，還會(huì)出現(xiàn)似然函數(shù)發(fā)散的情況。

本文提出了一種基于KF聯(lián)合EKF的參數(shù)辨識(shí)短時(shí)4D航跡預(yù)測(cè)方法，該方法在飛行器運(yùn)動(dòng)模型的基礎(chǔ)上，建立狀態(tài)方程和觀測(cè)方程，運(yùn)用實(shí)時(shí)航跡位置和時(shí)間數(shù)據(jù)依次采用KF和EKF濾波方法進(jìn)行狀態(tài)更新[12]，并結(jié)合等角航跡模型預(yù)測(cè)過(guò)點(diǎn)時(shí)間。

1　航跡預(yù)測(cè)模型及其實(shí)現(xiàn)

等角航線是地球表面與經(jīng)線相交成相同角度的曲線。地球表面除經(jīng)線和緯線以外的等角航線都是以極點(diǎn)為漸近點(diǎn)的螺旋曲線。在實(shí)地觀測(cè)中，航線與經(jīng)線的夾角易于測(cè)得；等角航線在航行圖（采用墨卡托投影）上表現(xiàn)為直線。因此，航空器在短距離飛行時(shí)多采用等角航跡的計(jì)算模式，建立等角航跡的運(yùn)動(dòng)模型，在該運(yùn)動(dòng)模型的基礎(chǔ)上選擇合適的濾波器狀態(tài)方程和觀測(cè)方程，然后根據(jù)觀測(cè)點(diǎn)的狀態(tài)信息（經(jīng)度、緯度、時(shí)間等）估計(jì)辨識(shí)參數(shù)，根據(jù)運(yùn)動(dòng)模型預(yù)測(cè)觀測(cè)點(diǎn)的過(guò)點(diǎn)時(shí)間，實(shí)現(xiàn)短期航跡推測(cè)。

1.1等角航跡飛行模型

等角航跡航線角和距離的計(jì)算公式由積分方法推導(dǎo)。假設(shè)一段航跡：α（φ1，λ1）為起始點(diǎn)，β（φ2，λ2）為終止點(diǎn)，其中φ為緯度，λ為經(jīng)度。設(shè)航線與經(jīng)度線的夾角為θ，距離為s，如果在該段航線內(nèi)任取一小段航段ds，則該航段與相應(yīng)的緯度、經(jīng)度（dφ，dλcosφ）之間構(gòu)成的三角形可看作平面直角三角形，如圖1所示。由此可得

圖1　等角航跡模型Fig.1　Isometric track model

或

航空器等校正空速巡航過(guò)程中τ時(shí)間的等角航線距離為

如果已知航空器當(dāng)前位置A（φ0，λ0）點(diǎn)的經(jīng)緯度、速度等，則航空器τ時(shí)間后的位置為

或

計(jì)算過(guò)程中，各個(gè)變量都轉(zhuǎn)換成弧度制單位，等角航線距離（sm）與弧度制換算關(guān)系為1m=π（/180×111120）。

1.2參數(shù)辨識(shí)模型與航跡預(yù)測(cè)

通過(guò)雷達(dá)或ADS-B接收機(jī)可獲取航空器實(shí)時(shí)航跡信息，包括航班號(hào)、飛機(jī)類型、飛行狀態(tài)、經(jīng)緯度、高度、時(shí)間等。在等速巡航階段，航空器的地速變化主要受氣象因素（溫度、風(fēng)速等）影響，會(huì)產(chǎn)生一定范圍內(nèi)的波動(dòng)[13]。本文采用KF和EKF聯(lián)合的方法對(duì)航空器地速進(jìn)行參數(shù)辨識(shí)，估計(jì)飛行器短時(shí)間內(nèi)準(zhǔn)確的飛行狀態(tài)變量。

航空器地速為vGS，則航空器在k時(shí)刻的運(yùn)動(dòng)狀態(tài)量為vG（Sk），短時(shí)間內(nèi)航空器所在空域的氣象條件保持相對(duì)穩(wěn)定，對(duì)航空器地速的影響也是固定的，地速僅存在隨機(jī)擾動(dòng)。因此，可得到遞推方程，以及描述時(shí)間間隔為τ的下一時(shí)刻的運(yùn)動(dòng)狀態(tài)方程為

觀測(cè)狀態(tài)方程為

觀測(cè)矢量包含2個(gè)觀測(cè)變量，并且第1個(gè)觀測(cè)量和狀態(tài)變量為線性關(guān)系，第2個(gè)觀測(cè)量和狀態(tài)變量為非線性關(guān)系。本文采用緯度KF更新然后經(jīng)度EKF更新的聯(lián)合觀測(cè)更新方法。緯度觀測(cè)變量KF更新為

其中

經(jīng)度觀測(cè)變量EKF更新為

其中

每個(gè)觀測(cè)點(diǎn)可獲得2個(gè)狀態(tài)變量的真實(shí)值，k時(shí)刻的真實(shí)值為

設(shè)定在一個(gè)觀測(cè)位置點(diǎn)的這兩個(gè)真實(shí)值權(quán)重相同，則

最后，將辨識(shí)到的航空器地速結(jié)合等角航跡的飛行模型，可以推測(cè)未來(lái)時(shí)刻航空器的位置或航路上的過(guò)點(diǎn)時(shí)間。

2　實(shí)驗(yàn)仿真與分析

本文實(shí)驗(yàn)使用ADS-B接收機(jī)采集數(shù)據(jù)，采集時(shí)間為2013年10月25日，航班由石家莊正定國(guó)際機(jī)場(chǎng)飛往大連周水子機(jī)場(chǎng)，航班號(hào)為9C8873，機(jī)型為波音A320，選取巡航階段的部分航跡數(shù)據(jù)參與計(jì)算，巡航高度大約為8 860 km，升降率為0，經(jīng)緯度、速度等部分?jǐn)?shù)據(jù)如表1所示。

表1　觀測(cè)點(diǎn)數(shù)據(jù)Tab.1　Observation Data

取初始化地速為100 m/s，通過(guò)KF和EKF對(duì)地速進(jìn)行參數(shù)辨識(shí)，最后得到辨識(shí)地速，結(jié)果如圖2所示。取數(shù)據(jù)中的25個(gè)時(shí)間點(diǎn)為觀測(cè)點(diǎn)，觀測(cè)點(diǎn)的部分經(jīng)緯度位置數(shù)據(jù)如表2所示，用辨識(shí)得到的地速預(yù)測(cè)過(guò)點(diǎn)時(shí)間，預(yù)測(cè)時(shí)間如圖3所示。

圖2　地速辨識(shí)結(jié)果圖Fig.2　Results of identified ground speed

表2　預(yù)測(cè)點(diǎn)數(shù)據(jù)Tab.2 Predicted data

圖3　預(yù)測(cè)過(guò)點(diǎn)時(shí)間和實(shí)際過(guò)點(diǎn)時(shí)間對(duì)比Fig.3　Comparison between predicted time and actual time

由圖3可知，該方法可有效地預(yù)測(cè)等速巡航階段的過(guò)點(diǎn)時(shí)間，預(yù)測(cè)過(guò)點(diǎn)時(shí)間在實(shí)際過(guò)點(diǎn)時(shí)間上下波動(dòng)，比較貼合實(shí)際過(guò)點(diǎn)時(shí)間。預(yù)測(cè)誤差隨時(shí)間的推移會(huì)增大，5 min航程的預(yù)測(cè)誤差可控制在50 s以內(nèi)。因此，為了保證航跡預(yù)測(cè)的準(zhǔn)確性，可根據(jù)實(shí)時(shí)航跡和飛行狀態(tài)對(duì)航跡進(jìn)行預(yù)測(cè)和修正。本文算法和傳統(tǒng)KF算法的復(fù)乘次數(shù)對(duì)比如表3所示，其中N表示迭代更新的次數(shù)。本文方法避免了矩陣求逆等運(yùn)算，減少了復(fù)乘次數(shù)，降低了運(yùn)算復(fù)雜度。由于算例采用緯度經(jīng)度分別進(jìn)行觀測(cè)更新的方法，運(yùn)算復(fù)雜度低、模型簡(jiǎn)單、動(dòng)態(tài)適應(yīng)性強(qiáng)，使得該方法適應(yīng)性強(qiáng)、適用范圍廣。

表3　運(yùn)算量分析Tab.3　Com putation analysis

3　結(jié)語(yǔ)

本文研究了一種KF和EKF聯(lián)合的方法用于辨識(shí)航空器飛行地速，并結(jié)合等角航跡運(yùn)動(dòng)模型進(jìn)行短時(shí)4D航跡預(yù)測(cè)。該方法采用觀測(cè)量分別更新的方法，避免了矩陣運(yùn)算、減少了運(yùn)算復(fù)雜度。本文暫未考慮更多的氣象、管制因素等復(fù)雜狀況，后期將進(jìn)一步完善巡航階段航跡預(yù)測(cè)方法，使該模型可適應(yīng)更為復(fù)雜的巡航飛行；結(jié)合航空器進(jìn)離港的爬升和下降階段的特點(diǎn)建立預(yù)測(cè)模型，實(shí)現(xiàn)對(duì)整個(gè)飛行過(guò)程的航跡預(yù)測(cè)和修正。

[1]SOLOMOS G,URLASS S,BHADRA D,et al.Capacity Needs in the National Airspace System:An Analysis of Airport and Metropolitan Area Demand and Operational Capacity in the Future[R].The MITRE Center for Advanced Aviation System Development,2004.

[2]盧奕羽.初始飛行計(jì)劃4D航跡生成及仿真系統(tǒng)研究[D].南京:南京航空航天大學(xué),2014.

[3]伊群.美國(guó)新一代空中交通管理系統(tǒng)運(yùn)行概念[J].中國(guó)民用航空, 2007,80(8):27-31.

[4]呂小平.中國(guó)民航新一代空中交通管理系統(tǒng)發(fā)展總體框架[J].中國(guó)民用航空,2007,80(8):24-26.

[5]PREVOST C G,DESBIENS A,GAGNON E.Extended Kalman Filter for State Estimation and Trajectory Prediction of a Moving Object Detected by an Unmanned Aerial Vehicle[C]//AmericanControlConference, New York,IEEE,2007:1805-1810.

[6]LYMPEROPOULOS I,LYGEROS J.Sequential Monte Carlo methods for multi-aircraft trajectory prediction in air traffic management[J].International Journal of Adaptive Control and Signal Processing,2010,24 (10):830-849.

[7]吳鹍,潘薇.基于數(shù)據(jù)挖掘的四維飛行軌跡預(yù)測(cè)模型[J].計(jì)算機(jī)應(yīng)用,2007,27(11):2637-2639.

[8]SLATTERY R,ZHAO Y.Trajectory synthesis for air traffic automation[J].Journal of Guidance,Control and Dynamics,1997,20(2):232-238.

[9]王超,郭九霞,沈志鵬.基于基本飛行模型的4D航跡預(yù)測(cè)方法[J].西南交通大學(xué)學(xué)報(bào),2009,44(2):295-300.

[10]CHESTER G,WILLIAM N C.Using Flight Manual Data to Derive Aero-Propulsive Models for Predicting Aircraft Trajectories[C]//AIAA Aircraft Technology,Integration and Operations(ATIO)2002 Technical.California:AIAA,2002:1-7.

[11]徐琴,湯新民,韓松臣,等.基于參數(shù)辨識(shí)的短期4D航跡預(yù)測(cè)[J].信息與控制,2014,43(4):501-505.

[12]楊宏,李亞安,李國(guó)輝.一種改進(jìn)擴(kuò)展卡爾曼濾波新方法[J].計(jì)算機(jī)工程與應(yīng)用,2010,46(19):18-20.

[13]吳曉光,張軍峰,蔣海行.GRIB數(shù)據(jù)及其在航跡預(yù)測(cè)中的應(yīng)用[J].航空計(jì)算技術(shù),2013,43(6):59-62.

（責(zé)任編輯：黨亞茹）

Short-term 4D trajectory prediction based on KF joint EKF parameter identification

ZHANG Tao1,GAO Yang1,ZHANG Chengwei2,WU Renbiao1
（1.Intelligent Signal and Image Processing Key Lab of Tianjin,CAUC,Tianjin 300300,China; 2.Shenzhen Air Traffic Management Station,Shenzhen518128,China）

A method based on KF-EKF jointed algorithm is proposed to identify motion model parameters for 4D trajectory prediction.On the basis of isometric track flight model,KF-EKF jointed algorithm is employed to identify the ground speed of aircraft,then the arriving time of aircraft’s scheduled position is calculated.Latitude and longitude of aircraft are taken as observed variables and are updated respectively.Simulation result shows that it could reduce complexity of the computation and accurately predict the short-term flight track of constant speed cruise phase.

Kalman filter;extended Kalman filter;parameter identification;isometric track;4D trajectory prediction

TN959.4

A

1674-5590（2016）05-0001-04

2015-11-10；

2016-01-07基金項(xiàng)目：國(guó)家自然科學(xué)基金項(xiàng)目（F011206）

章濤（1989—），男，河北衡水人，碩士研究生，研究方向?yàn)?D航跡預(yù)測(cè).