劉建歌，慕德俊

(西北工業(yè)大學(xué) 自動(dòng)化學(xué)院，陜西 西安 710072)

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

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

1 實(shí)驗(yàn)數(shù)據(jù)

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

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

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

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

2 曲波變換和無邊界主動(dòng)輪廓

2.1 曲波變換

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

(1)

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

2.2 無邊界主動(dòng)輪廓模型

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

(2)

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

3 冰水邊緣線提取算法

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

圖1 算法流程圖

3.1 SAR影像預(yù)處理

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

3.2 邊緣信息圖

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

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

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

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

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

(4)

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

(5)

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

4 實(shí)驗(yàn)及結(jié)果

4.1 實(shí)驗(yàn)步驟

實(shí)驗(yàn)步驟如下：

步驟1 對(duì)SAR HH和HV進(jìn)行預(yù)處理．

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

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

4.2 實(shí)驗(yàn)結(jié)果

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

圖3 20110216影像及結(jié)果

圖4 20110220影像及結(jié)果

圖5 20110405影像及結(jié)果

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

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

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

由圖3～圖5可以看出，在低海冰密集度區(qū)域，利用被動(dòng)微波海冰密集度數(shù)據(jù)大面積低估了冰水邊緣線的位置，然而文中提出的算法能夠有效地提取冰水邊緣線的位置．由于結(jié)果為不規(guī)則曲線，為驗(yàn)證算法的準(zhǔn)確性，文中通過

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

(6)

分別計(jì)算被動(dòng)微波海冰密集度數(shù)據(jù)得到的冰水邊緣線以及文中算法得到的冰水邊緣線與標(biāo)注數(shù)據(jù)冰水邊緣線的平均距離，其中，A(Rgt-r)表示結(jié)果與標(biāo)注數(shù)據(jù)之間的面積，M(dgt)表示標(biāo)注數(shù)據(jù)之間的最大距離，結(jié)果如表1 所示．可以看出，文中所提出的算法與采用海冰密集度數(shù)據(jù)得到的冰水邊緣線相比，準(zhǔn)確率得到了明顯的提高，能夠?qū)⒄`差縮減到 16 km 以內(nèi)．

表1 海冰密集度數(shù)據(jù)及文中提出的冰水邊緣線與人工標(biāo)注冰水邊緣線的平均距離Dkm

5 結(jié) 束 語

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