李軍軍, 曹建農(nóng), 廖娟, 程貝貝
(長(zhǎng)安大學(xué)地球科學(xué)與資源學(xué)院,西安 710054)
隨著遙感影像的空間分辨率提高,地物內(nèi)部細(xì)節(jié)得到了更好表征,空間結(jié)構(gòu)信息更加豐富,地物形狀以及與相鄰地物之間的關(guān)系得到更加詳盡的反映。因此,如何利用計(jì)算機(jī)更有效地處理這些高分影像成為當(dāng)前一項(xiàng)熱門研究課題。高空間分辨率遙感影像的邊緣提取及其相關(guān)的應(yīng)用,在地區(qū)規(guī)劃和防災(zāi)救災(zāi)等過(guò)程中占據(jù)著重要位置。傳統(tǒng)的邊緣檢測(cè)方法如Sobel算子、Robert算子和Canny算子在許多應(yīng)用中都具有較好的效果[1-2]。其中Canny算子在邊緣檢測(cè)中能得到較為清晰連續(xù)的邊界而得到許多學(xué)者的研究與改進(jìn)[3-6]。但上述算子的本質(zhì)仍是一種基于局部梯度的算法,容易受到噪聲的影響,遙感影像的獲取過(guò)程中由于受到自然天氣條件等處界因素的影響通常帶有噪聲。而使用濾波方法對(duì)圖像去噪會(huì)平滑掉部分弱邊緣,使得一部分邊緣丟失或不連續(xù)。因此,這些方法存在著抑制噪聲與邊緣提取結(jié)果精確性間的矛盾[7]。
小波變換作為非平穩(wěn)信號(hào)分析的有力工具,其基函數(shù)具有局部性和衰減性,不僅能夠得到信號(hào)分解后的頻率,而且能夠知道頻率發(fā)生的位置[8]。圖像邊緣通常對(duì)應(yīng)高頻信息,使用小波變換能夠有效檢測(cè)圖像邊緣,并且具有較好的邊緣定位性質(zhì),所以,采用小波變換提取圖像邊緣成為當(dāng)前一個(gè)研究方向。文獻(xiàn)[9-11]中提出以小波變換為基礎(chǔ)的不同影像邊緣提取方法,提取結(jié)果相較于傳統(tǒng)方法有一定提高。但是,標(biāo)準(zhǔn)小波變換方向數(shù)較少,小波系數(shù)僅被分解到水平、垂直和對(duì)角3個(gè)方向上,且對(duì)角子帶系數(shù)一般由于含有較多噪聲而舍棄。而高分遙感影像中地物內(nèi)部空間結(jié)構(gòu)信息豐富,地物邊緣存在于各個(gè)方向,使得基于小波變換的邊緣提取存在結(jié)果不完整問(wèn)題。在此基礎(chǔ)上,研究人員發(fā)展了一些二維圖像的表示方法,如Curvelet,Ridgelet,Contourlet和Brushlet等[12-15],相比于小波變換能更加精確地表示二維圖像,因此被廣泛應(yīng)用于圖像壓縮、分割及去噪等領(lǐng)域。文獻(xiàn)[16-18]等使用上述方法提取不同目標(biāo)的邊緣信息,但這些方法也存在著明顯的缺點(diǎn),濾波器的結(jié)構(gòu)設(shè)計(jì)較為復(fù)雜,計(jì)算復(fù)雜度較高等。例如Contourlet變換是過(guò)完備的,即圖像變換后的大小約為原始圖像大小的4/3; Curvelet變換冗余度高,存在加窗效應(yīng),實(shí)現(xiàn)效率差; Ridgelet變換是單尺度的,因此這些方法在圖像邊緣分析中都存在一定的不足[19],尤其是對(duì)于數(shù)據(jù)量較大的遙感影像。
針對(duì)上述方法中的部分缺陷,本文設(shè)想提出一種理想的方法,通過(guò)采用具有多方向性基函數(shù)對(duì)圖像進(jìn)行不同方向的一維小波變換,能夠從多個(gè)方向檢測(cè)圖像的奇異點(diǎn),獲取影像邊緣信息。這樣既具有一維小波變換的計(jì)算簡(jiǎn)單性,又具有二維不可分離小波變換的多方向性。相較于文獻(xiàn)[16-18]等方法具有較低的計(jì)算復(fù)雜度,同時(shí)也具有其多方向性的優(yōu)勢(shì),在獲得更加完整的邊緣提取結(jié)果的同時(shí)減少時(shí)間消耗。
圖1(a)顯示了使用某一方向性基函數(shù)檢測(cè)與其方法正交或近似正交的奇異曲線邊緣,以獲取較大的邊緣系數(shù)。這種方法處理的好處在于,使用多個(gè)方向處理后選取較大邊緣系數(shù)能夠得到較為連續(xù)的邊緣,不會(huì)因?yàn)闄z測(cè)方向數(shù)較少而產(chǎn)生一些曲線部分系數(shù)與兩側(cè)差異過(guò)大從而導(dǎo)致最后的邊緣不連續(xù); 而且不具備線狀結(jié)構(gòu)的噪聲則仍會(huì)在后續(xù)形態(tài)學(xué)處理中被消除。本文從Directionlet變換[20]出發(fā),根據(jù)其基本理論與方法完成圖1(b)所示的過(guò)程。
(a) 方向性基函數(shù)原始高分影像(b) 垂直基函數(shù)與經(jīng)過(guò)調(diào)制的高分影像
圖1影像方向性小波變換
Fig.1DirectionalwavelettransformofHSRimage
利用具有不同方向特性的小波基函數(shù)對(duì)高分影像進(jìn)行變換,等價(jià)于對(duì)影像預(yù)先調(diào)制再利用水平或垂直方向的小波變換。
Directionlet變換是一種多尺度幾何變換工具,相比于其他常用的多尺度幾何變換工具除了具有多尺度、多方向和局部化的優(yōu)點(diǎn)以外還具有標(biāo)準(zhǔn)二維小波變換的計(jì)算簡(jiǎn)單性。其基函數(shù)被稱為斜各向異性小波基S-AWT(n1,n2),n1和n2分別表示沿不同方向的變換次數(shù),當(dāng)n1=n2時(shí)表明變換是各向同性的,否則為各向異性變換。
定義由2個(gè)線性無(wú)關(guān)的方向向量d1=[a1,b1],d2=[a2,b2](斜率分別為r1=b1/a1和r2=b2/a2,且r1和r2都為有理數(shù)),d1和d2線性組合構(gòu)成滿秩整數(shù)格Λ,其基本形式用MΛ表示為
(1)
任一整數(shù)格Λ均為立方整數(shù)格(cubic integer lattice)Ζ2的子格,即Λ?Ζ2。
由格理論[21]可知,完整的格空間被MΛ劃分為|det(MΛ)|個(gè)關(guān)于Λ的陪集(co-set,CS)。每一個(gè)陪集CS[Sk,1 ,Sk,2]可由格Λ通過(guò)偏移矢量Sk=[s1,s2](k=1,2,...,|det(MΛ)| ,s1,s2∈Ζ)得到。位于格Λ及其CS上的像素點(diǎn)的坐標(biāo)可以表示為(i,j),即
i=c1a1+c2a2+Sk,1,
(2)
j=c1b1+c2b2+Sk,2,
(3)
式中:c1和c2為任意整數(shù),c1,c2∈Ζ,具體取值視影像大小而定;Sk,1=s1;Sk,2=s2。當(dāng)固定任意c2,遍歷c1時(shí)可以得到嚴(yán)格位于d1方向上的一條陪線(co-line,CL)。同樣,由偏移矢量可以得到其他CS上的CL。反之可以得到d2方向上CL。分別沿d1和d2方向進(jìn)行n1和n2次一維小波變換則完成了單尺度下的Directionlet變換,詳細(xì)內(nèi)容請(qǐng)參考文獻(xiàn)[20]。在此不再詳細(xì)討論Directionlet變換的整個(gè)過(guò)程,而是根據(jù)其基本理論方法獲取多方向一維CL用于進(jìn)行多方向小波變換。
以d1=[-1,2]和d2=[1,1]所構(gòu)成的一個(gè)基本單元格為例。圖2(a)中黑色點(diǎn)表示了由MΛ生成的一個(gè)基本單元格,|det(MΛ)|=3表示平面可以被MΛ分解為3個(gè)子集,偏移矢量分別為S1=[0,0],S2=[0,1],S3=[0,2]。變換c1和c2的值完成整個(gè)分解過(guò)程得到CS[0,0],其余2個(gè)CS[0,1]和CS[0,2]上的點(diǎn)可以由CS[0,0]上的點(diǎn)經(jīng)過(guò)偏移矢量得到。圖2(b)表示3條CL,CL上的點(diǎn)都嚴(yán)格位于d2方向上,同樣其余的CLSk也可以通過(guò)偏移矢量獲得。將上述方法應(yīng)用于圖像,即圖2中圓圈視為像素點(diǎn),完成整個(gè)分解過(guò)程獲取任意有理斜率方向的像素點(diǎn)集合。
(a) 一個(gè)基本單元格(b) 陪線
圖2格與陪線
Fig.2Latticeandco-line
小波變換是對(duì)傅里葉變換的發(fā)展,是一種在空間域和頻率域同時(shí)擁有分辨性的多尺度分析方法。小波變換使用小波函數(shù)族及相應(yīng)的尺度函數(shù)將原信號(hào)分解成不同頻率子帶。對(duì)于連續(xù)函數(shù)f(x),定義其小波變換為
(4)
根據(jù)Directionlet變換的基本理論,對(duì)于給定的圖像可以獲取任意有理斜率方向上的CL集合 。例如分別采用[1,0],[2,1],[1,1],[1,2]4個(gè)方向(分別對(duì)應(yīng)0°,30°,45°,60°)對(duì)圖像進(jìn)行多方向性小波變換,并給出變換后的子帶系數(shù),其結(jié)果如圖3所示。
(a) 變換方向(b) 高分影像(c) [1,0]方向
(d) [2,1]方向(e) [1,1]方向(f) [1,2]方向
圖3高分影像不同方向小波變換結(jié)果
Fig.3WavelettransformresultsofHSRimageindifferentdirections
從圖3結(jié)果可以看出,在不同方向小波的作用下其變換結(jié)果不同。由于變換方向與紅色圈中道路的方向一致,在Directionlet變換方向消失矩(directional vanishing moments ,DVM)的作用下,邊緣在沿30°方向的子帶中基本消失,即小波系數(shù)較小,在與30°方向相鄰的水平和45°方向上的子帶中重現(xiàn)但數(shù)值仍不大,而在60°方向子帶上系數(shù)的響應(yīng)相對(duì)較大,且呈高亮狀態(tài)。
圖像經(jīng)變換后最重要的是需要選合適的系數(shù)進(jìn)行后處理。以標(biāo)準(zhǔn)小波變換為例,fd1和fd2分別是圖像f經(jīng)過(guò)水平與垂直方向小波變換后的系數(shù),對(duì)于圖像中的任意一點(diǎn)(x,y),fd1(x,y)和fd2(x,y)正好構(gòu)成了一組梯度矢量,由此可得f(x,y)小波變換后的梯度模值Mf(x,y)和幅角Af(x,y)分別為
(5)
(6)
由于圖像中突變點(diǎn)主要位于圖像邊緣,這些突變點(diǎn)對(duì)應(yīng)于模極大值,梯度方向指向梯度模極大值,這就意味著可以通過(guò)梯度方向?qū)ふ夷O大值來(lái)檢測(cè)圖像邊緣[22]得到初始結(jié)果。
傳統(tǒng)的方法中,梯度分量只有2個(gè),對(duì)其求模和幅角就可以確定梯度大小和方向。而在本方法中變換方向一般為8個(gè)方向或4個(gè)方向。若以4個(gè)方向?yàn)槔?,則圖像變換后得到4個(gè)梯度分量分別為fd1,fd2,fd3,fd4。文獻(xiàn)[23]對(duì)此處理方法為
(7)
但這種做法太過(guò)籠統(tǒng),且當(dāng)變換方向較多時(shí)計(jì)算量會(huì)很大。在實(shí)際情況中,多數(shù)邊緣提取實(shí)驗(yàn)中發(fā)現(xiàn),只有少數(shù)幾個(gè)梯度分量的模值較大。這符合人們的認(rèn)知,由于邊緣一般出現(xiàn)在像素灰度值差異較大的地方,通常情況下,f(x,y)點(diǎn)處的梯度分量最大值出現(xiàn)在該點(diǎn)法線方向,而其他方向梯度分量的模沿法線兩側(cè)遞減,最小為0,如圖4所示。
圖4 多方向梯度分量Fig.4 Multi-direction gradient components
(8)
非極大值抑制是進(jìn)行邊緣檢測(cè)的重要步驟,用于尋找局部梯度最大值,細(xì)化邊緣。其基本思想是: 如果圖像(x,y)位置上的系數(shù)模值大于沿梯度方向兩端的值,則認(rèn)為該像素位置為候選圖像邊緣,否則為非邊緣點(diǎn)。而通常梯度方向介于2個(gè)梯度分量方向之間的亞像素點(diǎn)位置上(如圖5(a)所示),因此需要對(duì)梯度方向兩邊的梯度值進(jìn)行插值,得到沿梯度方向上的梯度值。在本文中,選取4個(gè)變換方向{[1,0],[1,1],[0,1],[-1,1]},對(duì)圖像進(jìn)行實(shí)驗(yàn)并加以說(shuō)明。
(a) 2個(gè)最大梯度分量方向相鄰(b) 2個(gè)最大梯度分量方向不相鄰
圖52個(gè)最大梯度分量方向的位置關(guān)系
Fig.5Positionrelationofdirectionofthetwomaximumgradientcomponents
多數(shù)情況下,梯度值最大的方向與次大梯度方向是相鄰的,如圖5(a)所示。但是當(dāng)采用多方向小波變換時(shí),由于變換數(shù)量較多及影像本身的幾何結(jié)構(gòu)等因素,存在部分2個(gè)最大梯度方向不相鄰的情況,如圖5(b)所示。對(duì)于第一種情形,采用一般方法,對(duì)梯度分量f1和f2插值得到梯度方向上的2個(gè)梯度值M1和M2,即
(9)
(10)
比較Mf(x,y)與M1,M2值的大小確定當(dāng)前位置是否為候選邊緣點(diǎn)。對(duì)于第二種情況,由于系數(shù)最大方向與系數(shù)第二大方向的位置不相鄰,相關(guān)性弱,無(wú)法直接對(duì)其插值求梯度方向模值。因此,將f(x,y)處的值分別與f1和f2方向前后值比較,若為局部最大值則保留f(x,y)為候選邊緣。圖6顯示了一幅高分影像經(jīng)非極大值抑制后的結(jié)果。最后采用雙閾值法[24-25]對(duì)非極大值抑制后的邊緣圖像進(jìn)行弱邊緣去除,得到主要邊緣。由于圖像中通常含有噪聲,雖然噪聲形成的邊緣像素個(gè)數(shù)較少,但會(huì)在部分邊緣產(chǎn)生毛刺現(xiàn)象,因此需要采用數(shù)學(xué)形態(tài)學(xué)的方法對(duì)其進(jìn)行一定的后處理[26],得到最終的圖像邊緣。
(a) 高分影像(b) 非極大值抑制結(jié)果
圖6高分影像與非極大值抑制后的結(jié)果
Fig.6HSRimageandresultofnon-maximumsuppression
本文方法具體步驟如下: ①對(duì)每個(gè)給定的方向向量d=[a,b],將輸入圖像進(jìn)行基于格的分解,獲取位于給定方向上的一維線集合; ②對(duì)第①步中每一條線進(jìn)行小波變換,使用Haar小波高通濾波器,將變換后的系數(shù)恢復(fù)為圖像矩陣形式,便于計(jì)算,從而完成方向性小波變換,得到不同方向小波變換的高頻子帶; ③重新定義小波系數(shù)模值的方法計(jì)算每個(gè)像素位置模和梯度方向,用非極大值抑制確定候選邊緣點(diǎn),分別采用不同方法處理2個(gè)梯度方向是否相鄰時(shí)的情況; ④采用雙閾值法得到邊緣結(jié)果; ⑤使用數(shù)學(xué)形態(tài)學(xué)方法對(duì)細(xì)節(jié)進(jìn)行適當(dāng)后處理。多方向小波變換的高分影像邊緣提取具體流程如圖7所示。
圖7 多方向小波變換的高分影像邊緣提取流程Fig.7 HSR image edge extraction flow chart based on multi-direction wavelet transform
本實(shí)驗(yàn)在Intel Core i5 2.7 GHz,4 GB內(nèi)存的計(jì)算機(jī)上進(jìn)行,編程環(huán)境為Matlab2014a。數(shù)據(jù)選擇3個(gè)地區(qū)空間分辨率為0.5 m的彩色航空影像,由于光譜特征對(duì)圖像邊緣提取影響不大,因此將影像預(yù)先轉(zhuǎn)為灰度影像,選擇地物主體類型為建筑物和道路。大小均為256像素×256像素。選擇{[1,0],[1,1],[0,1],[-1,1]}4個(gè)方向進(jìn)行方向性小波變換。將與基于Canny算子的邊緣提取結(jié)果和基于標(biāo)準(zhǔn)小波變換的邊緣提取結(jié)果進(jìn)行對(duì)比分析。
3幅實(shí)驗(yàn)數(shù)據(jù)的Canny算子邊緣提取結(jié)果、小波變換邊緣提取結(jié)果及本文方法邊緣提取結(jié)果如圖8所示。3幅實(shí)驗(yàn)數(shù)據(jù)Canny算子邊緣提取的閾值參數(shù)分別選擇0.2,0.23和0.23。
(a) 實(shí)驗(yàn)1高分影像(b) 實(shí)驗(yàn)1Canny算子結(jié)果(c) 實(shí)驗(yàn)1小波變換結(jié)果(d) 實(shí)驗(yàn)1本文方法結(jié)果
(e) 實(shí)驗(yàn)2高分影像(f) 實(shí)驗(yàn)2Canny算子結(jié)果(g) 實(shí)驗(yàn)2小波變換結(jié)果(h) 實(shí)驗(yàn)2本文方法結(jié)果
(i) 實(shí)驗(yàn)3高分影像(j) 實(shí)驗(yàn)3Canny算子結(jié)果(k) 實(shí)驗(yàn)3小波變換結(jié)果(l) 實(shí)驗(yàn)3本文方法結(jié)果
圖8實(shí)驗(yàn)結(jié)果
Fig.8Resultsofexperiment
由圖8(a)—(h)可以看出,由于這2類影像的幾何結(jié)構(gòu)較為簡(jiǎn)單,主要為規(guī)則線狀結(jié)構(gòu),邊緣與內(nèi)容部分灰度值差異明顯,3種方法均能識(shí)別主體邊緣,但前2種方法結(jié)果在某些較弱邊緣線處都發(fā)生斷裂現(xiàn)象,而本文方法則有效地避免了這種情況。而對(duì)于結(jié)構(gòu)較為復(fù)雜的高分圖像,如圖8(i),本文方法明顯要優(yōu)于Canny算子,從圖8(l)中可以看出,本文方法獲得的結(jié)果包含更多更完整的線狀邊緣。雖然傳統(tǒng)方法中選擇較小的閾值時(shí)可以獲得更多更完整的邊緣,但是由草地和樹(shù)冠等形成的紋理區(qū)域也包含大量的高頻信息,容易形成偽邊緣,因此結(jié)果很難同時(shí)兼顧邊緣完整性與去除偽邊緣。本文方法在控制由于光照和紋理區(qū)的噪聲和偽邊緣的同時(shí),能夠獲得更加完整的邊緣,這是因?yàn)楸疚姆椒ǜ玫乩昧藞D像本身的幾何結(jié)構(gòu),有針對(duì)性地檢測(cè)圖像的奇異信息,而噪聲和光照不均勻造成的偽邊緣通常不具備邊緣的幾何結(jié)構(gòu)。
對(duì)比采用小波變換的邊緣提取結(jié)果可以明顯看出,本文方法優(yōu)于標(biāo)準(zhǔn)小波方法,因?yàn)闃?biāo)準(zhǔn)小波變換只能檢測(cè)水平和垂直方向奇異點(diǎn),對(duì)于非垂直與水平的弱邊緣小波系數(shù)被分解到垂直與水平子帶,沒(méi)有充分利用圖像的幾何結(jié)構(gòu)信息。而基于多方向性小波變換可以檢測(cè)任意方向的奇異點(diǎn),從中選取較大系數(shù),有利于獲取豐富連續(xù)的圖像邊緣。
但是在多次的實(shí)驗(yàn)結(jié)果對(duì)比分析中,本文方法也存在一些缺陷。例如對(duì)于某些較弱的邊緣或者較小地物的邊緣提取結(jié)果中發(fā)現(xiàn)仍有部分像素點(diǎn)斷裂的情況??赡艽嬖诘脑蛟谟诓糠珠撝翟O(shè)置不恰當(dāng),或未考慮多尺度的方法。某種程度上地物的邊緣也可以視為是多尺度的統(tǒng)一,使用單一或區(qū)間閾值時(shí)能夠得到某個(gè)尺度下較好的邊緣結(jié)果,但在多尺度下提取結(jié)果會(huì)不可避免地出現(xiàn)問(wèn)題。
本文提出一種圖像邊緣提取方法,對(duì)高空間分辨率遙感影像進(jìn)行邊緣信息提取。首先分析了傳統(tǒng)邊緣提取方法和小波邊緣提取方法存在的不足之處,再根據(jù)高分遙感影像的自身幾何結(jié)構(gòu)特性提出多方向性小波變換邊緣提取方法,兼顧了一維小波變換計(jì)算簡(jiǎn)單性和二維多尺度幾何分析多方向性。對(duì)比實(shí)驗(yàn)結(jié)果表明,本文方法相對(duì)于傳統(tǒng)方法能夠更好地提取高分影像邊緣,具有應(yīng)用有效性。
小波變換是一種多尺度的表達(dá)工具,高分影像空間分辨率較高,地物內(nèi)部細(xì)節(jié)更加詳細(xì),使用單一尺度的邊緣提取過(guò)于籠統(tǒng)。在應(yīng)用中,有時(shí)需要忽略較小地物目標(biāo)而更多關(guān)注整體的輪廓信息,有時(shí)則不然。因此,充分利用高分影像的優(yōu)勢(shì),在不同尺度空間下提取圖像邊緣信息的更有效方法是今后仍需繼續(xù)研究的問(wèn)題。