趙慶平,于力剛,王 耀,張馨月
(1.淮北師范大學 a.物理與電子信息學院,b.實驗室與設備管理處,安徽 淮北 235000;2.塔斯馬尼亞大學 信息技術與系統學院,澳大利亞 塔斯馬尼亞 7005)
隨著全球氣候變暖,對于極地海冰的研究得到了越來越多的關注。海冰監測對氣候變化研究和航行安全保障具有重要意義。合成孔徑雷達(Synthetic Aperture Radar SAR)是海冰監測的有效工具,它是一種能產生高分辨率遙感圖像的主動式微波傳感器,具有全天時、全天候、多波段、穿透力強等特點[1],現已經得到了廣泛的應用。海冰SAR圖像的分割是圖像處理(解譯)的重要環節,加拿大冰署(CIS)的海冰分析員曾對Radarsat-1/2所提供的大量的海冰SAR圖像進行處理,制作冰況圖,以反映海冰的密集度、類型和尺寸等信息,但由于當時缺少有效的業務化分割工具,只能由海冰分析員對圖像進行人工分割,而人工分割費時耗力,且其精確度和分辨率均是有限的[2]。因此,實現海冰SAR圖像的自動分割變得尤為重要。
為了滿足海冰SAR圖像分割的實時性要求,計算機視覺技術在此領域得到了重要應用。目前已提出的海冰圖像分割方法中普遍采用了提取SAR圖像的紋理特征描述海冰的方法[3-5],紋理特征提取方法主要包括灰度共生矩陣、Gabor變換、小波變換和馬爾科夫隨機場模型等。和強度特征相比,紋理特征通過描述像素的位置分布特性能夠更好地實現海冰SAR圖像分割。
基于Radarsat-1數據,本文提出了一種集成Gibbs采樣標記步驟的海冰SAR圖像分割算法,該方法經由像素到區域再到大尺度區域這一途徑,把區域化標識、Gibbs采樣標記以及區域合并等操作組合起來,以有效提高分割算法的準確度。
假設把圖像分為n(冰)類,S是柵格上的一組位置,s∈S是柵格中的某個位置。
圖像分割是把原始圖像轉換成符號,即
R:{Is|s∈S}→{xs|s∈S}。
(1)
通過這種方式圖像空間被劃分成n個分割段Ω1,…,Ωn,得到
a)Ωi={s|xs=i,s∈S},
c)?i≠j:Ωi∩Ωj=φ。
(2)
在上面的定義中,任一Ωi可以是單個區域,也可以是一個不相鄰區域集合(區域集群)。每個區域中的像素具有相同的類標記。
圖1為算法流程圖,圖中虛線給出了本文所提方法的主要步驟,包括:預處理、區域化表示、計算區域特征、Gibbs采樣標記和區域合并。實線代表了未包括Gibbs采樣標記步驟的分割流程,稱之為基于區域的MRF(Markov Random Field,MRF)圖像分割算法。

圖1 算法流程圖
由于SAR圖像受相干斑噪聲的影響,如果直接進行分水嶺分割,過分割現象比較嚴重,影響后續的分割與分類的準確性。因此,在分水嶺分割之前需進行SARD濾波[6],在保持地物邊界的同時能夠有效地抑制過分割現象,減少分割時間。SARD濾波是將各向異性擴散引入傳統的自適應相干濾波器,構造偏微分方程來消除SAR圖像的相干斑噪聲。區域化表示的目的是消除雷達圖像的噪聲影響,并識別顯著的同質區域。具體步驟為:
對濾波結果進行分水嶺分割[7],獲得大量同質區域,每個區域v由一組屬于該區域的位置Sv構成,通過把各位置的特征矢量{ys|s∈Sv}平均到每一個區域的特征矢量yv上,可以有效降低斑點噪聲,此過程稱為結構去噪。接下來構建區域鄰接圖(RAG)[8],在RAG中,每個結點代表一個區域,結點之間的弧表示相鄰區域之間的共同邊界。分水嶺分割及RAG表示的示意圖如圖2所示,其中,U(x)表示區域的一元屬性,U(x,y)表示鄰接區域之間的二元屬性。

圖2 區域鄰接圖構建示意圖
假設特征模型服從高斯分布,則矢量ys的分布是n個高斯函數的混合,這里n是類的數量。通過期望最大化(EM)算法求解高斯混合模型的參數[2]。
根據給定的觀察數據和所有參數的當前估計,假設πi代表類Ωi的先驗,ωsi表示位置s屬于類Ωi的概率。EM算法通過計算E步驟進行迭代:
(3)
迭代M步驟:
(4)
(5)
(6)
經由上述EM算法獲得的參數就是初始使用的參數。在每一次迭代中進一步更新參數,使用M步驟計算參數。在本文中對大量不同初始值進行了隨機選擇,在少量迭代后,選擇具有最小能量的值,并經由更多次迭代進行細化,以獲得最終解。
Gibbs采樣標記的目的是要找出每個分水嶺區域的最優標記,這是通過找出標記的配置來完成的。最小化的成本函數表示如下:

(7)
其中,Ri表示第i個區域;ls表示區域Ri的類別標號;μls和σls表示區域Ri內所有像素的均值和標準差;C表示雙位團鄰域系統;g(st)表示邊界強度函數,即邊界懲罰函數;點s和點t組成屬于不同類別的鄰接區域的共享邊界的點對。
在每次迭代的同時進行區域合并的操作,減少了區域的數目,降低了優化過程最優解的搜索空間,提高了運算速度。鄰接區域可以合并的條件如下:

(8)
其中,Ri和Rj是兩個鄰接區域,同時Ri和Rj的類別標號不同;Rk=Ri∪Rj;NRi表示區域Ri的邊界點的集合;ηs表示s點的鄰域系統[2]。當δE為負數時,即區域合并前后能量減少,表示區域Ri和Rj可以合并,否則不能合并。
評價海冰圖像分割效果常用的兩個度量指標是分割精度和Kappa系數[9]。分割精度是指正確劃分像素的百分比。Kappa系數定義如下:
(9)
式中,P(A)為觀測一致率,P(E)為期望一致率。Kappa系數的范圍為[-1,1]。
為驗證本文提出的算法的有效性,選擇兩幅不同的海冰SAR圖像,應用基于區域的MRF分割算法和本文提出的算法進行分割測試,結果如圖3,圖4和表1所示。兩幅測試圖像中第一幅是由Radarsat-1模式于2003年2月18日采集的S. Laurence Bay圖像,包括初期冰和平滑冰兩類,如圖3(a)所示;第二幅由Radarsat-1模式于2006年11月27日在Beaufort海上空采集,包含初期冰和平滑冰兩類,如圖4(a)所示。
對比圖3中(b)(c)和圖4中(b)(c)可知,基于區域的MRF分割算法出現了海冰類型誤分的情況,即原屬于平滑冰的區域被錯誤地劃分為初期冰,同樣,原屬于初期冰的區域也被錯誤地劃分為平滑冰。而且定性地看,對于Beaufort海冰SAR圖像,基于區域的MRF算法的分割結果(圖4(b))顯示出了更多的歧義性,特別是在中、右區段,而本文提出的分割算法的分割結果(圖4(c))具有更好的區域連通性,在同樣的區域能夠比較好地兼顧邊界位置定位和區域一致性,圖像中孤立區域明顯減少。

(a)原始海冰SAR圖像;(b)基于區域的MRF算法分割結果;(c)本文算法分割結果;(d)地面實況圖

(a)原始SAR海冰圖像;(b)基于區域的MRF算法分割結果;(c)本文算法分割結果;(d)地面實況圖

方法S. Laurence Bay海冰SAR圖像分割精度(%)Kappa系數分割數所用時間/sBeaufort海海冰SAR圖像分割精度(%)Kappa系數分割數所用時間/s基于區域的MRF算法76.890.748 629 1532380.640.776 435 48525本文算法81.630.796 514 6861287.750.857 316 45718
由表1可知,對于S. Laurence Bay海冰SAR圖像,基于區域的MRF分割算法的Kappa系數是0.748 6,精度為76.89%,而本文提出的分割算法Kappa系數是0.796 5,精度為81.63%;對于Beaufort海海冰SAR圖像,基于區域的MRF分割算法的Kappa系數是0.776 4,精度為80.64%,而本文提出的分割算法Kappa系數是0.857 3,精度為87.75%。數據表明本文算法更有效地實現了海冰SAR圖像的準確分割。
本文算法利用Gibbs采樣標記實現了海冰SAR圖像的準確分割,克服了分割特征魯棒性差的問題,并從特征的魯棒性和分割過程的自適應性兩方面提高了對SAR海冰圖像結構的認知能力。另外,從表1還可以看出,本文提出的分割算法與基于區域的MRF分割算法相比,大幅度減少了分割區數目,從而降低了計算量,節省了計算時間。
本文針對Radarsat-1獲取的海冰SAR圖像提出了一種新的分割算法,該算法經由像素到區域再到大尺度區域這一途徑,把區域化、Gibbs采樣標記以及區域合并等操作組合起來實現對海冰SAR圖像的分割。將本文算法運用到由Radarsat-1所獲取的兩幅真實海冰SAR圖像進行測試,結果表明,在精確定位邊界和生成較大均勻區域方面該算法具有優越性,同時證明了該算法的可行性和準確性,這為實現海冰SAR圖像的準確分割提供了一條新途徑。