一种结合形态滤波和标记分水岭变换的遥感图像分割方法

江怡，梅小明，邓敏，陈杰，陈铁桥

（中南大学地球科学与信息物理学院，湖南 长沙 410083）

0 引言

遥感图像分割是遥感图像目标识别、分析和处理领域中一项极为关键的技术，分水岭变换作为一种经典而有效的图像分割工具，已经成为遥感图像分割的重要手段。在分水岭变换的发展过程中，国内外学者对其做了深入研究［1－5］，最具代表性的是Vincent等提出的基于浸水模型的分水岭变换快速算法（简称V－S算法）［4］以及 Mortensen等提出的基于降水模型的分水岭变换快速算法［5］。其中，V－S算法是一种基于地形学思想和区域增长理论的图像分割方法，它通过对梯度图像进行分水岭变换，从而得到单像素宽、定位精确且连续封闭的边缘轮廓。然而，由于图像背景噪声等因素影响，这种直接对梯度图像进行分水岭变换的方法容易产生过分割现象，导致物体边缘信息淹没在过度分割产生的复杂信息中，因而不能获得理想的分割结果。

针对V－S算法中的过分割问题，已出现较多改进方法［6－12］，大体分为3类：1）对原始图像或者梯度图像进行滤波处理［6，7］，将图像背景噪声、地物内部细密纹理消除，再实施分水岭变换分割。2）在分水岭分割之前对梯度图像添加标记以指导区域增长［8，9］，将由非标记的极小点形成的分割区域合并到其他由标记点形成的区域当中。该方法中阈值通常是根据待分割图像人为设定，缺乏严格的衡量标准。3）在分水岭分割之后，根据区域之间的相似性进行区域合并，进一步减少分割结果中区域数量［10－12］。这类方法虽然能够有效完成分割任务，但合并过程运算复杂、计算量大，而且终止合并时最优阈值的选取极为困难。为此，本文在V－S算法的基础上提出一种结合形态滤波和标记分水岭变换的遥感图像分割算法：采用非线性混合开闭重构滤波器对原始图像进行去噪处理，在滤除噪声的同时较好地保持了目标地物边缘；通过对梯度图像进行形态修正和极值标记，进一步滤除伪极小区域并限制极小区域数目；在此基础上实施分水岭变换，可以有效抑制分水岭变换的过分割现象。通过调整滤波尺寸和深度阈值参数，在一定程度上实现快速有效地分割，以满足不同遥感图像分割要求。

1 结合形态滤波和标记分水岭变换的分割方法

1.1 形态学混合开闭重构滤波

分水岭变换直接用于梯度图像时，噪声和梯度的局部不规则性通常会导致过分割，严重干扰了对真实目标的提取。因此，在进行梯度运算之前必须对原始图像加以平滑，以有效消除图像中的噪声。一般情况下采用线性滤波器（如均值滤波器、高斯滤波器等）可以滤除图像的高频分量，达到平滑图像、消除噪声的目的，但不能有效抑制纹理细节对后续图像分割带来的过分割影响。相比于传统的线性滤波方法，数学形态学的非线性滤波方式更适合图像分割及与图像几何特性相关方面的应用［13］。

形态学混合开闭重构滤波器是一种基于非线性区域连通算子的滤波器［14］，其在简化图像的同时，可以较好地保持图像中剩余连续区域的边缘，而且后续分割时不会有新的轮廓边缘产生［15］，从而得到更有效的分割结果。形态学混合开闭重构建立在测地膨胀和测地腐蚀的基础上，对于灰度图像f和参考图像r（取f－1），形态学开、闭重构定义如下［16］：

式中：f为原始图像；r为参考图像；b为结构元素；n为结构元素尺寸；rec表示重构运算；○和·分别为形态学开和闭运算；和分别表示测地膨胀和腐蚀运算收敛时的重构图像。形态学混合开闭重构运算在尺度n时定义为开－闭和闭－开重构运算的平均，表达为［16］：

通过分析发现，传统的形态学开闭运算只能去除图像中的部分高灰度和低灰度细节，而混合开闭重构运算在图像的平滑过程中，完全剔除了比当前尺度小的高灰度和低灰度区域细节。随着尺度n的递增，只会消除图像中原有的局部极小值，不会产生新的区域极小值，重构图像中保留的区域轮廓也不会发生位置偏移。

图1a为待处理的Canberra地区（Australia）分辨率为2.5m的SPOT5子图像（2003年2月拍摄），图1b和图1c分别表示利用线性滤波器（高斯低通滤波器）和形态学非线性滤波器（混合开闭重构滤波器）对图像滤波后的效果。从图1b可以发现，经过高斯低通滤波处理，并没有达到理想的滤波效果，部分地物边缘粗糙、模糊，地物内部仍存在大量纹理细节。图1c采用形态学混合开闭重构滤波器滤波后，图像得到了更大的简化，同时图像的边缘轮廓保持良好；但受未完全滤除的噪声和地物内部细密纹理的影响，仍存在由目标边缘附近灰度值均匀变化的过渡区域产生的伪极小区域。因此，要有效减少造成过分割现象的伪极小区域的数目，并且得到一个比较准确的分割结果，对梯度图像进行形态修正处理十分必要［15］。

图1 滤波效果图对比Fig.1 Filtering results of original image with Gaussian and morphological method

滤波器中结构元素形状和尺寸的选取对图像滤波效果具有重要影响。根据不同目的对遥感图像进行分割时，滤波应以不损失最小待分割地物为准则，结构元素尺寸设置应小于所有待分割地物尺寸的最小值，结构元素形状则选取与所需提取特定地物相同形状的结构元素。常用结构元素中，圆盘形结构元素具有各向同性，不会造成图像特征值的畸变，对各种不同形状的地物处理效果比较均衡。本文滤波器结构元素选取方法如下：根据分割要求确定图像中最小待分割地物，滤波器中结构元素尺寸r（0≤r＜Rmin，r取整数）设置为小于该最小待分割地物的外接圆半径Rmin的最大整数，若不规则地物无法求取外接圆，则Rmin选取能包含该地物的最小圆半径为结构元素尺寸；若r取值大于Rmin时，将会有部分待分割地物被当做噪声滤除，无法得到所需分割效果。

1.2 形态学标记

原始图像和梯度图像进行形态滤波之后，尽管局部极小值点和极大值点数目大幅减少，图像背景噪声和地物内部细密纹理得到有效抑制，但在此基础上直接进行分水岭变换，所得分割结果中仍会产生许多零碎区域。针对这一问题，本文引入扩展极小变换H－minima运算的形态标记方法对梯度图像进行极值标记，限制局部极小值点的数目。扩展极小变换的本质是形态学阈值算子，利用深度阈值提取显著的局部极小值，以此为依据对梯度图像进行重构，将大多数的微小无关区域标记为0，从而有效避免出现过分割问题。梯度▽f图像经过深度阈值为h的扩展极小变换运算，即［17］：

式中：▽f代表滤波后的梯度图像；▽fmark代表经过标记的梯度图像；Hmin代表形态学的 H－minima变换通过深度阈值h的设定，消除积水盆地低于给定阈值的局部极小值。

深度阈值满足尺度因果性，其取值直接影响到分割结果区域的数目，即深度阈值愈大，被标记的极小值点数目愈少，最终分割的区域数目愈少。因此可以根据具体的分割对象和分割要求设置相应的深度阈值，控制和指导合理分割结果的生成。

极值标记后，用形态学强制最小运算修正梯度图像，使得图像的局部极小区域仅出现在被标记位置。修正后的梯度图像用▽fws表示，公式如下［17］：

式中：IMmin代表形态极小值标定操作；▽fmark代表经过标记的梯度图像。

1.3 结合形态滤波和标记分水岭变换的分割方法

本文结合形态滤波和标记分水岭变换的新型分割算法具体步骤如下（图2）：1）计算最小待分割地物外接圆半径，设置相应尺寸r1的圆盘形结构元素的滤波器，对输入的待分割图像实施混合开闭重构滤波，消除噪声、平滑目标地物内部细密纹理；2）采用sobel梯度算子对滤波图像进行梯度运算，获取梯度图像；3）对梯度图像进行形态修正，以消除目标边缘附近灰度值均匀变化的过渡区域所形成的梯度极小值，结构元素尺寸r2选取法则同步骤1；4）根据文献［17］，设置深度阈值h，采用扩展极小变换和强制最小运算对梯度修正图像进行形态标记，获取标记的梯度极值图像；5）对经过标记的梯度图像实施分水岭变换运算得到最终的分割图像。

图2 本文分割算法流程Fig.2 The flowchart of the proposed segmentation approach

2 实验与分析

选取Reno地区（Nevada）的SPOT5影像作为测试图像（影像获取时间为2003年9月14日，空间分辨率为5m），在原图像中截取大小为256×256像元的具有典型地物的两部分子图像进行实验：图3a代表城区图像，主要地物有工厂、房屋、道路等；图4a代表郊区图像，主要地物有农田、水体等。实验平台为 Win7操作系统，CPU主频2.27GHz，内存2G，用Matlab编程。

对图3a原始待分割城区图像实施形态滤波后结果如图3b所示。经过混合开闭重构滤波后，在一定程度上有效消除了图像背景噪声，平滑了目标地物内部细密纹理，突出了目标地物边缘轮廓，很好地保持了目标地物的几何特性。图3c为经过形态修正后的梯度图像，目标地物边缘轮廓清晰可见且定位精确，基本避免了目标地物边缘模糊和轮廓偏移的现象。对图3c实施扩展极小变换运算并进行形态标记，梯度极值标记图像如图3d所示，图像中极小值点的数目进一步减少，原梯度图像中一部分极小值点消失，另一部分极小值点与周边更大地物合并，形成与地物大小、边缘轮廓位置相吻合的极小值区域。在此标记的极值梯度图像上进行分水岭变换，分割结果如图3e所示。图3f为将最终分割结果叠加到原始待分割图像上所得的叠加效果图，工厂、房屋、道路边缘轮廓清晰且与原待分割地物吻合度较高。同样采用本文分割算法对郊区图像进行实验（图4），典型地物如农田、水体的边缘轮廓与原始图像均有较高吻合度。

上述实验表明，采用本文算法能较为合理的分割出城区图像（图3）中的道路、工厂、房屋和郊区图像（图4）中的农田、水体等地物，待分割目标地物的边缘轮廓连续、封闭且位置准确。但对于建筑用地中面积较小的地块提取的精确性还需改进，主要原因是面积较小的地块像素的灰度值与周围地物的灰度值对比度不高，导致利用分水岭分割时无法将其从周围地物中准确分割出来。

本文提出的算法涉及3个参数的选取：1）对原始图像进行混合开闭重构滤波时结构元素尺寸r1。由图像分析，本文城区图像中最小待分割地物外接圆半径Rmin为5，故滤波器结构元素尺寸r1取4。当r1取值大于Rmin时，出现部分待分割地物被当做噪声滤除的现象。此参数用于消除所有小于该结构元素尺寸的图像背景噪声，使得整体目标地物内部的过分割现象得到有效改善，但部分目标地物边缘附近由于灰度单调变化的过渡区域造成的过分割现象仍然存在。2）对梯度图像进行处理的结构元素尺寸r2，用其可以消去梯度图像中尺寸小于该结构元素的梯度极小值区域（由目标地物边缘附近的过渡区域产生）。这项处理在改善边缘附近过分割问题的同时，确保边缘轮廓位置不发生偏移，从而得到单像素宽、定位精确且连续、封闭的边缘轮廓。根据与r1类似的选取原则，确定本文城区图像的r2取2。3）对梯度图像进行扩展极小变换运算的深度阈值h。根据文献［17］对本文待分割城区影像进行计算，最终确定深度阈值h取31，仅标记所有与其邻域像素亮度之差大于31的极小值点（块），从而使得梯度图像中极小值点（块）数目减少。

本文在多尺度下进行了实验结果对比（以城区图像为例），实验结果如图5所示。若只对原始图像进行平滑而不同时对梯度图像进行形态学滤波处理，过分割现象虽有减轻，但并未完全消除。随着本文算法中5个步骤的逐步进行，分割区域数目逐级减少，满足分割区域数对尺度的因果性关系。

为了更充分地验证本文分割算法，引入经典分割算法进行实验对比。分别利用Canny算子对图3a和图4a进行边缘检测（图6），由图6可知，Canny算子能够准确地提取出大部分边缘像素点；但由于图像背景噪声和地物内部细密纹理的影响，导致检测到的轮廓边缘局部断裂成若干连续的小段甚至孤立点，有些地方存在遗漏真正边缘像素点的现象。而本文提出的分割算法不但能够生成单像素宽、定位准确且连续、封闭的边缘轮廓，而且无需在分割之后再进行边缘追踪与连接过程，可直接得到所需的分割结果。

图5 城区图像在本文算法3个参数下的分水岭变换分割结果Fig.5 The segmentation results of urban image with the three parameters

图6 Canny算子的检测结果Fig.6 Results of Canny′s edge－based technique

表1记录了对两幅图像进行分割实验所选取的参数值大小、分割消耗的时间（t）及分割区域数（n）。对于城区图像，若采用V－S算法，最终所得的分割区域数为5 170块，所需时间为0.257s；而采用本文改进分割算法，实验中选取的各个参数分别为r1＝4、r2＝2、h＝31，最终的分割区域数为122块，所需时间为0.357s。对于郊区图像，若采用V－S算法，最终的分割区域数为6 100块，所需时间为0.271s；而采用本文改进分割算法，实验中选取的各个参数分别为r1＝4、r2＝3、h＝10，最终的分割区域数为69块，所需时间为0.405s。可见，采用本文的算法进行分割后所得最终分割区域数大大减少，有效克服了分水岭变换的过分割现象。

表1 城区、郊区图像分割参数及分割区域数对比Table 1 Values of the parameters and segmentation region numbers of urban image and rural image

3 结语

本文针对分水岭变换过分割问题，引入数学形态学理论对分水岭变换算法进行改进，提出了一种结合形态滤波和标记分水岭变换的遥感图像分割算法，有效解决了传统算法边缘定位不准、对弱边提取困难等问题，能够获得符合人类视觉特点、具有实际意义且与实际地物目标相一致的分割区域。整个分割过程不需进行分割后的区域合并处理，有效降低了算法的时间复杂度。针对特定的图像和图像中的特定目标，通过调整算法中3个参数的大小，有效控制分割处理过程，从而为遥感分类和信息提取奠定了良好的基础。从整个分割过程看，图像的最终分割结果取决于所选择的结构元素的尺寸，而尺寸的选择则根据对图像的具体分析需要而定。因此，本文算法具有较强的适应性，能根据不同的需求对遥感图像进行分割，是一种有效而实用的遥感图像分割方法。

［1］ BEUCHER S，LANTUEJOUL C.Use of watersheds in contour detection［A］.Proc.International Workshop on Image Processing［C］.1979.17－21.

［2］ MAYER F，BEUCHER S.Morphological segmentation［J］.Visual Communication Image Representation，1990，1（1）：21－46.

［3］ NAJMAMN L，SCHMITT M.Geodesic saliency of watershed con－tours and hierarchical segmentation［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1996，18（12）：1163－1173.

［4］ VINCENT L，SOILLE P.Watersheds in digital spaces：An efficient algorithm based on immersion simulations［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1991，13（6）：583－598.

［5］ MORTENSEN E N，BARRETT W A.Toboggan－based intelligent scissors with a four－parameter edge model［A］.Proc.IEEE Conf.Computer Vision and Pattern Recognition［C］.1999.452－458.

［6］ ANGULO J，VELASCO F S，CHANUSSOT J.Multi－scale stochastic watershed for unsupervised hyper－spectral image segmentation［A］.IEEE International Geoscience and Remote Sensing Symposium［C］.2009.93－96.

［7］ 陈景广，佘江峰，黄海涛.基于形态学的多尺度遥感图像分割方法［J］.地理与地理信息科学，2012，28（4）：22－24.

［8］ SUN Y，HE G J.Segmentation of high－resolution remote sensing image based on marker－based watershed algorithm［A］.The 5th International Conference on Fuzzy Systems and Knowledge Discovery［C］.2008.271－276.

［9］ PESARESI M，BENEDIKTSSON J A.A new approach for the morphological segmentation of high－resolution satellite imagery［J］.IEEE Transactions on Geoscience and Remote Sensing，2001，39（2）：309－320.

［10］ TARABALKA Y，CHANUSSOT J，BENEDIKTSSON J A，et al.Segmentation and classification of hyper－spectral data using watershed［J］.Pattern Recognition，2010，43（7）：2367－2379.

［11］ HARIS K，EFSTRATIADIS S N，MAGLAVERAS N，et al.Hybrid image segmentation using watersheds and fast region merging［J］.IEEE Transactions on Image Processing，1998，7（12）：1684－1699.

［12］ CHEN J，DENG M，XIAO P F，et al.Multi－scale watershed segmentation of high－resolution multi－spectral remote sensing image using wavelet transform［J］.Remote Sensing，2011，15（5）：908－926.

［13］ 罗玲，解梅，陈杉.基于多尺度形态滤波的分水岭图像分割方法［J］.计算机辅助设计与图形学学报，2004，16（2）：168－173.

［14］ SALEMBIER P，SERRA J.Flat zones filtering，connected operator，and filters by reconstruction［J］.IEEE Transactions on Image Processing，1995，4（8）：1153－1160.

［15］ 蒋志勇，陈晓玲，秦前清.一种基于重构算子的分水岭变换算法［J］.中国图象图形学报，2009，14（12）：2527－2533.

［16］ SALEMBIER P.Morphological multi－scale segmentation for image coding［J］.Signal Processing，1994，38（3）：359－386.

［17］ SOILLE P.形态学图像分析：原理与应用（第二版）［M］.北京：清华大学出版社，2008.83－150.