基于不同岸线类型的象山港岸线提取及变化研究

2022-07-08徐一凯李阳东朱宇航

海洋技术学报 2022年3期

徐一凯，胡松，李阳东，朱宇航

（上海海洋大学海洋科学学院，上海 201306）

海岸带是沿海地区地形变化最快的地区之一[1-2]。海岸线变化原因主要分为自然因素和人为因素[3]，获得岸线信息对于研究自然与人类环境尤为重要。提取海岸线可通过实地调查获得[4]，该方法精度较高，但进行大规模海岸线测量非常昂贵和耗时。此外，受到地理条件限制，有时很难对无法进入的地区进行测量。遥感（Remote Sensing）技术则可提供实地调查无法获得的即时、大规模的海岸线信息[5-7]。

利用遥感技术提取海岸线的方法主要有目视解译和自动提取。目视解译依据专业人员的主观判断，通过目视观察提取岸线，操作简单、提取岸线精度高，但工作效率低、工作量大，不适合快速提取大范围遥感影像数据。与目视解译相比，自动提取具有提取速度快、效率高的优势，但存在同物异谱和同谱异物的问题，因此，需要设计科学的提取算法。常用的自动提取方法有：基于边缘检测的方法[8]、基于指数分析的方法[9]、基于阈值分割的方法[10]、基于区域生长的方法[11]、基于神经网络的方法[12]和亚像素方法[13]。然而，在进行自动提取时，往往利用单一算法对岸线进行提取[14]，未考虑单一算法对不同类型岸线提取的适用性，所以本文针对不同海岸类型评测不同方法对其的适用性，采取适用方法进行该海岸类型的岸线提取，从而提高岸线提取精度，并对提取出的岸线进行分析，总结研究区近年来岸线变化的原因。

1 研究区域及数据方法

1.1 研究区域概况

研究区域象山港位于宁波市东南部（29˚24′N—29˚46′ N，121˚25′ E—122˚00′ E），地处穿山半岛和象山半岛之间，沿东北—西南走向为一个狭长的半封闭式海湾。全港长达60 多千米，水深一般10~15 m，海岸线蜿蜒曲折。港内有西沪港、黄墩港、铁港3 个分港。

1.2 数据获取

1.2.1 岸线提取数据

鉴于象山港区域近15 年间岸线变化较为剧烈，研究其岸线变迁需要进行较长时间跨度的对比分析，故本文从地理空间数据云（http：//www.gscloud.cn）选取了2002 年1 月3 日的Landsat 5 TM、2010年11 月1 日和2017 年3 月9 日Landsat 7 ETM 影像数据。其中，Landsat TM 影像包含7 个波段，Landsat ETM+影像数据包括8 个波段（表1）。

表1 Landsat TM 波段和Landsat ETM+波段

1.2.2 岸线分类数据

为比较不同的岸线提取方法在提取不同岸线时的精确程度，本文采用2017 年9 月的谷歌地球（Google Earth）全球高清影像进行岸线分类。通过目视解译方法，以离岸线500 m 以内区域的主要用地类型作为该段岸线的类型进行手动划分。

1.3 数据方法

不同类型的海岸有不同的地面特征，单一的算法通常难以保证岸线提取的准确性，故本文采用了不同的阈值分割法对象山港区域影像进行了图像二值化处理，再对二值化后的图像利用Canny 边缘检测算子进行了边缘线提取。具体技术流程如图1所示。

图1 岸线提取技术流程

1.3.1 图像二值化

图像二值化就是利用阈值将原始图像分割成前景和背景两幅图像。因此，图像二值化的关键是选择最优阈值，在取阈值时，背景与前景的差异应最大。

（1）大津法（OTSU 法）

最优阈值是衡量差值的标准，而在OTSU 算法中，衡量差值的标准是最大类间方差。

OTSU 算法的原理是将图像中每个像素点的灰度值与假定的阈值T进行比较，将像素灰度值小于阈值T的归类为C1，大于阈值T的归类为C2。假设归类为C1和C2的均值分别为m1、m2，图像全局均值为mG，而此时像素被分为C1和C2类的概率分别为p1、p2，计算公式如下所示。

类间方差σ2表达式如下所示。

把上式化简，可得

然后，遍历所有灰度级，求出使σ2最大的灰度值k就是OTSU 阈值。

（2）改进的归一化差异水体指数法（MNDWI 法）

遥感图像对于不同类型的地物反射率不同，从而呈现出的亮度也不同。如砂质海岸的干燥滩面反射率较高，在遥感图像上表现为较亮区域；与之相反，潮湿滩面反射率较低，在遥感图像上表现为较暗区域；由于海水反射率最低，在遥感图像上表现为最暗区域。遥感图像瞬时水边线应为潮湿滩面和干燥滩面的分界线，由于此分界线并不明显，若直接提取，则精度不高[15]。

针对此问题，本文采用MNDWI 来增大潮湿滩面和干燥滩面的差异，并进行岸线提取。其公式如下。

式中，P（Green）和P（MIR）分别代表绿波段（Landsat TM/ETM 中的波段2）和中红外波段（Landsat TM/ETM 中的波段5）的亮度值。

利用此公式将水体凸显出来，再设置阈值来提取水体，进行水陆分离，实现图像二值化。

1.3.2 Canny 算子提取水边线

Canny 边缘检测算子是一种多级边缘检测算法。尽管Canny 算子的操作较为复杂，但相比Sobel 边缘检测算子，其能处理包含噪声较多的图像，且具有Robert 算子定位精确的优点。Canny 算子具有提取位置精度高，包含信息多的优点，近年来已被广泛用于边缘提取等图形检测中[16]。Canny 算子原理包括以下4 个步骤：①用高斯滤波器平滑处理原图像；②梯度幅值和方向的计算；③消除冗余窗口；④双阈值算法检测和连接边缘。

2 数据处理与步骤

2.1 岸线分类步骤

本文根据海岸地理环境特点，将海岸类型分为人工海岸和自然海岸两个一级类别，其中人工海岸分为养殖区、农田、码头、建设用地、裸地5 个二级类别，自然海岸分为淤泥质海岸、砂质海岸、基岩海岸3 个二类级别，具体如表2 所示。

表2 海岸分类

参考文献[11]，根据象山港特点，本文对于淤泥质海岸和砂质海岸，非海水部分在涨潮作用下可以被海水覆盖，因此上述海岸的非海水部分应与海水一起标记为“可涨潮区域”，另一方面，基岩海岸在涨潮作用下不被海水覆盖，应标为“不可涨潮区域”。人工海岸是海岸的一种改良形式，以养殖区、农田、码头、施工围堰等为特征，在涨潮时不被海水覆盖，应标为“不可涨潮区域”，对于浮阀养殖区进行水下潜水或在涨潮时可被海水覆盖的区域，应标为“可涨潮区域”[17]。

在进行岸线划分时，Google Earth 高清影像部分区域颜色相近，如农田和围塘养殖（图2 框示区域），无法直接进行准确划分，所以对象山港整个海岸线进行了实地考察，并最终得出了较为准确的岸线分布情况。

本文利用了Google Earth 全球真彩色影像，结合上述实地考察，对岸线类型进行目视解译并提取岸线，将其导入ArcGIS 进行要素的整合，获得岸线类型及分布如图2 所示。

图2 象山港岸线类型及分布

2.2 岸线提取步骤

2.2.1 图像预处理

为获得象山港区域的完整图像，将同时期的两幅Landsat 影像并进行拼接处理，并按照象山港的经纬度范围进行裁剪。由于2003 年Landsat7 ETM+机载扫描行校正器故障，导致其之后图像会出现条带，所以本文使用了条带插件将2010 年及2017 年的影像进行修复。

由于需要对不同时间、不同传感器获取的图像进行比较，对去条带后的遥感图像进行了辐射定标，将图像的亮度灰度值转换为绝对的辐射亮度。

然后，对辐射定标之后的图像进行了大气校正，以消除大气影响所造成的辐射误差，将辐射定标值反演为地表真实信息，从而保证地物波谱信息的准确性。

最后，对各时期图像进行几何校正，消除及减弱由各种内外因素造成的遥感图像的几何畸变。

2.2.2 OTSU 结果

由于基岩海岸和人工海岸在遥感图像中都具有较明显的水陆分界，利用边缘检测便可提取水边线，其因不考虑海岸线的背景差异，因此获得的海岸线位置一般较为准确。

然而，在遥感图像采集过程中，周期性传感器偏移或电磁干扰引起的噪声会降低图像质量，因此在进行边缘检测之前，需要对图像进行去噪处理。因为中红外波段对于水陆分离情况较好，所以将经过大气校正后的图像中的MIR 灰度图像进行中值滤波以减少噪声对结果的影响，然后采用OTSU 法对该图像进行二值化处理。

2.2.3 MNDWI 结果

本文将影像导入ENVI 中，利用公式（5）将水体凸显出来，再尝试并最终设置了合适的阈值来提取水体，进行水陆分离，实现图像二值化。

2.2.4 Canny 算子提取水边线结果

基于在上文中提及的Canny 算子步骤，本文在OTSU 算法所得二值图像的基础上对其进行形态学操作，然后通过Matlab 编程利用Canny 算子分别对两种二值化方法得到的结果图进行边缘检测并消除孔洞，以保留主要岸线及岛屿轮廓。以2017 年为例，结果如图3 所示。

如图3 所示，利用Canny 算子提取出的水边线连续、清晰且较为平滑，虚假边缘被很好地抑制，受噪声干扰小，提取效果较好。将Canny 算子计算过后的图像导入ArcGIS 中并进行栅格数据转矢量数据处理得到矢量海岸线，并将不同年份的岸线进行叠合以进行后续的观察和分析。

图3 2017 年二值化后Canny 算子运算结果

3 结果分析

3.1 算法比较

如上所述，本文分别采用了OTSU 法与MNDWI 法对图像进行了二值化处理，使用MNDWI法提取出的岸线比使用OTSU 法提取的岸线更向海推进，尤其是箭头指向区域。将真彩色影像进行叠加分析，发现相差部分为淤泥质和砂质海岸（图4区域1 和区域2），从而验证了MNDWI 法相比OTSU 法更适合提取砂质和淤泥质海岸。而在提取基岩岸线及人工岸线时，OTSU 法更加精细及平滑（图4 区域3）。

基于上述情况，本文将砂质岸线、淤泥质岸线及浮筏养殖区按照MNDWI 法提取，将人工岸线和基岩岸线按照OTSU 法提取，将两者进行矢量裁剪并拼接，得到最终的象山港岸线。

3.2 岸线类型变化分析

岸线类型与经济发展、航行安全、生态保护有密切联系，不同的岸线类型往往对应了不同的用海功能，港口的开发往往选择基岩岸线，而淤泥质岸线地形可以用作天然的盐场。象山港地区在进行经济开发的同时，也秉着生态保护的理念，保证经济发展的同时，打造多元化的海岸区域[18]。

为研究各岸线类型的占比变化，本文利用ArcGIS 分别对2002 年、2010 年及2017 年各岸线类型进行了长度计算，结果如表3 至表5 及图5 所示。

图5 2002—2017 年各岸线占比变化

表3 2002 年各类岸线统计

表5 2017 年各类岸线占比

研究结果显示，2002—2017 年，各类自然岸线占比均呈现持续下降趋势，总占比由2002 年的45.83%下降至2010 年的38.30%，在2017 年仅剩31.94%；而各类人工岸线占比却在持续上升，总占比由2002 年的54.17%上升至2010 年的61.70%，在2017 年达到了68.06%。其中，自然岸线中的基岩岸线变化较小，而淤泥砂质岸线由2002 年至2017 同比下降44.79%；而人工岸线中的养殖区由2002 年至2017 年同比上升43.05%。

3.3 岸线变迁变化分析

表4 2010 年各类岸线占比

为具体分析其变化原因，本文按照行政区划图，把象山港分为郭巨—大嵩、大嵩—桐照、铁港、黄墩港、峡山—乌沙、西沪港、西泽—钱仓七个区域，然后根据上文所做出的2002 年、2010 年、2017 年岸线变迁图（图6），分析15 年间象山港岸线的变化情况。可以看到，2002—2017 年间，郭巨—大嵩、大嵩—桐照末尾、铁港、黄墩港、西沪港变化较大，其他地区变化较小。本文将变化较大区域用箭头标明，并对照高清卫星图对该区域变化原因进行分析。

图6 2002—2017 年岸线变迁图

主要变化如下：郭巨—大嵩：2002—2010 年大嵩附近岸线出现较大的扩张（区域2），2010—2017年出现巨大的扩张变化（区域1），查阅资料得知，2002—2010 年的岸线变化是由于围海造陆所致，而2010—2017 年的变化是由于凸出的部分是原有岛屿（梅山岛），在2010—2017 年间梅山水库的建设将岛屿和大陆连接，从而致使岸线出现了巨大变化（图7（a））。

大嵩—桐照末尾：2002—2010 年桐照西侧岸线出现较大扩张（区域4），2010—2017 年桐照东侧岸线出现较大扩张（区域3），2002—2010 年桐照西侧区域进行了围海造陆，而桐照东侧出现岸线变化原因如区域1 相似，即将原有岛屿规划进陆地（图7（b））。