CN106169086A - 导航数据辅助下的高分辨率光学影像损毁道路提取方法 - Google Patents

Abstract

本发明公开了一种导航数据辅助下的高分辨率光学影像损毁道路提取方法，以道路导航矢量数据作为辅助数据，将道路的多种特征结合学习‑检测法在高分辨率光学影像上进行毁损道路提取。由于灾前的道路导航矢量数据相对容易获取，本发明依托灾前道路导航矢量先验信息，可更高效、更准确地从高分辨率光学影像中提取出疑似损毁路段。

Description

导航数据辅助下的高分辨率光学影像损毁道路提取方法

技术领域

本发明属于遥感图像信息提取技术领域，具体涉及一种导航数据辅助下的高分辨率光学影像损毁道路提取方法。

背景技术

随着近年来卫星和航空遥感技术的发展，高分辨率影像在各方面得到了越来越多的应用。道路作为一类非常重要的军事和民用目标，当灾害或其它突发事件发生时，很可能发生损毁和阻塞，这给救灾人员和物资的输送造成了极大的阻碍。遥感影像快速获取和信息准确丰富的特征，使之成为监测和提取灾害信息的核心手段。

利用遥感影像提取道路损毁的经典方法主要为目视解译的方式，如基于灾后单时相的道路损毁提取方法。秦军等(2010)对无人机真彩色图像进行纠正之后利用目视判读的方式提取出灾害信息，然后利用道路损毁类别、损毁尺度、损毁百分比3个因子作为评价指标对道路的损毁级别评定分类。目视判读方法虽然精度高，但劳动强度大，也比较耗时。刘明众等(2014)采用面向对象法分割影像，根据道路分割对象的光谱信息、几何信息和空间信息建立知识库，然后依据知识库中的法则来获取影像中的道路。虽然分类的方法提高了效率，但是受灾后影像场景的复杂性，有时候精度不高，同时此方法常常会借助于易康等面向对象的分类软件进行，应用不方便。基于多时相的道路损毁提取主要利用灾前、灾后数据，采用变化检测的方式进行损毁提取。

导航数据是当前应用最为广泛的地图数据，其更新周期非常快，其是网上地图协作计划的产物，包括道路、水域轮廓、建筑物轮廓等矢量数据。考虑到道路是一种相对稳定的地物目标，这使得利用导航数据来辅助进行道路检测和提取具有很高的可行性和实用性。相对于传统的损毁检测方法大多缺乏对损毁信息的验证，导航数据能够避免无关信息的干扰，提高了道路损毁检测和提取的准确性。

发明内容

针对现有技术存在的不足，本发明提供了一种导航数据辅助下的高分辨率光学影像损毁道路提取方法。

本发明以矢量数据作为辅助数据，将道路的多种特征结合Learning-Detection(学习-检测)法，对受灾区域进行检测，并提取疑似损毁路段，基于道路上下文信息对疑似损毁路段进行验证，剔除虚警。

为解决上述技术问题，本发明采用如下技术方案：

一种导航数据辅助下的高分辨率光学影像损毁道路提取方法，包括：

S1获取受灾区域的高分辨率光学影像和道路导航矢量数据，高分辨率光学影像和道路导航矢量数据分别简记为影像和矢量数据；

S2对影像和矢量数据进行配准；

S3采用学习-检测法从影像中提取疑似损毁路段，本步骤进一步包括子步骤：

3.1构建搜索窗和匹配窗，搜索窗在影像中沿矢量道路搜索，匹配窗在搜索窗内搜索；

3.2在当前搜索窗内采集匹配窗图元，以匹配窗图元和基准图元为输入，以匹配窗图元和基准图元的像素灰度值方差为特征值，依次采用方差分类器、随机森林分类器对匹配窗图元进行分类，获得损毁区域图元；所述的基准图元即损毁模板库中图元；

3.3根据子步骤3.2所得损毁区域图元更新当前道路模板；

3.4沿矢量道路移动搜索窗，重复执行子步骤3.2～3.3，直至遍历完影像道路，输出损毁区域图元集；

3.5采用最邻近相关性分类器对子步骤3.4输出的损毁区域图元集和当前道路模板进行匹配检测，获得最终的损毁区域图元集，即疑似损毁路段；根据最终的损毁区域图元集更新当前道路模板；

S4检测原始影像中植被、阴影和车辆，获得植被、阴影和车辆的二值图像，对二值图像分别求并集和交集，并集和交集之差即上下文验证信息集；从疑似损毁路段中剔除疑似损毁路段与上下文验证信息集的交集，得损毁路段。

步骤S2具体为：

直接采取点选控制点进行配准，这里控制点即同名点。

步骤S2具体为：

选择若干沿道路方向且相互不平行的线段，以线段或线段延长线的交点作为同名点。

作为优选，匹配窗为多尺度匹配窗。

子步骤3.1中，搜索窗的高度H_s＝mw₀，匹配窗的高度H_p＝nw₀，搜索窗的宽度W_s＝m'w₀，匹配窗的宽度W_p＝n'w₀，其中，w₀为初始道路模板的宽度，m和n在[1,3]范围取值，且m＞n；m'和n'在[0.8,1.5]范围取值，且m'＞n'；

所述的初始道路模板的宽度采用如下方法获得：

在矢量数据中矢量道路上选定若干点C_i，以各点C_i分别为基准点建立道路足迹辐条，获得道路足迹辐条与影像中道路和非道路相连处的交点E_i，计算点C_i和对应交点E_i的距离；将矢量道路上所有点C_i和对应交点E_i的距离的平均值作为初始道路模板的宽度w₀。

子步骤3.5中，采用最邻近相关性分类器对子步骤3.4输出的损毁区域图元集和当前道路模板进行匹配检测，具体为：

Step1：采用最邻近相关性分类器分别计算损毁区域图元集中各损毁区域图元与当前道路模板在法线方向的距离，过滤掉距离小于距离阈值T_d的损毁区域图元，执行Step2；

Step2：对Step1筛选后的各损毁区域图元，计算损毁区域图元与当前道路模板的最邻近相对相似度，过滤掉最邻近相对相似度大于相似度阈值的损毁区域图元。

步骤S4中所述的检测原始影像中植被、阴影和车辆，获得植被、阴影和车辆的二值图像，具体为：

将原始影像从RGB色彩模型转换到HSI色彩模型获得强度影像，基于强度影像检测植被区域，并获得植被的二值图像；

将原始影像从RGB色彩空间转换至HSI色彩空间，采用阴影色彩不变量指数检测阴影区域，并获得阴影的二值图像；

对原始影像分别进行形态学开重构和闭重构运算，得顶帽变换影像和底帽变换影像，分别对顶帽变换影像和底帽变换影像进行阈值分割，获得车辆的二值图像。

和现有技术相比，本发明具有如下优点：

由于灾前的道路导航矢量数据相对容易获取，本发明依托灾前道路导航矢量先验信息，可更高效、更准确地从高分辨率光学影像中提取出疑似损毁路段。

附图说明

图1为实施例中影像和矢量数据的配准示意图，其中，图(a)为单点配准示意图，图(b)为线选配准示意图；

图2为提取疑似损毁路段的具体流程图；

图3为实施例中道路足迹辐条建立示意图；

图4为搜索框和匹配框示意图；

图5为三级级联分类器的具体流程图；

具体实施方式

下面将结合附图和具体实施方式进一步说明本发明技术方案。

步骤1：获取受灾区域的高分辨率光学影像，“高分辨率光学影像”后文简记为“影像”。

步骤2：搜索是否有受灾区域的道路导航矢量数据，道路导航矢量数据一般为shp文件，“道路导航矢量数据”后文简记为“矢量数据”。

若存在本地的矢量数据，则根据受灾区域进行矢量文件裁剪，以得到受灾区域的矢量数据。若不存在本地的矢量数据，则下载对应的矢量数据。

步骤3：矢量数据和影像的配准。

影像与矢量数据的配准主要是利用空间同名点，寻找矢量数据中点位坐标与影像上相应点位坐标间的空间几何变换关系。

(1)大场景下，矢量数据中道路的拐点和交叉点在影像上相应位置明显，因此，可直接采取点选控制点，所选控制点即同名点，见图1(a)。

(2)小场景、影像中道路较稀少或影像中道路弯曲度较大时，将导致道路的拐点和/或交叉点难以找到，此时，可使用线配准方式。选择若干沿道路方向且相互不平行的线段，以线段或线段延长线的交点确定配准同名点，见图1(b)。

步骤4：采用Learning-Detection法从影像中提取疑似损毁路段。

(1)根据道路足迹辐条(Jiuxiang Hu等，2007)计算初始道路模板的宽度。

在矢量数据中矢量道路上选定若干点C_i，对各点C_i分别进行如下操作：

以点C_i为基准点建立法线，以预设角度张量α建立道路足迹辐条，道路足迹辐条与法线夹角为α，见图3。由于影像道路局部具有灰度均质性，由点C_i发射的道路足迹辐条会在影像中道路与非道路相连处有交点E_i，即灰度特征突变点。计算点C_i和交点E_i的距离。

将矢量道路上所有点C_i和交点E_i的距离的平均值作为初始道路模板的宽度。

在确保灰度特征突变点在检测框内的前提下，角度张量越大越好；具体实施时，出于计算效率考虑，可适当减小角度张量。

(2)采用检测窗扫描影像，检测窗所包围图像区域称为图元。计算基准图元和图元的特征值，便于后续的匹配检测、识别损毁，基准图元即损毁模板库中图元。

特征值即图元的像素灰度值方差，如下：

$s^2=\frac{1}{n}\[{(g_1-\stackrel{\‾}{g})}^2+{(g_2-\stackrel{\‾}{g})}^2+{(g_3-\stackrel{\‾}{g})}^2+...{(g_n-\stackrel{\‾}{g})}^2\]---\left(1\right)$

式(1)中，s²为图元的像素灰度值方差；n为图元中像素数；g_i为图元中第i个像素的灰度值，i＝1,2,...n；为图元中像素灰度值的平均值。

本发明中，检测窗包括搜索窗和匹配窗两级模板。搜索窗定义为沿影像中沿矢量道路进行搜索的缓冲区，由于受配准精度和初始道路模板计算误差的影响，一般搜索窗尺寸可以设置的大些。匹配窗是在搜索窗内搜索的模板，需要上一步的特征检测器进行特征的计算寻找最优匹配的检索范围，见图4。

搜索窗和匹配窗尺寸的确定跟初始道路模板尺寸相关。记初始道路模板的高度H_org＝d_h，d_h取子步骤(1)所获得的初始道路模板的宽度。顾及匹配误差，搜索窗的高度H_s＝md_h，匹配窗的高度H_p＝nd_h，m和n为放大系数。本具体实施方式中，取m＝2，匹配窗放大系数n＝{0.8,1.0,1.2,1.5}，采用多尺度。

记初始道路模板的宽度W_org＝w，由于水平方向也需要跟踪，因此搜索窗和匹配窗的宽度可按照固定参数设置，本具体实施方式中，搜索窗的宽度W_s＝1.2W_org，匹配窗的宽度W_p＝0.8W_org。

上述，搜索窗和匹配窗的高度为垂直于矢量道路的边的长度，搜索窗和匹配窗的宽度为平行于矢量道路的边的长度。

见图5，对匹配窗中图元采用三级级联分类器进行匹配，检测损毁图元。三级级联分类器包括依次级联的方差分类器、随机森林分类器和最邻近相关性分类器。方差分类器是将匹配窗图元和基准图元的像素灰度值方差比较，首先进行粗筛选，即将像素灰度值方差大于方差阈值T的匹配窗图元过滤。为方便后续提高计算效率，再将与基准图元像素灰度值方差差值大于差值阈值的匹配窗图元过滤。随机森林分类器用来对方差分类器的输出结果进行进一步过滤。

随机分类器基于决策树实现，其中，节点代表对象，对象具有属性，分叉即对象可能的属性，叶子节点代表从根节点遍历到叶子节点所经过的对象的属性值。决策树仅有单一输出，是一种简单但广泛应用的分类器。最近邻相关分类器用来判断候选待匹配的图元与基准图元相关性的分类器，通过计算相关系数寻找最优匹配。

基于上述理论检测疑似损毁路段，具体流程见图2，包括：

从初始道路模板开始，沿着矢量道路方向迭代进行检测跟踪。迭代过程如下：

(1)采用道路足迹辐条法初始化道路模板；

(2)构建搜索窗和多尺度匹配窗，采用三级级联分类器检测匹配窗图元，获得损毁区域图元；

(3)采用损毁区域图元更新当前道路模板，即将损毁区域图元添加到当前道路模板的负模板库，非损毁区域图元添加到当前道路模板的正模板库；

(4)沿着矢量道路方向进行迭代追踪，直至遍历完影像道路，输出损毁区域图元集。

根据输出结果进行道路损毁检测，主要分为以下两步：

Step1：采用最邻近相关性分类器分别计算损毁区域图元集中各损毁区域图元与当前道路模板在法线方向的距离，过滤掉距离小于距离阈值T_d的损毁区域图元，然后执行Step2。

Step2：对Step1筛选后的各损毁区域图元，计算损毁区域图元与当前道路模板的最邻近相对相似度，过滤掉最邻近相对相似度大于相似度阈值的损毁区域图元，输出的损毁区域图元即疑似损毁路段。

距离阈值T_d和相似度阈值为通过重复试验确定的经验值。

相对相似度计算方法如下：

把道路模板(即正模板)和非道路模板(即负模板)的集合记为目标模型M，目标模型M可表示为p⁺和p^-分别代表正模板和负模板，下标代表模板序号。正模板根据其加入到道路模板库中的时间顺序排列，即代表首先被添加到道路模板库的正模板，代表最后加入道路模板库的模板。

图元p_i、p_j的相似度S(p_i,p_j)定义为：

S(p_i,p_j)＝0.5(NCC(p_i,p_j)+1) (2)

式(2)中，NCC(p_i,p_j)为图元p_i、p_j的归一化互相关系数。

将图元p_i、p_j归一化后，采用灰度进行计算图元p_i、p_j的归一化互相关系数NCC(p_i,p_j)，如下：

$NCC(x,y)=\frac{\underset{i=1}{\overset{n}{\Σ}}{(x_i-\stackrel{\‾}{x})}^2{(y_i-\stackrel{\‾}{y})}^2}{\sqrt{\underset{i=1}{\overset{n}{\Σ}}{(x_i-\stackrel{\‾}{x})}^2.\underset{i=1}{\overset{n}{\Σ}}{(y_i-\stackrel{\‾}{y})}^2}}---\left(3\right)$

式(3)中，x、y分别表示图元p_i、p_j；x_i、y_i分别表示图元p_i、p_j中第i个像素的灰度值；分别表示图元p_i、p_j的像素灰度值的平均值；n是图元中像素数。

正目标最近邻相似性S⁺(p,M)，也叫正最近邻相似性：

$S^{+}(p,M)={\max }_{p_i^{+}\∈M}S(p,p_i^{+})---\left(4\right)$

负目标最近邻相似度S^-(p,M)，也叫负最近邻相似性：

$S^{-}(p,M)={\max }_{p_i^{-}\∈M}S(p,p_i^{-})---\left(5\right)$

相对相似度S^r，取值区间在[0,1]，取值越大表示近邻性越强：

$S^r=\frac{s^{+}}{s^{+}+s^{-}}---\left(6\right)$

设相对相似度集S＝{S^r ₁,S^r ₂,S^r ₃...S^r _n}，S^r _i是相对相似度，S^r _max代表最邻近(最大)相对相似度。当最邻近相对相似度S^r _max大于阈值T_r，则当前匹配结果有效，判断为道路；否则认为当前点匹配结果无效，认为此处为道路损毁区域。

进行下次迭代前，更新当前道路模板。

当道路损毁检测成功时，表明此处为道路损毁区域。进行下次迭代时，当前道路模板应该使用上一次迭代的道路模板进行继续匹配，沿着道路矢量直到找到下一处完好道路时再进行更新。

当道路损毁检测失败时，即没有检测到损毁时，说明此处为完好道路，再进行下一次迭代时，当前道路模板更新为最邻近相关性分类器输出的当前最邻近匹配窗图元。

道路检测存在出错的可能，因此对模板库进行更新：

将检测到的道路图元放人正模板库；利用结构约束规则对分类器产生的道路图元和非道路图元进行判断，查找错误分类的图元；如果分类器将道路图元分类错误，则将更正后的正负目标样本放入相应的模板库进行更新；下次迭代检测时学习训练。

步骤5：道路损毁检测。

通过步骤4可获取疑似损毁路段，在疑似损毁路段缓冲区内，对植被、阴影和车辆进行检测，具体如下：

步骤5.1：将原始影像从RGB色彩模型转换到HSI色彩模型，H、S、I分别表示色调、饱和度、强度。

归一化强度影像I，采用阈值T_I分割强度影像I得到明暗分类标记影像；将根据色彩不变量分割后标记影像V(i,j)与明暗标记反转后影像～I(i,j)做逻辑与操作，得植被检测影像VI(i,j)，见式(7)：

VI(i,j)＝V(i,j)&(～I(i,j)) (7)

其中，(i,j)表示像素位置。

步骤5.2：将原始影像从RGB色彩空间转换至HSI色彩空间，由于被阴影遮挡像素将导致路面的亮度降低，有阴影路面的色调、饱和度则同路面近似，通过阈值化处理方式将阴影色彩不变量指数在一定范围内的像素判定为阴影，见式(8)：

${\mathrm{\ψ}}_s(i,j)=\frac{4}{\mathrm{\π}}\·arctan\left(\frac{I(i,j,r)-\sqrt{I{(i,j,r)}^2+I{(i,j,g)}^2+I{(i,j,b)}^2}}{I(i,j,r)+\sqrt{I{(i,j,r)}^2+I{(i,j,g)}^2+I{(i,j,b)}^2}}\right)---\left(8\right)$

其中，(i,j)表示像素位置，ψ_s(i,j)表示像素(i,j)的阴影色彩不变量指数，I(i,j,r)、I(i,j,g)、I(i,j,b)分别表示像素(i,j)在红光波段、绿光波段与蓝光波段的像素值。

步骤5.3：按照特定角度间隔定义一系列线状结构元素，基于线状结构元素分别对原始影像进行形态学开重构和闭重构运算得顶帽变换影像和底帽变换影像，开重构运算用于滤除相对于背景较亮的小尺寸对象，即长度小于线状结构元素的对象；闭重构则用于滤除较暗的对象。

将基于不同方向线状结构元素重构得到的影像融合，则保证对不同方向车辆对象兼顾。重构后影像中分别滤除了尺寸短于线状结构元素的偏亮或偏暗对象。将重构后影像与原始影像做差值，实现对被关注地物的增强。

增强后影像中亮对象与暗对象相对于背景的反差被增强，被关注地物得以凸显。分别对顶帽变换影像和底帽变换影像执行基于Otsu(OTSU,N,1979)的阈值分割处理，得到车辆的二值图像。

步骤5.4：对子步骤5.1～5.3生成的植被、阴影和车辆的二值图像分别进行求并集和求交集运算，将并集和交集求差，得上下文验证信息集。

步骤5.5将疑似损毁路段与上下文验证信息集求交集，从疑似损毁路段剔除疑似损毁路段与上下文验证信息集的交集，得损毁路段。

Claims (7)

1.一种导航数据辅助下的高分辨率光学影像损毁道路提取方法，其特征是，包括步骤：

S1获取受灾区域的高分辨率光学影像和道路导航矢量数据，高分辨率光学影像和道路导航矢量数据分别简记为影像和矢量数据；

S2对影像和矢量数据进行配准；

S3采用学习-检测法从影像中提取疑似损毁路段，本步骤进一步包括子步骤：

3.1构建搜索窗和匹配窗，搜索窗在影像中沿矢量道路搜索，匹配窗在搜索窗内搜索；

3.2在当前搜索窗内采集匹配窗图元，以匹配窗图元和基准图元为输入，以匹配窗图元和基准图元的像素灰度值方差为特征值，依次采用方差分类器、随机森林分类器对匹配窗图元进行分类，获得损毁区域图元；所述的基准图元即损毁模板库中图元；

3.3根据子步骤3.2所得损毁区域图元更新当前道路模板；

3.4沿矢量道路移动搜索窗，重复执行子步骤3.2～3.3，直至遍历完影像道路，输出损毁区域图元集；

3.5采用最邻近相关性分类器对子步骤3.4输出的损毁区域图元集和当前道路模板进行匹配检测，获得最终的损毁区域图元集，即疑似损毁路段；根据最终的损毁区域图元集更新当前道路模板；

S4检测原始影像中植被、阴影和车辆，获得植被、阴影和车辆的二值图像，对二值图像分别求并集和交集，并集和交集之差即上下文验证信息集；从疑似损毁路段中剔除疑似损毁路段与上下文验证信息集的交集，得损毁路段。

2.如权利要求1所述的导航数据辅助下的高分辨率光学影像道路损毁提取方法，其特征是：

步骤S2具体为：

直接采取点选控制点进行配准，这里控制点即同名点。

3.如权利要求1所述的导航数据辅助下的高分辨率光学影像道路损毁提取方法，其特征是：

步骤S2具体为：

选择若干沿道路方向且相互不平行的线段，以线段或线段延长线的交点作为同名点。

4.如权利要求1所述的导航数据辅助下的高分辨率光学影像道路损毁提取方法，其特征是：

所述的匹配窗为多尺度匹配窗。

5.如权利要求1所述的导航数据辅助下的高分辨率光学影像道路损毁提取方法，其特征是：

子步骤3.1中，搜索窗的高度H_s＝mw₀，匹配窗的高度H_p＝nw₀，搜索窗的宽度W_s＝m'w₀，匹配窗的宽度W_p＝n'w₀，其中，w₀为初始道路模板的宽度，m和n在[1,3]范围取值，且m＞n；m'和n'在[0.8,1.5]范围取值，且m'＞n'；

所述的初始道路模板的宽度采用如下方法获得：

在矢量数据中矢量道路上选定若干点C_i，以各点C_i分别为基准点建立道路足迹辐条，获得道路足迹辐条与影像中道路和非道路相连处的交点E_i，计算点C_i和对应交点E_i的距离；将矢量道路上所有点C_i和对应交点E_i的距离的平均值作为初始道路模板的宽度w₀。

6.如权利要求1所述的导航数据辅助下的高分辨率光学影像道路损毁提取方法，其特征是：

子步骤3.5中，采用最邻近相关性分类器对子步骤3.4输出的损毁区域图元集和当前道路模板进行匹配检测，具体为：

Step1：采用最邻近相关性分类器分别计算损毁区域图元集中各损毁区域图元与当前道路模板在法线方向的距离，过滤掉距离小于距离阈值T_d的损毁区域图元，执行Step2；

Step2：对Step1筛选后的各损毁区域图元，计算损毁区域图元与当前道路模板的最邻近相对相似度，过滤掉最邻近相对相似度大于相似度阈值的损毁区域图元。

7.如权利要求1所述的导航数据辅助下的高分辨率光学影像道路损毁提取方法，其特征是：

步骤S4中所述的检测原始影像中植被、阴影和车辆，获得植被、阴影和车辆的二值图像，具体为：

将原始影像从RGB色彩模型转换到HSI色彩模型获得强度影像，基于强度影像检测植被区域，并获得植被的二值图像；

将原始影像从RGB色彩空间转换至HSI色彩空间，采用阴影色彩不变量指数检测阴影区域，并获得阴影的二值图像；

对原始影像分别进行形态学开重构和闭重构运算，得顶帽变换影像和底帽变换影像，分别对顶帽变换影像和底帽变换影像进行阈值分割，获得车辆的二值图像。