王法景，李俊锋

(1.杨凌职业技术学院，陕西 杨凌712100；2.陕西测绘地理信息局，陕西 西安710054)

自2007年开始，我国黄海海域每年都发生浒苔为主的绿潮(green tide)灾害。绿潮虽无毒，过量增殖却会阻塞海上航道等，造成巨大的经济损失，尤其是2008年夏天青岛黄海海域暴发的绿潮，覆盖面积达到600 km2，严重影响了第29届奥帆赛的顺利进行[1]。已有专家学者把绿潮、水母和赤潮称之为“海洋三大灾害”。由于绿潮具有漂移速度快、覆盖面积广、生长周期短等特点，卫星遥感观测具有探测波长广、覆盖范围宽和时间分辨率高等优势，已成为绿潮宏观监测的重要手段[2-5]。2006年Gower等[6]首次使用中分辨率成像光谱仪(medium resolution imaging spectrometer，MERIS)数据，采用最大叶绿素指数(maxinum chlorophyll index，MCI)算法对墨西哥湾马尾藻进行了监测。随后，中分辨率成像光谱仪(moderate-resolution imaging spectroradiometer，MODIS)、环境一号(huanjing-1，HJ-1)、陆地资源卫星、SAR和北京一号等卫星数据也逐步应用到绿潮的监测中。Son等[5]利用韩国海洋水色卫星(geostationary ocean color imager，GOCI)与现场的遥感数据分析绿潮的光谱特征，并提出了基于GOCI波段设置的“GOCI绿潮漂浮指数”(index of floating green algae for GOCI，IGAG)，结果表明GOCI绿潮漂浮算法更稳定；Shanmugam等[7]提出了海面藻华提取指数(ocean surface algal bloom technique，OSABT)算法，并将该算法应用到阿拉伯海、墨西哥湾的绿潮灾害暴发区域；陈莹等[8]基于GOCI遥感数据设计了一种绿潮提取算法，并应用此算法和2017年的GOCI数据对黄海的绿潮发生区域进行了时空变化分析；施英妮等[9]基于HJ-1A/1B数据利用归一化植被指数(normalized difference vegetation index，NDVI)和绿潮的光谱特性对东海、黄海的绿潮进行了监测提取；徐福祥等[10]基于无人机遥感影像数据对烟台海域的绿潮进行了监测，以验证不同指数对绿潮信息监测的适宜性；郭宇龙等[11]基于线性混合分解，利用野外的现场光谱数据为端元，构建了一种新的叶绿素a光谱指数；张海龙等[12]以近黄海绿潮漂移暴发区域作为研究区，设计了一种基于HJ-CCD和高分一号卫星数据的多光谱绿潮指数(multispectralgreentideindex，MGTI)-多波段差值耦合算法。绿潮灾害业务化监测所利用的数据源仍以光学影像为主，但光学影像受云雾影响严重，无法满足绿潮日变化监测的精度要求。Carkeer[13-18]等也从不同方面对绿潮进行了研究。

1 材料及方法

1.1 RadarSat-2数据

RadarSat-2是于2007年12月14日由加拿大太空署发射的一颗搭载C波段雷达传感器的高分辨率商用卫星，卫星寿命7～12年。全极化SAR影像采用RadarSat-2卫星全极化精细模式SAR影像，相比于RadarSat-1，RadarSat-2更加灵活，重访周期更短，立体成像能力增加，具有高空间分辨率等优势，主要应用于绿潮、溢油等监测、防灾减灾、冰川监测等领域(图1)。

图1 RadarSat-2数据影像

1.2 研究区域

根据从光学影像上目视解译出的中国黄海海域绿潮的分布范围，确定研究区域为33 °N～37 °N、119 °E～123 °E(图2方框区域)。

图2 绿潮提取研究区域

1.3 研究方法

光学遥感受云雾影响严重，难以满足绿潮业务化监测需求；而传统的SAR绿潮提取算法以目视解译结合阈值法对绿潮进行分类提取，但是阈值法有其难以克服的缺点：人为因素影响重、耗时费力，无法对影像进行批处理等，而且由于绿潮分布的不规则性，用统一的一个阈值对绿潮进行分类提取会造成以绿潮聚集区域提取的绿潮覆盖面积较真实值偏大、以零散分布区域提取的绿潮会比真实覆盖范围偏小。因此，发展新的绿潮提取算法是SAR绿潮监测必须克服的瓶颈。

1.3.1 DRLSE模型 DRLSE(Dr. Li Chunming’s improved algorithm for local region segmentation)模型[12]是一种基于边缘信息的变分水平集分割模型，主要为克服水平集需不断重新初始化为符号距离函数的缺点，在水平集模型中添加了一个距离保持的正则化项，该项的引入提高了分割的演化效率和速度，能量函数公式如下：

ε(φ)=uRP(φ)+εext(φ)

(1)

式(1)中：RP(φ)是水平集函数φ的内部能量泛函(Ω为图像区域)，用来定量地纠正符号距离函数与水平集函数之间的偏差：

(2)

式(2)中：p(s)为能量密度函数。

εext(φ)是外部能量函数，依赖图像特征，用来驱使零水平集向目标边界演化，并在零水平集到达边界时获得最小值。其定义为公式(3)：

(3)

式(3)中：u、λ、α分别为水平集正则化项、长度项、边缘停止项的系数；φ为水平集轮廓线，δ(φ)为一维Dirac函数；g为边缘指示函数；H(-φ)为一维Heaviside函数。

I(x，y)为图像的灰度函数，则g(I)为公式(4)：

(4)

公式(4)中：Gσ表示标准方差为σ的二维高斯滤波器，作用是去除噪声；表示空间梯度算子。当边缘指示函数坐标变量(x，y)位于图像中目标的边缘时，图像梯度的值|GσI(x，y)|相对很大，因而g(x，y)的值会趋近0；相反，当远离图像中目标边缘时，图像梯度的值|GσI(x，y)|相对很小，因而g(x，y)的值会趋近1。由此，可以根据函数值g(x，y)是接近0还是接近1来判断对应像素I(x，y)是否位于图像的目标边界。

利用偏微分方程求解曲线演化过程，通过极小化能量泛函δε(φ)可得控制演化的偏微分方程(5)：

(5)

DRLSE模型法用来分割边界清晰图像，但由于SAR图像斑点噪声严重等特点，当其运用在SAR图像分割时存在如下缺点：一旦目标与水平集初始轮廓曲线交叉，DRLSE模型法则无法正确分割，而SAR图像弱边缘时存在边缘不封闭问题。

1.3.2 阈值法和水平集模型相结合的分割算法 阈值法是一种常用的有效图像分割算法。灰度图像的阈值分割首先确定一个处于图像灰度范围中的阈值，然后将影像中各像素的灰度值与这阈值相比较，根据比较结果将相应的像素分为两类，从而达到图像分割的目的。最优阈值的确定是阈值分割的最关键步骤，现有的大部分算法都集中在如何确定阈值的研究上。本文采用的方法首先使用阈值分割的方法确定阈值，然后将这个阈值代入水平集DRLSE模型的初始轮廓函数，由此确定水平集DRLSE模型的初始轮廓，最终准确分割出绿潮灾害区域。

2 结果与分析

依据改进的水平集绿潮提取算法，对RadarSat-2有绿潮覆盖的区域进行分块，使分块区域既有聚集区域也有零散分布区域，分别用本文发展的绿潮提取算法进行绿潮提取，以验证算法的有效性(图3、图4)。

a为原始影像；b为水平集初始轮廓线。图3 水平集算法

a为绿潮零散分布原始影像；b为对应的提取影像；c为绿潮聚集分布原始影像；d为对应的提取影像。图4 水平集算法绿潮提取

2.1 偏差验证

为了更好地验证水平集算法提取的绿潮精度，本文采用传统的单波段阈值法作为绿潮提取，“真值”的参考，对本研究提出的改进水平集算法进行精度验证，将本研究提出的方法提取的绿潮覆盖面积与单波段阈值法提取的绿潮覆盖面积之间的偏差值E作为该算法精度的评价指标：

(6)

公式(6)中Sestimate为本文算法提取的绿潮面积，Strue为单波段阈值法提取的绿潮覆盖面积，为参考真值。分别利用单波段阈值法和改进水平集算法对绿潮进行提取，图3展示了RadarSat-2灰度图像基于水平集算法的绿潮提取图。从中我们可以看出，本研究应用的算法提取的绿潮分布轮廓与目视解译高度一致，与目视解译绿潮的分布范围基本对应。

为防止试验结果的偶然性，本文采用不同的研究区域进行对比，对比结果如表1和图5。

图5 用两种算法提高取的面积关系

表1 绿潮覆盖面积水平集算法与单波段阈值法偏差对比

本研究所用水平集算法提取的绿潮平均覆盖面积为46.0 km2，单波段阈值结合目视解译法提取的平均覆盖面积为45.1 km2。根据偏差计算公式，E=1.12%，与单波段阈值法提取的绿潮覆盖面积基本一致。传统阈值法受人为因素影响严重，由于人为主观因素容易导致误差，且计算速率较慢，无法满足应急业务化监测的应急要求；同一景影像若既有绿潮聚集区又有分散区则很难将绿潮很好地分离开来，造成的误差较大。本文发展的算法较好地解决了上述问题，无论在绿潮零星分散区域还是聚集区域都能很好地将绿潮分离开来，这也是传统单波段阈值法同一阈值提取无法做到的。

2.2 基于误差矩阵的精度评价

随机生成750个检查点，通过目视解译结合专家先验知识，在这750个检查点中有152个是绿潮像元，598个是水体像元；本文所采用的改进水平集算法对绿潮的分类后影像中，有154个绿潮像元，596个水体像元。根据各项精度指标的计算公式可以得出如下结果(表2～表3)：

表3 基于改进水平集算法的绿潮分类精度评价(Kappa系数)

表2 基于改进水平集算法的绿潮分类精度评价(总体精度)

(1)绿潮信息的制图精度为100%，表明在实际为绿潮的范围内所有的绿潮被划分成正确的一类；海水的制图精度为99.7%，表明99.7%海水被正确划分到海水类中。

(2)绿潮的用户精度为98.7%，表明在分类影像中绿潮像元98.7%为实际绿潮的覆盖像元；海水用户精度为100%，表明在分类后的影像中水体的像元中所有的海水像元对应实际的海水覆盖像元。

(3)Kappa系数为0.992，结果表明原始影像与多项式逻辑回归拟合算法分类提取的影像吻合度为99.2%，分类结果极佳。

通过试验我们可以得出：发展基于改进水平集的绿潮覆盖面积提取算法与传统单波段阈值法有很好的一致性，契合度达到99.73%，而且无论是总体精度，还是制图精度、用户精度、Kappa系数等都能很好满足绿潮业务化监测精度的需要，是一个理想的绿潮覆盖面积提取算法；用本文发展的基于改进水平集算法提取绿潮覆盖面积的精度高并且稳定性好，相比于传统的绿潮提取算法，人为因素干预更少，在时间效率上有较大的优势，可较好地满足绿潮应急监测的需要。

2.3 时间效率比较

应急业务化绿潮监测算法优劣的另一个评价指标为算法提取绿潮所需要的时间，因此本文将发展的改进水平集算法跟传统阈值法提取同一区域绿潮所用的时间进行比较，以验证算法时间上的优势(表4)。

表4 时间效率对比

通过对比我们可以发现：对于相同的6个研究区域，基于传统的波段阈值法提取绿潮耗时最长的是区域5，所用时间达到6.8 s；耗时时间最短的是区域6，所用时间达到4.7 s。而本文发展的改进水平集绿潮提取算法最长耗时(区域5)仅为2.12 s，最短耗时(区域2)1.32 s，跟传统阈值法相比，具有更好的时间优势，更适于绿潮业务化监测。

3 讨论

基于雷达影像的绿潮信息提取，传统的方法仍以波段阈值法为主。基于阈值法的绿潮提取存在诸多方面的问题，首先，阈值法需要耗费大量的人力，在绿潮暴发期间，科研人员需每景影像都要人为地进行提取，效率偏低，同时受科研人员的人为因素影响较大。影像在提取前需要进行预处理，由于科研素质的差异，不同科研人员对影像的提取可能存在差异性，导致绿潮的业务化监测存在不稳定性。针对SAR绿潮提取算法相对单一、提取结果受人为因素影响大、不稳定、效率低等的问题，本文基于图像分割的改进水平集方法进行了绿潮信息提取研究。

通过目视判读可以发现，本文改进的水平集算法的初始轮廓线能够较好地将绿潮与海水进行分离并且能够很好地将绿潮信息进行提取，无论是绿潮零星分布的区域还是绿潮聚集区域，该算法均能很好地将绿潮像元提取出来。

为了验证该算法的提取精度，本文通过利用发展的算法分不同区域进行绿潮提取，从与波段阈值法提取的偏差、精度评价、时间效率等方面与传统方法进行比较分析。实验结果表明，主流的波段阈值法提取的绿潮覆盖面积与本文改进的水平集算法具有很好的相关性，偏差仅为1.12%。但是阈值法在零星分布、聚集区域需要人为设置不同的阈值，导致其提取效率低下，因此本文发展的算法在绿潮业务化监测过程中具有阈值法无法比拟的效率优势。通过总体精度、Kappa系数等进行比较，本文发展的基于改进水平集的绿潮覆盖面积提取算法与传统单波段阈值法有很好的一致性，契合度达到99.73%，而且无论是总体精度，还是制图精度、用户精度、Kappa系数等都能很好地满足绿潮业务化监测精度的需要，是一个理想的绿潮覆盖面积提取算法；最后通过提取时间分区域进行比较，可以发现该算法比传统阈值法对同一区域进行提取所用的时间更少，具有很好的效率优势。用本文发展的基于改进水平集算法提取绿潮覆盖面积的精度高并且稳定性好，为绿潮的监测提取提供了新的方法。