CN103778413B - 一种遥感影像欠分割对象自动识别方法 - Google Patents

Abstract

本发明公开了一种遥感影像欠分割对象识别方法，包括如下步骤：对图像进行数据降维获取分割对象；将所述分割对象进行聚类；根据聚类计算所述分割对象的混杂度指标；根据所述混杂度指标在所述分割对象中识别欠分割对象。本发明的判定方法结合空间纹理和光谱维信息，所得到的图像分割结果更符合真实地物的分布情况，提高了识别欠分割对象的速度和识别结果的准确率。

Description

一种遥感影像欠分割对象自动识别方法

技术领域

本发明属于遥感数据处理领域，涉及一种高分辨率遥感影像欠分割对象自动识别方法。

背景技术

图像分割是计算机视觉研究中的一个经典难题，是实现图像对象信息提取的重要环节，通过图像分割生成图像对象后，进一步的图像分析(例如面向对象的特征提取和目标识别与图像分类)与理解才能展开。所以，图像分割的好坏将直接影响到后续的图像分析精度。如何高效准确地使分割效果与实际地物相匹配已成为图像理解领域关注的一个热点。

目前，图像分割技术主要可分为包括基于区域的分割技术、基于边缘的分割技术和两者结合在内的传统图像分割方法，以及与特定理论结合的图像分割方法，包括基于数学形态学、模糊理论集、神经网络、支持向量机、小波分析和多尺度图像分割等等分割方法。虽然众多新方法已经在原基础上不断改进以避免图像分割的欠分割和过分割带来的矛盾，但仍没有一个通用的分割方法可以解决图像局部欠分割的问题。由于分割对象是图像后续信息提取和目标识别的最小单位，若存在欠分割现象将给后续图像处理带来不可逆转的误差。

以往对于欠分割现象的改善方法大多是在图像分割方法本身做改进以减少欠分割现象，或通过以增加图像计算量为代价把图像做过分割处理以保证最小对象的纯度。这些方法在减少欠分割现象的同时，又增添了其他矛盾(例如增加了进行图像分割的时间成本和后续建立图像分类规则的复杂化)。也有从分割对象角度为研究出发点，通过光谱值来聚类，适合于聚成两类的判定为欠分割对象，否则为单一对象。但该方法无法解释对象纹理形似斑马线状或明暗交织的灌木丛实为单一对象却被判为欠分割对象的情况。

发明内容

本发明克服现有技术中图像分割计算量大、分割过程中易产生其他矛盾以及产生错误的欠分割对象等缺陷，提出了一种基于对象纹理混杂度的高分辨率遥感影像欠分割对象识别方法。

本发明提出了一种遥感影像欠分割对象识别方法，包括如下步骤：步骤a.对图像进行数据降维，获取分割对象；步骤b.将所述分割对象中的像元进行聚类；步骤c.根据聚类计算所述分割对象的混杂度指标(MD，Mixed Degree)；步骤d.根据所述混杂度指标在所述分割对象中识别欠分割对象。

本发明提出的遥感影像欠分割对象识别方法中，所述步骤a包括：a1.通过主成分变换压缩所述图像得到待处理的波段与分割对象信息；a2.根据所述分割对象信息分割所述图像；a3.筛选分割后的图像得到分割对象。

本发明提出的遥感影像欠分割对象识别方法中，所述步骤b包括：b1.确定所述分割对象的像元值的初始聚类中心C1,C2；b2.对所述初始聚类中心进行Kmeans聚类，得到类别U1与类别U2。

本发明提出的遥感影像欠分割对象识别方法中，所述初始聚类中心如下式表示：

其中，b＝1,2,…k；b＝1,2,…k；

式中，C1与C2表示初始聚类中心，k表示主成分层数，n为分割对象的像元数，xbi表示第b主成分层的第i个像元的光谱值，x_b表示第b主成分层所有像元的平均光谱值，s_b表示第b主成分层的像元值的标准差。

本发明提出的遥感影像欠分割对象识别方法中，所述Kmeans聚类包括如下步骤：b21.分别计算像元与所述初始聚类中心C1,C2之间的距离，若所述像元与初始聚类中心C1的距离ri1小于像元与初始聚类中心C2的距离ri2，则该像元归入类别U1，否则该像元归入类别U2；b22.利用均值方法重新计算类别U1和类别U2的中心，若计算前后的中心保持不变则执行下一步骤，否则重新开始执行步骤b21。

本发明提出的遥感影像欠分割对象识别方法中，所述步骤b中进一步包括：b3.对所述分割对象进行去除噪声处理。

本发明提出的遥感影像欠分割对象识别方法中，所述步骤c包括：c1.在所述图像上叠加网格图，所述网格图的尺寸为所述图像的一半；c2.统计每个网格内原始像元数信息；所述原始像元数信息包括网格内背景值像元数T0、网格内类别U1的像元数T1和网格内类别U2的像元数T2；c3.根据所述原始像元数信息计算纯像元网格数和混合像元网格数；c4.根据所述纯像元网格数和混合像元网格数计算混杂度指标。

本发明提出的遥感影像欠分割对象识别方法中，在所述步骤c3中，当网格内T1＝0，T2≠0或T1≠0，T2＝0时，该网格是纯网格，其纯像元网格数Pure_i如下式表示：，式中，i表示第i个纯网格，N为网格内像元数，T0为网格内背景值像元数，Pure_i表示纯像元网格数；当网格内T1≠0，且T2≠0时，该网格是混合网格，其混合像元网格数Mix_j,如下式表示：式中，j表示第j个混合网格，N为网格内像元数，T1表示网格内类别U1的像元数，T2表示网格内类别U2的像元数，Mix_j表示混合像元网格数。

本发明提出的遥感影像欠分割对象识别方法中，所述混杂度指标如以下公式所示：

式中，MD表示混杂度指标,n为纯网格数，m为混合网格数。

本发明提出的遥感影像欠分割对象识别方法中，所述步骤d包括：d1.按所述混杂度指标对所述分割对象进行升序排序；d2.设定阈值，将序号小于阈值的分割对象作为欠分割对象。

本发明利用主成分变换原理，将更多的信息重组在更少的数据量中。通过度量待判定分割对象被聚为两类后，根据混杂程度进行欠分割对象的判定。本发明通过叠加更粗精度的网格，计算每个网格内的纯像元网格数和混合像元网格数，求和整个对象所有网格内的纯像元网格数和混合像元网格数得到分割对象的混杂度指标，从而根据混杂度指标判定是否为欠分割对象。本发明提出的混杂度指标MD值小的，说明混杂度小，可分性好，应判定为欠分割对象。

本发明对高分辨率遥感影像的欠分割对象检测具有更高的效率，现代高精度的高分遥感影像至少具有4个波段(有些如WorldView2达8个甚至更多波段)，计算量大，采用传统分割方法需要对每个图像层进行计算，而本方法仅需3个主成分层的计算量，效率明显提高；且利用混杂度指标判定欠分割对象精度更高，误判率和漏判率明显降低。

附图说明

图1是本发明遥感影像欠分割对象识别方法的流程图

图2是实施例中一个欠分割对象在聚类后计算得到的混杂度指标示意图。

图3是实施例中一个非欠分割对象在聚类后计算得到的混杂度指标示意图。

图4是欠分割对象识别方法在WorldView2高分辨率遥感图像分割中的应用实例。

图5是非欠分割对象识别方法在WorldView2高分辨率遥感图像分割中的应用实例。

具体实施方式

结合以下具体实施例和附图，对本发明作进一步的详细说明。实施本发明的过程、条件、实验方法等，除以下专门提及的内容之外，均为本领域的普遍知识和公知常识，本发明没有特别限制内容。

本发明遥感影像欠分割对象识别方法首先确定参与图像分割的光谱层，得到分割后的分割对象层，判断每个分割对象的面积(像元数)，小于一定值则可忽略不计。若对象面积大于一定阈值，则通过K均值聚类方法将对象分为两类，然后计算对象的混杂度指标(MD)作为判断欠分割的依据。对象的混杂度越低其对象被分为的两类的聚合程度越高，该对象越趋向于欠分割。

参阅图1，本发明遥感影像欠分割对象识别方法的具体步骤如下：

步骤a：数据降维

a1.通过主成分变换压缩图像得到待处理的波段与分割对象信息。主成分变换是设法将原来众多具有一定相关性的变量，重新组合成一组新的互相无关的及格综合变量，同时根据实际需求从中可以取出较少的综合变量尽可能多地反映原来变量的信息。主成分分析是一种多元统计分析方法，也是数学上用来降维的一种方法。利用主成分变换后的综合光谱图像(第一、二、三主成分)作为图像分割的编辑层，并以各主成分的特征值百分比作为权重值。该方法将遥感影像的多波段图像信息压缩到比原波段更有效的少数几个转换波段中,仅利用前三个主成分就可以解释95％以上的图像信息，不仅有利于减少后续处理和分析的数据量，同时主成分变换得到的百分比(即主成分变换层中所涵盖信息量的比重)可作为图像分割时进行权重设置(w_c)的客观指标。

a2.根据分割对象信息分割图像。本发明中使用前三个主成分变换层作为图像分割的编辑层，进行有权重设置的图像分割。

a3.筛选分割后的图像得到分割对象。本发明中通过排除整幅分割图像中面积较小的对象以减少计算成本，面积阈值的确定视遥感影像空间分辨率、图像分割的要素和分割精度而定。例如本实施例中采用空间分辨率为0.5m的WorldView2高分遥感影像，分割的要素是城市地物，那么当对象面积小于500像元时(相当于一个羽毛球片场面积)，该对象不参与运算。

步骤b：将分割对象中的像元进行聚类。

b1.确定分割对象的初始聚类中心C1,C2。

本实施例采用空间分辨率为0.5m的WorldView2高分遥感影像，其图像信息是混合的地物，利用对象的光谱信息在多光谱特征空间中做固定初始聚类中心C1，C2为中心的两个聚类集类别U1和类别U2。确定初始聚类中心C1和C2，初始聚类中心的计算公式如下：

其中，b＝1,2,…k；b＝1,2,…k

式中，C1与C2表示初始聚类中心，k表示主成分层数，n为对象像元数(即面积)，x_bi表示第b主成分层的第i个像元的光谱值，表示第b主成分层所有像元的平均光谱值，s_b表示第b主成分层的所有像元值的标准差，根据原始数据的光谱数可以定义k维多光谱特征空间。

b2.对初始聚类中心进行Kmeans聚类，得到类别U1与类别U2。计算每个像元Pi到C1和C2的距离r_i1和r_i2，其中，如果r_i1<r_i2，Pi点归入U1类，否则归入U2类，并利用均值方法重新计算U1和U2的中心，如此反复迭代更新后，直到新老类别中心的值保持不变，则迭代结束，否则继续迭代。

b3.去除椒盐噪声的干扰。对聚类后的类别U1作数学形态学的开启闭合操作，去除椒盐点噪声的干扰。本实施例中使用尺寸大小为2的SE对分割对象做开启闭合操作，消除图像中的凸点和凹陷。

步骤c.根据聚类计算分割对象的混杂度指标；

c1.在图像上叠加网格图，网格图的尺寸为图像的一半。若原始对象的外接矩形为m*n，则在原始对象上叠加一张的网格图，若m或n为奇数，则加1后除以2，以保证能被2整除。

c2.统计每个网格内原始像元数信息；原始像元数信息包括网格内背景值像元数T0、网格内类别U1的像元数T1和网格内类别U2的像元数T2。

c3.根据原始像元数信息计算纯像元网格数和混合像元网格数。先计算单个网格内的纯像元网格数Pure_i和混合像元网格数Mix_j：

若T1＝0，T2≠0或T1≠0，T2＝0，则说明该网格内只含有类别U2或类别U1，该网格是纯网格，计算纯像元网格数Pure_i，N为网格内像元数，此处N＝4。

若T1≠0，且T2≠0，则说明该网格是混合的，即U1类和U2类的交界处，计算混合像元网格数Mixj,

以上两个公式中皆排除了外接矩形的背景值，即排除了分割对象的边界长度和对象形状的干扰，只考虑分割对象内部被强制分为两类时混杂网格数与总的分割对象网格数的比值情况。

c4.根据纯像元网格数和混合像元网格数计算混杂度指标。

累计求和计算上述步骤c3中整个分割对象内所有网格的纯像元网格数Pure和混合像元网格数Mix，然后计算该分割对象的混杂度指标(MD)：

式中，MD表示混杂度指标,n为纯网格数的序号上限，m为混合网格数的序号上限。

步骤d.根据混杂度指标在分割对象中识别欠分割对象。

d1.按混杂度指标对分割对象进行升序排序。计算得到整幅影像所有参与计算的分割对象的混杂度指标(MD)，该MD值是一个相对值，若MD_x<MD_y，说明相较于x对象，y对象聚类后的两类在空间分布上没有明显界限，更趋向为同一种地物，而x对象两类界限所占整个分割对象的比重较低，说明对象本身的光谱值在空间分布上分异性好于分割对象y，根据实际经验可判断分割对象x比分割对象y更有可能是欠分割对象。

d2.设定阈值，将序号小于阈值的分割对象作为欠分割对象。一般设定阈值为10-25％，取MD值从小到大排列中的前10％～25％分割对象为欠分割对象进行再分割。该MD值提供了筛选欠分割对象的有效依据，以便后续处理中对欠分割对象进行分割。

以下结合示意图对本发明实施混杂度指标MD的过程进行详细描述。

图2和图3为分割对象1(欠分割对象)和分割对象2(非欠分割对象)经过步骤a和步骤b后的聚类示意图。对象1和对象2都接有20*20的外接矩形，在此基础上叠加一张10*10的网格图(即步骤c1)，网格图如图2与图3中的粗黑边。

欠分割对象往往是材质相接近的不同地物，比如水泥制房屋与房屋周围水泥道路的混合对象。欠分割对象的像元被强制分为两类后，两类的分界线比较明确，判为欠分割对象的概率较高。

非欠分割对象往往是由于灌木丛的高低差形成阴影，从而体现在图像像素值上的差异，或类似斑马线纹路的地物，容易分为两类，实则一类地物。非欠分割对象的像元被分为两类后，两类比较混杂，彼此交融混合度较高，界限不明显，判为非欠分割对象的概率较高。

根据步骤c2与c3计算得到每个网格内的T0、T1和T2值，并计算得到每个网格的纯像元网格数(淡色字体)和混合像元网格数(黑色字体)。

根据步骤c4对图2的Pure_i和Mix_j累计求和，得到Mix＝10.25，Pure＝65.75，MD＝0.13486；同理，对于图3，Mix＝37.5，Pure＝27.75，MD＝0.57471。

当图2和图3对象被聚类分为两类U1与U2后，判断出图2中U1和U2的分异性很高，混杂度很低，其为将不同属性的两种地物被欠分割为同一个分割对象(比如，遥感影像中的房屋和其阴影)，而图3中U1和U2的混杂度很高，因此强制分成两类是不合理的，分割区域更趋向是单一地物(比如，遥感影像中的灌木丛等)。从定量的角度去筛选欠分割对象的话，就是计算整幅影像中所有分割对象的MD指标，MD越小，被判为欠分割对象的可能越大。图2中的MD远小于图3中的MD值，相较图3的对象而言，图2对象为欠分割对象，需要进一步分割。

1.本发明的欠分割对象识别方法在WorldView2高分辨率遥感图像分割中的应用实例。

选取市区一复杂地块(图4)中的三块区域作为研究区，包含居民楼房、单位楼房、绿地、道路等多种地面目标，以检验本文提出方法的适用性。将WorldView2影像的8波段进行主成分变换后(表1)，取前三个主成分层(PC1的累积特征百分比70.58％，PC2的累积特征百分比96.80％，PC3的累积特征百分比99.11％)在Definien8.7中进行多尺度分割，并分别将主成分的特征值百分比70、26和2作为三个图像编辑层的权重值参与到图像分割的运算中。

表1研究区主成分变换后的统计属性

获得分割尺度为88时的最优分割图像，含有534个分割对象，排除其中面积(即像元数量)小于500个像元的对象，一共得到400个有效分割对象。将每个分割对象计算得到的MD值从小到大排列，取前20％的对象(80个对象)在图4中进行放大并以高亮轮廓显示。有图4可见，房屋与水泥地面混杂和树木与水泥地面等混杂在一起的欠分割对象都被筛选出来。

反之，当选取用本发明所得的所有MD值中较大的20％对象(80个对象)，可以发现筛选出来的多是灌木丛、屋顶和草地等仅包含单一地物的对象。如图5放大区域内高亮轮廓显示的对象。

通过对比图4和图5可发现，依据本发明所得的MD值排序，可准确地筛选出绝大多数欠分割对象，且排除仅包含单一地物的纯对象。

本发明的保护内容不局限于以上实施例。在不背离发明构思的精神和范围下，本领域技术人员能够想到的变化和优点都被包括在本发明中，并且以所附的权利要求书为保护范围。

Claims (7)

1.一种遥感影像欠分割对象识别方法，其特征在于，包括如下步骤：

步骤a.利用主成分变换对图像进行数据降维获取分割对象，用主成分变换后的综合光谱图像的各主成分的特征值百分比作为权重值；

步骤b.将所述分割对象的像元值进行聚类，得到类别U1与类别U2；

步骤c.根据聚类计算所述分割对象的混杂度指标；其中，包括如下子步骤：

c1.在再分割为两类的图像上叠加一张以长为a，宽为b像元的网格图作为单元网格；

c2.统计每个网格内原始像元数信息；所述原始像元数信息包括网格内背景值像元数T0、网格内类别U1的像元数T1和网格内类别U2的像元数T2；当网格内T1＝0，T2≠0或T1≠0，T2＝0时，该网格是纯网格，其纯像元网格数Pure_i如下式表示：

${\mathrm{Pure}}_i=1-\frac{T0}{a*b};$

式中，i表示第i个纯网格，a为单元网格的长，b为单元网格的宽，T0为单元网格内背景值像元数，Pure_i表示纯像元网格数；

当网格内T1≠0，且T2≠0时，该网格是混合网格，其混合像元网格数Mix_i,如下式表示：

${\mathrm{Mix}}_j=\frac{T1+T2}{a*b};$

式中，j表示第j个混合网格，a为单元网格的长，b为单元网格的宽，T1表示单元网格内类别U1的像元数，T2表示单元网格内类别U2的像元数，Mix_i表示混合像元网格数；

c3.根据所述原始像元数信息计算纯像元网格数和混合像元网格数；

c4.根据所述纯像元网格数和混合像元网格数计算混杂度指标；所述混杂度指标如以下公式所示：

$MD=\frac{{\&Sigma;}_{j=1}^m{\mathrm{Mix}}_j}{{\&Sigma;}_{j=1}^m{\mathrm{Mix}}_j+{\&Sigma;}_{i=1}^n{\mathrm{Pure}}_i};$

式中，MD表示混杂度指标,n为纯网格数的序号上限，m为混合网格数的序号上限；

步骤d.根据所述混杂度指标在所述分割对象中识别欠分割对象。

2.如权利要求1所述的遥感影像欠分割对象识别方法，其特征在于，所述步骤a包括：

a1.通过主成分变换压缩所述图像得到待处理的波段与分割对象信息；

a2.根据所述分割对象信息分割所述图像；

a3.筛选分割后的图像得到分割对象。

3.如权利要求1所述的遥感影像欠分割对象识别方法，其特征在于，所述步骤b包括：

b1.确定所述分割对象的像元值的初始聚类中心C1,C2；

b2.对所述初始聚类中心进行Kmeans聚类，得到类别U1与类别U2。

4.如权利要求3所述的遥感影像欠分割对象识别方法，其特征在于，所述初始聚类中心如下式表示：

$C1=({\stackrel{\&OverBar;}{x}}_1+s_1,{\stackrel{\&OverBar;}{x}}_2+s_2,{\stackrel{\&OverBar;}{x}}_3+s_3,...),C2=({\stackrel{\&OverBar;}{x}}_1-s_1,{\stackrel{\&OverBar;}{x}}_2-s_2,{\stackrel{\&OverBar;}{x}}_3-s_3,...);$

其中，b＝1,2,…k；b＝1,2,…k；

式中，C1与C2表示初始聚类中心，k表示主成分层数，n为分割对象的像元数，x_bi表示第b主成分层的第i个像元的光谱值，表示第b主成分层所有像元的平均光谱值，s_b表示第b主成分层的像元值的标准差。

5.如权利要求3所述的遥感影像欠分割对象识别方法，其特征在于，所述Kmeans聚类包括如下步骤：

b21.分别计算像元与所述初始聚类中心C1,C2之间的距离，若所述像元与初始聚类中心C1的距离ri1小于像元与初始聚类中心C2的距离ri2，则该像元归入类别U1，否则该像元归入类别U2；

b22.利用均值方法重新计算类别U1和类别U2的中心，若计算前后的中心发生变化则重新开始执行步骤b21。

6.如权利要求3所述的遥感影像欠分割对象识别方法，其特征在于，所述步骤b2之后进一步包括：b3.对聚类后的所述类别U1作数学形态学的开启闭合操作，去除椒盐点噪声的干扰。

7.如权利要求1所述的遥感影像欠分割对象识别方法，其特征在于，所述步骤d包括：

d1.按所述混杂度指标对所述分割对象进行升序排序；

d2.设定阈值，将升序排序后的序号小于阈值的分割对象作为欠分割对象。