CN105160355B - 一种基于区域相关和视觉单词的遥感图像变化检测方法 - Google Patents

Abstract

本发明公开了一种基于区域相关和视觉单词的遥感图像变化检测方法，具体如下：首先选取参考图像和待检测图像；对参考图像和待检测图像分别进行直方图均衡化处理；对处理后的参考图像和待检测图像分别进行形态学重建的分水岭分割；对分割后得到的两幅图像进行融合处理，获得融合图像；融合图像中共包括k个子区域，对其中的每一个子区域j，计算加权相关系数C_j；设定加权相关系数阈值Δ，若第j个子区域的加权相关系数C_j大于Δ，则第j个子区域即为不变区域，直接将该不变区域提取出来；若第j个子区域的加权相关系数C_j小于Δ，则通过基于视觉单词理论的方法进行提取。使用该方法能够在满足检测任务的基础上提高检测的实时性和鲁棒性。

Description

一种基于区域相关和视觉单词的遥感图像变化检测方法

技术领域

本发明属于遥感图像处理技术领域，具体涉及一种基于区域相关和视觉单词的变化检测方法。

背景技术

多时相遥感图像变化检测是遥感技术应用的一个重要方向。针对同一区域在不同时间获取的遥感影像，自动提取两幅图像中的不同区域并标记出来。这项技术在民用和军用领域都发挥着重要的作用，如防灾减灾、地震预警和战场毁伤评估等。

目前在变化检测方面已经开展了许多研究工作，提出了多种变化检测方法。常见的变化检测方法可以分为三类：第一类是直接比较法，如影像差值法，相关系数法，变化向量分析法等；第二类是分类法，如分类后比较法，多时相直接分类法等；第三类是模型法，如预测模型法，统计模型法等。

直接比较法主要是应用两幅图像的代数运算进行变化检测的方法。影像差值法是传统的变化检测方法，通过设定差分图像灰度阈值得到变化区域，但通常效果不理想；相关系数法利用变化区域相关性小的特点提取变化信息，复杂度低，不足之处是需要选择合理的计算区域；而分类法可以获取变化的类别信息，但是分类精度直接影响变化检测的结果，并且通常监督分类器需要准确的变化样本，这种样本较难提取，影响了此种方法的实时性。模型法的优点是可以有效克服噪声和光照影响，简化复杂问题。其缺点是建模困难，并且一般的模型如马尔科夫模型等都比较复杂，计算量非常大，不能满足快速实时的处理要求。

综上，变化检测的自动快速高精度处理成为现今亟待解决的问题。

发明内容

有鉴于此，本发明提供了一种基于区域相关和视觉单词的变化检测方法，通过分割选择出合理的相关计算区域，使用加权相关系数提出可靠度高的不变区域，并且引入视觉单词理论提取变化区域，可以在满足检测任务的基础上提高检测的实时性和鲁棒性。

为了达到上述目的，本发明的技术方案包括如下步骤：

步骤一：选取参考图像和待检测图像；其中待检测图像与参考图像中场景相同。

步骤二：对参考图像和待检测图像分别进行直方图均衡化处理。

步骤三：对处理后的参考图像和待检测图像分别进行形态学重建的分水岭分割。

步骤四：对分割后得到的两幅图像进行融合处理，获得融合图像。

步骤五：融合图像中共包括k个子区域，对其中的每一个子区域j，计算加权相关系数C_j，

其中，是第j个子区域相关系数的权重，N_j是第j个子区域的像素个数，N是所有子区域像素个数的总和，R_1jm为参考图像中第j个子区域中第m个像素的灰度值，R_2jm为待检测图像中第j个子区域中第m个像素的灰度值，是参考图像中第j个子区域的灰度均值；是待检测图像中第j个子区域的灰度均值。

设定加权相关系数阈值Δ，若第j个子区域的加权相关系数C_j大于Δ，则第j个子区域即为不变区域，直接将该不变区域提取出来。

若第j个子区域的加权相关系数C_j小于Δ，则通过基于视觉单词理论的方法进行提取。

进一步地，以参考图像或待检测图像为原图像，步骤二中直方图均衡化处理的过程具体为：

步骤201)统计原图像中各灰度级出现的次数，计算每个灰度级出现的概率，其中第l个灰度级出现的概率为P_l(r_l)＝N_l/N，l＝0,1,2,…L-1

其中，N_l为第l个灰度级像素个数，N为图像像素总个数，L为图像中出现的灰度级总数；

步骤202)计算图像的累积分布函数S₀～S_L：l＝0,1,2,…L-1

步骤203)计算直方图均衡化处理后的结果S₀'～S_L'：S_l'＝round(S_l*L+0.5)；其中，round为取整操作。

进一步地，以步骤二中直方图均衡化处理后的参考图像和待检测图像作为本步骤的待处理图像，步骤三中形态学重建的分水岭分割具体为：

步骤301)用形态学算子对待处理图像进行腐蚀操作，用腐蚀前的图对腐蚀后的图进行重建；

步骤302)用形态学算子对步骤301)所得的重建图进行膨胀操作，接着用步骤301)所得的重建图对膨胀操作所得的图进行重建；

步骤303)用梯度算子计算待处理图像的梯度图像；

步骤304)用步骤302)所得的重建图对步骤303)得到的梯度图像进行极小值标记，最后用分水岭算法对标记后的梯度图像进行分割处理；

进一步地，以步骤三中分割后的参考图像和待检测图像为待融合图像，步骤四中的融合处理过程具体为

步骤401)对待融合图像进行二值化处理，二值化规则为将分水岭上的像素置为1，非分水岭上的像素置为0；

步骤402)对二值化后的两幅图像进行逻辑或操作，即得到分割区域融合后的结果。

进一步地，形态学算子选择圆盘型的结构。

进一步地，基于视觉单词理论的方法包括局部特征提取、字典生成、视觉单词词频统计以及使用Bhattacharyya距离提取变化区域。

有益效果：

1、本发明通过计算分割融合子区域的加权相关系数和视觉单词直方图相似度来提取变化区域，结合了两种方法特点，能带来更高的检测精度。加权相关系数的提出，衡量了相关系数的可靠程度，增强了合理性。

2、本发明对参考图像和待检测图像进行直方图均衡化处理，减少了光照因素对检测结果的影响。

3、选择基于形态学重建的分水岭算法对参考图像和待测试图像进行分割，对梯度图像的极小值进行抑制，解决了分水岭算法的过分割问题，得到目标级的分割轮廓，使融合子区域包含了目标级的信息；同时可以调整形态学算子的大小以适应不同分辨率的图像，提高了检测的鲁棒性。

4、本发明在融合分割区域之前，首先对分割结果进行了二值化处理，保留了变化的有效信息的同时大大降低了融合的运算量。

5、本发明提供了针对不规则区域的视觉单词特征描述的提取方法，具有一定的应用价值。

附图说明

图1是变化检测流程示意图；

图2是子区域融合示意图；

图3本发明中的词频统计决策机制。

具体实施方式

下面结合附图并举实施例，对本发明进行详细描述。

如图1所示，基于区域相关和视觉单词的变化检测方法，具体步骤为：

步骤一、选取参考图像和待检测图像。待检测图像必须和参考图像是同一个场景，并且待检测图像和参考图像的配准误差要在2个像素以内。

步骤二，分别对参考图像和待检测图像进行直方图均衡化处理。

下面对直方图均衡化技术进行说明：

直方图均衡化技术是一种图像变换方法，可以使图像的灰度值均匀分布在所有灰度级上，增强图像对比度，同时对参考图像和待测图像进行处理可以均衡两者的光照强度。具体实施过程如下：

步骤201)统计原图像中各灰度级出现的次数，计算每个灰度级出现的概率P_l(r_l)＝N_l/N，l＝0,1,2,…L-1

其中，表示第l个灰度级出现的概率，N_l为第l个灰度级像素个数，N为图像像素总个数，L为图像中出现的灰度级总数；

步骤202)计算图像的累积分布函数：l＝0,1,2,…L-1

步骤203)S_l'＝round(S_l*L+0.5)；其中，round为取整操作，则S_l'为直方图均衡化处理结果。

步骤三、分别对处理后的参考图像和待检测图像进行形态学重建的分水岭分割。下面对形态学重建的分水岭分割算法进行说明：

以步骤二中直方图均衡化处理后的参考图像和待检测图像作为本步骤的待处理图像，步骤三中形态学重建的分水岭分割具体为，：

步骤301)用形态学算子对待处理图像进行腐蚀操作，用腐蚀前的图对腐蚀后的图进行重建；

步骤302)用形态学算子对步骤301)所得的重建图进行膨胀操作，接着用步骤301)所得的重建图对膨胀操作所得的图进行重建；

步骤303)用梯度算子计算待处理图像的梯度图像；

步骤304)用步骤302)所得的重建图对步骤303)得到的梯度图像进行极小值标记，最后用分水岭算法对标记后的梯度图像进行分割处理。

在本实例中，需要分别对直方图均衡化后的参考图像和待检测图像进行两次重建操作，第一次是腐蚀重建，第二次是膨胀重建。本实例使用圆盘型，大小为4的腐蚀和膨胀算子。

2)用梯度算子计算出参考图像和待检测图像的梯度图像。本实例使用sobel算子进行计算。

3)使用1)中重建后的图像对参考图像和待检测图像进行局部极小值标记，对标记后的图像使用分水岭分割算法。分水岭算法的原理是每一个局部极小值及其影响区域可以看作集水盆，而集水盆的边界则形成分水岭。它的一个的经典实现过程是由L.Vincent提出的模拟浸润过程，它包括了两个步骤：首先对每个像素的灰度级进行从低到高排序，然后在从低到高实现淹没过程中，对每一个局部极小值在h阶高度的影响域采用先进先出结构进行判断及标注。这样分割得到的各个子区域都有着封闭的轮廓。

步骤四、以步骤三中分割后的参考图像和待检测图像为待融合图像，步骤四中的融合处理过程具体为

步骤401)对待融合图像进行二值化处理，二值化规则为将分水岭上的像素置为1，非分水岭上的像素置为0；

步骤402)对二值化后的两幅图像进行逻辑或操作，即得到分割区域融合后的结果。

图2显示了子区域融合的示意图，假定2(a)为参考图像分割结果，2(b)为待测图像分割结果，那么区域融合结果为2(c)。

步骤五、鉴于分割融合后每个子区域的形状不规则，并且大小一般也不相同，传统的相关系数法的合理性大大降低，这里提出一种针对本发明的加权相关系数法。设融合结果T中一共有k个子区域，对应于T中的k个子区域，对参考图像和待检测图像做相关处理，计算加权相关系数，加权相关系数表达式如下：

其中，是第j个子区域相关系数的权重，衡量了该区域相关系数的可靠程度；由于增加了权重这一项，计算出来的绝对值较大的加权相关系数相对于原始的相关系数值的可靠程度更高，亦即根据加权相关系数判定区域为不变区域的可信度更高(基本判定规则为加权相关系数大于阈值的为不变区域)。N_j是第j个子区域的像素个数，N是所有子区域像素个数的总和，C_j是第j个子区域的加权相关系数，R_ijm(i＝1为参考图像，i＝2为待检测图像)是第j个子区域中第m个像素的灰度值，是第j个子区域的灰度均值。

计算出各子区域加权相关系数后，设定加权相关系数阈值Δ，大于Δ的区域判定为不变区域，由于这里加入了权重系数w_j，因此加权相关系数C_j小的区域并不意味着该区域变化的可能性大(还可能是因为权重小)，也就是说通过加权相关系数法可以更可靠地判定出不变的区域，这里并不使用加权相关系数法的提取变化的区域，原因如前所述。接下来本发明通过视觉单词理论提取位于小于Δ的区域中的变化区域。

其中，对于加权相关系数小于Δ的区域，使用视觉单词理论提取变化区域，具体步骤如下：

1)局部特征提取：特征提取是将图像小块(即单词)转变为特征空间多维向量的过程，本发明采用16*16的窗口对整幅参考图像和待测图像进行均匀分块，使用SIFT特征作为图像小块的特征描述，得到特征集合X＝{x₁,x₂,...,x_N}，其中N为图像小块个数，x_N是第N个图像小块的SIFT描述子，是128维的向量。

2)字典生成：字典生成过程是将所有图像小块用视觉单词表示的过程，一个视觉单词就是许多相似的图像小块的聚类中心。本发明采用无监督K-means聚类算法对直方图均衡化之后的参考图像和待测图像的图像小块特征，也就是集合X的元素，进行聚类，得到K个聚类中心：V＝[v₁,v₂,...,v_K]，即K个视觉单词，V被称为词典。

3)视觉单词词频统计：分别计算参考图像和待测图像对应于加权相关系数小于Δ的区域内的局部特征，每一个子区域一般包含多个16*16窗口，因此存在多个局部特征描述，形成一个特征描述集：

D_i＝(x₁,x₂,...,x_M)，其中M是第D_i个子区域的特征个数。接着对于子区域中每一个特征x_m，找到词典中与其距离最近的视觉单词，即这样，每一个子区域都对应于M个视觉单词。接下来统计M个单词中v₁到v_k各自出现的次数，形成一个词频直方图。特别地，步骤1)和步骤2)与经典的视觉单词理论一致，但在词频统计这一步骤，由于加权相关系数小于Δ的子区域的形状并不是规则的矩形，因此无法直接用标准16*16矩形窗口来统计词频。这里给出一种针对本发明的决策机制：如图3所示，在对每一个加权相关系数小于Δ的子区域A进行局部特征描述时，先求得该子区域的最大内接矩形，统计最大内接矩形的特征描述集，接着在余下的区域里寻找最大内接矩形，统计特征描述集，不停的寻找，直到矩形小于16*16不能统计为止，合并所有内接矩形的特征描述集作为该子区域特征描述集。此决策机制会损失一小部分的特征描述信息，但提供了一种针对不规则形状区域的词频统计解决方案。

4)Bhattacharyya距离提取变化区域：对应于加权相关系数小于Δ之间的每一个子区域，计算参考图像和待测图像的词频直方图，使用Bhattacharyya距离描述直方图之间的相似性。设定Bhattacharyya距离阈值δ，大于δ的说明直方图相似性很高，亦即此区域变化可能小，判定为不变区域，小于δ的区域判定为变化区域，这样加权相关系数小于Δ的子区域的变化区域就被提取出来了。

综上，以上仅为本发明的较佳实施例而已，并非用于限定本发明的保护范围。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (5)

1.一种基于区域相关和视觉单词的遥感图像变化检测方法，其特征在于，该方法包括如下步骤：

步骤一：选取参考图像和待检测图像；其中待检测图像与参考图像中场景相同；

步骤二：对参考图像和待检测图像分别进行直方图均衡化处理；

步骤三：对处理后的参考图像和待检测图像分别进行形态学重建的分水岭分割；

步骤四：对分割后得到的两幅图像进行融合处理，获得融合图像；

步骤五：所述融合图像中共包括k个子区域，对其中的每一个子区域j，计算加权相关系数C_j，

其中，是第j个子区域相关系数的权重，N_j是第j个子区域的像素个数，N是所有子区域像素个数的总和，R_1jm为参考图像中第j个子区域中第m个像素的灰度值，R_2jm为待检测图像中第j个子区域中第m个像素的灰度值，是参考图像中第j个子区域的灰度均值；是待检测图像中第j个子区域的灰度均值；

设定加权相关系数阈值Δ，若第j个子区域的加权相关系数C_j大于Δ，则第j个子区域即为不变区域，直接将该不变区域提取出来；

若第j个子区域的加权相关系数C_j小于Δ，则通过基于视觉单词理论的方法进行提取；

所述基于视觉单词理论的方法包括局部特征提取、字典生成、视觉单词词频统计以及使用Bhattacharyya距离提取变化区域。

2.如权利要求1所述的一种基于区域相关和视觉单词的遥感图像变化检测方法，其特征在于，以参考图像或待检测图像为原图像，所述步骤二中直方图均衡化处理的过程具体为：

步骤201)统计原图像中各灰度级出现的次数，计算每个灰度级出现的概率，其中第l个灰度级出现的概率为P_l(r_l)＝N_l/N，l＝0,1,2,…L-1 ；

其中，N_l为第l个灰度级像素个数，N为图像像素总个数，L为图像中出现的灰度级总数；

步骤202)计算图像的累积分布函数S₀～S_L-1：

步骤203)计算直方图均衡化处理后的结果S₀'～S_L-1'：S_l'＝round(S_l*L+0.5)；其中，round为取整操作。

3.如权利要求1或2所述的一种基于区域相关和视觉单词的遥感图像变化检测方法，其特征在于，以步骤二中直方图均衡化处理后的参考图像和待检测图像作为本步骤的待处理图像，所述步骤三中形态学重建的分水岭分割具体为：

步骤301)用形态学算子对待处理图像进行腐蚀操作，用腐蚀前的图对腐蚀后的图进行重建；

步骤302)用形态学算子对步骤301)所得的重建图进行膨胀操作，接着用步骤301)所得的重建图对膨胀操作所得的图进行重建；

步骤303)用梯度算子计算待处理图像的梯度图像；

步骤304)用步骤302)所得的重建图对步骤303)得到的梯度图像进行极小值标记，最后用分水岭算法对标记后的梯度图像进行分割处理。

4.如权利要求3所述的一种基于区域相关和视觉单词的遥感图像变化检测方法，其特征在于，以所述步骤三中分割后的参考图像和待检测图像为待融合图像，所述步骤四中的融合处理过程具体为：

步骤401)对待融合图像进行二值化处理，二值化规则为将分水岭上的像素置为1，非分水岭上的像素置为0；

步骤402)对二值化后的两幅图像进行逻辑或操作，即得到分割区域融合后的结果。

5.如权利要求3所述的一种基于区域相关和视觉单词的遥感图像变化检测方法，其特征在于，所述形态学算子选择圆盘型的结构。