当航班遭遇恶劣天气、空域容量下降等情况时，初期的战术级空中流量管理如地面等待，会造成大面积航班延误影响航班运行。因此，通过动态变更航空器的飞行路径来提高航班运行效率变得越来越重要，即以少量的空中飞行换取长时间的地面和空中等待。国外在改航路径规划研究上已经取得了一些具有代表性的成果：Krozel等^[1]提出了自由飞行改航路径规划方法；Mukherjee等^[2]研究了动态改航路径规划的模型；Taylor等^[3]研究了基于模拟退火算法的改航路径规划模型；Chen等^[4]提出一种基于网络的空中交通流模型，并将其用于实时的航路飞行流量预测，这对空中交通的改航路径规划具有相当大的借鉴意义；Sadovsky等^[5]提出了一种针对天气不确定性的Risk-Hedged改航方法；Agustin等^[6]提出了0–1模型，设计相关风险函数针对航班取消、改航策略进行研究以提高空中交通运行效率；Sheth等^[7]在专利中也针对空域单元相关限制提出了行之有效的改航策略与方法。Souza等^[8]针对交通拥堵问题进行研究，提出了相关预警和改变交通流向的方法以提高城市地面交通的容量和效率。国内的宋柯^[9]、李雄等^[10]、王飞等^[11]等也对改航策略、改航模型以及算法做了大量具体研究，提出了基于几何和Maklink图的改航算法。

上述大多数方法只考虑某些单一因素进行改航研究，未能建立合理全面的约束条件。此外，在改航路径选取上没有考虑关键航段的历史运行情况，不能很好的衡量选取改航路径的合理性。因此，作者在总结已有研究的基础上，提出了一种基于改航约束和挖掘历史雷达数据的运行特征的改航算法。算法首先根据经济性、安全性、航路特征和容量影响等因素粗选改航区间的航路通路以得到满足条件的关键点序列集合。随后从雷达数据中分析频繁飞行航段（利用率）确定最终改航关键点序列。算法考虑了航班真实运行情况，频繁飞行航段（利用率高的航段）意味着有较高安全性、较少飞行冲突以及更加完善的飞行保障条件。最后通过机器学习方法挖掘经过该航段的航班运动特征（运动状态向量），采用运动学方程预测出完整的改航路径。

1 改航约束条件

1）经济性约束

主要考虑当前飞机燃油量是否满足改航造成的额外燃油消耗，描述如下：

${R_{\rm f}} = w(\sum\limits_{i = 1}^N {{l_i}} - {l_{{\rm{src}}}}) \le 0.6{R_{{\rm{surplus}}}}$

(1)

式中， $w$ 为飞机飞行单位距离的燃油消耗量， $N$ 为改航规划路径的航路段数量， ${l_i}$ 为第 $i$ 段的飞行距离， ${l_{{\rm{src}}}}$ 为改航区间在改航前的飞行距离， ${R_{{\rm{surplus}}}}$ 为航班按原路飞行的理论剩余油量。

2）安全性约束

主要考虑改航路径与受限区（危险区、限制区、天气影响区域、军事活动区和流量管制区域等）之间的距离间隔是否满足最低要求，将整个改航路径分段计算，因此只需要航段起始和终止点与受限区的距离满足如下条件：

${R_{\rm s}} = \min ({S_{\!\!mn}}) \ge {S_{{\rm{thresold}}}}$

(2)

式中： ${S_{\!\!mn}}$ 为第 $m$ 个受限区域与第 $n$ 个改航航段之间的最小距离， $1 \le m \le M$ ， $1 \le n \le N$ ， $M$ 是受限区数量， $N$ 是改航航段数量； ${R_{\rm s}}$ 为指改航路径与受限区之间的最小间隔； ${S_{{\!\!\rm{thresold}}}}$ 为安全距离门限值。

3）航路特征约束

主要考虑改航路径的关键点数量、航段长度以及在各个关键点的转弯角度，包含进入改航规划路径和回到原航路的角度，约束描述如下：

改航路径关键点数量： $N \le 6$ ；

改航路径的航段长度： ${l_i} \ge 2.5{R_w},i \in [1,N]$ ，其中 ${R_w}$ 为飞机性能参数中的最小转弯半径；

飞行进入和离开改航路径的转弯角度： ${30^ \circ } \le {w_i} \le {60^ \circ },i \in [1,N]$ 。

4）对空中交通流量的影响

改航后的路径不能造成与该路径相关的空域单元（关键点、航段、空域等）出现流量告警或者预警问题。约束条件为 ${V_{\rm d}} < \alpha {V_{\rm s}}$ ，其中， ${V_{\rm s}}$ 为空域单元的静态理论飞行容量值， ${V_{\rm d}}$ 为当前时刻的飞行流量值。 $\alpha $ 为告警阈值参数。

2 算法流程

Step1　根据改航问题描述的区域获取改航起始关键点和最小改航区间终止关键点；

Step2　根据空管基础数据获取改航区间的关键点以及航段连通性，建立搜索树得到所有从起始点到终止点的改航路径集；

Step3　根据约束条件对改航路径集的所有路径进行判别，去掉不满足所提出的改航约束条件的路径；

Step4　若改航路径集中不存在满足条件路径，将终止关键点设置为计划航路的下一关键点，重复Step2和Step3，否则转到Step5；

Step5　基于历史雷达数据进行航段利用率分析，选取最优改航路径关键点序列；

Step6　根据改航路径关键点序列信息，基于历史雷达数据规划改航路径。

3 改航关键点选择 3.1 变量描述

1）历史航班信息集合 $F = \{ {f_1},{f_2},\cdots,{f_m}\} $ ，每个航班用 ${f_i}$ 表示，包含如下信息 ${f_i} = \{ id,P,R\} $ ，其中： $id$ 为航班唯一标识，航班计划航路关键点序列为 $P = \langle p{r_1},p{r_2},\cdots\!,p{r_m}\rangle $ ， $p{r_j}$ 为空管基础数据中的航路关键点。该航班实际飞行过程中的雷达监视数据轨迹点序列为 $R = \langle { r}{{ t}_{ 1}},{ r}{{ t}_{ 2}},\cdots\!,{ r}{{ t}_{ n}}\rangle $ ， $r{t_j}$ 为该航班的实际轨迹点（本文被分析对象），数据中的经纬度位置已转换在同一直角坐标下，表示为 ${r}{{t}_i} = [{x_i},{y_i},{{\textit{z}}_i},{v_x},{v_y},{v_{\textit{z}}}]$ 。

2）针对某空域单元 $au$ ， $F$ 中根据雷达数据实际经过的航班集合表示为 ${F_r}(au)$ ，根据飞行计划中的计划航路应该经过的航班集合表示为 ${F_p}(au)$ 。

3） $KP = \{ k{p_1},k{p_2},\cdots,k{p_n}\} $ ：空域中关键点集合。每个关键点 $k{p_j}$ 包含位置信息和航班执行情况 ${F_r}(k{p_j})$ 和 ${F_r}(k{p_j})$ ，表示为 ${k}{{p}_j} = [x,y,{\textit{z}},{F_r}(k{p_j}),{F_r}(k{p_j})]$ 。

4） $AR = \{ a{r_1},a{r_2},\cdots,a{r_n}\} $ ：空域中航路段集合，每个航路段为关键点序列。

5） $d$ ：距离门限值，若 $l(r{t_i},k{p_j}) \le d$ 认为航班经过该关键点/航段。

6） ${\delta _i}$ ：航段 $a{r_i}$ 利用率，在实际飞行过程中，实际经过航段 $a{r_i}$ 的航班数量与飞行计划中预计经过航段的航班数量的比值。

3.2 算法实现

算法核心思想为：遍历历史航迹轨迹点 $r{t_i}$ ，通过航班位置信息计算其与航段集合 $AR$ 中各元素的距离，当满足距离门限 $d$ 时，记录经过了某航段数量 ${N_{{\rm{track}}}}$ ，并对比飞行计划中该航班的计划航路信息，得到预计经过该航段的航班数量 ${N_{{\rm{plan}}}}$ ，最后计算各航段的利用率。

3.2.1 距离

$S[{f_i}(r{t_j}),a{r_k}]$ 表示航班 ${f_i}$ 的轨迹点 $r{t_j}$ 与航段 $a{r_k}$ 之间的最小欧式距离，因改航多发生在巡航阶段，因此仅考虑平面距离。当 $S[{f_i}(r{t_j}),a{r_k}] \le d$ 时，认为该轨迹点经过了航段；当经过航段的轨迹点数量 $N \ge {N_{{\rm{std}}}}$ 时，认为航班经过了航段。

3.2.2 航段利用率

令 $pn(a{r_i},{F_r}(a{r_i}))$ 和 $pn(a{r_i},{F_p}(a{r_i}))$ 分别表示历史航迹数据中经过航段 $a{r_i}$ 的航班数量以及飞行计划中计划航路经过 $a{r_i}$ 的航班数量，则本文定义的航段利用率为：

${\delta _i} = \frac{{pn(a{r_i},{F_r}(a{r_i}))}}{{pn(a{r_i},{F_p}(a{r_i}))}}$

(3)

航段利用率定义主要考虑在航段实际经过航班数量与计划经过航班数量的比值。一般来说，若无特殊情况航班将按照计划航路执行整个飞行过程，其经过的航段利用率应该为1。若某一航段的利用率小于1，则说明该航段可能因为某些情况导致了航班改航，即该航段存在改航的潜在可能（如经常进行飞行管制、导航设施覆盖等问题），在确定改航路径时必须避免这些具有潜在风险的路径。同时，航段利用率大除与正常飞行次数相关外，还可以衡量其被选择为改航路径的可能性。若一条航路多次被选择为航段路径，则说明该航段满足本文提出的改航约束条件（如经济、安全等），这也是在确定改航路径应该重点考虑的因素。总之，利用率高的航段更加符合改航路径的要求。

4 路径规划

经过历史航迹数据分析得到改航路径关键点之后，进一步预测完整的改航飞行路径。本文路径规划算法主要思想是以历史航班飞行的速度矢量 ${{v}} = [{v_x},{v_y},{v_{\textit{z}}}]$ 和运动学方程为基础进行航班位置预测。首先建立包含未知参数的速度矢量分布函数，并以历史雷达数据中经过改航关键点序列的航班轨迹集为训练数据学习和优化模型未知参数。由于机型对飞机运动性能有较大影响，因此机器学习时仅选取本次改航航班机型对应的轨迹数据。随后根据已学习的速度分布模型采样速度矢量序列^[12]，并以雷达数据更新周期为时间间隔采用匀加速运动方程预测航班位置，最终形成改航路径。

4.1 速度矢量分布

由于各关键点航段之间的飞行条件差异，因此对经过相邻关键点组成的航段区间的轨迹速度矢量单独建模。采用高斯分布对航段速度矢量分布进行描述^[13]，其分布函数表示为：

$p({ V}|\psi ) = \frac{1}{{\sqrt {\det 2{\text{π}} {{\varSigma}}} }}{{\rm e}^{ - \frac{{{{(V - {{\mu }})}^T}{{{\varSigma}} ^{ - 1}}(V - {{\mu }})}}{2}}}$

(4)

式中， $\psi $ 为高斯分布中待学习参数，包含均值向量 ${{\mu }} $ 和方差矩阵 ${{\varSigma}} $ 。由于本文的速度矢量包含3维空间分量，因此 ${{\mu }} $ 为3维向量， ${{\varSigma}} $ 为3 $ \times $ 3的正定矩阵。

4.2 模型参数学习

采用贝叶斯估计方法对高斯分布对未知参数进行学习^[14]，贝叶斯估计方法认为未知参数也服从某一分布，假设其分布为 $\psi = N({{{\mu }} _0},{{{\varSigma}} _0})$ ，且所有样本数据是独立同分布（iid）采样的结果。其联合分布和贝叶斯规则下的后验概率分布分别为：

$p({{V}}|\psi ) = \mathop \prod \limits_{i = 1}^N p({{{v}}_i}|\psi )$

(5)

$p(\psi |{{V}}) = \frac{{p({{V}}|\psi )p(\psi )}}{{\int_{\psi \in \Omega } {p({{V}}|\psi )p(\psi ){\rm d}\psi } }}$

(6)

在平方误差损失下的最优贝叶斯估计量为：

$\overline \psi = d\int_{\psi \in \Omega } {\psi p(} \psi |{{V}}){\rm d}\psi $

(7)

式中： $V$ 为训练样本的速度矢量序列，长度为 $N$ ；积分 $\int_{\psi \in \Omega } {p({{V}}|\psi )p(\psi ){\rm d}\psi } $ 与待学习参数无关。

4.3 速度矢量序列预测

在已学习到各航段区间速度矢量分布函数之后，需要通过采样得到预测的速度矢量序列，本文采用了一种效率较高的MCMC（Markov chain Monte Carlo）采样方法^[15]，其基本思想如下：

1）初始化速度矢量为上一阶段的最终值，假设为 ${{{v}}^0} = [{{v}}_x^0,{{v}}_y^0,{{v}}_{\textit{z}}^0]$ ；

2）在 $t = 1,2,\cdots,n$ 时刻采用如下方法进行速度矢量采样， $n$ 为预测周期数：

$\left\{\begin{aligned}& v_x^{t + 1} \sim p({v_x}|v_y^t,v_{\textit{z}}^t),\\& v_y^{t + 1} \sim p({v_y}|v_x^{t + 1},v_{\textit{z}}^t),\\& v_{\textit{z}}^{t + 1} \sim p({v_{\textit{z}}}|v_x^{t + 1},v_y^{t + 1})\end{aligned}\right.$

(8)

在对速度矢量序列采样之后，采用匀加速运动方程计算航班各预测周期的位置形成最终的改航路径。设航班在 $t$ 时刻的位置和速度矢量分别为 ${{{p}}_t} = [{x_t},{y_t},{{\textit{z}}_t}]$ 和 ${{{v}}^t} = [v_x^t,v_y^t,v_{\textit{z}}^t]$ ，根据匀加速运动学方程，航班在 $t + 1$ 时刻位置可表示为 ${{{p}}_{t + 1}} = [{x_{t + 1}},{y_{t + 1}},{{\textit{z}}_{t + 1}}]$ ，其中（ $T$ 为轨迹更新周期）：

$\left\{\begin{aligned}& {x_{t + 1}} = v_x^tT + 0.5(v_x^{t + 1} - v_x^t)T,\\& {y_{t + 1}} = v_y^tT + 0.5(v_y^{t + 1} - v_y^t)T,\\& {{\textit{z}}_{t + 1}} = v_{\textit{z}}^tT + 0.5(v_{\textit{z}}^{t + 1} - v_{\textit{z}}^t)T\end{aligned}\right.$

(9)

通过预测航班在改航飞行阶段不同的预测周期的位置，我们可以得出航班的改航航迹。

5 应用实践

本文提出的改航算法自研究开始已在某大型流量管理系统中得到应用并逐步完善，通过事后轨迹分析验证算法的正确性和有效性。下面简要说明本文提出算法的仿真过程，仿真数据为2014年7月14日某已改航航班，该航班在SUBUL-P248航段区间受天气影响需改航飞行，系统历史数据包含2013年1月1日至2014年7月1日的所有航班信息，航班数据包含航班计划信息和雷达航迹信息，评价指标为算法改航路径的航班位置信息与雷达采集的航班改航区间内位置信息的误差。本文算法的“改航关键点选择”部分由于与现有系统相关联采用C++实现，并为后续实现提取出所需数据；“路径规划部分”在Matlab中实现并仿真。

5.1 关键点选择

首先从最小改航区间开始进行航路连通性判别，得到满足条件的改航关键点序列有4个候选，其航路点序列分别为：

1）〈SUBUL, TOREG, P248〉；

2）〈SUBUL, TOREG, GAO, VENON〉；

3）〈SUBUL, TOREG, OGOMO, JTG〉；

4）〈SUBUL, TOREG, WFX, OMGUS〉。

随后，按照约束条件对关键点集合进行判别，结果如表1所示，其中经济性以改航造成的额外飞行距离表示，km；安全性以关键点与天气影响区域的最小距离表示，km。

然后，对表1中满足约束条件的关键点序列进行利用率分析可知，由第3和4个序列点组成的航段的利用率相对较高，约为93%左右。但第4个序列的额外飞行距离较短（经济性更高），因此选择关键点序列4作为改航航段，如图1所示，其中，实线为计划航路，虚线为改航关键点的连线。

5.2 路径规划

根据算法步骤，在历史数据中查找所有经过“SUBUL-TOREG-WFX-OMGUS”关键点序列的雷达数据，并提取飞行各时刻在三维空间内的速度矢量，将其按飞行日期和航班号保存在文件中。在Matlab中实现GMM参数学习和速度矢量预测功能，利用匀加速运动方程进行航班位置点计算，并与历史数据中该航班在改航区间采集的位置点进行比较。本文的改航发生在巡航阶段，因此不考虑高度变化（ ${\textit{z}}$ ）的影响。

图2和3分别描述了本文算法的改航路径与真实采集轨迹的误差曲线图。通过图2可知，由本文算法选择的最优改航关键点序列（路径）与实际运行过程中采取的改航策略一致。进而从图3可知基于机器学习的改航路径规划算法与历史数据采集的轨迹点具有较高的相似度，其距离误差较小且稳定。通过更加准确的改航路径轨迹预测，改航路径相关的空域单元（机场、空域等）能够提前预测局部空中交通态势并制定相关管理措施。

6 总结与展望

针对航班运行过程中的改航问题进行研究，综合考虑经济性、安全性、航路特征等约束以及真实历史数据的运行特点，应用数据挖掘和机器学习的思想解决改航关键点选择以及路径规划问题。通过实践证明了算法的正确性，且具有很高的实用价值。本文的主要工作如下：

1）建立了改航约束条件，为改航区间的粗选航段提供了一种有效的方法。采用的约束条件全面的考虑了改航过程的经济性、安全性、航路特征以及局部流量影响等因素，能更加合理的选择改航路径；

2）提出了基于历史雷达数据分析关键点利用率的路径规划算法以获取改航路径的最优关键点。算法考虑关键点（航段）在历史数据中的运行情况，能够筛选出更加安全的改航路径。

3）建立速度矢量分布函数，通过机器学习算法学习模型参数并预测速度矢量序列，最终利用运动学方程预测规划路径。机器学习模型能够更加准确的预测航班的飞行路径，也为改航路径相关空域单元的流量管理提供数据支撑。

同时，在运行实践中遇到的问题也指明了本课题下一步研究的重点：建立更加合理高效的历史雷达数据搜索和分析模型，并基于已选关键点进行更加准确的路径规划等。