基于改进标记分水岭的高分辨率遥感影像海岸水边线提取方法

栾奎峰，刘 帅，潘与佳，朱卫东*，李丕学,裘 诚，邱振戈，沈 蔚，3，王 洁，王振华

(1.上海海洋大学 海洋科学学院，上海 201306；2.上海市海洋监测预报中心，上海 200062；3.上海河口海洋测绘工程技术研究中心，上海 201306)

0 引言

海岸线是海洋与陆地的分界线[1]。如何快速、准确提取海岸线具有重要意义。遥感技术具有全天候、全天时、覆盖面积广和数据来源丰富等优点，已成为海岸线提取的重要手段[2-3]。现有基于遥感手段的海岸线提取方法主要有：阈值分割法[4]、边缘检测法[5]、活动轮廓模型法[6]、面向对象法[7]、区域生长法[8]和元胞自动机法[9]等。分水岭算法是一种基于形态学的图像分割算法，定位准确，适用于分割粘连的目标，得到单个像素宽度和连续封闭分割边缘，适用于海岸线的提取。VINCENT et al[10]最早提出一种比较经典的分水岭检测算法，以梯度图像作为输入进行分割；GAO et al[11]先采用数学形态学中的面积算子对图像进行简化，然后再进行分水岭分割；NING et al[12]在Mean-Shift聚类分割后进行分水岭分割，并通过人工干预进行区域合并得到分割结果；林振荣 等[13]先采用形态学滤波进行降噪处理，然后在分水岭分割后进行区域合并；PUISSANT et al[14]基于形态学工具“击中击不中变换”定义模板，通过双重阈值对结果进行组合，再进行标记的分水岭变换，得出准确的海岸线。随着国产高分辨率多光谱卫星数据的发展，影像上像元的纯净度提高，地物更加清晰，但细节也被放大，噪声变得更多，同物异谱和异物同谱现象大量存在，使影像中出现许多局部灰度极值区域，原有的标记分水岭算法存在过分割现象，产生许多无效数据。

针对分水岭方法在高分辨率多光谱数据提取岸线中存在的过分割和局部极值问题，本文提出了基于扩展极大值和极小值变换的标记分水岭水边线提取方法，采用形态学重建和扩展极值变换等方法抑制噪声点的影响，以达到较好的水边线提取效果。

1 改进的扩展极值变换标记分水岭水边线提取方法

改进的扩展极值变换标记分水岭水边线提取方法主要步骤如图1所示：首先对待分割的高分辨率多光谱影像进行辐射定标和大气校正等预处理，其次对预处理后的影像进行归一化水体指数(Normalized Difference Water Index，NDWI)计算和梯度图像计算，并对NDWI计算后的影像进行形态学重建。然后，使用本文提出的改进的扩展极值变换方法对图像进行变换，获得前景和背景标记，并对原梯度图像进行强制最小值标定处理，修正梯度图像。最后，对修正后的梯度图像进行分水岭变换，完成水边线的提取。

图1 算法流程图Fig.1 Algorithm flow chart

1.1 形态学重建滤波

分水岭变换是利用图像灰度梯度幅值作为分割函数，对每个像素点的灰度级进行排序，然后按照由低到高的顺序进行处理[15]。为突出水体信息、增强水体与陆地部分像元灰度值的差异，采用NDWI计算方法，对高分辨率多光谱影像数据进行波段运算。之后，利用形态学重建滤波方法对处理后的影像进行平滑，去除影像中的噪声。形态学基本运算包括腐蚀、膨胀、开运算和闭运算，开运算是先腐蚀后膨胀，闭运算是先膨胀后腐蚀。膨胀和腐蚀运算定义如下

(g⊕l)(x,y)=max{g(x+i,y+j)-
l(i,j)|(x+i,y+j)∈Dg,(i,j)∈Dl}

(1)

(gΘl)(x,y)=min{g(x+i,y+j)-
l(i,j)|(x+i,y+j)∈Dg,(i,j)∈Dl}

(2)

式中：g为待处理图像，l为结构元素，(x,y)为待处理图像中运算的原点，(i,j)表示运算过程中结构元素范围内像元的位置，Dg表示待处理图像的范围，Dl表示结构元素的范围，⊕、Θ分别表示膨胀和腐蚀运算。

(3)

(4)

结构元素l的选取影响着影像形态学重建后的结果。在常用的结构元素中，圆盘型结构元素对各种地物的处理效果都比较均衡，具有各向同性，不会对图像特征值造成畸变[17]。结构元素大小m应大于目标区域内部单体地物，同时小于目标区域的最大内切圆半径，目的是使目标区域内部更为平滑，同时保持其整体的外形轮廓。通过先进行开重建运算再进行闭重建运算的二次重建滤波，可有效平滑图像和削弱噪声。

1.2 基于扩展极值变换的形态学标记

经过形态学重建处理的多光谱影像数据，虽然大部分局部噪声被去除，但仍存在一些对影像梯度修正的准确性造成干扰的噪声。扩展极值变换[18]是一种形态学阈值算子，是通过设定阈值，限制灰度局部极值点的数量，提取出显著的灰度极值区域，将大多数无关区域标记为0。一般的研究是在影像分割过程中仅采用扩展极小值变换对影像中的灰度极小值区域进行限制，然后进行标记提取[17,19-20]。但由于高分辨多光谱卫星影像细节丰富，同物异谱和异物同谱等因素，影像中不仅灰度极小值区域会对分割结果造成影响，过多的灰度极大值点和区域也会对分割结果造成一定程度的干扰。因此，本文提出了基于扩展极小值和极大值变换的标记分水岭方法，以此对形态学重建后的图像进行重构，从而解决高分辨多光谱数据的过分割问题。

形态学重建滤波处理后的图像，陆地范围的像元灰度值低于水体像元的灰度值，则采用扩展极小值变换标记陆地，其公式为

If=Hmin(M|h1)

(5)

式中：If为扩展极小值变换后的图像，Hmin为形态学扩展极小值变换运算，M为形态学重建后的图像，h1表示扩展极小值变换运算的灰度阈值。扩展极小值变换后的图像中陆地区域标记为1，其余部分标记为0，消除陆地部分低于阈值的局部灰度极小值区域，并与已标记为陆地的区域合并。

为去除水体区域中存在的局部噪声点，本文采用扩展极大值变换的方法，对图像进行扩展极值变换，其公式为

Ib=Hmax(M|h2)

(6)

式中：Ib为扩展极大值变换后的图像，Hmax表示形态学的扩展极大值变换运算，h2为扩展极大值变换运算的灰度阈值。扩展极大值变换后的图像中水体区域标记为1，其余部分标记为0，消除水域部分高于阈值的局部灰度极大值区域，并与已标记为水域部分的区域进行合并。

扩展极小值和极大值变换过程中，h1和h2的设定是扩展极值变换的关键，阈值过小则处理后的图像仍会存在过分割现象，阈值过大会使得最终的结果出现错误分割。本文采用基于最大类间方差法[21]设定上述两个阈值的初值，其公式为

h1=t-Imin

(7)

h2=Imax-t

(8)

式中：t为采用最大类间方差法计算的图像分割阈值，Imin和Imax分别为全图像中像元灰度的最小值和最大值。由于高分辨率多光谱数据的差异，h1和h2的值在实际计算中需要根据标记计算的效果再进行细微调整。

本文对扩展极值变换后的图像进行形态学开、闭运算和腐蚀运算的处理，基于前景和背景标记的结果，采用形态学强制最小值标定的方法修正原始梯度图像。修正后的梯度图像f1为[20]

f1=IMmin(f|I)

(9)

式中：IMmin表示形态学强制最小值标定运算，f表示原始梯度图像，I表示前景和背景标记图像。f1对应于标记图像I中为1的点，均被标定为0，其余像元点的灰度值与原始梯度图像保持不变。修正后的梯度图像中保留边缘轮廓信息，其余区域均标记为0，完成局部噪声的去除。在此基础上再开展分水岭变换的处理，得到单像素宽的准确水边线[22]。

2 实验与结果分析

2.1 研究区概况及数据源

本文选取南海海域的南薰岛和鸿庥岛作为2种岸线类型的研究区(图2)。南薰岛是吹沙填海而来的人工岛，岛屿面积约为0.18 km2，岸线类型主要为人工岸线。鸿庥岛为自然形成的小岛屿，面积约为0.08 km2，外观呈椭圆形，岸线类型主要为沙质岸线。两个岛屿海岸线类型特征明显。

图2 研究区地理位置示意图Fig.2 Geographical position of the study area

遥感数据选用GF-2卫星遥感影像，2个岛屿的卫星影像获取时间分别为2017年5月29日和2017年3月21日。高分二号全色影像空间分辨率为1 m，多光谱影像空间分辨率为4 m。其中全色影像波段为450～900 nm；多光谱影像分为4个波段，分别是450～520 nm(蓝波段)、520～590 nm(绿波段)、630～690 nm(红波段)、770～890 nm(近红外波段)[23]。本文使用其多光谱影像数据进行实验。

2.2 实验结果与分析

2.2.1 数据处理

实验首先对南薰岛原始遥感图像进行辐射定标、大气校正等预处理(图3a)，再进行NDWI运算，并建立Sobel梯度图像(图3b)，图像边缘信息得到了加强和突出，但由于存在过多的噪声点，会出现严重的过分割现象(图3c)。在NDWI运算后进行形态学重建滤波(图3d)，图像中海陆像元灰度值对比差异显著，噪声在一定程度上得到消除，同时边缘轮廓信息也能得到较好的保留。其次，针对陆地部分的噪声，在形态学重建后进行扩展极小值变换运算，为保证能得到分割结果，再进行形态学开、闭运算和腐蚀运算来收缩边缘，得到陆地标记结果(图3e)，图像中陆地部分的灰度极小值点被消除，但水体部分由于存在过多的灰度极大值点，会对分割结果造成干扰。针对水体部分存在过多的灰度极大值点，采用本文提出的扩展极大值变换运算方法，在形态学重建后进行扩展极大值变换运算，同时收缩边缘，得到水体标记结果(图3f)，可发现图像中水体部分的灰度极大值点和区域得到有效的去除，水体部分过分割现象基本消除。结合前景(陆地)、背景(水体)标记，采用形态学强制最小值标定方法对图3b进行修正，得到修正后梯度图像(图3g)，并在此基础上进行分水岭变换，图3h为分割结果，中间单像素宽的分割线即为水边线，叠加在原图上如图3i所示。

(a)预处理后图像

(d)形态学重建图像

(g)修正后梯度图像

对比本文方法提取的水边线和图像中的水边线，轮廓封闭，形状较为规则，对于伸向海洋的码头边缘定位也较准确，两者基本重合，说明本文所用算法对南薰岛的陆地与海水部分能够较为准确地分割与提取。

同样，用本文算法对以沙质岸线为主的鸿庥岛进行实验，实验结果如图4所示。在对预处理后图像进行NDWI运算和形态学重建后，图像上的大部分噪声点已经消除，但依旧存在部分灰度极值区域对影像梯度修正的准确性造成干扰。采用扩展极小值变换对陆地部分的灰度极小值区域进行限制，同时采用扩展极大值变换对水体部分的灰度极大值区域进行限制，完成了前景和背景标记，图像中陆地和水体内部的灰度极值区域分别得到有效的去除。然后，采用形态学强制最小值标定方法和分水岭变换后得到单像素宽的水边线。从结果中发现，沙质岸线同样可以得到较为准确的海陆分割和水边线提取结果，对沙滩能够进行准确识别，水边线提取结果与图像中的水边线基本吻合。

(a)预处理后图像

(c)过分割图像

(e)扩展极小值变换后陆地标记

(g)修正后梯度图像

(h)分割结果叠加图像

扩展极值变换运算中，对于h1和h2的设定是变换的关键所在。根据公式(7)和(8)的阈值确定方法，南薰岛影像的h1和h2值分别为55和60。在上述阈值基础上，分别设定过小和过大的阈值进行实验，实验结果如图5所示。从结果中发现，当h1和h2均设定过小时(图5a)，图像中仍会存在过多的灰度极值区域，最终结果会出现过分割现象；当h1设定过小、h2设定过大(图5b)或h1设定过大、h2设定过小(图5c)时，会使得图像中的灰度极值区域被过分消除或合并，最终均会出现错误分割的结果；当h1和h2均设定过大时，会使得图像中的灰度极值区域被完全消除，最终无法得到分割结果，使整幅影像呈现单一灰度值。由此说明扩展极值变换阈值设定得过大或过小均会使得最终结果出现过分割、错误分割现象，甚至无法得到分割结果。

在以上实验过程中，由于部分水陆边缘光谱信息较复杂，水陆差异对比不显著，导致分水岭变换时部分海陆边界定位不够准确，对水边线的提取精度存在一定影响。

2.2.2 精度评价

为定量评价本文方法提取水边线的精度和可靠性，以目视解译提取的水边线作为参考，采用缓冲区分析法[24]开展精度验证。该方法中，先对参考水边线以一定的像素宽度建立缓冲区，提取出的水边线如果落入到该缓冲区内，则认为该段水边线提取结果与参考水边线相匹配，其长度记为TP1；如若不在缓冲区内，则认为不相匹配，长度记为FP。对采用本文方法提取出的水边线也建立相同像素宽度的缓冲区，如果参考水边线落入到缓冲区内，则认为二者相匹配，记其长度为TP2；如果不在缓冲区内，记其长度为FN。基于上述分析，建立完整度(Complete)、正确度(Correct)和质量(Quality)这三个评价指标：

Complete=TP2/(TP2+FN)

(10)

Correct=TP1/(TP1+FP)

(11)

(12)

长度相对误差也常用于评价海岸线提取精度，其计算方法为：(提取的海岸线长度—参考海岸线长度)/参考海岸线长度。该值越低，说明长度误差越小，精度越高。

通过以上4个评价指标验证本方法提取水边线的精度，其结果如表1和表2所示。

表1 南薰岛提取结果精度分析Tab.1 Accuracy analysis of Nanxun Dao extraction results

表2 鸿庥岛提取结果精度分析Tab.2 Accuracy analysis of Hongxiu Dao extraction results

从上述结果可以发现，南薰岛的长度相对误差为0.57%，鸿庥岛的长度相对误差为-1.09%，说明两个岛屿水边线提取结果的长度精度均较高。对于南薰岛，当缓冲半径为0.5个像素时，完整度为92.08%，正确度为91.55%，质量为84.87%，提取效果良好；当缓冲半径为1个像素时，完整度为96.67%，正确度为96.00%，质量为92.93%，提取的准确性较好，定位精度较高，说明该方法对于人工岸线的提取精度能够达到1个像素。对于鸿庥岛，当缓冲半径为0.5个像素时，质量仅有69.09%，提取效果较差；当缓冲半径为1个像素时，正确度达到90%以上，说明提取结果较为准确；当缓冲半径为1.5个像素时，完整度、正确度和质量均在90%以上，此时提取效果较好，说明该方法提取沙质岸线的精度能达到1.5个像素。

3 结论

针对分水岭变换方法在高分辨率多光谱卫星影像数据应用中存在同物异谱和异物同谱等因素引起过分割的问题，本文提出了一种基于扩展极值变换标记分水岭的海岛水边线提取算法，并采用高分二号卫星多光谱影像，在南海海域的南薰岛和鸿庥岛进行了水边线提取的方法验证，在现有的标记分水岭算法基础上，增加了扩展极大值变换的影像数据背景标记，再经过强制最小值标定和分水岭变换后提取海岛水边线信息。实验中，人工岸线为主的海岛水边线提取精度在1个像素之内高于90%，沙质岸线为主的海岛水边线提取精度在1.5个像素之内高于90%。实验结果表明：本文所用的改进算法能够较好抑制过分割现象，对于人工和沙质岸线提取的精度均较高，对于面积较小的岛屿也能准确提取出其水边线信息。部分区域受自然因素和影像光谱信息影响，使得分水岭变换后定位不准确，存在水边线提取结果出现细微偏差的现象。同时，本研究缺少海岛的多期观测数据和潮汐数据，未对其提取的结果进行潮汐改正计算，仅为海岛水边线。