基于改进分水岭算法的空域实时流量监测*

朱代武 刘 豪 朱瑀然

（1.中国民用航空飞行学院 广汉 618307）（2.中国民航西南地区空管局 成都 610000）

1 引言

随着航空器数量日益增多，空域趋于饱和或超负荷使用，使得空中交通流量管理压力激增，为解决该问题需对空域容量和实时流量进行监测和告警。传统交通监测方法从战略和战术角度对描述对象的数量进行统计，其中战略级交通监测在地面交通领域使用较为广泛，主要包括基于线性系统理论的方法［1～2］、基于知识发现的智能模型监测方法［3］、基于非线性系统理论的方法［4～5］和基于组合的监测方法［6］等，并逐渐推广到空中交通流量监测［7～8］；战术级交通监测主要应用在空中交通领域，根据航迹监测不同时段内的航空器数量［9～10］。但是上述方法都无法对重叠航空器（同一空域不同高度层配备）进行有效监控，使得航空器监测统计上会有缺失，进而对航行安全及飞行数据管理造成影响。为解决监测统计缺失问题，在对分水岭算法研究基础上结合结合航路航线、尾流间隔及移交间隔等进行改进，提出基于改进分水岭算法对空域实时流量进行监测告警。

2 分水岭算法

分水岭算法依据图像拓扑结构进行形态切割，将在空间位置上相近且灰色度相近的像素点相互连接形成封闭轮廓［11］。

2.1 分水岭算法模型构建

设定分水岭算法的具体步骤为N1、N2、N3、…、Nk为图像（fx，y）的局部特征极值点，A1、A2、A3、…、Ak为局部极值点对应的汇聚点集合，为在n 阶段汇水盆地Ak中的点集，An表示n 阶段所有极值特征点的集合，fmax和fmin表示图像灰度的最大值和最小值，得式（1），同时可得汇水盆地集合为式（2）：

若Tn表示图像中低于平均水平面的集合点即Tn=｛（i，j）|f（i，j）＜n｝，则有

通过式（3）可构造单二值化图像，若原图的点f（i，j）ϵTn，则设该点为0，反之为1。通过对0～1 变量计数可知在n 阶段处于设定水位下的点数目［12］，具体方法如下。

设定初始汇水盆地的合集Afmin+1=Tfmin+1，然后进入递归阶段。将水位从n=fmin+1以整数值递增至n=fmax+1，记录相对阶段被覆盖的点集Tn。设第n阶段Tn中有P 个连通向量，其中p≤P，设第n 阶段已经完成An的构造，对于p 值比较连通变量与An-1，两者需满足同等关系执行同类操作，反之亦然。

通常情况下图像在采集过程中因偶然误差或必然误差存在噪声，分水岭算法在对目标进行处理时会将噪声误以为目标信息，会使得图像处理质量变差。因此需要对图像进行降噪处理，其中定义理想图像为u（x，y），含噪图像为u0（x，y），n（x，y）为噪声，则两者关系为

为去除图像中存在的噪声问题，可建立偏微分去噪模型，即

式（8）需满足如下条件：

基于拉格朗日方程和快速下降法可得去噪图像I（x，y），其中x，y属于噪声位置点坐标：

2.2 分水岭算法改进

因分水岭算法对边缘依赖性较强，易受到干扰因素的影响产生过分切割现象，因此需对分水岭分割算法进行改进。因梯度可以增强图像中灰度变化比较明显的边缘，可突出区域轮廓，因此可通过梯度对分水岭算法进行改进［12］。

设去噪图像为I（x，y），结构元素为s（i，j）；图像定义域设置为RI，结构元素定义域设定为Rs，则数学形态学的膨胀运算和腐蚀运算如下。

1）膨胀

2）腐蚀

则形态学梯度可表示为

针对不同类型噪声可选取多元结构元素进行滤波保证，同时为加强梯度效果，本文拟选两种结构元素进行梯度精确取值。因信息熵可对图像承载的丰富信息进行描述，即信息熵值与图像丰富程度成正比，可引入不同元素的信息熵进行表示：

式（15）中qi（j）为梯度图像内灰度值为j 概率，梯度信息权重可表示为

式（16）中m表示元素数量，不同结构元素梯度信息可通过权重进行融合，可描述为

因拟选两种结构元素进行梯度精确取值，则n为元素组别编号；通过式（18）可得彩色图像各分类梯度，选取最大值作为彩色图像梯度取值，则梯度为

3 空域条件约束

本文拟选取西南贵州高空空域作为研究对象，在巡航过程中因前航空器涡流会对后航空器运行安全造成影响，需设定尾流间隔进行安全保证。根据《空中航行服务程序-空中交通管理》（PANS-ATM，Doc4444），雷达间隔下的尾流间隔最低标准为表1所示［13］。

表1 达间隔下的尾流间隔最低标准（单位：km）

如表1 所示，可知尾流间隔大小排序为（6，8，10，12）（单位：km），选取尾流间隔最大值对应元素作为全平面结构半径R，则上述尾流间隔标准转化为自适应图像分割率C，即

其中Ti表示每张图像标称航空器总数，Ni表示每张图像中航空器实际总数，pi表示每张图像的有效分割率，R0为设定的平面圆盘形的结构元素标准半径，R1为平面圆盘形的结构元素自适应全局半径。

4 仿真分析

为测试改进分水岭算法在空域实时流量监测的分割性能，选取硬件运行环境为Intel Core i7 处理器，16GB 内存容量，软件运行平台为MATLAB2018a。令平面圆盘形的结构元素标准半径R=3，最大外接矩形取监测航空器边缘取值最大处，缺省值为1，N 最大取值≤255，随机从flighttradar24截取部分西南贵州高空空域，图片像素尺寸为1000×800，如图1所示。

图1 部分西南贵州空域航空器分布

对图1 进行灰色处理，并根据模型及所设参数对图1 利用imbinarize 函数进行自适应阈值和反运算处理，如图2，图3所示。在反运算处理基础上对其进行腐蚀后开运算［14］，其中选取平面圆盘形的结构元素指定半径为3，如图4所示。

图2 部分西南空域图灰色处理

图3 图像自适应阈值与发运算处理

图4 图像腐蚀处理与开运算

通过对图像进行灰色处理、自适应阈值二值化处理与反运算和腐蚀开运算，可得出基于传统分水岭法图像分割的基本轮廓，并按照式（18）～（20）进行连通变量计算和边缘点标记，每个连通域的面积保存在数组中，见图5。

图5 传统分水岭算法标记的航空器数量

已知该空域内航空器实有架数为31 架，传统分水岭算法标记为29 架，准确率为93.5%，分析可知处于重叠的两架或多架航空器（同一空域不同高度层）无法进行精确标记，如图5 黑色方框所示。为此引入ICAO 尾流间隔对自适应结构半径R 进行调整，可对同一空域内遵循最小尾流间隔标准的航空器进行自适应标记［15］，以减小监测统计误差。通过imshowpair 函数将蒙版叠加在原始图像，改进后的分水岭算法可对欧氏距离进行变换，使得仅在所需位置具有最小值［16］，并重复循环直到所有航空器在遵循尾流间隔条件下皆可标记，改进后自适应阈值二值化处理图像及分水岭分割图像如图6、图7所示。

图6 改进后的图像二值化处理

图7 为欧氏距离变换及采用改进分水岭算法的图像丰富程度，通过对信息熵进行采集标记如图8，可知标记航空器数量为31 架，已知该空域内航空器实有架数为31 架，准确率为100.0%。对分水岭算法改进前后准确率进行对比，准确率如表2 所示。对比可知利用改进分水岭算法进行二值图像的距离变化，并在理想情况下梯度分割，使得边缘点、中心点与最大外接矩形进行最优化处理，使得区域总数能够被精确标记。

图7 距离变换及改进分水岭变换

图8 改进分水岭算法标记的航空器数量

表2 改进前后空域流量准确率表

因航空器在高空空域巡航时主要依靠地面通信导航监视设备对其进行信号传输，且航空器自身载有监视通信导航设备，在一定程度上可作为独立个体存在，因此引入几何分布对改进前后分水岭算法监测的准确率进行检验分析［17］。已知pi表示每张图像的有效分割率，无效分割率为（1-pi），m为重复循环次数，设定m次有效分割概率的倒数Lm为误差值：

通过将图像所含元素数据进行处理并带入算法模型，分水岭算法改进前后的误差值Lm为图9 所示，基于传统分水岭算法的原始误差为1800，改进算法的原始误差为1760。随着图像分割循环次数的增加，图像有效分割误差减少，且改进后图像的误差皆小于传统算法误差。如当循环次数为4 时，传统分水岭算法下的实时误差为1580，改进算法下的实时误差为1500，间接反映了改进分水岭算法在实时流量预测方面更具有合理性。

图9 分水岭算法改进前后误差值

5 结语

为解决噪声对图像的干扰问题及标记元素重叠问题，在遵循图形分割原则上提出一种结合空域尾流间隔约束的空域流量实时监测方法。该算法能够有效对重叠粘连的元素进行自适应处理分割与识别，使得监测精度结果可达到100%。目前国内主要从容量评估、流量预测和飞行程序优化三方面提高航行安全和效率，较少利用图像处理技术对高空实时流量进行监测告警，因此本文在一定程度上为基于图像处理的自由航行（TBO/4D 航迹）提供理论依据，也为缓解现实中的管制负荷提供了新途径。