从遥感影像上自动提取道路是地理信息系统(GIS)快速获取和更新道路数据的一种有效途径。随着高分辨率遥感技术的快速发展，在高分辨率遥感影像上提取道路成为近年的研究热点^{[1, 2, 3, 4, 5, 6]}。在高分辨遥感影像上，道路两侧的树木阴影、路面上行驶的汽车等都会对道路提取产生不利影响，如何滤除这些不利影响是高分辨率影像道路提取的关键。

目前遥感影像道路提取方法按自动化程度可分为全自动提取和半自动提取^[7]。全自动提取虽然经过多年的研究，但仍未取得突破性的进展^[8]。现有的全自动算法稳健性差，提取结果需要大量人工后处理，效果不理想，因此采用人机交互方式的半自动提取是当前更好的选择^[9]。在各种半自动提取方法中，基于模板匹配的方法被认为是更实用的一类方法^[10]。但传统的模板匹配方法主要采用相关系数作为相似性测度^[11]，因此对遮挡敏感，稳健性差，当路面上出现树荫、车辆时会产生较大的匹配误差，不适用于高分辨率遥感影像道路提取。针对这一问题，国内外的学者提出了各种方法。这些方法主要通过引入几何约束条件^[12]、运用贝叶斯滤波^{[13, 14]}或神经网络^[15]等方式来优化匹配结果，提高道路追踪的稳健性，取得了一定的效果，但这些方法并没有从本质上解决模板匹配对遮挡敏感的问题。

通过分析传统模板匹配对遮挡敏感的原因，本文设计了一种基于均值漂移的道路中心点匹配算法。

1 基于均值漂移的道路中心点匹配

本文设计了一种基于均值漂移的道路中心点匹配算法，该算法采用一种稳健的概率相似性测度，能够在目标区域内根据道路中心点模板寻找相似度最大点，作为道路中心点，对遮挡不敏感。

1.1 核密度估计

核密度估计^[16]是在概率论中用来估计未知的概率密度函数，是一种非参数估计方法。假定n个数据点x_i,i=1,2，…,n在d维空间R^d，核密度估计采用如下形式

式中，(x)为x处的概率密度估值；h为带宽；K(x)为核函数。本文中使用的核函数为高斯核函数

1.2 概率相似性测度

相关系数是匹配算法中常用的相似性测度，但对粗差敏感，这导致很多道路追踪算法对车辆、树荫敏感，稳健性差。本文采用概率相似性测度^[17]，这种相似性测度对粗差具有稳健性，具体计算方式如下：

设样本点集I_x={x_i}_i=1ⁿ、I_y={y_j}_j=1^m分别为描述对象x和y的点集，n、m为点集中的样本数，那么x与y的相似度定义为

式中，K(x)为高斯核函数。相似性测度公式(3)的实质是利用式(1)计算对象y的每一个样本点在对象x中的概率密度估值，再将它们的平均值作为对象y与对象x的相似度，可以直观地理解为计算对象y是对象x的概率。相似性测度公式(3)由于使用了核函数，因而能够有效地限制粗差的影响。

1.3 均值漂移算法

均值漂移算法是一种步长自适应的梯度上升算法^{[18, 19, 20]}，能够快速地收敛于d维空间内一系列离散点对应的概率密度函数的局部极值。均值漂移向量的标准形式为

式中，m(x)为x处的均值漂移向量；K(x)为核函数；W(x)为权函数；S为样本点集，s为样本点。

1.4 道路中心点匹配的计算过程

道路中心点匹配的基本过程为：首先，给定一个道路中心点模板和一个匹配初始点；然后，算法以匹配初始点为中心建立目标区域；最后，利用均值漂移算法寻找目标区域内相似度最大的点，作为道路中心点。计算过程如下所示。

设x₀为道路中心点模板的中心点坐标，y₀为匹配初始点的坐标，I_x={x_i,u_i}_i=1ⁿ为道路中心点模板，I_y={y_j,v_j}_j=1^m为待匹配的目标区域，其中x_i、u_i为模板内样本点的坐标和颜色信息，y_j、v_j为目标区域内样本点的坐标和颜色信息，n和m均表示样本数。由式(3)知，初始点与道路中心点的相似度为

式中，K(x)和W(x)都为高斯核函数。

令初始点与道路中心点的相似度函数为

相似度函数l(y₀)的极大值点₀由梯度上升算法迭代求解，具体如下

而

式中，Δy_j=y_j-y₀；Δx_i=x_i-x₀。则

在梯度∇l(y₀)的基础上乘以正数h²不会改变梯度的方向，因而用∇L(y₀)取代∇l(y₀)不会改变梯度上升算法的性质。迭代求解₀的过程为

式(7)、式(8)中，∇L(y_0j)便是式(4)定义的均值漂移向量。文献[20]已证明，当核函数K(x)具有凸性并且它的剖面函数k(x)严格单调递减时，均值漂移算法能够快速地收敛于概率密度函数的局部极大值点。本文使用高斯核函数(2)，满足收敛条件，因而能够保证匹配初始点收敛于相似度函数l(y₀)(局部)极大值点。

基于均值漂移的道路中心点匹配算法根据公式(8)移动匹配初始点y₀，最终使y₀收敛于相似性函数l(y₀)的(局部)极大值点₀，准确匹配出道路中心点，如图 1(a)和1(b)所示。图 1(b)所示的曲面表示目标区域内每一点与道路中心点的相似度大小。

在图 1(a)所示的匹配过程中，匹配初始点y₀移动的步长与次数的关系如图 1(c)所示。从图 1(c)可以看出，初始点距离道路中心点较远时(移动次数小于16)，移动步长较大，随着初始点越来越接近道路中心点(移动次数超过16次)，移动步长迅速递减并趋于零，初始点y₀最终收敛于道路中心点。

2 结合卡尔曼滤波的道路中心线追踪

基于均值漂移的道路中心点匹配算法充分利用道路中心点模板信息和目标区域内的样本信息来寻找道路中心点，对车辆和树荫等遮挡具有稳健性，但没有考虑道路中心点的先验信息。因此，为了进一步提高道路中心线追踪的稳健性，本文在道路中心点匹配算法的基础上结合卡尔曼滤波，利用卡尔曼滤波的预测过程给道路中心点匹配算法提供匹配初始点，而匹配得到的道路中心点作为卡尔曼滤波的观测值，进而迭代追踪道路中心线。具体流程如图 2所示。

卡尔曼滤波^[21]是一种高效率的递归滤波器，它能够从一系列不完全且包含噪声的测量中，估计动态系统的状态。本文的卡尔曼滤波状态参数由道路中心点坐标(x_t,y_t)、中心处的道路方向Φ_t组成。在卡尔曼滤波预测过程中，t+dt处的状态向量预测值由t处的状态向量估值通过状态方程来得到，X _t+dt表示t+dt处的预测值，_t表示t处的状态估值。状态方程为

由式(9)可知，状态预测过程不属于线性过程，需要使用扩展卡尔曼滤波^[21]。预测状态向量的协方差矩阵为

式中，Φ_t+dt表示状态方程线性化后的系数矩阵；Q_t+dt表示系统噪声的协方差矩阵。在预测过程之后，则利用道路中心点匹配算法获取观测值，即

式中，Φ_t+dt表示t+dt处的道路方向观测值，可由道路中心点坐标观测值和上一状态的中心点坐标估值计算得到。则观测方程为

式中，z_t+dt表示t+dt处的状态向量。那么t+dt处的状态向量估值为

式中，K_t+dt为卡尔曼滤波增益，计算公式为

式中，R_t+dt为观测值的协方差矩阵。

由式(14)可知，系统误差的协方差Q和观测值的协方差R在卡尔曼滤波中发挥着非常重要的作用。在本文中，系统误差主要由实际道路的曲率决定，因而系统误差的协方差Q应根据实际道路曲率而定。对于观测值协方差R的确定，采用如下自适应策略：在每匹配一个道路中心点之后，算法将利用式(5)计算匹配点与道路中心点的相似度，根据相似性实时确定R，同时还结合道路中心线应为光滑曲线这一特征，即如果观测值中的方向与上一个道路方向估值的差值超过阈值，则增大R。

最后，利用式(9)—式(14)迭代追踪高分辨率遥感影像上的道路中心线。

3 试验与分析

本文设计了4组试验：第1组试验验证概率相似性测度的稳健性；第2组试验测试本文的道路中心点匹配算法的性能；第3组试验验证卡尔曼滤波的优化作用；第4组试验为高分辨率遥感影像道路中心线追踪试验。试验数据采用QuickBird影像。

3.1 相似性测度的稳健性测试

试验将本文所使用的相似性测度与相关系数对比，验证本文的相似性测度的稳健性，具体如下。

(1) 试验1。对目标旋转的稳健性测试。图 3(a)为目标模板，图 3(b)、(c)为待匹配图像。由图 3(b)、图 3(c)可知，待匹配图像相对目标模板发生了25°的旋转。试验运用两种相似性测度逐像素的计算待匹配图像内每一个位置与目标模板的相似度，然后选取相似度最大的点作为最佳匹配点。结果如图 3所示。

(2) 试验2。对遮挡的稳健性测试。图 4(a)为目标模板，图 4(b)、4(c)为待匹配图像，其中待匹配图像中的白色区域为遮挡，灰度值为255。过程与试验1相同，结果如图 4所示。

由图 3、4可知，当目标旋转或被遮挡时，相关系数无法准确度量出目标的真实位置，而本文所使用的相似度不受目标旋转、遮挡的影响，准确度量出目标的位置。因此，由本组试验可知，本文所使用的相似性测度比相关系数更稳健。

3.2 道路中心点匹配试验

(1) 试验3。测试基于均值漂移的道路中心点匹配算法对车辆、树荫等遮挡的稳健性。试验3在各种存在遮挡的道路场景中进行单个道路中心点匹配试验。试验中，先由人工给定道路中心点模板和一个匹配初始点，然后算法在以初始点为中心的目标区域内匹配道路中心点，结果如图 5所示。

从图 5(a)—图 5(e)可以看出，基于均值漂移的道路中心点匹配算法对车辆遮挡具有很强的稳健性。而当道路存在小面积树荫遮挡时，本文的道路中心点匹配算法能准确匹配出道路中心点，如图 5(f)-图 5(h)所示。但当道路存在大面积树荫遮挡时，道路中心点匹配算法出现误匹配，产生较大的匹配误差，如图 5(l)、图 5(m)所示。针对这些误匹配点，本文将在后续部分使用卡尔曼滤波对其进行修正。

(2) 试验4。将基于均值漂移的道路中心点匹配算法与经典的相关系数模板匹配进行对比。试验中，使用同一模板，分别利用两种匹配算法对一幅存在汽车遮挡的道路影像进行中心线追踪，具体过程为:首先，以道路中心线上一点为中心定义模板窗口，并将这一点作为道路中心线追踪的起始点；然后，根据起始点和道路方向预测下个道路中心点的初始位置；接着，在道路中心点初始位置运用匹配算法匹配道路中心点，并利用匹配得到的道路中心点继续预测下个道路中心点的初始位置；最后，重复前述过程，提取道路中心线。框架为道路中心点模板，试验结果如图 6所示。

由图 6(b)可知，利用本文的道路中心点匹配算法能够准确提取道路中心线，对车辆遮挡具有稳健性。对比图 6(a)可知，两种匹配算法在没有车辆遮挡的路段上都能准确匹配出道路中心点，但当路面出现车辆时，相关系数匹配受车辆的影响会产生较大的匹配误差。而基于均值漂移的道路中心点匹配算法能够排除车辆的影响，准确匹配道路中心点。

3.3 卡尔曼滤波优化试验

试验5，目的是为了验证卡尔曼滤波对道路中心点匹配算法的优化作用。试验过程为：①单独利用基于均值漂移的道路中心点匹配算法对两幅存在树荫遮挡的道路影像进行中心线追踪；②将基于均值漂移的道路中心点匹配算法与卡尔曼滤波结合对同样两幅影像进行道路中心线追踪。本文使用中心线提取的位置偏差作为精度衡量指标^[3]，道路中心线提取和精度统计结果分别如图 7和表 1所示。

由图 7(a)可知，仅使用道路中心点匹配算法提取的道路中心线在树荫遮挡路段明显偏离了实际道路中心线，存在较大匹配误差；而在图 7(b)中，在原有的道路中心点匹配算法基础上结合卡尔曼滤波后，修正了在树荫遮挡路段的误匹配。在图 7(c)中，仅使用道路中心点匹配算法提取的道路中心线受大面积树荫的影响出现一定的波动；而在图 7(d)中，在结合卡尔曼滤波后算法有效地排除了树荫的影响。由表 1可知，在结合卡尔曼滤波后道路中心线的提取精度得到显著提高。

由本试验可知，卡尔曼滤波能够有效地修正由道路中心点匹配算法产生的误匹配，提高道路中心线追踪的稳健性和精度。

3.4 道路中心线追踪试验

试验6，选取存在车辆、树荫遮挡的道路影像对算法进行测试。试验道路分为高速公路和乡村道路，其中高速公路上存在大量的车辆和树荫遮挡，在乡村道路上存在树枝遮挡。试验中，需由人工选定道路中心点模板和指定道路初始方向，试验结果如图 8所示。

从图 8可知，无论是高速公路还是乡村道路，本文的道路中心线追踪方法都能够排除车辆、树荫等遮挡的影响，提取道路中心线。因此，本文所提出的道路中心线追踪方法是一种稳健的道路追踪算法。

为了进一步测试算法的性能，试验7选取尺寸为2428×2504像素的影像进行试验。影像中有两条相交的道路，总长约为3 km，路面上存在较多车辆、树荫遮挡。试验采用人机交互的方式，由人工输入一系列种子点，算法以第一个种子点为起始点，根据当前种子点和下一个种子点计算道路初始方向，追踪至下一个种子点时更换道路中心点模板，如图 9所示。图 9(c)中的黑色框架为道路中心点模板。

试验7中，人工总共输入6个种子点，算法共提取1283个道路中心点，提取的道路中心点的最大位置偏差为4像素，平均位置偏差为1.0像素，位置偏差的中误差为0.8像素，如图 9所示。由图 9(b)可知，本文提出的道路追踪算法能够排除车辆、树荫的干扰，准确提取道路中心线。当道路一侧出现分岔口时，算法能够不受分岔口影响，继续按原先道路方向追踪，而当道路A追踪完毕后，再由人工重新输入种子点，更换模板，追踪道路B，如图 9(c)所示。本文算法运用卡尔曼滤波，结合道路先验信息，因此能够适应道路宽度发生短暂变化的情况，如图 9(d)所示。

4 结 论

本文提出了一种基于均值漂移和卡尔曼滤波的高分辨率遥感影像道路中心线追踪方法。该方法是一种对车辆、树荫等遮挡稳健的模板匹配方法，首先利用概率相似性测度和均值漂移算法设计了一种稳健的道路中心点匹配算法，取代了传统追踪算法中常用的相关系数匹配；然后结合卡尔曼滤波实现道路中心线的追踪。试验表明，该方法能够排除车辆、树荫遮挡的影响并能适应道路一侧存在分岔口的情况，对高分辨率道路影像取得较好的效果。但该方法只适用于道路存在明显边缘的道路影像，对道路形态的一致性有一定的要求，同时目前仍需人工提供道路中心点模板，且每次只能提取模板所在的道路。因此为了提高方法的自动化程度，下一步将研究模板的自动选取和更新。