航空器軌跡變點檢測與判識技術(shù)研究

蘇志剛，李志強，郝敬堂

（中國民航大學(xué)a.天津市智能信號與圖像處理重點實驗室；b.中歐航空工程師學(xué)院，天津 300300）

航空器軌跡變點檢測與判識技術(shù)研究

蘇志剛a，b，李志強a，郝敬堂b

（中國民航大學(xué)a.天津市智能信號與圖像處理重點實驗室；b.中歐航空工程師學(xué)院，天津 300300）

航空器飛行階段的有效劃分是支撐航空排放主動監(jiān)測的重要技術(shù)手段，采用變點檢測與判識的方法可有效地實現(xiàn)對航空器飛行階段的劃分。本研究根據(jù)多變量T2統(tǒng)計量的特點，在數(shù)據(jù)整體參數(shù)一致性假設(shè)下，提出了一種采用順序雙滑窗的航空器軌跡變點檢測與判識的方法。仿真與實驗結(jié)果表明，該方法對窗口長度較為敏感，在適宜的窗口長度下可有效地解決軌跡的變點檢測問題。根據(jù)航空器軌跡的特點，選擇適當?shù)拇翱陂L度前提下，本方法是一種有效的變點檢測與判識方法。

航空排放；排放監(jiān)測；航跡分割；變點檢測；T2統(tǒng)計量

國際民航組織（ICAO，international civil aviation organization）預(yù)測2050年航空排放將比2010年增長300%左右[1]，因此，在國際氣候變化談判中，民用航空的節(jié)能減排成為其前沿問題和各國斗爭的焦點之一。中國民航2016年運輸總周轉(zhuǎn)量比上年同期增長12.8%，隨著運輸量的增加，中國民航也面臨嚴峻的航空節(jié)能減排壓力，推進民航節(jié)能減排成為近年來中國民航管理部門的重點工作。民航節(jié)能減排需要有效的主動核查與監(jiān)測手段。根據(jù)航空器不同飛行階段的能耗模型，結(jié)合航空器軌跡信息，可以實現(xiàn)對航空器排放的主動監(jiān)測。因此，基于航空器軌跡信息的飛行階段劃分技術(shù)是實現(xiàn)航空排放主動監(jiān)測任務(wù)的核心技術(shù)之一。

在民航監(jiān)視領(lǐng)域，二次雷達系統(tǒng)和廣播式自動相關(guān)監(jiān)視（ADS-B，automatic dependent surveillance-broadcast）系統(tǒng)可獲得航空器呼號、經(jīng)度、緯度、高度、地速、航向和升降率等信息。利用與航空器軌跡相關(guān)的高度、地速和升降率等信息可進行飛行階段的劃分。目前航空器飛行階段的劃分主要涉及對特定飛行階段的提取和對全航程飛行階段的劃分兩類。

對特定飛行階段的提取主要是用于研究在特定飛行階段，如爬升、巡航和下降等，航空器所涉及的飛行動作和姿態(tài)的變化、燃油消耗[2]和噪聲[3]等問題，相應(yīng)飛行階段的提取通常采用變點分析方法[4]、遺傳算法[5]和隱馬爾可夫模型[6]等方法實現(xiàn)。對全航程飛行階段的劃分主要采用支持向量機方法[7-8]，該方法依賴人工獲取的標簽，僅是對歷史數(shù)據(jù)的回顧性分析。單窗口Hotelling’s T2統(tǒng)計量方法可對航跡中的變點進行檢測[9-10]，現(xiàn)有方法重點在變點檢測上，一旦檢測出變點就停止，無法實現(xiàn)連續(xù)地變點檢測。雙窗口距離測度方法通過比較兩個窗口內(nèi)數(shù)據(jù)的差異性來完成變點檢測[11]，但該方法未能有效解決數(shù)據(jù)更新問題，也無法實現(xiàn)對實時航跡的變點檢測。因此，現(xiàn)有對特定飛行階段的提取和對全航程飛行階段劃分的研究成果均不具備對航跡數(shù)據(jù)實時處理的能力，無法為航空排放的主動監(jiān)測提供實時的階段劃分。

本文為解決航空器飛行階段實時劃分問題，根據(jù)飛行參數(shù)特點，提出軌跡的變點檢測方法。該方法利用民航監(jiān)視系統(tǒng)提供的航空器高度、地速及升降率等信息，采用順序雙滑窗的方式對航跡數(shù)據(jù)進行切割，在一致性假設(shè)條件下計算雙窗口數(shù)據(jù)的多變量Hotelling’s T2統(tǒng)計量，通過 Hotelling’s T2統(tǒng)計量是否超過門限值來檢測變點，并采用最大似然估計（MLE，maximum likelihood estimation）方法確定變點的位置。最后對仿真數(shù)據(jù)和B737-800航線數(shù)據(jù)進行了實驗分析。

1 問題描述

航空器航跡信息中的高度、地速和升降率可用于描述飛行階段，設(shè)航空器航跡在第n個位置點的測量矢量為

其中：hn、vn和rn分別表示航空器在第n個位置點的高度、地速與升降率的測量值。由于測量矢量受到諸多因素的影響，可以假設(shè)其服從多維高斯分布。序列{xn′}反映了航空器的航跡，可利用其進行相應(yīng)的變點檢測。

由于序列{xn′}具有時間相關(guān)性，對于采用獨立樣本進行變點檢測的方法不利，因此需采取差分或?qū)?shù)化的方式消除樣本間的時間相關(guān)性，先對原始數(shù)據(jù)一階差分處理，即

本研究將根據(jù)差分序列{xn}進行變點檢測與變點位置估計。

航空器航跡變點是順序兩個飛行階段的轉(zhuǎn)換點，即變點前后的樣本在統(tǒng)計特征上存在差異，利用此差異能夠進行變點存在與否的檢測。

基于上述考慮，可以順序地建立兩個數(shù)據(jù)窗口，相應(yīng)的數(shù)據(jù)矩陣分別為

其中：n1和n2分別為兩個窗口的數(shù)據(jù)起始位置；l1和l2分別為兩個窗口的數(shù)據(jù)長度。由于建立的是順序滑動窗口，所以，n2=n1+l1。設(shè)總體數(shù)據(jù)的協(xié)方差矩陣為C，數(shù)據(jù)矩陣W1的均值矢量與協(xié)方差矩陣分別為μ1和C1，數(shù)據(jù)矩陣W2的均值矢量與協(xié)方差矩陣分別為μ2和C2，相應(yīng)的估計分別為

航空器航跡變點檢測的核心問題就是判斷兩窗口內(nèi)數(shù)據(jù)W1和W2的均值矢量是否發(fā)生顯著變化。

2 變點檢測模型

根據(jù)兩個窗口數(shù)據(jù)W1和W2，在一定顯著度水平α下構(gòu)建假設(shè)檢驗

式（9）假設(shè)檢驗問題的原假設(shè)即為所有變量的均值在兩個窗口內(nèi)都相等，備擇假設(shè)為至少有一個變量的均值不相等。在原假設(shè)H0條件下，相應(yīng)的Hotelling’s T2統(tǒng)計量為

其中，符號Cov（·）表示協(xié)方差運算。在樣本獨立條件下

若兩窗口數(shù)據(jù)W1和W2協(xié)方差矩陣和總體的協(xié)方差矩陣相等，即C1=C2=C，則相應(yīng)的合并協(xié)方差矩陣為

顯然，式（13）只是在窗口數(shù)據(jù)W1和W2具有相同協(xié)方差矩陣條件下的特例。

式（10）中的Hotelling’s T2統(tǒng)計量的分布為[12]

其中：d=3是式（1）中數(shù)據(jù)矢量的維數(shù)。

在原假設(shè)H0條件下，由顯著度水平α確定的門限為c，則當T2＞c時，拒絕原假設(shè)H0，判定備擇假設(shè)H1成立，即至少有一個參數(shù)的均值在窗口數(shù)據(jù)W1和W2中發(fā)生了顯著變化，則判決有變點發(fā)生，其中窗口數(shù)據(jù)W2中可能存在變點。

3 變點位置估計

當備擇假設(shè)H1成立時，窗口數(shù)據(jù)W1和W2屬于兩個不同均值的總體，即數(shù)據(jù)W2中可能存在變點。假設(shè)變點發(fā)生在 W2中的第 k（k=1，…，l2）個數(shù)據(jù)點，則

其中：μ21、C21和 μ22、C22分別為變點前后的序列均值矢量與協(xié)方差矩陣。

考慮到變點有可能發(fā)生在k=1位置，因此在數(shù)據(jù)W2前補充2d個數(shù)據(jù)，即[xn2-2d，…，xn2-1]。參考式（10）可得

當k處于真實變點位置時，k點前后數(shù)據(jù)間的差異性最大，即式（16）中的T2（k）取最大值，因此

4 數(shù)據(jù)窗口長度選擇

兩窗口數(shù)據(jù)長度l1與l2的大小會影響Hotelling’s T2統(tǒng)計量的變點檢測性能。較長窗口有利于保持數(shù)據(jù)統(tǒng)計特征的穩(wěn)定性，但對于相距較近的變點檢測造成困難，可能會導(dǎo)致漏檢。較短窗口能較好地刻畫航跡參數(shù)變化特征，但易受到隨機因素的影響，虛警率較高。在一般的變點檢測問題中，為了滿足樣本有效性，數(shù)據(jù)窗口的長度至少要滿足l1＞2d，l2＞2d，其中d是參數(shù)的維度。在航空器航跡變點檢測中，窗口長度的選擇需充分考慮航空器航跡變化的特點，盡可能地選擇較長的窗口，本文主要以航空器階段持續(xù)時間和ADS-B接收機數(shù)據(jù)更新時間來選擇窗口長度。

航空器的起、降階段主要在相對地面高度（AGL，altitude above ground level）3 000 ft（1 ft=0.304 8 m，下同）以下空域。在此空域內(nèi)，航空器主要涉及慢車、起飛、爬升和進近4個狀態(tài)，ICAO對各狀態(tài)持續(xù)時間的建議參數(shù)如表1所示。由表1可見，在AGL 3 000 ft以下空域內(nèi)，狀態(tài)持續(xù)時間最短的起飛過程也需0.7 min，即42 s。

表1 起降階段各狀態(tài)持續(xù)時間Tab.1 Duration of each stage of taking-off and landing

在AGL 3 000 ft以上空域，航空器逐漸爬升進入巡航階段。在巡航階段時，航空器會因沖突、天氣等原因進行高度層的調(diào)整。根據(jù)中國民航高度層劃分規(guī)則，同方向高度層間的最小垂直間隔為2 000 ft，以在30 000 ft高空巡航的總重50 t的B737為例，當進行高度層遷徙時，升降率為2 000 ft/min，則該航空器至少需要1 min完成鄰近高度層間的變遷。所以，航空器的狀態(tài)持續(xù)時間均不小于0.5 min。

ADS-B地面站所提供的航跡信息為每0.5 s更新一次，即在0.5 min內(nèi)可提供60組數(shù)據(jù)，因此窗口的長度設(shè)置為

在每次變點檢測后，后續(xù)的檢測將以本次檢測變點位置為起始重新建立窗口。

5 實驗結(jié)果

本部分將分別利用仿真數(shù)據(jù)及實際記錄的航跡數(shù)據(jù)驗證本文提出方法的有效性。考慮仿真數(shù)據(jù)具有相同協(xié)方差矩陣

生成三維正態(tài)數(shù)據(jù)矢量1 000條，其中數(shù)據(jù)矢量的3個分量AD_hn、AD_vn和AD_rn分別模擬航空器爬升過程的飛行高度、地速和升降率等參數(shù)的一階差分。假設(shè)在仿真數(shù)據(jù)中設(shè)置6個變點，相應(yīng)的位置分別為第 301、501、601、801、831 和 881 時刻，其中變點最小間隔為30。因此，相應(yīng)的變點將仿真數(shù)據(jù)分割為7段，每段數(shù)據(jù)的均值情況如表2所示。

表2 7段仿真數(shù)據(jù)的起始間隔及其對應(yīng)仿真變量Tab.2 Mean vectors of 7 segments simulation data

首先研究顯著度水平對變點檢測性能的影響。設(shè)置兩個順序滑動窗口長度l1=l2=30，仿真數(shù)據(jù)的Hotelling’s T2統(tǒng)計量曲線如圖 1（a）所示，統(tǒng)計量的幾個峰值位置準確地標識出變點的位置。將圖1（a）中虛線框內(nèi)區(qū)域放大后，如圖1（b）所示，隨著顯著度水平α的減小，Hotelling’s T2統(tǒng)計量的檢測門限也相應(yīng)提高，從而降低了虛警概率；當α＞10-3時，能準確檢測出所有變點，但存在虛警；當α≤10-6時，可以有效地抑制虛警的發(fā)生，但可能存在漏檢。因此顯著度水平取10-5≤α≤10-3。

其次研究窗口長度對變點檢測性能的影響。設(shè)置等長的兩滑動窗口長度（l1=l2）為不同值時，Hotelling’s T2統(tǒng)計量曲線如圖2所示。由圖2可見，隨著窗口長度的增大，Hotelling’s T2統(tǒng)計量在數(shù)值上增大，且峰的寬度加大，使得變點的分辨能力降低，可能造成變點的漏檢。窗口長度減小，Hotelling’s T2統(tǒng)計量在數(shù)值上減小，且峰的寬度變窄，分辨變點的能力提高，但由于在峰值數(shù)值上的減小，可能會導(dǎo)致無法超過門限，也會產(chǎn)生漏檢。因此，窗口長度的選擇具有依賴于實際數(shù)據(jù)的特點。

對仿真數(shù)據(jù)變點判識的結(jié)果如圖3所示。由于設(shè)計仿真數(shù)據(jù)的變點最小間隔為30，因此，設(shè)置窗口長度為l1=l2=30，顯著度水平α=10-3，判識出的變點位置由圖中虛線表示。由圖3可見，本文所提出的方法準確地檢測出6個變點，并估計出相應(yīng)的位置。

圖1 檢測門限和仿真數(shù)據(jù)統(tǒng)計量曲線Fig.1 Thresholds and statistic curves of simulation data

圖2 不同窗口長度的統(tǒng)計量曲線Fig.2 Statistic curves for different window sizes

下面利用由ADS-B接收機獲取的真實航跡數(shù)據(jù)進行變點檢測與估計性能分析。該數(shù)據(jù)描述航空器B737-800由天津—昆明的航跡信息，選取巡航階段中約25 min的航跡數(shù)據(jù)，相應(yīng)數(shù)據(jù)曲線如圖4所示。圖4（a）是高度、地速與升降率等參數(shù)的實際變化原始數(shù)據(jù)，航空器高度從約22 000 ft經(jīng)過兩次調(diào)整爬升到32 000 ft，在此期間地速保持在420 kts左右，升降率存在兩次明顯的波動。圖4（b）是相應(yīng)的差分結(jié)果。由圖4（b）可見，高度差分數(shù)據(jù)有明顯的規(guī)律，能直觀地反映出變點的位置。當有變點發(fā)生時，地速差分數(shù)據(jù)波動的頻率顯著增高，升降率差分數(shù)據(jù)波動的幅度明顯增大。

圖3 仿真數(shù)據(jù)的變點判識Fig.3 Change points recognition of simulation data

圖4 巡航段部分數(shù)據(jù)的參數(shù)變化和一階差分曲線Fig.4 Curves of original data and first-order differential data of cruise phase

順序首尾相接兩窗口沿航跡數(shù)據(jù)滑動，得到Hotelling’s T2統(tǒng)計量的變化曲線，如圖 5（a）所示。Hotelling’s T2統(tǒng)計量曲線的4個峰值說明有4個顯著的變點發(fā)生，表明了數(shù)據(jù)整體參數(shù)存在突變。由圖5（a）可見，當有變點發(fā)生時，統(tǒng)計量曲線便迅速攀升。而在240和1 300附近，出現(xiàn)了較低的峰值，相應(yīng)的峰值反映了地速和升降率的突變。圖5（a）中Hotelling’s T2統(tǒng)計量的最大峰值超過了1 000，由圖4（b）中的一階差分數(shù)據(jù)可以看出，高度數(shù)據(jù)的均值跳變幅度最大約為40 ft，地速的約為2 kts，升降率的約為5 ft/s，由于高度數(shù)據(jù)的均值存在較大的跳變幅度，這就導(dǎo)致了Hotelling’s T2統(tǒng)計量的顯著增大。當顯著度水平設(shè)置為α=10-3時，變點檢測結(jié)果如圖5（b）所示。圖中共檢測出7個變點，其中4個變點反映了整體參數(shù)的突變，3個變點反映了地速或升降率的變化。

圖5 巡航段航空數(shù)據(jù)統(tǒng)計量曲線與變點判識（l1=l2=60）Fig.5 Statistic curve and segmentations of cruise phase（l1=l2=60）

當兩窗口的長度設(shè)置為l1=l2=60時，顯著度水平與變點判識情況如表3所示。由于檢測到變點后，在數(shù)據(jù)W2前補充2d個數(shù)據(jù)，所以變點的位置可能大于60。由表3的統(tǒng)計結(jié)果可見，當顯著度水平α較大時，門限值較低，檢測出的變點數(shù)量較多，發(fā)生虛警的概率較大；當α≤10-3時，當前數(shù)據(jù)段變點的數(shù)量不再發(fā)生變化，但估計的變點位置會有差別，這是由于檢測門限升高，檢測到變點的時刻發(fā)生延時引起的；當α≤10-5時，檢測到的變點位置不再發(fā)生大的變化。因此，對于航跡數(shù)據(jù)的變點檢測的顯著度水平取α=10-5。

表3 不同顯著度水平α下檢測到的變點位置分布Tab.3 Distribution of detected change points with different significant levels of α

最后對全航程的各個飛行階段進行變點檢測與判識。采用上述數(shù)據(jù)集的全部數(shù)據(jù)，共11 686條，兩窗口長度設(shè)置為l1=l2=60，顯著度水平α=10-5時，共檢測出30個變點，如圖6所示。由圖6可見，檢測出的變點主要反映了航空器整體參數(shù)的顯著變化，但也存在虛警，如第1 140點處的結(jié)果。

圖6 航空器一次完整飛行的變點判識Fig.6 Change point detection and segmentation of a complete flight of an aircraft

6 結(jié)語

根據(jù)航空器軌跡參數(shù)變化特點，將航空器軌跡的實時劃分問題轉(zhuǎn)化為實時變點檢測問題，利用順序雙滑動窗口截取航跡數(shù)據(jù)，在整體參數(shù)協(xié)方差一致性的假設(shè)前提下，引入多變量Hotelling’s T2統(tǒng)計量，通過判識統(tǒng)計量的大小實現(xiàn)相應(yīng)的變點檢測。實驗結(jié)果表明，結(jié)合航空器軌跡的具體特點，選擇適度的窗口長度和顯著度指標，可以有效地解決航空器軌跡的實時劃分問題。

[1]ICAO.ICAO Environ-Mental Report 2016[R].Montreal:Environment BrandoftheICAO,2016.

[2]陳靜杰,肖冠平.飛機油耗分析工具設(shè)計[J].計算機工程與設(shè)計,2014,35(11):4012-4016.

[3]張召悅,王 超,孟 娜,等.基于航跡分段的飛行程序噪聲評估方法[J].科學(xué)技術(shù)與工程,2014,14(3):116-120.

[4]ECKSTEIN A.Automated Flight Track Taxonomy for Measuring Benefits from Performance Based Navigation[C]//2009 Integrated Communications,Navigation and Surveillance Conference.Arlington:IEEE Computer Society Press,2009:1-12.

[5]尹文軍,倪世宏.一種基于遺傳算法的飛行動作識別方法[J].電腦開發(fā)與應(yīng)用,2006,19(4):50-51.

[6]王潔寧,孫 禾,趙元棣.面向終端區(qū)航空器飛行狀態(tài)識別的HMM方法[J].航空計算技術(shù),2013,43(2):1-5.

[7]LU H,MAO K F.Research on Stage Classification of Flight Parameter Based on PTSVM[C]//13th IEEE International Conference on Computational Science and Engineering.Hong Kong:IEEE,2010:55-63.

[8]李丹丹,路 輝,郎榮玲.基于主成分分析和支持向量機的飛參階段劃分研究[J].現(xiàn)代電子技術(shù),2010,33(8):134-137.

[9]HOTELLING H.Multivariate Quality Control—Illustrated by the Air Testing of Bombsights,in Techniques of Statistical Analysis[M].New York:McGraw-Hill,1947:110-122.

[10]ZAMBA K D,HAWKINS D M.A multivariate change-point model for statistical process control[J].Technometrics,2006,48(4):539-549.

[11]KUNCHEVA L I.Change detection in streaming multivariate data using likelihood detectors[J].IEEE Transactions on Knowledge and Data Engineering,2013,25(5):1175-1180.

[12]JOHNSONRA,WICHERNDW.AppliedMultivariateStatisticalAnalysis[M].6th ed.Upper Saddle River,NJ:Prentice Hall,2007:285-286.

（責(zé)任編輯：楊媛媛）

Research on change point detection and recognition of aircraft trajectory

SU Zhiganga,b,LI Zhiqianga,HAO Jingtangb
（a.Intelligent Signal and Image Processing Key Lab of Tianjin;b.Sino-European Institute of Aviation Engineering,CAUC,Tianjin 300300,China）

Effective division of aircraft flight phase is an important technical measure of supporting active monitoring of aviation emission,which can be realized via change point detection and recognition.According to the characteristics of multivariate T2statistics,under the assumption of variance consistency of overall data parameters,a method with sequential double sliding window is proposed for detecting and estimating the change points.Simulation and experimental results show that the proposed method is more sensitive to the size of the window,change point detecting problem of trajectory can be effectively solved with well selected window length.According to the characteristics of aircraft trajectory,the proposed method is validated effective under the premise of appropriate window length.

aviation emission;emission monitoring;trajectory segment;change point detection;T2statistics

V557；TB114.2

：A

：1674－5590（2017）04－0001－06

2017-03-28；

：2017-04-10

：中央高校基本科研業(yè)務(wù)費專項（3122017111）

蘇志剛（1972—），男，黑龍江尚志人，教授，工學(xué)博士，研究方向為陣列信號處理、信息融合及空管新技術(shù).