在气候变暖、海平面逐渐上升导致的全球岸线普遍后退的背景下, 国内很多沿海地区大规模的围垦导致岸线向海推进^[1]。这对快速准确获取海岸线提取出要求, 传统人工获取岸线的方法已无法满足需求。用面向对象的方法提取海岸线, 一定程度上可去除影像分类的“椒盐效应”^[2]。面向对象处理的是影像对象而不是像元, 因此可以充分利用对象的空间特征和光谱特征, 可以有效地克服基于像元层次分类的不足^[3]。近年来国内外很多学者在基于面向对象分割的基础上结合不同的指数确定阈值和分类方法对海岸线进行提取。贾明明等^[4]利用多源遥感数据在面向对象分割的基础上对Landsat卫星数据选取SWIR波段进行阈值分割提取杭州湾海岸线, 精度较高。巢子豪等^[5]基于改进化水体指数和潮位校正模型, 使用数字岸线分析系统提取了不同年份的海岸线, 杨红等^[6]以多源遥感影像为数据源, 使用基于像素的方法对大连新港溢油区域进行溢油面积提取, 获得较好结果。吴小娟等^[7]基于多尺度分割的基础上使用“高分二号”卫星数据根据在近红外波段水体具有强吸收的特性, 选取近红外波段合适的阈值提取了海岸线。鉴于目前, 很少有研究使用LBV变换中的B分量作为阈值分割条件提取海岸线。B分量值相当于水体特征图像, B值反映的地面水分的多少B值越大, 表明地面水分含量越高^[17]。选取合适的B分量阈值可以较好地实现水陆分离。本文在根据相关研究基础上, 推导出适合Sentinel-2A的B分量公式对如东海岸的淤泥质海岸进行海岸线提取。为验证该方法的适用性, 同时对交通围堤海岸和养殖围堤海岸进行岸线提取及验证, 由于本文使用的图像获取时刻处于高潮时, 实际提取的海陆分界线为瞬时水边线。

1 研究区域与数据 1.1 区域概况

研究区域位于江苏省如东县(32°12~32°36′N, 120°42′~121°22′E), 地处长江三角洲北翼, 位于江苏省东部, 东北方向濒临黄海, 境内海岸线全长102.59 km, 海域面积大于6 000 km²。

1.2 遥感数据源

欧洲航天局的哥白尼计划于2015年6月23日发射的Sentinel-2A, 重访周期为10 d, 携带的多光谱仪(MSI)可测量从可见光到短波红外(SWIR)的13个光谱波段, 其中4个波段(Blue, Green, Red和NIR)空间分辨率为10 m, 6个波段(包括SWIR)空间分辨率为20 m^[8]。Sentinel-2A影像数据可免费从ESA(欧空局)的Open Access Hub(https://scihub.copernicus.eu/)分发系统中下载。本文采用2017年11月15日Sentinel-2A的Level-1C级遥感影像, 天气状况晴好, 云量小于5%。是经正射校正和亚像元几何精校正的大气表观反射率产品, 几何纠正误差控制在一个像元以内。

2 基于面向对象方法的海岸线提取

由于影像在获取的过程中受到各种因素的影响会产生几何和辐射畸变^[9]。为正确的评价地物的辐射特征和几何特征, 需要对影像进行预处理操作, 预处理的主要步骤包括辐射定标、大气校正、影像剪裁和影像融合。面向对象提取岸线的过程中使用e-Cognition软件, 提取流程主要为多尺度分割、光谱差异性分割、B分量阈值选取、分类、分类后处理。具体技术路线见图 1。

2.1 遥感数据预处理

提取岸线所用的指数是依据影像的光谱值计算得到, 利用欧空局提供的Sentinel-2 Toolbox软件对具有大气表光反射率的产品影像进行辐射定标和大气校正。为减少运算成本, 在进行融合图像之前先用ENVI软件对图像进行剪裁。在四个10 m波段中, NIR波段和PAN相似波段具有最好的相关性系数^[10-11], 选择NIR波段作为PAN相似波段与20 m多光谱波段通过ENVI中的Gram-Schmidt算法进行融合, 通过Gram-Schmidt方法融合得到的影像与原始影像多光谱信息比较接近, 色彩失真较小, 并且可以清晰地辨别地物的细节信息。通过ENVI的Layer Stacking功能对10 m的四个波段和融合后20 m的六个波段进行处理, 最后得到10 m的研究区多光谱影像。对空间分辨率为10 m的波段进行标准假彩色合成(波段8, 4, 3标准假彩色合成), 影像结果见图 2。

2.2 B分量公式推导

用预处理后Sentinel-2A的多光谱影像作为B分量推导的基础影像, 通过对不同的波段组合做对比试验, 选取提取水体信息效果较好并且精度较高的Blue、Red、Vegetation Red Edge、SWIR波段组合, 其对应中心波长分别为λ₁=0.490 μm、λ₂=0.665 μm、λ₃=0.842 μm、λ₄=1.610 μm, 对应的灰度值记为D_N。推导出适用于Sentinel-2A多光谱影像水体信息的B分量公式, 用一次线性回归方程(1)即可

$ {D_\lambda } = a + b\lambda . $

(1)

根据文献[12-13]和公式(1)求B分量公式:

第一步, 根据公式(1)和各中心波长计算出方程的斜率。

$ b = {\rm{ }}-0.5731{D_1}{-_{}}0.3329{D_{2}}-0.0583{D_3} + 0.9643{D_4}. $

(2)

第二步, 根据文献[13-14]可用b的相反数表示B₀,

$ {B_0} = {\rm{ }}0.5731{D_1}{ + _{}}0.3329{D_2} + {\rm{ }}0.0583{D_3}-0.9643{D_4}{._{}} $

(3)

第三步, 将4个波段的灰度值分别拉伸到0~255范围, 对拉伸后的波段应用公式(3)做波段运算得到B₀图像, 对B₀图像做统计分析得到标准差(σ)为40.9507, 根据文献[13]令κ=25/σ, κ与B₀相乘得到B₁。

$ {B_1} = {\rm{ }}0.3525{D_1}{ + _{}}0.2032{D_2} + {\rm{ }}0.0356{D_3}-{\rm{ }}0.5887{D_4}{._{}} $

(4)

第四步, 对B₁图像做统计分析得到B₁图像的灰度平均值为–2.8138, 根据文献[13]令A=128–m, A与B₁相加得到研究区Sentinel-2A多光谱影像的B分量公式。

$ B = {\rm{ }}0.3525{D_1}{+}0.2032{D_2} + {\rm{ }}0.0356{D_3}- \\ 0.5887{D_4} + {\rm{ }}130.8138. $

(5)

不同区域Sentinel-2A多光谱影像的B分量公式是有差异的, 根据需要目标区域融合后的Sentinel-2A多光谱影像完成第三步到第五步的操作, 即可以得到相应区域的B分量公式。Sentinel-2A影像波段经过B分量运算后, 水体的细节信息变得更加突出, 见图 3。

2.3 海岸线提取 2.3.1 影像分割

本文采用多尺度分割与光谱差异分割相结合方式对影像进行分割, 多尺度分割是一种自下而上的分割算法, 通过合并相邻的像素或小的分割对象, 在保证对象与对象之间平均异质性最小、对象内部像元之间同质性最大的前提下, 基于区域合并技术实现分割^[15]。光谱差异分割是在多尺度分割的基础上进行分割的, 多尺度分割可以将尺度参数设置的较小, 实现对影像的过分分割。光谱差异分割对多尺度分割的结果进行优化, 减少分割数量, 这样的组合方式在很大程度上减少了手动选取阈值的次数。根据SWIR波段对水体具体强吸收性, 淤泥具有强反射性^[16]。在分割之前, 先对SWIR(Band11)波段进行Canny边缘检测, 再把边缘检测的结果作为一个新的波段参与到多尺度分割当中去^[2]。经过不断实验, 得到多尺度分割的最适分割尺度为10, 形状指数为0.6, 紧凑度为0.5, 各波段权重均为1。光谱差异分割的最适分割尺度为100。影像分割的局部结果见图 4。

2.3.2 B阈值选取及分类后处理

B反映的地面水分的多少, B值越大, 表明地面水分含量越高。当地面完全被水占据时, 具有最高的一些B值^[17]。在Sentinel-2A的B分量公式中, SWIR波段前系数的绝对值权重最大, 这是由于Sentinel-2A中的SWIR波段相对于Sentinel-2A中的其他波段, SWIR波段对土壤湿度的变化更加敏感, 可以更准确地反映淤泥海岸地表水含量^[18]。将影像中B大于138的分为海水, B小于138的分为陆地; 对陆地使用查找封闭类运算进行二值化处理, 不对海水使用查找封闭类运算, 是因为部分碎边缘在合并运算时会被直接归类为海水; 使用合并区域算法对陆地和海水进行合并运算; 研究区域内有部分河流直接与海水相连, 划分的原则为:保留大型河口的港湾特征, 反映河口处沙咀、沙洲和瀉湖地貌, 表现出河口的喇叭形状^[7]。

考虑到碎边缘的保留问题, 在进行形态学封闭运算的过程中选取较小的Mask值, 对保留河口的喇叭形状有利, 选取较大的Mask值, 对封闭碎边缘有利。经过多次试验, 选取最适的Mask值为11, 最终得到淤泥质海岸的水边线提取结果。为验证该方法的适用性, 通过调整B分量阈值分别对交通围堤海岸和养殖围堤海岸进行水边线提取。

3 精度验证

淤泥质岸线受涨落潮的影响较大, 需经过潮汐校正位置高程后才能得到较为准确的岸线位置, 由于缺乏验潮站数据, 本文不对提取的淤泥质水边线做定量精度评价。将提取的淤泥质水边线与原始影像叠加见图 5中a和b, 可以看出提取的淤泥质水边线较为理想。

交通围堤海岸和养殖围堤海岸受涨落潮影响较小, 为方便验证选取部分的交通围堤水边线和养殖围堤海岸水边线, 以同期的亚米级Google Earth卫星影像上研究区域内相应海岸类型的水边线作为基准线, 在ArcGIS里将提取的水边线以50 m为间隔离散成点做近邻分析。对分析的结果做最小值、最大值、平均值和标准差做统计分析, 验证结果见表 1。

从表 1统计的结果可以看出, 提取的交通围堤海岸和养殖围堤海岸的水边线平均值相近。对于提取的养殖围堤海岸类型的水边线的最大值和平均值都大于交通围堤海岸类型的水边线, 其最小值小于交通围堤海岸类型的水边线。这表明提取的养殖围堤海岸类型的水边线整体效果不如交通围堤海岸类型的水边线, 主要原因是养殖围堤海岸类型与淤泥质潮滩的植被接壤, 植被生长的边界较为模糊, 这对提取养殖围堤海岸类型的水边线带来较大误差。对于与海水直接接壤的交通围堤海岸类型的水边线, 公路的B分量值与海水的B分量值具有较大的差异, 并且公路具有清晰的轮廓, 所以提取的水边线结果相对准确。将提取的水边线与具有相应海岸类型的影像叠加可以看出提取的结果较为准确(见图 5中c和d)。

4 结论与讨论

1) 通过对Sentinel-2A多光谱影像波段进行组合实验, 提出了适用于该传感器的B分量公式。并且使用B分量作为分割阈值可以有效地进行水陆分离。

2) 通过提取的水边线与影像叠加, 同时对提取的受潮汐影响较小的水边线做精度验证, 证明提取的结果较为准确。

3) 对于淤泥质海岸类型提取的水边线验证, 需要进行潮位高程校正后精度验证的结果才较为客观和准确。对于其他类型的海岸水边线的提取有待实验。