(中國(guó)電子科技集團(tuán)公司第三十八研究所,安徽合肥230088)

0 引言

變化檢測(cè)是通過分析同一地區(qū)不同時(shí)刻獲取的圖像來識(shí)別一個(gè)目標(biāo)或者現(xiàn)象的狀態(tài)差異的過程[1]。近年來,由于SAR圖像相比光學(xué)圖像具有不受氣候和太陽光照條件影響的特點(diǎn),越來越多的學(xué)者通過分析SAR圖像來識(shí)別相關(guān)的變化。雖然SAR圖像數(shù)據(jù)的內(nèi)在復(fù)雜性,既要求一個(gè)有效的降斑預(yù)處理,又要求一個(gè)有效的能處理該斑點(diǎn)噪聲的分析技術(shù),使得通過SAR圖像識(shí)別變化區(qū)域的難度有所增加,但仍有大量的相關(guān)技術(shù)形成,并用于地震災(zāi)區(qū)的標(biāo)識(shí)[2]、水災(zāi)監(jiān)控[3]和栽種面積的評(píng)估[4]等方面。

一般地,SAR圖像變化檢測(cè)由預(yù)處理、差異圖構(gòu)造和變化檢測(cè)結(jié)果生成[5]三步完成。預(yù)處理主要是對(duì)原始兩時(shí)相SAR圖像進(jìn)行諸如匹配、輻射校正、幾何校正和噪聲去除等處理。對(duì)已配準(zhǔn)、已校正的兩時(shí)相SAR圖像,預(yù)處理的主要目的是降低斑點(diǎn)噪聲。現(xiàn)有變化檢測(cè)方法大多采用Lee或Gamma MAP等濾波器進(jìn)行降斑。差異圖構(gòu)造則是為生成可后續(xù)處理的特征圖,常用的方法有比值法、對(duì)數(shù)比值法和區(qū)域均值法等。變化檢測(cè)結(jié)果生成,根據(jù)其所利用的信息,可分為上下文不敏感和上下文敏感兩大類。上下文不敏感方法只利用像素自身的信息,如文獻(xiàn)[5]方法,通過假設(shè)類統(tǒng)計(jì)分布模型和推導(dǎo)出相應(yīng)的閾值選擇標(biāo)準(zhǔn),閾值分割差異圖生成變化檢測(cè)結(jié)果。上下文敏感方法利用像素的一個(gè)鄰域內(nèi)信息。由于上下文敏感方法考慮了相鄰像素間具有的相關(guān)性,其變化檢測(cè)性能通常優(yōu)于上下文不敏感方法,諸如通過馬爾可夫模型(MRF,Markov random field)模型[6]、Hopfield神經(jīng)網(wǎng)絡(luò)(HNN,Hopfield neural network)以能量的形式將空間關(guān)系引入變化檢測(cè)結(jié)果生成階段。結(jié)果顯示,由于鄰域像素信息的利用,變化檢測(cè)性能都有所提高,但該鄰域的大小和形狀都被固定,導(dǎo)致變化檢測(cè)方法對(duì)幾何細(xì)節(jié)不敏感,造成變化區(qū)域的邊緣有過光滑現(xiàn)象。

為提高變化檢測(cè)性能和改善變化檢測(cè)效果,本文提出一種新的差異圖構(gòu)造方法,并基于主動(dòng)輪廓模型分割技術(shù)完成了兩時(shí)相SAR圖像的水災(zāi)變化檢測(cè)。所提方法避免了上下文不敏感技術(shù)[5-6]的變化檢測(cè)性能受影響于統(tǒng)計(jì)分布模型的合理性的缺陷,突破了上下文敏感技術(shù)只利用固定形狀和大小鄰域內(nèi)像素信息的限制。本文方法因無需假設(shè)統(tǒng)計(jì)分布模型和能利用非緊湊鄰域內(nèi)像素信息對(duì)當(dāng)前像素進(jìn)行決策,其變化檢測(cè)性能優(yōu)于其他相關(guān)算法。兩組真實(shí)SAR圖像數(shù)據(jù)的變化檢測(cè)結(jié)果,驗(yàn)證了本文方法的有效性。

1 差異圖構(gòu)造及變化檢測(cè)結(jié)果生成

1.1 模糊貼近度差異圖

模糊貼近度原本為度量?jī)蓚€(gè)模糊子集間相似性測(cè)度。Weken等人[7]成功將其應(yīng)用于圖像比較,并為改善效果,將基于像素的相似性測(cè)度推廣到基于鄰域的相似性測(cè)度。貼近度度量公式的形式有多種,但本文選用如式(1)相似性度量公式構(gòu)造差異圖,是在于其采用的比值形式和能利用窗口內(nèi)像素信息增強(qiáng)差異圖的類內(nèi)一致性,這有利于克服SAR圖像的乘性斑點(diǎn)噪聲和后續(xù)變化檢測(cè)結(jié)果生成。

考慮到窗口在降低斑點(diǎn)噪聲和保持細(xì)節(jié)信息方面的兩難兼顧,對(duì)已配準(zhǔn)、已校正的兩時(shí)相SAR圖像I1和I2,選取3×3的窗口N(x)構(gòu)造差異圖DI,其中N(x)為與x在位置上相鄰的位置元素構(gòu)成的集合,諸如二階鄰域。對(duì)位置x上的像素,取其空間上的一個(gè)鄰域,把由鄰域元素構(gòu)成的集合看作是模糊向量,用式(1)度量?jī)赡：蛄块g的相似度。式(1)值越大,位置x上的對(duì)應(yīng)像素間的相似度越高,則該點(diǎn)屬于非變化類的可能性更大;反之,該點(diǎn)屬于變化類的可能性更大。為方便后面陳述,將依據(jù)相似度公式構(gòu)造的差異圖稱為模糊貼近度差異圖。

1.2 變化檢測(cè)結(jié)果生成

主動(dòng)輪廓模型是一種基于能量的分割方法,其基本思想是進(jìn)化一條或多條受限于給定圖像的曲線,使之向目標(biāo)的邊緣逐漸靠近,以達(dá)到分割的目的。對(duì)已配準(zhǔn)、已校正的SAR圖像構(gòu)造差異圖DI后,本文以主動(dòng)輪廓模型為原型的分割技術(shù)對(duì)差異圖分析生成變化檢測(cè)圖,即把差異圖分成兩類：變化類wc和非變化類wn。若把DI∶Ω→R看成一個(gè)包含兩區(qū)域的連續(xù)的標(biāo)量場(chǎng),并用Ω表示其圖像域,C0表示DI中兩類的邊界曲線,即wc的邊緣;則對(duì)在圖像域Ω中的任意一條封閉曲線C,我們希望曲線C進(jìn)化到邊界曲線C0時(shí)停止,此時(shí)曲線C包圍的區(qū)域?qū)?yīng)wc,C之外而Ω內(nèi)的區(qū)域?qū)?yīng)wn,也就完成了變化檢測(cè)。

以主動(dòng)輪廓模型為原型的分割技術(shù),其本質(zhì)為把分割的問題映射為能量的問題,通過求能量的最小化問題完成圖像分割,能量表達(dá)式為

式中,C為任意假設(shè)的在Ω中的封閉曲線,若假設(shè)C為Ω中一個(gè)真子集w的邊界(即w?Ω,C=?w),則inside(C)為區(qū)域w,outside(C)為區(qū)域Ω/ˉw;Per(inside(C))是給曲線C施加光滑約束的正則項(xiàng);c1和c2是依賴于C的常數(shù),其中c1表示曲線C包含的像素灰度均值,c2表示落在曲線C外像素的灰度均值。不難看出,當(dāng)C收縮在wc的邊緣C0時(shí),式(2)達(dá)到最小值,即以主動(dòng)輪廓模型為原型的分割技術(shù)對(duì)差異圖分析的問題轉(zhuǎn)化為求式(2)最小值的問題。

直接求解式(2)的最小值存在易陷于局部最小值的缺點(diǎn),因此利用加權(quán)的總差分范數(shù)(TV,total variation)將經(jīng)典的主動(dòng)輪廓模型和無邊界的主動(dòng)輪廓模型(ACWE,active contour without edges)模型結(jié)合,變換如式(2)所示的能量函數(shù)為

式中,r(x,c1,c2)=(c1-DI(x))2-(c2-DI(x))2;表示以經(jīng)典主動(dòng)輪廓模型中的邊緣指示函數(shù)g為權(quán)重的總差分范數(shù)。

求解如式(3)所示的能量表達(dá)式的最小值時(shí),能確保收斂到全局最小,但仍存在收斂速度比較慢的問題。為此,用對(duì)偶規(guī)劃將式(3)變?yōu)槭?4)[8]：

式中,參數(shù)θ為大于零的小數(shù)值。求解如式(4)表示的能量最小化問題時(shí),已有的數(shù)值實(shí)現(xiàn)能確保較快的收斂到全局最小;該策略通過迭代的解決如下兩個(gè)能量最小值問題：

(1)固定v不變,求u

(2)固定u不變,求v

來求如式(4)表示的能量最小化問題,而式(5)的解為

式(6)的解為

更多的信息可參見文獻(xiàn)[8]。一定的迭代次數(shù)完成后,u即為所求的變化檢測(cè)結(jié)果。從函數(shù)形式轉(zhuǎn)變的整個(gè)過程來看,如式(4)所示的能量形式使曲線的進(jìn)化趨勢(shì)由經(jīng)典主動(dòng)輪廓模型中的指示函數(shù)g與被ACWE模型利用的區(qū)域同質(zhì)性特點(diǎn)共同決定,以致即使差異圖中兩類的對(duì)比度不高,該模型也能使曲線逐漸進(jìn)化收斂于變化區(qū)域的邊緣輪廓。

為使模糊貼近度差異圖近似滿足分片恒值(piecewise-constant)的條件,可先采用雙邊濾波器處理模糊貼近度差異圖。

2 實(shí)驗(yàn)分析

為了證明本文方法的有效性,我們分別證明了如下兩個(gè)方面：(1)差異圖構(gòu)造方法的有效性;(2)變化檢測(cè)性能改善的有效性。為定量地比較變化檢測(cè)效果,依據(jù)參考圖,給出了不同變化檢測(cè)結(jié)果的漏檢數(shù)(變化的像素被判為未變化的像素個(gè)數(shù))、虛警(未變化的像素被判為變化的像素個(gè)數(shù))、總錯(cuò)誤檢測(cè)數(shù)(漏檢數(shù)+虛警數(shù))和Kappa系數(shù)等四個(gè)性能評(píng)價(jià)指標(biāo),

式中,N為差異圖像素單元總數(shù),E c為變化區(qū)域內(nèi)像素單元被正確檢測(cè)的個(gè)數(shù),T c為非變化區(qū)域內(nèi)像素單元被正確檢測(cè)的個(gè)數(shù)。漏檢及虛警的像素單元個(gè)數(shù)越少,Kappa系數(shù)越大,變化檢測(cè)性能越好。

2.1 實(shí)驗(yàn)數(shù)據(jù)描述

圖1給出了用于驗(yàn)證變化檢測(cè)性能的第一組真實(shí)SAR圖像數(shù)據(jù)集。該數(shù)據(jù)集由ERS-2分別于1999年4月和5月獲取的瑞士Bern城市水災(zāi)發(fā)生前后兩時(shí)相SAR影像構(gòu)成。圖1(a)為洪水發(fā)生前SAR影像,圖1(b)為受水災(zāi)后SAR影像。SAR影像尺寸均為301×301。圖1(c)為實(shí)際變化參考圖,圖中的白色像素表示兩時(shí)相間發(fā)生變化的區(qū)域,像素個(gè)數(shù)為1 155。

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

第二組真實(shí)的SAR圖像數(shù)據(jù)集原始圖像及參考變化圖如圖2所示。其中數(shù)據(jù)集的原始圖像分別為加拿大Ottawa地區(qū)的1997年5月和1997年8月Radarsat SAR影像,如圖2(a)和(b)所示。兩幅圖像的大小均為350×290,灰度級(jí)為256,配準(zhǔn)誤差為1.5個(gè)像素左右。前后時(shí)相的變化在于雨季(5月)期間部分陸地區(qū)域被淹沒。圖2(c)為參考變化檢測(cè)圖,圖中的白色像素表示兩時(shí)相間發(fā)生變化的區(qū)域,像素個(gè)數(shù)為16 049。

圖2 Ottawa實(shí)驗(yàn)數(shù)據(jù)

2.2 差異圖構(gòu)造方法有效性驗(yàn)證

依據(jù)參考圖,采用一種有監(jiān)督的方式,計(jì)算差異圖中變化類和未變化類的方差。方差越小,類內(nèi)一致性越好。圖3展示了比值法、對(duì)數(shù)比值法及區(qū)域均值比值法等經(jīng)典構(gòu)造方法的差異圖以及本文的模糊貼近度差異圖,其中構(gòu)造區(qū)域均值比法差異圖時(shí),同樣采用3×3的窗口。從效果圖上看,比值法差異圖的變化區(qū)域很明顯,但未變化區(qū)域存在易產(chǎn)生錯(cuò)誤檢測(cè)數(shù)的像素,因其灰度值與變化區(qū)域內(nèi)像素相當(dāng);對(duì)數(shù)比值法差異圖的未變化區(qū)域的一致性較好,但是變化區(qū)域信息有所丟失,易造成變化區(qū)域像素單元的漏檢較多;區(qū)域均值比差異圖,相比于本文差異圖,變化類與非變化類的對(duì)比度較強(qiáng),但未變化類的一致性略差于本文差異圖。

圖3 兩組實(shí)驗(yàn)數(shù)據(jù)各方法構(gòu)造的差異圖

表1 兩組實(shí)驗(yàn)數(shù)據(jù)差異圖中兩類的方差

分析表1數(shù)據(jù)可知：對(duì)Bern實(shí)驗(yàn)數(shù)據(jù),本文方法構(gòu)造的差異圖,變化類和未變化類的方差均優(yōu)于經(jīng)典方法構(gòu)造的差異圖;對(duì)Ottawa實(shí)驗(yàn)數(shù)據(jù)的結(jié)果亦如此。這從定量的角度驗(yàn)證了本文模糊貼近度差異圖構(gòu)造方法的有效性。

2.3 變化檢測(cè)性能改善有效性驗(yàn)證

本文所提變化檢測(cè)方法具有不涉及統(tǒng)計(jì)分布模型假設(shè)和上下文敏感兩個(gè)特點(diǎn)。為了驗(yàn)證不涉及統(tǒng)計(jì)分布模型假設(shè)的有效性,給出了涉及統(tǒng)計(jì)分布模型變化檢測(cè)方法的檢測(cè)結(jié)果圖,為陳述方便,假設(shè)類分布符合高斯模型的分析技術(shù)記為GKIT,假設(shè)類分布為廣義高斯模型的分析技術(shù)記為GGKIT;假設(shè)類分布符合對(duì)數(shù)正態(tài)分布的分析技術(shù)記為L(zhǎng)N-GKIT;假設(shè)類分布符合“Weibull-Ratio”模型的分析技術(shù)記為 WR-GKIT。為了體現(xiàn)本文方法利用上下文信息對(duì)檢測(cè)效果的貢獻(xiàn),給出了兩種已有上下文敏感變化檢測(cè)方法的結(jié)果,分別記為EM+MRF,HNN。

圖4給出了Bern實(shí)驗(yàn)數(shù)據(jù)集的變化檢測(cè)結(jié)果;圖5給出了Ottawa實(shí)驗(yàn)數(shù)據(jù)集的變化檢測(cè)結(jié)果。表2給出了對(duì)應(yīng)的定量分析評(píng)價(jià)指標(biāo)。

表2 兩組實(shí)驗(yàn)數(shù)據(jù)變化檢測(cè)結(jié)果的定量分析指標(biāo)

圖4 Bern實(shí)驗(yàn)數(shù)據(jù)變化檢測(cè)結(jié)果

圖5 Ottawa實(shí)驗(yàn)數(shù)據(jù)變化檢測(cè)結(jié)果

分析表2的數(shù)據(jù),本文所提變化檢測(cè)方法對(duì)兩組實(shí)驗(yàn)數(shù)據(jù)變化檢測(cè)結(jié)果的Kappa系數(shù)都最高。對(duì)Bern實(shí)驗(yàn)數(shù)據(jù),本文方法的Kappa系數(shù)比GGKIT法高0.0172,比LN-GKIT方法高0.0367,比EM+MRF法高0.1731,比HNN法高0.2439;對(duì)Ottawa實(shí)驗(yàn)數(shù)據(jù),本文方法的Kappa系數(shù)比GKIT法高0.063,比 WR-GKIT法高0.064 1,比EM+MRF法高0.046 2,比HNN法高0.089 8。這是因?yàn)楸疚姆椒ū苊饬祟惤y(tǒng)計(jì)分布選擇的合理性對(duì)檢測(cè)結(jié)果的影響,同時(shí),通過非緊湊鄰域像素的信息增強(qiáng)了對(duì)幾何細(xì)節(jié)信息的敏感性。

兩組實(shí)驗(yàn)數(shù)據(jù)的分析,從定量的角度驗(yàn)證了本文方法的有效性。從定性的角度看,本文方法檢測(cè)的變化區(qū)域的邊緣更接近于真實(shí)變化區(qū)域的輪廓。

3 結(jié)束語

基于模糊數(shù)學(xué)概念,本文構(gòu)造了一種新的差異圖,并基于主動(dòng)輪廓模型,突破了只利用固定形狀和大小鄰域內(nèi)像素信息的限制,利用了非緊湊局部鄰域內(nèi)的像素信息,形成一種非緊湊鄰域的上下文敏感變化檢測(cè)方法。真實(shí)SAR圖像數(shù)據(jù)的變化檢測(cè)結(jié)果表明,本文所提變化檢測(cè)方法既提高了變化檢測(cè)性能,又改善了視覺效果。

[1]王剛,劉智,王番,等.基于KL距離的SAR影像變化檢測(cè)[J].雷達(dá)科學(xué)與技術(shù),2012,10(1)：59-63.WANG Gang,LIU Zhi,WANG Fan,et al.A Change Detection Method of SAR Image Based on Kullback Leibler Divergence and Morohology[J].Radar Science and Technology,2012,10(1)：59-63.(in Chinese)

[2]BRETT P T B,GUIDA R.Earthquake Damage Detection in Urban Areas Using Curvilinear Features[J].IEEE Trans on Geoscience and Remote Sensing,2013,51(9)：4877-4884.

[3]GIUSTARINI L,HOSTACHE R.A Change Detection Approach to Flood Mapping in Urban Areas Using TerraSAR-X[J].IEEE Trans on Geoscience and Remote Sensing,2013,51(4)：2417-2430.

[4]SAKAMOTO T,WARDLOW B D,GITELSON A A.Detecting Spatiotemporal Changes of Corn Development Stages in the U.S.Corn Belt Using MODIS WDKVI Data[J].IEEE Trans on Geoscience and Remote Sensing,2011,49(6)：1926-1936.

[5]MOSER G,SERPICO S B.Generalized Minimum-Error Thresholding for Unsupervised Change Detection from SAR Amplitude Imagery[J].IEEE Trans on Geoscience and Remote Sensing,2006,44(10)：2972-2983.

[6]BRUZZONE L,PRIETO D F.Automatic Analysis of the Difference Image for Unsupervised Change Detection[J].IEEE Trans on Geoscience and Remote Sensing,2000,38(3)：1171-1182.

[7]VAN DER WEKEN D,NACHTEGAEL M,KERRE E E.Using Similarity Measures and Homogeneity for the Comparision of Images[J].Image and Vision Computing,2004,22(9)：695-702.

[8]BRESSON X,ESEDOGLU S,VANDERGHEYNST P,et al.Fast Global Minimization of the Active Contour/Snake Model[J].Journal of Mathematical Imaging and Vision(JMIV),2007,28(2)：151-167.