一種SAR影像冰水邊緣線提取算法

劉建歌，慕德俊

(西北工業大學 自動化學院，陜西 西安 710072)

星載被動微波傳感器以及合成孔徑雷達(Synthetic Aperture Radar，SAR)不受天氣和光照的影響，能夠全天候大范圍獲取海洋信息，成為分析、跟蹤和監測海洋冰的常用手段．冰水邊緣線[1]作為海冰區域和海水區域的分界線，與全球氣候及生態系統存在著密不可分的關系[2]，其位置的確定能夠為研究海冰提供重要的參數，同時能夠降低人類活動的風險．遙感方法確定冰水邊緣線最常用的手段是將被動微波獲取的海冰密集度值為15%的界線定義為冰水邊緣線[2]．全球氣候變暖導致從固結冰到開闊水域過渡區域的海冰邊緣區(Marginal Ice Zone，MIZ)[3]出現的時間越來越長，覆蓋范圍越來越廣．在這種條件下，冰水邊緣線所在區域的海冰密集度較低，采用被動微波數據獲取的海冰密集度值往往會低估其真實值，導致以海冰密集度值為15%定義的冰水邊緣線位置的不準確[2，4]．隨著SAR能夠提供多極化高分辨率的雷達影像，在低海冰密集度區域依然能夠捕捉海冰特征，提取冰水邊緣線．另外，海冰邊緣區的動態性造成其位置變化快，人工解譯不能滿足實時性的需求，因此基于SAR影像的冰水邊緣線自動提取算法成為需求[5]．SAR影像自動提取冰水邊緣線的難點主要有兩點:SAR自身系統問題造成的影像邊緣提取及分割的難點[6]，包括相干斑噪聲和對于入射角非常敏感等[7]，造成海冰種類內和種類間的后向散射值強度變化大；由于海冰邊緣區海冰種類多，以及海冰表面粗糙度、鹽度、濕度等不同，且動態性強，會導致SAR后向散射值的強度變化大，海冰形態各異，難以提取合適的紋理特征．

SAR影像中的海冰邊緣區包含大量的曲線特征，但是曲線特征強度和對比度變化大，直接利用曲線特征會造成對比度小的區域被忽略，導致真實的冰水邊緣線被低估．采用邊緣檢測的方法，將對比度不同的曲線特征統一為邊緣信息則可解決該問題．根據上述分析，文中采用RADARSAT-2(RS-2)HH(horizontal polarization of transmitted and received)和HV(horizontal polarization of transmitted and vertical received)雙極化影像，利用Canny算子對影像進行邊緣提取，得到海冰邊緣區的邊緣信息圖，根據邊緣信息豐富程度的差異提取冰水邊緣線．由于得到的邊緣信息為長短不一、方向不同的離散曲線，且冰水邊緣線通常為不規則曲線，因此為達到平滑影像的同時保留邊緣位置及形態，文中采用曲波域多尺度主動輪廓模型，將邊緣信息圖轉換到曲波域，經過從粗糙尺度到精細尺度的迭代得到冰水邊緣線的位置．

1 實驗數據

實驗采用的影像位于加拿大拉普拉多海岸，所在區域為北緯65°至55°，西經64°至56°，獲取時間分別為2011年2月16日(20110216)、2月20日(20110220)和4月5日(20110405)．通過分析比較SAR影像和海冰分析圖表，可以看出，影像中的冰水邊緣線存在于低海冰密集度的海冰邊緣區，并且使用被動微波海冰密集度數據在影像中明顯低估了冰水邊緣線的位置．

實驗使用的SAR影像是由加拿大冰面監控部(Canadian Ice Service，CIS)提供的加拿大RS-2衛星寬幅掃描模式C波段的雙極化HH和HV影像．每幅SAR影像約為 10 000× 10 000 個像素，空間分辨率約為 50 m．

被動微波獲得的海冰密集度數據采用高級微波掃描輻射計-地球觀測系統(Advanced Microwave Scanning Radiometer-Earth observing system，AMSR-E)設備得到，根據ARTIST海冰算法實現的，空間分辨率約為 6.25 km．

海冰分析圖表作為標注數據，用來比較和驗證結果的準確性[8]．由CIS海冰分析專家根據RS-2 SAR影像人工估計的海冰密集度值，通過與SAR影像比較，將插值后的人工分析圖表海冰密集度為零的邊界線作為冰水邊緣線的標注數據．由于人工海冰分析圖表的分辨率較低，作為冰水邊緣線的標注數據存在一定的誤差，該誤差通過SAR影像的分辨率與海冰分析圖表分辨率的比值粗略估計δ為 3.5 km．

2 曲波變換和無邊界主動輪廓

2.1 曲波變換

文獻[9]提出的曲波變換是能夠將圖像分解為多尺度多方向信息的多尺度幾何分析方法之一．曲波變換具有各向異性，能夠很好地表達圖像的邊緣、曲線等高維奇異性，克服了小波函數只能表達點奇異性的局限．2004年提出的第2代曲波變換及兩種快速離散算法，包括不等間距快速傅里葉變換和基于卷繞的變換，促進了曲波變換在圖像處理各個方面的應用[10]．基于卷繞的方法計算速度快，冗余小，文中的曲波變換采用基于卷繞的方法，假設大小為n×n的圖像I(x，y)，曲波變換母函數為φj,l,K，則曲波系數矩陣Cj,l,K為

(1)

其中，j為尺度參數(j≤ns,ns表示可分解的尺度個數)，l和K分別為方向參數和位移參數．

2.2 無邊界主動輪廓模型

主動輪廓模型廣泛地應用于圖像分割．文中采用基于區域的Chan-Vese無邊界主動輪廓模型(C-V模型)[11]對圖像進行分割，模型不依靠圖中的梯度信息而是滿足全局能量最優化，通過最小化能量函數E驅動輪廓的移動．因此，當邊緣為模糊或者離散的狀態時，仍然可以有效地進行分割．假設模型的外力場為F(x，y)，通過邊界R將圖像劃分為區域A和B，則能量函數為

(2)

其中，UA和UB分別表示兩個區域的均值；L(R)和S(R)分別為閉合曲線R的長度及所圍區域的面積；u>0，v>0，λ1>0，λ2>0分別為各能量項的權值．

3 冰水邊緣線提取算法

SAR影像中冰水邊緣線為具有大量不同尺度、不同方向曲線特征的海冰邊緣區和相對平滑的開闊水域的邊界線，然而海冰邊緣區曲線特征的強度及對比度變化較大，直接利用曲線特征無法得到準確的結果．文中首先通過Canny邊緣檢測的方法將曲線特征統一為邊緣信息，再利用海冰邊緣區與開闊水域邊緣信息豐富程度的差異，采用C-V模型提取冰水邊緣線．如圖1所示，算法主要分為3個步驟:SAR影像預處理;對影像進行曲波域中尺度濾波，采用Canny算子對濾波后的影像進行邊緣檢測，得到邊緣信息圖;基于離散二值邊緣信息圖，采用曲波域多尺度C-V模型提取邊緣信息豐富區域的輪廓，得到冰水邊緣線．

圖1 算法流程圖

3.1 SAR影像預處理

將RS-2 HH和HV數據分別進行預處理，包括降采樣和采用鏡像的方式來去除影像中的陸地部分．降采樣采用影像4×4大小非重疊塊的平均值，達到降噪和減少數據計算量的目的，降采樣后的大小為原圖的 1/4，空間分辨率為 200 m．根據海陸邊界的標識，將SAR影像中陸地區域對應的像素，通過水平鏡像的方式，改為對應的海洋區域的像素值．通常情況下，海洋像素值和陸地像素值具有較大的差異，因此海陸邊界存在對比度較大的曲線特征，在曲波域表現為較大的曲波系數，對邊緣檢測產生干擾．通過鏡像的方式能夠有效去除海陸邊界的曲線特征．

3.2 邊緣信息圖

其中，jl和jh表示保留的中尺度范圍．該方法利用曲波系數的內部聯系，能夠達到銳化中曲線特征同時去除噪音的效果．

由于有些海冰特征只能在HH或HV影像中被捕捉到，因此為了能夠利用HH和HV之間信息的互補性，將HH和HV邊緣檢測結果采用“或”運算進行合并，得到完整的邊緣信息圖B=BHH‖BHV．

3.3 曲波域多尺度C-V模型

由3.2節得到的邊緣信息為二值離散圖像，采用該模型的原因有:海冰邊緣區邊緣信息的密集程度不同，因而分析一定尺度下的邊緣信息密度的方法很難找到合適的閾值進行準確的分割結果．采用多尺度分析的方法，利用從粗糙尺度到精細尺度的迭代在保證全局一致性的基礎上，逐步細化邊緣區域；冰水邊緣線通常為不規則的曲線，曲波變換不僅可以分解為多尺度信息和保留冰水邊緣線的幾何特征，而且其位置參數能夠準確定位曲線位置；C-V模型能夠滿足全局能量的最優化，而不依靠影像中的梯度信息．因此，當邊緣為模糊或者離散的狀態時，可以有效地進行分割．因此，文中采用曲波域多尺度C-V模型對得到的邊緣信息圖進行分割．

圖2 曲波域多尺度C-V模型

(4)

將各個方向上系數矩陣的平均值作為該尺度下的外力場Fj(x,y)，L表示尺度j下的方向個數．將Fj(x,y)代入式(2)，得到

(5)

通過最小化Ej搜索該尺度下目標輪廓Rj，將該結果作為下一尺度j+1 的初始值．其中，最粗糙尺度j=1 的初始值R0采用被動微波海冰密集度數據為15%的邊界線．由于最精細尺度通常對應影像的孤立點以及噪聲，曲波域的C-V模型迭代到尺度ns-1．

4 實驗及結果

4.1 實驗步驟

實驗步驟如下：

步驟1 對SAR HH和HV進行預處理．

步驟2 將預處理后的HH和HV影像分別進行曲波變換，觀察海冰邊緣區中的海冰物理特征約小于 3.2 km，表現在預處理后的SAR影像中約為小于24個像素．因此中尺度濾波范圍選擇jl=ns-3，jh=ns-1．利用式(3)得到新的曲波系數，將逆變換后的重構影像分別采用Canny算子進行邊緣檢測得到BHH及BHV，則邊緣信息圖B=BHH‖BHV．

步驟3 對B進行曲波變換，按照圖2所示的流程以及式(4)～式(5)，其中初始值R0為采用被動微波海冰密集度為15%的輪廓線，在曲波域從尺度j=1 到j=ns-1，將Rj作為尺度j+1 的初始值，采用C-V模型迭代得到冰水邊緣線位置．C-V模型參數設置: 通過實驗可知，每個尺度下的C-V模型在迭代 1 000 次后趨于穩定，因此這里將迭代次數設置為 1 000，反映輪廓長度的權值u= 0.01N，N表示影像像素點的總個數;v反映輪廓內部區域的面積，通常不考慮面積，因此v=0；另外，λ1=λ2= 1．

4.2 實驗結果

圖3～圖5分別為影像20110216、20110220以及20110405的RS-2 HH、HV影像以及冰水邊緣線的提取結果．

圖3 20110216影像及結果

圖4 20110220影像及結果

圖5 20110405影像及結果

圖3(a)、圖4(a)、圖5(a)和圖3(b)、圖4(b)、圖5(b)分別為降采樣后的HH和HV SAR影像(白線表示陸地掩膜)；圖3(c)、圖4(c)、圖5(c)為邊緣信息圖E，并且陸地部分標記為零；圖3(d)、圖4(d)、圖5(d)為冰水邊緣線提取結果(白星線為文中提出的方法，黑點線為海冰密集度數據為15%的邊界線，黑星線為標記數據)．

圖4(d)中黑框部分為存在稀疏海冰的區域，在這種情況下，采用C-V模型由于需要考慮區域的一致性，因此該區域往往被分割為海水，導致結果不正確．

圖5(d)中黑框部分表示由于對比度明顯小于其他海冰區域的情況下，不能準確得到邊緣信息，使得該區域被誤認為海水，導致結果的不正確．

由圖3～圖5可以看出，在低海冰密集度區域，利用被動微波海冰密集度數據大面積低估了冰水邊緣線的位置，然而文中提出的算法能夠有效地提取冰水邊緣線的位置．由于結果為不規則曲線，為驗證算法的準確性，文中通過

D=A(Rgt-r)/M(dgt)-σ,

(6)

分別計算被動微波海冰密集度數據得到的冰水邊緣線以及文中算法得到的冰水邊緣線與標注數據冰水邊緣線的平均距離，其中，A(Rgt-r)表示結果與標注數據之間的面積，M(dgt)表示標注數據之間的最大距離，結果如表1 所示．可以看出，文中所提出的算法與采用海冰密集度數據得到的冰水邊緣線相比，準確率得到了明顯的提高，能夠將誤差縮減到 16 km 以內．

表1 海冰密集度數據及文中提出的冰水邊緣線與人工標注冰水邊緣線的平均距離Dkm

5 結 束 語

文中基于SAR影像設計了一種冰水邊緣線提取算法，該算法能夠在低海冰密集度區域準確提取冰水邊緣線的位置，解決采用被動微波數據獲取的海冰密集度數據不能準確提取冰水邊緣線位置的問題．通過實驗可知，該算法能夠有效地提取冰水邊緣線的位置，與被動微波海冰密集度數據得到的冰水邊緣線相比，準確率有大幅提升．冰水邊緣線的準確提取能夠為研究、跟蹤和預測海冰邊緣區范圍等方面提供有效參數，對研究、分析海冰及環境變化以及保證人類活動安全等方面起著重要的作用．