基于遥感和OpenCV 的海域使用疑点疑区自动提取与变化检测方法研究

付弘涛，张昊睿，秦 平，王丽琳

（自然资源部东海预报中心，上海 200136）

自2006 年以来，原国家海洋局在全国11 个沿海省（自治区、直辖市）启动国家海域动态监视监测管理系统建设，并逐步开展海域动态监视监测工作。该系统利用“岸—海—空—天”一体化监测手段，实现对我国海域开发利用活动的立体和动态监视监测，全面掌握海域资源、海域使用和海域管理等状况[1]。作为海域动态监视监测工作的重点之一，海域使用疑点疑区监测对全国疑似违法违规用海行为进行监测，是贯彻落实《国务院关于加强滨海湿地保护严格管控围填海的通知》（国发〔2018〕24 号）重要精神的具体手段。开展基于卫星遥感的海域使用疑点疑区全覆盖业务化监测，可迅速、准确和及时地掌握用海情况，保护海洋生态环境，避免国有资产流失，更好地践行“绿水青山就是金山银山”的理念。

海域使用疑点疑区是指通过前后两期卫星遥感影像的对比分析，海岸线向海一侧某区域未经批准（依据海域法审批）由海水变为陆地或露出水面的其他形态的区域。海域使用疑点疑区监管手段主要依靠卫星遥感，数据源以国产高分系列卫星为主。国内外针对遥感影像变化识别已提出了大量方法，目前较为主流的方法即基于地物光谱特性的比较法和基于像元特征的比较法。前者把不同时相影像先分类再比较，简单且应用广泛，但耗时耗力，地面验证工作难；后者以像素作为检测单位直接进行变化检测，再进行判别分析，利用像素级多特征融合，使检测结果更加可靠，但变化检测的指标规则、变化阈值难判定[2]。

基于遥感的自动变化检测作为最重要的关键技术之一，越来越多的应用于国土资源管理、地物变化的监测等领域[3]。但多数变化检测算法主要是针对特定的条件（区域和影像质量） 具有较好的效果，还没有较为通用性的算法，且现有算法在效率、精度、智能性等方面存在的问题，大多算法解决的问题及理论相对分散。本文旨在落实自然资源部提出的利用卫星遥感等多种手段实现对我国管辖海域“15 天一覆盖”的巡查总体要求，研究探索一套适用于海域使用疑点疑区监管工作的自动提取与动态变化检测方法，同时探究影像像元“椒盐现象”，潮间带潮位不一致对自动提取精度的影响等关键问题，实现大尺度、快速且较为准确地掌握显著用海情况，遏制违法违规用海行为，保护海洋生态环境。

1 区域概况

金山城市沙滩，位于杭州湾北岸，上海金山区南端，全长约2 km，由东北向西南徐徐延伸，是一条集海景、海趣、海乐、海韵为一体的城市景观海岸线，地理范围在30˚42′N—30˚43′N、121˚20′E—121˚22′E 之间。作为长三角最具有海派风格的城市海岸景观，围水面积达1.5 km2，该区域集聚填海、围海、构筑物、房屋建筑、码头、潮间带滩涂等海岸带典型解译标志，围填海工程边界清晰、轮廓成熟。近两年新增简易码头、构筑物修缮时有发生，基本覆盖了准业务化工作中目视解译提取的面状要素类型。

2 研究方法

2.1 数据选择

遥感影像变化检测首先需要选择数据源。目前变化检测所利用的数据源有多种，包括低空无人机航拍影像、SAR 图像、遥感卫星影像（如高分辨率彩色图像、DSM（Digital Surface Model）、Landsat、Pleiades 等）、红外影像、GIS 及多影像融合数据等；在数据的时相上，有单时相影像、多时相影像；在光谱的波段上，有单波段影像、多波段影像；在数据的类型上，有影像间变化检测、以影像为主源辅助其他数据的变化检测、非影像数据与影像数据综合考虑的变化检测。在数据选取时，考虑到检测区域的自然特征、光照及光谱信息等因素，最好选择同一传感器的数据，在时空上尽可能使数据减少物候和太阳高度角等因素引起的误差[4]。

本文选取的高分辨率卫星影像具有经济性、稳定性、实时性，影像清晰度高，可提供彩色、纹理等信息，较好地反映结构信息，而且应用范围广；缺点在于数据量大、处理数据难度大且耗时、卫星影像视场角一般在15˚~20˚，范围较小；该类星源广泛应用于各个领域，且通常与其他数据融合。

2.2 数据来源与预处理

本研究采用高分一号（GF-1） 卫星作为数据源，其所搭载的两台PMS 传感器2 m 分辨率全色，8 m 分辨率多光谱。常用的遥感影像处理软件有Envi、Erdas、Arcgis 等。基于Envi 软件平台，影像经辐射校正、几何粗校正，并经波段合成与影像融合、标准差光谱增强、感兴趣区影像裁剪及几何精校正等处理后，作为变化检测的本底数据。

2.3 疑点疑区遥感解译

为保证提取精度，对2019 年12 月10 日和2021 年11 月2 日获取的影像进行了融合处理及粗裁剪，选取10 个像控点进行GPS 控制点测量，完成几何精校正。选取单期影像进行栅格数据矢量化水边线提取岸线及感兴趣区划定，剔除岸线向陆一侧对整体提取效率的影响，同时利用矢量数据裁切出研究区域，形成一组相同区域不同时相的正射影像集[4]。GPS 踏勘选取的控制点如图1 所示。

图1 控制点站位

按照《我国近海海洋综合调查与评价专项海岛海岸带卫星遥感调查技术规程》要求，基于Arcgis软件平台，采用目视解译方法，完成研究区域海域使用疑点疑区遥感判读，再结合现场踏勘情况，比对解译标志特征，最终确定了2 处新增海域使用疑点疑区。实测点位坐标如表1 所示。

表1 经现场验证的遥感解译疑点疑区

2.4 计算方法

2.4.1 AdaptiveThreshold 函数

在影像中灰度是不均匀的，因此不同区域的阈值也是不一样的。基于此原因，需要一种方法根据影像不同区域亮度或灰度的分布，计算其局部阈值。这种方法就是自适应阈值化图像处理，也称为局部阈值法[5-6]。

AdaptiveThreshold（Input Array src，Output Array dst，double maxValue，int adaptive Method，int thresholdType，int blockSize，double C）

其中，inputArray src 为输入图像，outputArray dst 为输出图像，double maxValue 为满足条件的像素非零值（本文取值255），int adaptiveMethod 为自适应阈值算法，分别是adaptive_thresh_means_c 和adaptive_thresh_guassian_c[7-9]。

均值滤波即对邻域内的像素取平均值后减去相同加权值；高斯平均是在符合正态分布的情况下，每个像素点的权重离中心点越近的，其权值占比也就越重，即中心点权值最大，越往外所占越小，比均值滤波的好处是可以突出重点。本研究为追求对图像的全幅二值化，选取均值滤波更为适合[7]。

针对阈值类型参数int thresholdType，选取以高电平模式开启的形式作为最大值，低电平为最小值零的THRESH_BINARY，计算公式如下。

int blockSize 是自适应阈值函数需要计算的领域块的大小，即领域块内的像素点数量（本文取值11）。double C 是偏移调整量，当计算得到阈值后，用阈值再减去这个偏移调整量得到最终阈值，该偏移调整量的取值正负皆可（本文中取值-2）。对于对比度和明暗变化有较大差异的影像，前期使用普通非自适应的Threshold 函数计算二值图，其输出轮廓效果并不理想，输出的二值图质量差异较大，影响了接下来变化检测的效果[7-9]。

2.4.2 FindContours 函数

经自适应阈值化处理后得到的二值图，作为FindContours 函数的输入图像进行处理。数字图像以坐标系方式确定像素点坐标，因此函数将每个像素点以“点向量”方式存储，即vector

Void findContours（InputArray image，OutputArrayOfArrays contours， OutputArray hierarchy， int mode，int method，Point offset=Point（））

其中，InputArray image 为输入图像；OutputArrayOfArrays contours 为找到的轮廓图；OutputArray hierarchy 为层级关系，定义为vector

2.4.3 DrawContours 函数

在屏幕上绘制检测到的轮廓。绘制轮廓算法不涉及调节参数，所以对算法的细节可以忽略，实现效果很好。

Void draw Contours（InputOutputArray image，Input Array Of Arrays contours，int contourIdx，const Scalar& color，int thickness=1，int lineType=8，InputArray hierarchy =noArray（）， int maxLevel =INT_MAX，Point offset=Point（））

其中，InputOutputArray image 表示目标图像，Input Array Of Arrays contours 表示输入的轮廓组，每一组轮廓由点vector 构成；int contourIdx 指明绘制轮廓的数量，本文设定该参数为负值（-1），画全部轮廓，const Scalar&color 为轮廓的颜色，三通道色组，int thickness为轮廓的线宽，如果为负值或CV_FILLED表示填充轮廓内部，int lineType 为线条类型，InputArray hierarchy 和int maxLevel 共同作用可以输出轮廓结构信息。Point offset 为偏移参数[7-9]。

3 结果与分析

3.1 变化检测结果

利用OpenCV 普通二值化函数（Threshold）对影像进行基本二值化的阈值方法得到的轮廓效果并不理想，后采用自适应阈值化函数（Adaptive Threshold）对每一个像素点单独计算阈值，根据像素的邻域块的像素值分布来确定该像素位置上的二值化阈值，对于不关注的提取要素，运用霍夫变换、二次线段剔除等特征提取方法进行噪点消除[7]。金山城市沙滩轮廓提取效果如图2 所示。

图2 金山城市沙滩自适应阈值图（2019 年12 月10 日）

在阈值图基础上，使用FindContours 函数，对基期和末期分别进行轮廓线提取，使用DrawContours函数对轮廓上色，并比较前后像素点色度差异，在末期图的轮廓线中，标记变化区域，如图3 所示。为验证算法对疑点疑区目标检测的性能，需结合实地验证后的遥感目视解译结果进行一次验证。

图3 金山城市沙滩自动变化检测结果（2021 年11 月2 日）

3.2 结果分析

通过上述自动提取结果可以看到，首先，利用OpenCV 自适应阈值化函数得到的轮廓结果能够将围填海、构筑物、房屋建筑、码头、潮间带滩涂等海岸带典型解译标志提取出来，因实际业务运行获取到的影像分辨率较高，虽然提取精度得到了一定保障，但需滤除船舶、海面垃圾、管道、浮阀等不关注的信息，以上要素在卫星影像中基本以简单的点状呈现，运用“霍夫变换”、二次线段剔除等特征提取方法进行尽可能的噪点消除后，非关注信息得到了有效净化，岸线向海一侧的主要地物得到保留。

其次，实际工作中会遇到两期卫星影像潮位不一致的情况，受限于影像质量差异，自动判读出的潮位变化区域基本为大面积不规则形状，边缘侦测效果并不理想，故自动提取出来的连续横向不规则条状变化信息结合人工判读被舍弃。

最后，提取出的轮廓边界较为准确，变化检测结果表现良好，但仅试验性开展方法研究并不能证明其整体有效性，而实际工作中监测范围更大，影像覆盖面更广，如何有效且快速的提取出变化区域仍需开展更多的案例研究，同时进行大尺度的精度验证分析。

4 结论及建议

本文提出了一种适应实际业务化工作需求的海域使用疑点疑区自动提取与变化检测方法，并以上海金山城市沙滩为研究区，选用研究区2019 年、2021 年两景高分一号（GF-1）卫星遥感影像作为数据源，进行了相关案例实验与精度评价分析。研究结果表明，所提融合遥感和OpenCV 方法自动提取后的图像空间结构边界清晰、信息增强效果好，且变化检测结果较为准确，该方法较适合于监测区域地物类型较少、时效性要求高，但调查面积较大的信息提取。而且不同的阈值函数对同一影像所得提取效果不同，利用OpenCV 普通二值化函数（Threshold）对影像进行基本二值化的阈值方法得到的轮廓效果并不理想，采用自适应阈值化函数（AdaptiveThreshold） 对每一个像素点单独计算阈值，根据像素的邻域块的像素值分布来确定该像素位置上的二值化阈值，得到的轮廓效果较好。

对研究区的海域使用疑点疑区动态分析结果表明，由融合遥感和OpenCV 方法所得出的疑点疑区变化信息，其精度达到了业务化开展海域使用疑点疑区变化检测的要求，且形成了一套基于自动变化检测方法的海域使用疑点疑区遥感解译业务化试运行基本操作流程，如图4 所示。

图4 基于人机交互的海域使用疑点疑区业务化试运行流程

OpenCV 是一种先进的计算机视觉库，用于数字图像及视频处理操作，能够较为轻松地实现算法目标和应用。OpenCV 函数库核心代码由C/C++代码编写而成，具有轻量化且高效的特点，同时提供Python、Ruby、MATLAB 等语言的接口程序，实现了图像处理和计算机视觉方面的很多通用算法，有很强的扩展性、通用性和可移植性。

针对海域海岛此类大尺度的目标变化检测，利用OpenCV 可高效完成应用目标。相比传统人工识别，OpenCV 对基于精校正的基末期遥感影像，进行轮廓提取并准确定位坐标，比较前后坐标差异，判定变化区域，可减少人工识别带来的误差。对于变化识别过程中，类似滩涂这样的干扰区域，使用OpenCV 中的滤波算法和霍夫变化算法进行特定区域的轮廓消除，同时配合光谱反射率数据，能够准确判定区域地物特征类型，避免人工误判。

相较于Envi 或ArcGIS 这样的专业遥感处理软件，OpenCV 作为计算机视觉库，在处理多光谱影像数据时，没有针对性的函数快速提取数据，需自行研究并编写对应的代码进行提取；基于不同分辨率影像底图比较中，使用OpenCV 进行分辨率调整后，基末期影像轮廓比较存在一定的误差。未来在差异变化区域的精度方面，需要进一步探究。