陈宏昆,察 豪,刘立国,孟 薇
(1.海军工程大学 电子工程学院,武汉 430033; 2.陆军武汉军事代表局,武汉 430000)(*通信作者电子邮箱chenhk93@163.com)
船舶在海上航行时通常遵循特定的航道,利用丰富的海上船舶航行信息提取出这些航道,判断船舶的航行轨迹是否遵循特定的航道作为该船舶的航行行为是否正常和是否进一步进行专家分析的依据,对分析海上船舶行为、异常检测和海洋监视具有重要意义。船舶自动识别系统(Automatic Identification System, AIS)在船舶、陆基站点和卫星之间发送与接收包含船舶静态信息(海上移动业务标识码、船名、船舶类型等)和动态信息(经纬度、航速、航向等)的AIS报文来进行船舶的识别和定位[1-3],是船舶航行数据的一种重要来源。但海上船舶数量巨大,导致AIS数据量太大而无法直接用于人工分析,对数据的存储和计算也提出了挑战,因此要求研究快速有效的智能算法对海量AIS数据进行处理并获取知识。
目前,国内外利用AIS数据获取航道信息主要是通过轨迹聚类和机器学习等方法来实现: 文献[4-5]提出了为异常检测提取航道(Traffic Route Extraction for Anomaly Detection, TREAD)的算法,该方法基于AIS信息,获取航道的起点和终点以及中途的停驻点,利用基于密度的空间聚类(Density-Based Spatial Clustering of Applications with Noise, DBSCAN)方法对这些点进行聚类,将连接类的矢量线段视为航道,分析航道的时空分布,并在此基础上对船舶行为进行预测;文献[6]提出用机器学习的聚类、分类和离群点检测等方法对轨迹相似性进行度量;文献[7]提出了基于船舶航向、航速和船舶间距离的相似性,利用k均值聚类方法获取航道内船舶的行为特征的算法。文献[8-11]提出建立轨迹结构距离模型,优化轨迹间相似性度量的算法,并结合经典DBSCAN聚类算法实现船舶轨迹聚类。但是上述算法均要求连续接收AIS信息以形成船舶轨迹,对数据质量要求高,而远海海域的AIS信息只能通过卫星平台接收,目前星载AIS接收设备不多,无法保证AIS信号的全程接收,因此不适用于远海海域的航道提取。文献[12]提出将船舶密度分布视为图像,运用Hough变换提取线性航道,并用高斯拟合法分析航道宽度,有效地提取了大面积海域的线性航道,但该方法提取航道的精度不高,且对噪声敏感。
针对上述情况,本文提出一种基于历史AIS数据,利用图像处理方法进行降噪,利用Hough变换提取线性航道中心线,并用核密度估计法估计航道宽度的航道提取算法,成功地实现对大面积远海海域的航道提取。最后采用真实的历史AIS数据进行实验验证,实验结果证明了所提方法的有效性,且提取航道的精度较文献[12]有了提升。
Hough变换是图像处理中检测直线的技术中最受欢迎的技术之一[13-14]。在远海开阔海域,航道上的船舶密度远高于其他区域,可将船舶密度分布视为Hough变换中的图像空间,分布在航道上的船舶构成图像中的直线,其他区域的船舶视为噪声,据此利用Hough变换提取线性航道。算法的总体流程如图1所示。
图1 航道提取算法流程
船舶密度是指单位面积内船舶数量。对目标区域划分网格,将单个网格的大小视为单位面积,计算每个网格内船舶的数量,即可得到该区域的船舶密度分布。受AIS数据量限制,且AIS数据经常存在误报、数据遭恶意篡改[15]的现象,船舶密度分布没有回归统计学规律,航道提取准确度较低, 因此利用图像处理中的中值滤波和形态学滤波方法对密度分布进行优化、降低噪声干扰、提高航道提取准确度。
从AIS数据库中提取目标区域的船舶位置信息,存入集合L中:
L= { (lon_i,lat_i);i= 1,2,…,N}
(1)
式中:lon表示经度,lat表示纬度,N为目标区域内AIS报文总数。将目标区域划分为n×n的网格,网格的宽度为:
(2)
式中:lon_min、lon_max和lat_min、lat_max分别表示目标区域的经纬度边界。计算每个网格内的AIS报文数量,获得密度矩阵M:
Mr,c=count(L{ (lon,lat)|lon_min+ (c-1)Δlon≤lon (3) 式中:r,c=1, 2, …,n分别表示矩阵的行和列,count(·)表示计数函数。 下面对密度矩阵进行降噪处理。首先进行中值滤波: Mr,c′=med{Mr-i,c-j;i,j∈W} (4) 式中:med(·)表示取中值,W为2×2的方形窗口。 然后利用形态学滤波中的开运算对M′进行形态学滤波:对M′作二值化处理,得到二值化矩阵C: (5) 利用Matlab软件中的开运算函数对C作开运算得到C′=imopen(C,W),其中W取为2×2的方形窗口。矩阵M′的形态学滤波的输出结果为: M″ =M′⊙C′ (6) 式中:⊙表示相同维数的矩阵对应元素相乘。 最后对矩阵M″作标准化处理: (7) Hough变换是图像处理中检测直线或曲线的一种有效方法,其原理是将图像空间中的一条直线或曲线映射到Hough空间中的一点。将船舶密度分布视为Hough变换中的图像空间,图像空间中的直线可由原点到直线的距离d和该垂线与c轴(矩阵行)的夹角θ唯一确定,直线的参数形式为: d=ccosθ+rsinθ (8) 图像空间同一条直线上的点作Hough变换后在Hough空间(d,θ)相交于同一点,计算该点位置,就可以确定船舶航道中心线。 应用Hough变换检测线性航道的算法步骤如下: 步骤1 根据1.1节船舶密度分布分析方法计算标准化密度矩阵MN。 步骤2 设置密度门限ρth,将矩阵MN中超过ρth的点记录到点集P= {(c,r)}。 步骤3 将点集P中的点(c,r)按照式(8)映射至点集Q={(d,θ)},式中:θ∈(0,2π)。 步骤4 为点集Q中的任意一点(d,θ)设置窗口G(d,θ): G(d,θ)={(dk,θk)||d-dk|≤Δd,|θ-θk|≤Δθ,(dk,θk)∈Q} (9) 式中:Δd和Δθ为窗口的大小。 步骤5 计算Hough空间中每个窗口G(d,θ)内点的数量,选取包含点数最多的窗口Gmax。 步骤6 提取Gmax内所有元素对应的点集P中的点(c,r),这些点近似共线,将其视为航道内的点,作直线拟合,得到航道的中心线l。 同陆地上的车道类似,海上的船舶航道也有一定的宽度,但是受气象、水文和海底地形等因素影响,海上航道宽度通常不始终保持一致,因此按照密度矩阵的列逐段分析航道宽度。为避免目标区域内其他航道对密度分布的影响,首先确定航道宽度分析区间。密度矩阵第c列的宽度分析区间I的定义为:设航道中心线l与矩阵第c列交于(γ,c),以γ为中心,向两侧延伸至密度值首次低于门限Ith所对应的行的集合。取该列航道中心线处密度值的0.2倍作为门限值,即Ith= 0.2Mγ,c,如图2所示。 图2 宽度分析区间 通常宽度分析区间内的船舶密度近似服从以γ为中心的高斯分布,因此用高斯函数对密度分布的密度函数进行拟合。设密度矩阵第c列的宽度分析区间I=(r1,r2,…,rm),其中m为宽度分析区间的宽度,密度分布的均值μc为: (10) 方差σc为: (11) 将μc和σc代入高斯函数,得到密度函数f(x): (12) 其中: (13) 在实际情况中,用高斯函数对密度矩阵部分列的宽度分析区间内船舶密度分布的拟合效果不好,因此采用核密度估计法。利用核密度估计法得到的概率密度函数为: (14) 式中:xi为独立同分布的样本点,h为带宽,K(·)为核函数。核函数采用高斯核函数,具体表达式为: (15) 式中:xc为核函数中心,σ为函数的宽度参数。根据密度矩阵构造样本点xi: xi= (16) 式中:{rk;k= 1,2,…,m}为密度矩阵第c列的宽度分析区间,m为宽度分析区间的宽度。将构造的样本点{xi}代入式(16),得到宽度分析区间的密度函数f(x)。 下面根据船舶密度分布的密度函数f(x)估计航道边界。参照3 dB带宽规则,取密度函数f(x)最高值的1/2作为边界门限λth: λth=fmax(x)/2 (17) 令 xl=f-1(λth) (18) 设xl的取值集合为Xl={x1,x2,…,xn},且x1 选取印度洋海域坐标为(10°N, 95°E)、(20°S, 95°E)、(20°S, 85°E)和(10°N, 85°E)四点依次连接形成的矩形海域为目标海域,提取该海域2017- 03- 01— 03- 10的历史AIS数据对本文算法进行实验验证。AIS数据统计结果如表1所示,共筛选出1 163艘船舶88 509条AIS数据。通过划分网格并设置密度门限ρth,将图像空间的对象数量从88 509降低到1 193,大幅降低了Hough变换的运算量。目标海域在地图上的位置和船舶在该海域的分布如图3(a)所示。 运用Hough变换提取线性航道,并分别用文献[12]的高斯拟合法和本文核密度估计法分析航道宽度,共提取了5条航道,如图3(b)~(c)所示。图像背景表示划分的网格,用网格色调的冷暖表示船舶密度值高低;图中实线表示航道中心线,虚线表示航道边界。图3(d)为采用DBSCAN轨迹聚类法获得的轨迹簇。航道内船舶轨迹数量占区域轨迹总数的比重和航道面积在目标区域的占比如表2所示。 表2 三种方法航道内船舶数量和航道面积的占比 从图3和表2可以看出,用DBSCAN轨迹聚类法对该海域船舶轨迹进行聚类,只能在轨迹密度极高的区域形成轨迹簇,簇内轨迹数量仅占区域轨迹总数的29.81%,而在船舶密度较为稀疏的区域则没有形成轨迹簇, 这是因为星载AIS接收机没有连续接收到船舶密度分布稀疏的区域的AIS信号,导致这些区域无法形成船舶航行轨迹。结果证明该方法不适用于大面积远海海域的航道提取。用Hough变换法提取航道,航道内船舶轨迹占区域轨迹总数的比重超过95%,而航道面积占比不到40%,没有位于航道内的轨迹占比不到5%。结果证明利用Hough变换可以有效地提取远海大面积海域的线性航道。 针对文献[12]的高斯拟合法和本文提出的核密度估计法,对这两种航道宽度分析方法作进一步分析比较。选取其中一条航道的宽度分析区间,分别用高斯拟合法和核密度估计法获得的密度函数及宽度分析结果如图4所示。 图4中航道宽度1为通过高斯拟合法分析的结果,航道宽度2为核密度估计分析的结果。从图中可以看出,在该宽度分析区间内,密度分布没有显示出对称特点,因此用高斯函数拟合的密度函数与实际密度分布有明显的偏差,而核密度估计法得到的密度函数与实际密度分布拟合的效果好于高斯分布,航道宽度2内包含的船舶密度值高于航道宽度1。从表2也可以看出,利用核密度估计法分析航道宽度,得到的航道面积比高斯拟合法小,而航道内船舶轨迹数量更多,用核密度估计法分析的航道内船舶密度高于高斯拟合法。实验结果证明了用核密度估计法估计航道宽度的精度高于高斯拟合法。 大面积远海海域航道提取对全球海上态势评估和海洋监视具有重要意义。本文提出基于AIS数据,用划分网格的方法分析海上船舶密度分布,并通过中值滤波和形态学滤波优化密度分布,将船舶密度分布视为图像,利用Hough变换提取线性航道的算法,并用核密度估计法分析航道宽度。选取印度洋海域的真实历史AIS数据对该算法进行验证,结果表明,利用Hough变换可以有效地提取远海大面积海域的线性航道,用核密度估计法分析航道宽度的结果比高斯拟合法更加精确。利用Hough变换提取远海海域的线性航道可为海洋监视和异常检测提供支持,将航行轨迹不符合航道的船舶视为异常进行重点关注,可以大幅减少海上目标监视的工作量,提高工作效率。1.2 Hough变换检测航道
2 航道宽度估计
2.1 高斯拟合法
2.2 核密度估计法
2.3 航道边界和宽度估计
3 实验验证
3.1 航道提取
3.2 算法比较
4 结语