CN102831423B - 一种sar图像道路提取方法 - Google Patents

Abstract

一种SAR图像道路提取方法，基于均值比率的线目标检测算子提取SAR图像上道路的特征，获得道路线特征响应图和道路线特征方向图；在SAR图像的道路上选择一个点作为粒子滤波的起始点；基于起始点，采用粒子滤波方法跟踪出道路上的种子点；以道路线特征响应图作为Snake模型中的图像能量，以种子点作为Snake模型的初始值进行拟合，得到一条道路目标。本发明结合粒子滤波和Snake模型，考虑道路在SAR影像上的方向、后向散射和边缘等特征因素解决道路提取问题，并且改进基于均值比率的线目标检测算子，精度好，准确率高。

Description

一种SAR图像道路提取方法

技术领域

本发明涉及遥感影像处理技术领域，尤其是涉及一种结合粒子滤波和Snake模型的SAR图像道路提取方法。

背景技术

SAR(合成孔径雷达)传感器的全天时、全天候工作能力使其成为光学传感器的有效补充，SAR图像道路提取在战时、灾后等恶劣天气情况上发挥着重要作用。粒子滤波和Snake模型是两种比较好的道路目标提取算法，但两者都有各自的优缺点。粒子滤波受初始点影响不大，只需一个起始点即可跟踪出整个线目标，但其受粒子权重和跟踪步长影响较大；而Snake模型能实现准确定位，保持线条光滑，但其对初始位置敏感，需要在关键位置选取高精度控制点才行。并且，SAR传感器特殊的相干距离成像方式，使得SAR图像与光学图像具有明显的差异，通常现有的光学图像道路提取方法难以直接用于SAR图像。

发明内容

针对上述问题，本发明提出了一种结合粒子滤波和Snake模型的SAR图像道路提取方法，以拟合出平滑、连续的道路目标。

本发明的技术方案为一种SAR图像道路提取方法，包括以下步骤：

步骤1，基于均值比率的线目标检测算子提取SAR图像上道路的特征，所述特征包括响应值和方向，获得道路线特征响应图和方向图；

步骤2，在SAR图像的道路上选择一个点作为粒子滤波的起始点；

步骤3，基于起始点，采用粒子滤波方法跟踪出道路上的种子点；

步骤4，以步骤1得到的道路线特征响应图作为Snake模型中的图像能量，以离散的道路点作为Snake模型的初始值，采用Snake模型进行拟合，得到一条道路目标；所述离散的道路点包括步骤2所得起始点和步骤3所得种子点。

而且，步骤3采用粒子滤波方法时，采用以下公式得到任一粒子W的权重w_t，

w_t＝αw₁-βw₂-γw₃

其中，α,β，γ分别是权重因子w₁、w₂、w₃的系数，α,β，γ均大于0；

$w_1=\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\σ}}{e}^{\frac{-{(1-r)}^2}{2{\mathrm{\σ}}^2}}$

其中，r是以粒子W为检测点采用基于均值比率的线目标检测算子提取得到的响应值，σ为方差；

w₂＝|θ_t-1-θ_t|

其中，θ_t-1是前一时刻得到的道路种子点的方向，θ_t是当前时刻以粒子W为检测点采用基于均值比率的线目标检测算子提取得到的方向；

$w_3=|\stackrel{\‾}{g}-g_t|$

其中，g_t表示当前时刻粒子W对应的后向散射值，表示当前时刻前跟踪得到的所有道路种子点的后向散射值的均值。

而且，所述基于均值比率的线目标检测算子采用D1或D2，基于均值比率的线目标检测算子提取SAR图像上道路的特征的实现过程中，计算某检测点在某个方向的响应值前，先判断检测点所在区域的后向散射值均值是否均小于两边区域的后向散射值均值，如果是则继续计算该检测点在这个方向的响应值，如果不是则直接设该检测点在这个方向的响应值为0。

本发明将线检测算子得到的响应值作为Snake模型中的图像能量，使Snake模型更适合于SAR图像道路提取。而且，本发明还提出利用道路的后向散射特性，改进现有的基于均值比率的SAR图像线特征提取算法，降低亮线目标提取对道路的虚警率；在计算粒子权重时，综合考虑道路的方向、后向散射和边缘特征，提高粒子滤波跟踪的精度。

附图说明

图1为本发明实施例的流程图。

图2为本发明实施例在竖直方向的9×9检测窗口示意图。

图3为本发明实施例的线检测方向示意图。

图4为本发明实施例的待检测SAR图像。

图5为本发明实施例的道路线特征响应值。

图6为本发明实施例中提取出来的道路在图4上的叠加显示示意图。

具体实施方式

本发明首先利用改进的线检测算子来获得道路的特征（包括道路响应和道路方向），然后基于蒙特卡洛原理，采用粒子滤波方法跟踪出道路的种子点，最后利用这些点作为Snake模型的初始控制点，通过寻找最小能量使之收敛到SAR图像中道路目标的真实位置，来拟合出平滑、连续的道路目标。

本发明技术方案可采用计算机软件技术实现自动运行流程。以下结合附图和实施例详细说明本发明技术方案。

如图1，实施例包括以下步骤：

步骤1，基于均值比率的线目标检测算子提取SAR图像上道路的特征，所述特征包括响应值和方向，获得道路线特征响应图和道路线特征方向图。

基于均值比率的线目标检测算子为现有技术，法国Florence Tupin的研究小组专门针对线型目标提出了基于均值比率的局部线特征检测算子D1和D2，这里以D2算子为例来说明。

基于均值比率的线特征检测算法是利用检测窗口在待检测的SAR图像上进行卷积运算，获取线特征的响应值和方向。如图2是一个竖直方向9×9像素大小的检测窗口，将检测窗口应用在SAR图像上，取窗口中心像素点为待检测点（图中黑色实心点为中心像素点），沿着过该点某直线的方向（图2所示为竖直方向，即为图3所示的方向“4”），当窗口在图像上滑动时，将检测窗划分为不重叠的3个区域，首先计算标号为1，2，3各区域内所有点的灰度强度均值μ₁，μ₂，μ₃和各区域内所有点的灰度强度的标准差σ₁，σ₂，σ₃。

定义R_ij为检测器关于区域i和区域j的D2线特征检测算子的输出响应：

$R_{\mathrm{ij}}=\sqrt{(1/(1+(n_i+n_j)*(n_i*{v_i}^{*}*\frac{{{\mathrm{\μ}}_i}^2}{{{\mathrm{\μ}}_j}^2}+n_j*{v_j}^2)/(n_i*n_j*{(\frac{{\mathrm{\μ}}_i}{{\mathrm{\μ}}_j}-1)}^2)))}---\left(1\right)$

其中i,j分别表示两相邻区域的标号，n_i，n_j分别表示两个区域i,j内像素的个数，μ_i，μ_j分别表示两个区域i,j内所有点的灰度强度均值。v_i,v_j分别表示两个区域i,j的局域方差系数，即

$\left\{\begin{array}{c}{v}_i={\mathrm{\σ}}_i/{\mathrm{\μ}}_i\\ v_j={\mathrm{\σ}}_j/{\mathrm{\μ}}_j\end{array}\right.---\left(2\right)$

根据式(1)分别计算区域2与区域1，3这两个区域的线特征检测响应值R₁₂和R₂₃。取线特征检测响应值R₁₂和R₂₃之间的最小值作为该方向的响应值R，即

R＝min(R₁₂,R₂₃)                (3)

考虑到线特征的不同取向，按图3所示8个方向0、1、2、3、4、5、6、7，对每个检测点构建对应方向的检测窗口，计算每个方向的响应值R，保留8个R值中最大的结果作为该点的线特征响应值r，其对应的方向为该点的线特征方向。

在现有技术的基于均值比率的线目标检测算子基础上，本发明提出了改进。由于通常道路属于镜面反射，后向散射系数较小，在SAR图像上表现为暗线目标，本发明提出在计算某个方向的响应值前，根据该方向相应对检测窗划分所得区域的后向散射值均值进行判断：先判断该检测点所在区域的后向散射值均值是否均小于两边区域的后向散射值均值，如果不是均小于，则直接认为该检测点在这个方向上不是道路特征点，即该检测点在这个方向上的道路特征值为0，即：

r＝0 if μ₂≥μ₁ orμ₂≥μ₃            (4)

否则该检测点在这个方向上可能是道路特征点，则继续根据式（1）、（2）、（3）计算该检测点在这个方向上响应值。这样可以避免亮线目标检测得到高响应值而带来的高虚警率，同时降低计算量。

根据基于均值比率的线目标检测算子提取SAR图像上道路的特征所得各检测点的响应值和方向，即可得到道路线特征响应图和道路线特征方向图。

步骤2，在SAR图像的道路上选择一个点作为粒子滤波的起始点。

具体实施时，只需在SAR图像的任意一条道路上选择一个点作为粒子滤波的起始点，该点可以是道路的起始点或者终止点，也可以是道路上的任意一点。具体实施时，可以采用人工或依据一定的规则在图像上自动选择道路上的点。

步骤3，基于起始点，采用粒子滤波方法跟踪出道路上的种子点。

以步骤2所得起始点，给定一个步长，例如8～10个像素，基于粒子滤波方法，跟踪出道路上的种子点。粒子滤波方法为现有技术，粒子滤波利用了随机样本来描述目标状态的后验概率密度函数，这些样本被称为“粒子”，通过调节各粒子权值的大小和样本的位置，来近似实际概率分布，并以样本的均值作为目标的估计值。也就是说，道路的每个种子点都是通过所有粒子加权平均得到的。而粒子的权重描述了该粒子的重要程度，用来衡量每个粒子与真值的相似性。本发明进一步提出，在计算每个粒子权重时，可以考虑道路的三个特征：

1）道路点的线检测响应值大于非道路点的响应值。由于线检测响应是一个0到1之间的实数，并且r值越大，该点属于道路的可能性越大。然而，如果直接使用r值作为权值，则各个粒子之间的区分度较小，不能够对粒子进行有效地选择，因此实施例用方差参数为σ的高斯分布函数对线检测响应强度来定义第一个权重因子：

$w_1=\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\σ}}{e}^{\frac{-{(1-r)}^2}{2{\mathrm{\σ}}^2}}---\left(5\right)$

其中，r是以粒子W为检测点采用基于均值比率的线目标检测算子得到的响应值，从步骤1所得道路线特征响应图中即可取得；

2）道路在局部范围内不会出现较大的弯曲，其形状近似于直线，道路的方向信息基本一致，即使有一点弧度，其变化范围也很小。因此，第二个权重因子定义为前一时刻的方向与当前时刻方向的差值，即：

w₂＝|θ_t-1-θ_t|            (6)

式中，θ_t-1是前一时刻得到的道路种子点的方向，θ_t当前时刻以粒子W为检测点采用基于均值比率的线目标检测算子提取得到的方向，实施例就是当前时刻该粒子所在像素经D2算子检测处理得到的方向。这些方向信息从步骤1所得道路线特征方向图中即可取得。

3）道路内部后向散射比较均匀、变化较小。因此，第三个权重因子定义为当前时刻粒子对应的后向散射值与前面道路种子点的后向散射值均值的差值，即：

$w_3=|\stackrel{\‾}{g}-g_t|---\left(7\right)$

式中，g_t表示当前时刻该粒子对应的后向散射值，表示当前时刻前，跟踪得到的所有道路种子点的后向散射值均值。

将这三个特征各自构成的权重归一化，然后通过加权平均来获得每个粒子的权重。

w_t＝αw₁-βw₂-γw₃                (8)

式中，α,β，γ是不同权重因子的系数，α,β，γ均大于0。具体实施时，α,β，γ的取值可以由本领域技术人员根据具体情况自行设定。由于第一个特征是最重要的特征，因此建议设定对应的系数α最大，而β大于γ。

步骤4，以步骤1得到的道路线特征响应图作为Snake模型中的图像能量，以离散的道路点作为Snake模型的初始点，采用Snake模型进行拟合，得到一条道路目标。所述离散的道路点包括步骤2所得起始点和步骤3所得种子点，采用Snake模型进行拟合可以得到一条光滑、连续的道路目标。本步骤基于Snake模型的初始点，将步骤1通过改进算子得到的道路线特征响应图作为Snake模型中的图像能量，通过迭代使整个Snake能量函数最小化，即可提取出高精度的道路目标。Snake模型具体实现为现有技术，本发明不予赘述。

上述步骤是针对图像上只有一条道路的情况，如果SAR图像上有多条道路，可分别在每条道路上人工或依据一定的规则选择一个点作为该条道路的起始点，通过多个粒子滤波跟踪可以提取出图像上的多条道路。

为理解本发明效果参考起见，采用实施例所提供流程对图4所示待检测SAR图像提取道路特征，得到如图5所示道路线特征响应值，根据粒子滤波所得离散道路点拟合提取出来的道路如图6所示，白色实线为提取出来的道路。

本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (2)

1.一种SAR图像道路提取方法，其特征在于，包括以下步骤： 

步骤1，基于均值比率的线目标检测算子提取SAR图像上道路的特征，所述特征包括响应值和方向，获得道路线特征响应图和道路线特征方向图； 

步骤2，在SAR图像的道路上选择一个点作为粒子滤波的起始点； 

步骤3，基于起始点，采用粒子滤波方法跟踪出道路上的种子点；采用粒子滤波方法时，通过以下公式得到任一粒子W的权重w_t， 

w_t＝αw₁-βw₂-γw₃

其中，α,β,γ分别是权重因子w₁、w₂、w₃的系数，α,β,γ均大于0； 

其中，r是以粒子W为检测点采用基于均值比率的线目标检测算子提取得到的响应值，σ为方差； 

w₂＝|θ_t-1-θ_t| 

其中，θ_t-1是前一时刻得到的道路种子点的方向，θ_t是当前时刻以粒子W为检测点采用基于均值比率的线目标检测算子提取得到的方向； 

其中，g_t表示当前时刻粒子W对应的后向散射值，表示当前时刻前跟踪得到的所有道路种子点的后向散射值的均值；

步骤4，以步骤1得到的道路线特征响应图作为Snake模型中的图像能量，以离散的道路点作为Snake模型的初始值，采用Snake模型进行拟合，得到一条道路目标；所述离散的道路点包括步骤2所得起始点和步骤3所得种子点。 

2.如权利要求1所示SAR图像道路提取方法，其特征在于：所述基于均值比率的线目标检测算子采用D1或D2，基于均值比率的线目标检测算子提取SAR图像上道路的特征的实现过程中，计算某检测点在某个方向的响应值前，先判断检测点所在区域的后向散射值均值是否均小于两边区域的后向散射值均值，如果是则继续计算该检测点在这个方向的响应值，如果不是则直接设该检测点在这个方向的响应值为0。