王法景，李俊鋒

(1.楊凌職業(yè)技術學院，陜西 楊凌712100；2.陜西測繪地理信息局，陜西 西安710054)

自2007年開始，我國黃海海域每年都發(fā)生滸苔為主的綠潮(green tide)災害。綠潮雖無毒，過量增殖卻會阻塞海上航道等，造成巨大的經濟損失，尤其是2008年夏天青島黃海海域暴發(fā)的綠潮，覆蓋面積達到600 km2，嚴重影響了第29屆奧帆賽的順利進行[1]。已有專家學者把綠潮、水母和赤潮稱之為“海洋三大災害”。由于綠潮具有漂移速度快、覆蓋面積廣、生長周期短等特點，衛(wèi)星遙感觀測具有探測波長廣、覆蓋范圍寬和時間分辨率高等優(yōu)勢，已成為綠潮宏觀監(jiān)測的重要手段[2-5]。2006年Gower等[6]首次使用中分辨率成像光譜儀(medium resolution imaging spectrometer，MERIS)數據，采用最大葉綠素指數(maxinum chlorophyll index，MCI)算法對墨西哥灣馬尾藻進行了監(jiān)測。隨后，中分辨率成像光譜儀(moderate-resolution imaging spectroradiometer，MODIS)、環(huán)境一號(huanjing-1，HJ-1)、陸地資源衛(wèi)星、SAR和北京一號等衛(wèi)星數據也逐步應用到綠潮的監(jiān)測中。Son等[5]利用韓國海洋水色衛(wèi)星(geostationary ocean color imager，GOCI)與現場的遙感數據分析綠潮的光譜特征，并提出了基于GOCI波段設置的“GOCI綠潮漂浮指數”(index of floating green algae for GOCI，IGAG)，結果表明GOCI綠潮漂浮算法更穩(wěn)定；Shanmugam等[7]提出了海面藻華提取指數(ocean surface algal bloom technique，OSABT)算法，并將該算法應用到阿拉伯海、墨西哥灣的綠潮災害暴發(fā)區(qū)域；陳瑩等[8]基于GOCI遙感數據設計了一種綠潮提取算法，并應用此算法和2017年的GOCI數據對黃海的綠潮發(fā)生區(qū)域進行了時空變化分析；施英妮等[9]基于HJ-1A/1B數據利用歸一化植被指數(normalized difference vegetation index，NDVI)和綠潮的光譜特性對東海、黃海的綠潮進行了監(jiān)測提?。恍旄Ｏ榈萚10]基于無人機遙感影像數據對煙臺海域的綠潮進行了監(jiān)測，以驗證不同指數對綠潮信息監(jiān)測的適宜性；郭宇龍等[11]基于線性混合分解，利用野外的現場光譜數據為端元，構建了一種新的葉綠素a光譜指數；張海龍等[12]以近黃海綠潮漂移暴發(fā)區(qū)域作為研究區(qū)，設計了一種基于HJ-CCD和高分一號衛(wèi)星數據的多光譜綠潮指數(multispectralgreentideindex，MGTI)-多波段差值耦合算法。綠潮災害業(yè)務化監(jiān)測所利用的數據源仍以光學影像為主，但光學影像受云霧影響嚴重，無法滿足綠潮日變化監(jiān)測的精度要求。Carkeer[13-18]等也從不同方面對綠潮進行了研究。

1 材料及方法

1.1 RadarSat-2數據

RadarSat-2是于2007年12月14日由加拿大太空署發(fā)射的一顆搭載C波段雷達傳感器的高分辨率商用衛(wèi)星，衛(wèi)星壽命7～12年。全極化SAR影像采用RadarSat-2衛(wèi)星全極化精細模式SAR影像，相比于RadarSat-1，RadarSat-2更加靈活，重訪周期更短，立體成像能力增加，具有高空間分辨率等優(yōu)勢，主要應用于綠潮、溢油等監(jiān)測、防災減災、冰川監(jiān)測等領域(圖1)。

圖1 RadarSat-2數據影像

1.2 研究區(qū)域

根據從光學影像上目視解譯出的中國黃海海域綠潮的分布范圍，確定研究區(qū)域為33 °N～37 °N、119 °E～123 °E(圖2方框區(qū)域)。

圖2 綠潮提取研究區(qū)域

1.3 研究方法

光學遙感受云霧影響嚴重，難以滿足綠潮業(yè)務化監(jiān)測需求；而傳統(tǒng)的SAR綠潮提取算法以目視解譯結合閾值法對綠潮進行分類提取，但是閾值法有其難以克服的缺點：人為因素影響重、耗時費力，無法對影像進行批處理等，而且由于綠潮分布的不規(guī)則性，用統(tǒng)一的一個閾值對綠潮進行分類提取會造成以綠潮聚集區(qū)域提取的綠潮覆蓋面積較真實值偏大、以零散分布區(qū)域提取的綠潮會比真實覆蓋范圍偏小。因此，發(fā)展新的綠潮提取算法是SAR綠潮監(jiān)測必須克服的瓶頸。

1.3.1 DRLSE模型 DRLSE(Dr. Li Chunming’s improved algorithm for local region segmentation)模型[12]是一種基于邊緣信息的變分水平集分割模型，主要為克服水平集需不斷重新初始化為符號距離函數的缺點，在水平集模型中添加了一個距離保持的正則化項，該項的引入提高了分割的演化效率和速度，能量函數公式如下：

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

(1)

式(1)中：RP(φ)是水平集函數φ的內部能量泛函(Ω為圖像區(qū)域)，用來定量地糾正符號距離函數與水平集函數之間的偏差：

(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 閾值法和水平集模型相結合的分割算法 閾值法是一種常用的有效圖像分割算法。灰度圖像的閾值分割首先確定一個處于圖像灰度范圍中的閾值，然后將影像中各像素的灰度值與這閾值相比較，根據比較結果將相應的像素分為兩類，從而達到圖像分割的目的。最優(yōu)閾值的確定是閾值分割的最關鍵步驟，現有的大部分算法都集中在如何確定閾值的研究上。本文采用的方法首先使用閾值分割的方法確定閾值，然后將這個閾值代入水平集DRLSE模型的初始輪廓函數，由此確定水平集DRLSE模型的初始輪廓，最終準確分割出綠潮災害區(qū)域。

2 結果與分析

依據改進的水平集綠潮提取算法，對RadarSat-2有綠潮覆蓋的區(qū)域進行分塊，使分塊區(qū)域既有聚集區(qū)域也有零散分布區(qū)域，分別用本文發(fā)展的綠潮提取算法進行綠潮提取，以驗證算法的有效性(圖3、圖4)。

a為原始影像；b為水平集初始輪廓線。圖3 水平集算法

a為綠潮零散分布原始影像；b為對應的提取影像；c為綠潮聚集分布原始影像；d為對應的提取影像。圖4 水平集算法綠潮提取

2.1 偏差驗證

為了更好地驗證水平集算法提取的綠潮精度，本文采用傳統(tǒng)的單波段閾值法作為綠潮提取，“真值”的參考，對本研究提出的改進水平集算法進行精度驗證，將本研究提出的方法提取的綠潮覆蓋面積與單波段閾值法提取的綠潮覆蓋面積之間的偏差值E作為該算法精度的評價指標：

(6)

公式(6)中Sestimate為本文算法提取的綠潮面積，Strue為單波段閾值法提取的綠潮覆蓋面積，為參考真值。分別利用單波段閾值法和改進水平集算法對綠潮進行提取，圖3展示了RadarSat-2灰度圖像基于水平集算法的綠潮提取圖。從中我們可以看出，本研究應用的算法提取的綠潮分布輪廓與目視解譯高度一致，與目視解譯綠潮的分布范圍基本對應。

為防止試驗結果的偶然性，本文采用不同的研究區(qū)域進行對比，對比結果如表1和圖5。

圖5 用兩種算法提高取的面積關系

表1 綠潮覆蓋面積水平集算法與單波段閾值法偏差對比

本研究所用水平集算法提取的綠潮平均覆蓋面積為46.0 km2，單波段閾值結合目視解譯法提取的平均覆蓋面積為45.1 km2。根據偏差計算公式，E=1.12%，與單波段閾值法提取的綠潮覆蓋面積基本一致。傳統(tǒng)閾值法受人為因素影響嚴重，由于人為主觀因素容易導致誤差，且計算速率較慢，無法滿足應急業(yè)務化監(jiān)測的應急要求；同一景影像若既有綠潮聚集區(qū)又有分散區(qū)則很難將綠潮很好地分離開來，造成的誤差較大。本文發(fā)展的算法較好地解決了上述問題，無論在綠潮零星分散區(qū)域還是聚集區(qū)域都能很好地將綠潮分離開來，這也是傳統(tǒng)單波段閾值法同一閾值提取無法做到的。

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%，分類結果極佳。

通過試驗我們可以得出：發(fā)展基于改進水平集的綠潮覆蓋面積提取算法與傳統(tǒng)單波段閾值法有很好的一致性，契合度達到99.73%，而且無論是總體精度，還是制圖精度、用戶精度、Kappa系數等都能很好滿足綠潮業(yè)務化監(jiān)測精度的需要，是一個理想的綠潮覆蓋面積提取算法；用本文發(fā)展的基于改進水平集算法提取綠潮覆蓋面積的精度高并且穩(wěn)定性好，相比于傳統(tǒng)的綠潮提取算法，人為因素干預更少，在時間效率上有較大的優(yōu)勢，可較好地滿足綠潮應急監(jiān)測的需要。

2.3 時間效率比較

應急業(yè)務化綠潮監(jiān)測算法優(yōu)劣的另一個評價指標為算法提取綠潮所需要的時間，因此本文將發(fā)展的改進水平集算法跟傳統(tǒng)閾值法提取同一區(qū)域綠潮所用的時間進行比較，以驗證算法時間上的優(yōu)勢(表4)。

表4 時間效率對比

通過對比我們可以發(fā)現：對于相同的6個研究區(qū)域，基于傳統(tǒng)的波段閾值法提取綠潮耗時最長的是區(qū)域5，所用時間達到6.8 s；耗時時間最短的是區(qū)域6，所用時間達到4.7 s。而本文發(fā)展的改進水平集綠潮提取算法最長耗時(區(qū)域5)僅為2.12 s，最短耗時(區(qū)域2)1.32 s，跟傳統(tǒng)閾值法相比，具有更好的時間優(yōu)勢，更適于綠潮業(yè)務化監(jiān)測。

3 討論

基于雷達影像的綠潮信息提取，傳統(tǒng)的方法仍以波段閾值法為主。基于閾值法的綠潮提取存在諸多方面的問題，首先，閾值法需要耗費大量的人力，在綠潮暴發(fā)期間，科研人員需每景影像都要人為地進行提取，效率偏低，同時受科研人員的人為因素影響較大。影像在提取前需要進行預處理，由于科研素質的差異，不同科研人員對影像的提取可能存在差異性，導致綠潮的業(yè)務化監(jiān)測存在不穩(wěn)定性。針對SAR綠潮提取算法相對單一、提取結果受人為因素影響大、不穩(wěn)定、效率低等的問題，本文基于圖像分割的改進水平集方法進行了綠潮信息提取研究。

通過目視判讀可以發(fā)現，本文改進的水平集算法的初始輪廓線能夠較好地將綠潮與海水進行分離并且能夠很好地將綠潮信息進行提取，無論是綠潮零星分布的區(qū)域還是綠潮聚集區(qū)域，該算法均能很好地將綠潮像元提取出來。

為了驗證該算法的提取精度，本文通過利用發(fā)展的算法分不同區(qū)域進行綠潮提取，從與波段閾值法提取的偏差、精度評價、時間效率等方面與傳統(tǒng)方法進行比較分析。實驗結果表明，主流的波段閾值法提取的綠潮覆蓋面積與本文改進的水平集算法具有很好的相關性，偏差僅為1.12%。但是閾值法在零星分布、聚集區(qū)域需要人為設置不同的閾值，導致其提取效率低下，因此本文發(fā)展的算法在綠潮業(yè)務化監(jiān)測過程中具有閾值法無法比擬的效率優(yōu)勢。通過總體精度、Kappa系數等進行比較，本文發(fā)展的基于改進水平集的綠潮覆蓋面積提取算法與傳統(tǒng)單波段閾值法有很好的一致性，契合度達到99.73%，而且無論是總體精度，還是制圖精度、用戶精度、Kappa系數等都能很好地滿足綠潮業(yè)務化監(jiān)測精度的需要，是一個理想的綠潮覆蓋面積提取算法；最后通過提取時間分區(qū)域進行比較，可以發(fā)現該算法比傳統(tǒng)閾值法對同一區(qū)域進行提取所用的時間更少，具有很好的效率優(yōu)勢。用本文發(fā)展的基于改進水平集算法提取綠潮覆蓋面積的精度高并且穩(wěn)定性好，為綠潮的監(jiān)測提取提供了新的方法。