基于超像素分割变化特征提取的高铁环境变化监测

2019-08-01王凯

铁道勘察 2019年4期

王凯

(中国铁路设计集团有限公司，天津 300251)

随着我国高铁建设的快速发展，高铁里程和覆盖范围不断增加，预计到2020年，我国高铁的总里程将超过3万公里。高速铁路这种高速交通运输方式在惠及千家万户的同时，其运营安全也引起了社会的高度关注。如何有效监测高铁沿线环境并快速提取其环境的变化情况，成为国内众多学者新的研究方向。

高速铁路时速高，对运营环境要求比较苛刻。在极端恶劣天气(飓风、台风、强降雨、冰雹等)下，铁路周围不稳定的建筑体，如彩钢瓦、农作物大棚、高耸的电塔、线杆、深挖的取土坑、密集堆积的货场等，都可能对高铁的运行造成严重威胁，轻则造成高速列车紧急停车，重则造成重大安全事故。因此，密切关注高铁沿线运营环境的变化，随时掌握高铁沿线风险源的分布与位置对高铁运营安全至关重要。现行高铁沿线变化监测的方式主要为人工巡视，该方法受人视野和作业环境等因素的限制，作业效率低下，精度也难以保证，无法满足作业任务的需求。遥感技术具有大尺度、多时相等特点，使用遥感数据可以快速获取大区域空间信息。随着我国高分卫星数据平台的研发和应用，遥感数据的空间分辨率得到了显著的提高，

通过高分数据平台可以获取地表米级、亚米级的遥感影像数据。以下采用1A级高分二号遥感影像数据作为数据源(空间分辨率为1 m)，对京津城际高速铁路天津市武清段的环境变化进行研究。

遥感变化监测的主要方法有:影像直接比较法、分类后比较法和直接分类法三种[1]。直接比较法使用配准的两个时相单波段遥感影像直接进行比较(结果为比值或差值)，得到地表变化信息的分布图像，如黄阅智等使用DEM数据差值对冰川的变化情况进行研究[2]。Khalsa S J S应用比值法进行了地表变化的监测工作[3]。无论是比值法还是差值法都易受地表其他因素的干扰，且只能利用单个特征图层的信息，精度较低。分类后比较法使用配准前后时相的遥感影像分别进行分类，然后将分类后的矢量数据进行空间叠加，得出地表变化分布，如赵恒谦等利用ETM+和Sentinel-2A数据对北京通州土地利用和生态环境变化进行了监测和分析[4]，王琎等基于Landsat等数据源，采用最大似然和面向对象的分类方法，对珠口地区的土地利用情况进行提取，得到了该区域的空间格局变化趋势[5]。分类后比较法过程较为繁琐，需要对两个时相的地物进行分类，提取精度取决于分类的精度，通常情况下，无法避免因分类细节差异而产生的噪声，其结果多用于类别变化的统计分析。直接分类法将前后两个时相的信息进行配准和叠加，综合利用叠加的影像波段组进行变化特征提取，并利用变化特征图像直接进行变化区域、未变化区域的分类。直接分类法具有噪声少、提取精度高、分类简便等优势，如于珊珊等利用差值法获取变化特征，在大连等地使用面向像元的直接分类法进行提取实验(kappa系数达85%)，取得了较好的提取效果[6]。图像直接比较法、分类后比较法和直接分类法均基于像元单位进行数据操作，对中低分辨率遥感影像较为适用，对于高分率遥感影像，容易产生“椒盐噪声”现象，使得提取效率下降[8]。

根据以往研究情况，采用面向对象直接分类方法中的最小噪声分离算法(MNF)，综合利用两个时相的多个波段信息提取变化特征，基于SLIC超像素分割算法将像素单元转化为地物对象单元，有效避免了“椒盐噪声”的出现。

1 多时相遥感图像配准

图像配准又称图像匹配，是将不同传感器、不同时间或不同视角拍摄图像中的同名地物点匹配到相同的目标坐标下。图像自动匹配方法可分为灰度匹配、特征匹配及变换域匹配。通常情况下，灰度匹配计算量大，计算速度慢，易受外界光线环境变化影响。特征匹配具有计算量小、匹配效果好等特点。特征匹配又可分为基于点和线的匹配，角点的提取要好于线段的检测(这是由于光圈问题决定的)，通过一个点观测一条直线，只有垂直于线段的时候才能观测到，而角点的观测则不受角度的限制[9]。角点为图像领域内各个方向上灰度变化的局部最大值。角点检测有多种算法，使用较多的为Moravec、Harris和Forstner算法等。其中，Forstner算法运算速度较快，精度较高，在图像匹配领域得到广泛的应用[10]。其算法的主体思想为：逐像元计算Robert梯度以及以该像元为中心的滑动窗口内的灰度协方差矩阵，并以此为依据选取特征点。

本研究以Forstner算法为基础进行影像匹配。首先构建匹配图像与被匹配图像的影像金字塔，分别利用Forstner算法提取最粗尺度上金字塔图层中的特征角点，对特征点进行相似性度量得到匹配特征点，之后进入更为精细的下一层金字塔图层匹配。利用最终匹配特征点作为两景影像的连接点，具体步骤为：首先逐层计算各个像元的Robert梯度(见图1)，有

fx=f(x+1，y+1)-f(x，y)

(1)

fx=f(x+1，y)-f(x，y+1)

(2)

f(x,y)f(x,y+1)f(x+1,y)f(x+1,y+1)

图1 Robert算子邻域

①逐像元计算滑动窗口内灰度协方差矩阵

(3)

其中dx，dy分别为x方向和y方向的微分。

②计算确定角点的特征值

(4)

(5)

其中trN为协方差矩阵的迹，设定阈值TW、TU，遍历图像，获取同时满足W>TW和U>TU的像元点。

③设定窗口的大小，以最大特征值点作为特征角点，选择特定的滑动窗口，利用马氏距离进行图层间特征角点的相似性度量，选择相似性最高的点对作为连接点，进入下一级金字塔继续进行匹配。对最终匹配连接点的误差进行估算，去除误差较大的连接点，(见图2)。使用研究区遥感影像进行配准，得到配准后的影像数据，其影像配准前后的对比效果见图3。由图2、图3可知，房屋等具有明显角点的区域得到了较为精准的匹配，道路等线状地物角点特征不是十分明显，其匹配效果亦较建筑地物差，但总体上得到较大改善。本次实验生成280组匹配点，图像匹配总体均方根误差(RMS)为0.31个像元，符合变化检测要求。

图2 连接点提取

图3 影像配准前后对比

2 基于最小噪声分离变换的地表变化特征提取

最小噪声分离变换(MNF)本质上是两次层叠的主成分变换，第一步变换用于分离和调节噪声波段数据间的相关性，使得噪声的方差最小[11]。第二次变换的目的是对噪声白化数据作进一步处理，通过特征值判断数据的内在维度。最终，数据被划分为两大部分：第一部分是与最大特征值相关的特征图层，其物理意义为图像中的主要信息；第二部分为相近特征值对应波段和噪声占主体的波段图层。对两个时相遥感影像图层组进行最小噪声分离，其图层最大特征波段为影像的主要信息，即未发生变化的空间主体地物信息，而变化部分的波谱信息不一致，表现出噪声的特征。因此，不同时相间空间变化信息在噪声图层表现出较为明显的特征。将精确配准后的2015年和2017年的遥感影像(见图4)进行空间叠加并进行最小噪声分离，得到8个变换后图层(见图5)。

图4 研究区遥感影像

由图5可知，主要地物信息(如建筑、水体、道路等)在前4个特征图像上得到了清晰的体现，而噪声信息在后4个特征波段较为明显。通过对比观测前后两个时相的变化区域，其变化特征在MNF6中有明显的体现，故选择MNF6图层作为区域变化信息提取的特征图层(以下称之为变化特征图层)。对变化特征图层进行9邻域窗口中值滤波处理，去除图像异常值单点噪声的影响。

3 SLIC超像素分割与阈值法变化区域提取

SLIC算法可将RGB颜色空间图像转化为CIELAB空间和笛卡尔坐标系下的特征向量，然后利用多维特征向量进行特征距离计算，以距离为依据进行像元局部聚类[12]。SLIC算法具有运算速度较快、对地物的轮廓保持较好等特点[13]。以2015年、2017年基础遥感影像图层为数据源，利用SLIC超像素分割算法进行图像分割、地物边界线提取，将图像的像素单元转化为地物对象单元，并根据光谱的相似度对分割结果进行融合。

以分割对象为统计单元，分别对变化特征图层进行分区均值和标准差统计，得到对象级别的变化特征值。分别选取50个变化区域和非变化区域的对象并对其变化特征图层的均值和标准差进行统计，得到其均值和标准差的区间范围(见表1)。

表1 样本变化特征图层统计

由表1可知，变化区域在特征图层上的取值均为负值，未变化区域主体上为正值，部分区域受植被、水体状态等因素影响出现负值。变化区域特征值标准差较小，未变化区域较大。以变化区域的均值和标准差为依据设置分类阈值，其分类的规则见表2。

表2 分类阈值

其中，Feature_value为变化特征图层的特征值，changemax、changestdev分别为变化区域的最大值和标准差，unchangemin和unchangestdev为未变化区域的最小值和标准差。应用该阈值进行特征图层分类，将图层区域分为变化区域、未变化区域和混淆区域。基于2015年、2017年遥感影像，分别计算两个时相下的归一化植被指数(NDVI)、归一化水体指数(NDWI)和亮度值(bright)，并以此为特征选择样本，使用决策树分类器，对混淆域进行二次分类，最终将区域划分为变化和未变化两类[14]。新增的彩钢瓦房(见图6(a)、(b)、(e)、(f))以及拆除的房屋(见图6(c)、(d))、楼房建筑(见图6(g)、(h))得到较明显的识别。受数据获取时相、拍摄角度等影响，高耸建筑和房屋阴影无法做到完全的配准，导致部分区域存在错分的现象[15]。以变化区域和非变化区域为目标区域，在两个目标区域分别随机生成100个检验点，通过目视叠加分析的方式对变化监测精度进行检验，统计得到混淆矩阵(见表3)。由表3可知，变化监测的总体精度为91%，错分率为4%，漏分率为5%。