基于多算法水边线提取的潮滩DEM构建

李天祺，王建超，吴芳，赵政，张文凯

(中国自然资源航空物探遥感中心，北京 100083)

0 引言

淤泥质潮滩通常指平均大潮高潮线与平均大潮低潮线之间的地带，也称潮间带[1-2]。在我国约有四分之一的海岸属于淤泥质海岸[3]，获取淤泥质潮滩地形对研究岸线变迁、海岸带生态变化及沿海工程建设有着重要意义[4-5]。淤泥质潮滩受波浪、泥沙、沿岸流和地质等诸多因素的影响，具有面积宽广、水浅滩平、变化频繁的特点，使得该区域地面调查与地形测绘的难度较大、成本较高，且存在一定风险[6-7]。遥感技术具有高时效、大范围、高频率的特点，利用多期遥感水边线复合潮位数据构建潮滩数字高程模型(digital elevation model,DEM)，成为获取大范围淤泥质潮滩地形信息的有效途径。

遥感水边线是卫星过境时获取的瞬时水边线，准确提取水边线是潮滩DEM构建的关键。现有的水边线提取方法主要有： 边缘检测法[8]、阈值分割法[9]、区域生长法[10]、主动轮廓模型法[11]、面向对象分类[12]等。Mason等[13]基于ERS SAR影像，通过纹理分割提取多期水边线，应用水动力模型对水边线高程进行赋值，构建了英国东岸Humber/Wash区域潮滩DEM； 沈芳等[14]基于Landsat TM影像，对比了不同波段的水边线提取精度，利用阈值分割法提取水边线，通过潮汐预测数据插值对水边线赋值，构建了长江口九段沙DEM； 穆敬等[15]基于BJ-1影像，利用面向对象分类法提取水边线，通过潮汐网格数据对水边线赋值，构建了黄骅市潮滩DEM。目前国内外学者在构建潮滩DEM的研究中大多关注于准确模拟水边线的瞬时潮位和提升单一方法的水边线提取精度[8-15]。但是，瞬时水边线受卫星过境时的潮情、天气以及滩面等因素的影响，在影像中的光谱与纹理差异较大，基于单一方法提取水边线，难以保证水边线的提取精度，从而影响潮滩DEM反演结果。吴迪等[16]虽提出了基于多源多算法的水边线提取模型，但未对水边线进行系统划分，且没有将其应用于潮滩DEM构建研究。

针对上述问题，本文从边缘类型与滩面噪声2方面划分水边线，根据瞬时潮情对水边线进行分类提取，通过潮汐数据构建多项式计算瞬时潮位，提出基于多算法水边线提取的潮滩DEM构建方法，并通过野外实测数据对潮滩DEM反演精度进行验证。

1 研究区概况及其数据源

1.1 研究区概况

研究区位于N37°21.5′～37°32.0′，E118°53.5′～118°58.5′，位于莱州湾西侧，潍坊港验潮站以北，如图1所示。研究区主要为淤泥质潮滩，滩面广阔，坡度平缓，南北长约20 km，宽度为2～4 km，近岸区域有人工堤坝和大量盐田，滩面生长高耐盐红海滩植被，近海区域无植被生长。潮汐类型属于不规则半日潮混合潮，平均潮差为125 cm，最大潮差为259 cm。

图1 研究区地理位置Fig.1 Location map of the study area

1.2 研究数据源及预处理

本文选用GF-1 WFV数据，空间分辨率为16 m，重访周期为2 d，幅宽为800 km。数据获取时间为2017年2—7月，筛除受厚云、雪及海冰污染的数据，共选取了18期数据用于研究。数据预处理包括绝对辐射定标、几何精纠正及研究区裁剪，其中几何精纠正误差控制在0.5个像元以内。

由于受野外测量条件限制，缺少研究区的实测潮位数据，本文使用潍坊港验潮站的潮汐预测数据用于水边线的分类与卫星过境瞬时的水边线潮位计算。潮汐预测数据是在潮汐观测数据的基础上，利用调和分析法推算的预测潮位[14]。潍坊港验潮站所在海域为不规则半日混合潮，潮汐预测数据包含每日整时的潮位值，以及每日的潮位极值与对应时刻，潮高基准面为平均海平面以下120 cm。

2 多算法水边线提取

2.1 基于潮情的水边线分类

水边线作为遥感影像中一种重要的边缘信息，是水陆区域的分界。不同遥感水边线提取算法具有不同的适用范围[17]，其提取效果主要取决于目标边缘的类型与边缘邻域内的噪声水平。根据边缘类型与噪声水平对水边线分类，选择适用算法进行提取，有利于准确地获取水边线信息。

水边线两侧水陆区域的灰度变化梯度主要受涨落潮的影响： 涨潮时滩面干湿分明，水边线两侧区域灰度差异显著； 落潮时受滩面残余水和水中泥沙影响，水边线两侧灰度差异较小。因此按灰度变化梯度，边缘类型可分为强边缘与弱边缘。水边线的邻域噪声水平，主要由潮位高低决定： 由于滩面地物种类由近海向内陆逐渐增多，高潮位时水边线在潮沟、植被、建筑物分布较多的近岸区域，滩面破碎复杂，邻域内噪声较高； 低潮位时水边线近海，滩面相对均质单一。

综上，本文通过分析潮汐数据，确定各期水边线的瞬时潮情： 涨、落潮与高、低潮位，根据瞬时潮情将水边线划分为单一滩面弱边缘、复杂滩面弱边缘、单一滩面强边缘和复杂滩面强边缘4类，具体分类规则如图2所示。

Test and optimization on heat transfer performance of electric heater under lower air pressure

图2 分类规则

在本实验中，影像中的薄云或雾明显削弱灰度变化梯度，因此通过人工目视对影像质量进行判断，将包含薄云或雾影像的水边线分为弱边缘水边线； 通过分析潮汐数据对晴朗影像的水边线进行划分。若H(10),H(11)和H(12)3个时刻的潮位依次升高，则瞬时潮情(11—12时之间)为涨潮，滩面不存在残余水的影响，则水边线划分为强边缘水边线，若3个时刻的潮位依次下降(落潮)或10时的潮位高于11时潮位，则存在滩面残余水或水中泥沙的影响，划分为弱边缘水边线。

获得强、弱边缘水边线后，按潮位高低对水边线所在滩面特征分类： 若水边线获取当日的11时与12时的平均潮位高于多期平均值H，则该期水边线近岸，对应复杂滩面； 若低于多期平均值H则该期水边线近海，对应单一滩面。最后得到4类水边线分别为： 单一滩面弱边缘水边线、复杂滩面弱边缘水边线、单一滩面强边缘水边线与复杂滩面强边缘水边线。对本文所用的18期GF-1 WFV影像的水边线进行分类，分类结果如表1所示。

表1 水边线分类结果Tab.1 Waterline classification results

2.2 基于边缘检测水边线提取

边缘检测是基于图像灰度的一阶或二阶导数来提取图像中灰度急剧变化的区域边界[18]，具有模型参数少、运行速度快的优点，但对噪声敏感，在复杂滩面提取的水边线存在严重的毛刺或不连续。本文选择定位准确度高、单边响应良好的Canny边缘检测算法[19]，以20170429影像为例，对单一滩面强边缘水边线进行提取，结果如图3所示。

图3 基于Canny边缘检测的水边线提取结果Fig.3 Waterline extraction resultbased on Canny edge detection

2.3 基于阈值分割水边线提取

由于弱边缘水边线两侧灰度变化不明显，边缘检测不能有效提取边缘信息，阈值分割通过设定阈值，如归一化水体指数(normalized different water index,NDWI)[20]，可以快速有效地区分水陆区域，但对于高潮位的复杂滩面，受潮沟或其他积水区域的影响，水边线提取不准确，因此，本文通过阈值分割提取单一滩面弱边缘水边线。

以20170617影像为例，首先，计算NDWI，增强影像中的水体信息，计算公式为：

(1)

式中Green与NIR分别代表GF-1 WFV影像中的绿光波段与近红外波段的反射率。NDWI的取值范围为[-1, 1]； 然后，通过OTSU算法[21]计算阈值，对NDWI灰度影像进行二值化； 最后，对二值图做开运算处理，对边界进行平滑，得到水边线提取结果，如图4所示。

图4 基于NDWI阈值分割的水边线提取结果Fig.4 Waterline extraction result based onNDWI threshold segmentation

2.4 基于面向对象法水边线提取

面向对象法是通过影像中的光谱信息、空间信息和纹理信息对影像进行分割，然后以分割后的对象为基本单元进行分类，从而提取水边线的一种方法[22]。面向对象法具有较强的抗噪能力，有效抑制复杂滩面上其他地物的影响，可准确保留强边缘水边线的细节信息，但该方法对纹理变化过于敏感，提取的弱边缘水边线不平滑，因此，本文通过面向对象法提取复杂滩面强边缘水边线。

首先，通过多尺度分割将影像分割为多个对象，然后根据各个对象的光谱相似度进行合并，多次迭代后最终得到水边线提取结果，以20170609影像为例，初始分割结果和水边线提取结果如图5所示。

2.5 基于改进分水岭水边线提取

图6 改进分水岭算法流程Fig.6 Flow chart of improved watershed algorithm

以20170310影像为例，对影像NDWI低通滤波去除高频噪声后，潮滩上的潮沟、河口处的沙洲及较小的人工建筑物均已去除； 二值化处理后通过数学形态学处理，将盐田、河流等内陆水体完全去除，在边界处做缓冲区处理，生成待分割区域的标记信息； 利用原始影像B4(R)，B3(G)，B2(B)假彩色图像结合标记信息进行分水岭分割，提取结果如图7所示。

图7 基于改进分水岭的水边线提取结果Fig.7 Waterline extraction result basedon improved watershed algorithm

3 潮滩DEM构建与精度验证

3.1 潮滩DEM构建

利用潍坊港验潮站的潮汐预测数据对影像获取的瞬时潮位进行插值。本文所用GF-1 WFV数据的获取时间均在每日上午11—12时之间，利用10—13时的4个整时潮位值构建三次多项式，若在10—13时之间存在潮位极值，则用潮位极值替代距影像获取时间较远的整点潮位值，通过多项式插值计算水边线的瞬时潮位。在研究资料有限时，这是一种行之有效的方法[14]，水边线对应瞬时潮位如表2所示。

表2 水边线瞬时潮位Tab.2 Instantaneous tidal level of waterlines (cm)

受滩面上潮沟或河流入海口处泥沙的影响，提取的多期水边线在局部存在交叉或重叠的现象。针对上述问题，本文按照单一滩面强边缘、复杂滩面强边缘、单一滩面弱边缘和复杂滩面弱边缘的优先顺序选用水边线，结果如图8所示。

图8 水边线选用结果Fig.8 Selection result of waterlines

将筛选后的水边线重采样为16 m×16 m的水边点，利用不规则三角网进行空间插值，并将插值结果从潮高基准面转换为以WGS84为基准的高程数据，得到潮滩DEM反演结果，如图9所示。

图9 潮滩DEM反演结果Fig.9 Inversion result of tidal flat DEM

3.2 精度验证

为验证基于多算法水边线提取的潮滩DEM构建方法有效性，本文利用野外实测高程数据对潮滩DEM反演结果进行精度评价。实测高程数据获取时间为2017年12月，采用Trimble R5 GPS接收机测量，静态观测不少于20分钟。实测高程数据通过精密星历解算为WGS84坐标数据，18个验证点分布如图10所示。

图10 验证点分布Fig.10 Distribution of verification points

对实测高程与反演结果进行相关分析，实测高程与对应DEM反演结果的相关性较高，R2为0.864 9，如图11所示，误差分布在0.31～0.78 m之间，中误差为0.173 4 m。

图11 相关性分析Fig.11 Correlation analysis

4 结论

本文针对遥感水边线法潮滩DEM构建中，应用单一方法提取多期水边线精度不高的问题，提出了基于多算法水边线提取的潮滩DEM构建方法，并以莱州湾西侧潮滩为例验证了方法的有效性。主要结论如下：

1)分析了不同潮情下水边线的边缘类型与滩面噪声的差异，并基于潮汐数据对水边线进行分类，将研究数据中的多期水边线分为： 单一滩面弱边缘水边线、复杂滩面弱边缘水边线、单一滩面强边缘水边线与复杂滩面强边缘水边线4种类型，该结果可进一步用于水边线提取与潮滩DEM构建。

2)基于多算法提取水边线，通过潮汐预测数据插值获得多期水边线瞬时潮位，并根据水边线类别进一步筛选，构建了潮滩DEM。研究区DEM反演结果精度较高，可准确反映潮滩的近似地形，方法可有效用于潮滩地形信息的获取。

但本文所选研究区与研究数据的代表性有限，可进一步完善水边线分类规则，以适用于其他区域。通过潮汐预测数据插值计算瞬时潮位与水边线实际潮位存在一定偏差，使用实测潮位数据可进一步提高研究区潮滩DEM的反演精度。