李瑩瑩，吳 昊，俞 雷，丁一帆，鮑 莉，李 豐

(北京市遙感信息研究所，北京 100192)

一、引 言

近年來可見光和SAR遙感已成為獲取全球空間信息的主要手段，二者分別反映了目標(biāo)區(qū)域的不同層次特性，互補(bǔ)性較強(qiáng)，因此其信息融合應(yīng)用也受到越來越多的關(guān)注，而圖像融合處理的基礎(chǔ)就是利用圖像間同名點(diǎn)的對(duì)應(yīng)關(guān)系進(jìn)行圖像配準(zhǔn)。

多源圖像自動(dòng)匹配技術(shù)[1-2]發(fā)展尚不完善，對(duì)于SAR和光學(xué)異源圖像之間的匹配，由于成像機(jī)理不同及SAR圖像斑點(diǎn)噪聲的影響，通常很難達(dá)到像素級(jí)的配準(zhǔn)精度，因此研究一種高精度、高效率的同名點(diǎn)自動(dòng)匹配技術(shù)顯得格外重要。目前自動(dòng)匹配算法歸納起來主要分為：①基于灰度的方法[3-4]直接利用圖像的灰度信息，對(duì)輻射特性比較敏感，不太適合輻射失真較大的情況，更不適合SAR和光學(xué)兩種不同傳感器獲得的圖像之間的配準(zhǔn)。②基于特征的方法[5-6]提取圖像的幾何特征作為配準(zhǔn)基元，克服了灰度方法的缺點(diǎn)，已成為主要發(fā)展方向，但是計(jì)算量較大且特征提取結(jié)果與圖像質(zhì)量密切相關(guān)。近年來還出現(xiàn)了一些復(fù)合的圖像配準(zhǔn)方法，主要有兩種類型：一是結(jié)合灰度和特征的圖像配準(zhǔn)方法，二是結(jié)合多種特征的圖像配準(zhǔn)方法。Zhang等[7]提出了一種組合邊緣特征和光流估計(jì)的配準(zhǔn)方法；文獻(xiàn)[8]利用標(biāo)記特征和灰度信息配準(zhǔn)醫(yī)學(xué)圖像；Dare[9]利用點(diǎn)、區(qū)域、邊緣三種特征和由粗到精的策略配準(zhǔn)SPOT 和SAR 圖像；彭文[10]分別利用點(diǎn)和輪廓特征、點(diǎn)和曲線特征進(jìn)行圖像配準(zhǔn)。

以上文獻(xiàn)中涉及的匹配算法往往較傳統(tǒng)方法耗時(shí)長,且算法對(duì)如何高效獲取高精度的同名點(diǎn)對(duì)等流程也未能作詳細(xì)論述。鑒于此，本文給出了一種基于不同邊緣檢測(cè)算子的光學(xué)和SAR圖像自動(dòng)匹配方法，系統(tǒng)地闡述了仿射變換粗配準(zhǔn)、特征點(diǎn)提取及同名點(diǎn)精匹配的全過程，并進(jìn)行了匹配精度的目視和定量檢驗(yàn)。試驗(yàn)結(jié)果證明，該算法不僅處理速度快，而且提取的同名點(diǎn)對(duì)精度高，分布均勻。

二、算法流程

本文采用從粗到精的多級(jí)匹配策略，首先利用仿射變換進(jìn)行光學(xué)和SAR圖像的粗配準(zhǔn)，解決兩幅圖像之間的旋轉(zhuǎn)和尺度不一致問題；然后采用Harris算子提取光學(xué)圖像的特征點(diǎn)，以這些特征點(diǎn)的區(qū)域子圖像作為匹配模板，提取其Sobel特征圖像，在SAR粗匹配圖像對(duì)應(yīng)的搜索范圍內(nèi)用ROA邊緣檢測(cè)算法提取SAR特征圖像；最后聯(lián)合兩幅圖像基于歸一化互相關(guān)原理進(jìn)行精匹配，確定出SAR圖像的同名點(diǎn)。具體實(shí)現(xiàn)如圖1所示。

1. 利用仿射變換進(jìn)行原始圖像粗匹配

由于光學(xué)和SAR圖像存在較大的分辨率和視角差異，在進(jìn)行同名點(diǎn)精匹配之前，需要解決原始圖像間的尺度和旋轉(zhuǎn)問題，具體過程如圖2所示。

1) 根據(jù)光學(xué)圖像的位置，依賴四角位置的定位精度并考慮粗匹配的精度，確定出需要裁切區(qū)域的對(duì)應(yīng)4個(gè)角點(diǎn)位置A1、A2、A3和A4。

2) SAR原始圖像的4個(gè)角點(diǎn)N1、N2、N3、N4位置已知，根據(jù)衛(wèi)星軌道參數(shù)和均勻分布的少量控制點(diǎn)，就可以利用仿射變換計(jì)算出A1～A44個(gè)點(diǎn)對(duì)應(yīng)于SAR圖像上的4個(gè)位置B1、B2、B3、B4。

3) 采用三次卷積灰度插值對(duì)B1B2B3B4區(qū)域的SAR圖像按照光學(xué)圖像尺度(分辨率)進(jìn)行重采樣。

圖1 同名點(diǎn)自動(dòng)匹配步驟

圖2 基于仿射變換的圖像粗匹配

最終得到SAR圖像的粗匹配圖像，在裁切時(shí)確保了該區(qū)域包括參考光學(xué)圖像區(qū)域。使用粗匹配系數(shù)記錄上述從SAR原始圖像到匹配圖像的點(diǎn)對(duì)點(diǎn)變換流程，用于精匹配步驟中通過反變換確定同名點(diǎn)在SAR原始圖像上的候選區(qū)域。

2. 特征點(diǎn)自動(dòng)提取

提取圖像特征點(diǎn)的方法有不少，如Moravec、Forstner等算子。Harris算子是由C．Harris提出的一種點(diǎn)特征提取算子，它采用了信號(hào)處理中的自相關(guān)函數(shù)矩陣，其特征值是自相關(guān)函數(shù)的一階曲率，如果兩個(gè)曲率值都高，那么就認(rèn)為該點(diǎn)是特征點(diǎn)。其基本思想是從局部小窗口中觀測(cè)圖像特征，因?yàn)榻屈c(diǎn)定義為包括灰度值變化明顯的點(diǎn)，所以窗口的移動(dòng)會(huì)導(dǎo)致灰度值出現(xiàn)明顯的變化。

首先對(duì)灰度圖像的每個(gè)點(diǎn)，分別計(jì)算橫向和縱向的一階導(dǎo)數(shù)及二者的乘積，然后計(jì)算其興趣值。Harris算法要提取的特征點(diǎn)就是局部范圍內(nèi)的極大興趣值所對(duì)應(yīng)的像素點(diǎn)。Harris算子是最穩(wěn)定的點(diǎn)特征提取算子之一，它的計(jì)算公式中只涉及灰度的一階導(dǎo)數(shù)，因此不需設(shè)置閾值；提取的點(diǎn)特征均勻而且合理，因?yàn)樗鼘?duì)圖像的每個(gè)點(diǎn)都計(jì)算興趣值，所以在鄰域中選擇最優(yōu)。

3. 同名點(diǎn)自動(dòng)精匹配

欲對(duì)兩幅圖像進(jìn)行匹配，首先要對(duì)兩幅圖像進(jìn)行特征描述。SAR圖像與可見光圖像對(duì)同一地面特征的像素灰度不具有高度的相關(guān)性，其主要區(qū)別是，在SAR圖像中道路顯示為黑色，而在可見光圖像上則呈白色；SAR圖像上對(duì)建筑物顯示非常清楚，而在可見光圖像上則可能有陰影效應(yīng)；SAR圖像有明顯斑點(diǎn)噪聲，而可見光圖像有豐富的表面細(xì)節(jié)特征等。只有道路、建筑物、植被區(qū)等景物的邊緣輪廓特征是兩種圖像的主要相似點(diǎn)。

(1) 光學(xué)Sobel梯度圖像

圖像的梯度強(qiáng)度正反映了邊緣和輪廓的情況，把圖像看成二維離散函數(shù)，梯度其實(shí)就是對(duì)這個(gè)二維離散函數(shù)的求導(dǎo)。即通過計(jì)算圖像鄰域的微分獲得灰度變化率，相當(dāng)于對(duì)圖像進(jìn)行銳化，從而增強(qiáng)圖像的灰度反差。同時(shí)由于梯度強(qiáng)度有一定的相似性，因此在匹配測(cè)度上采用互相關(guān)優(yōu)勢(shì)，最終形成了基于圖像梯度和歸一化互相關(guān)的同名點(diǎn)精匹配方法。

梯度定義為

振幅為

計(jì)算梯度的關(guān)鍵問題是如何數(shù)字化地估計(jì)Gx和Gy，有幾種經(jīng)典的估計(jì)方法，如Sobel、Roberts、Prewitt。本文采用Sobel方法提供了兩個(gè)3×3矩陣

Gx=(z7+2z8+z9)-(z1+2z2+z3)

Gy=(z3+2z6+z9)-(z1+2z4+z7)

其中,z1～z9為圖像上特征點(diǎn)周圍3×3大小鄰域內(nèi)的灰度值。將上述定義的兩個(gè)矩陣Gx和Gy代入振幅計(jì)算公式近似計(jì)算出光學(xué)圖像梯度。

(2) SAR ROA特征圖像

由于SAR圖像具有典型的乘性斑點(diǎn)噪聲，常規(guī)的梯度算子很難反映圖像的邊緣輪廓信息,因?yàn)榇祟惙椒ɑ诰植科骄抵頪6]，它們的虛警率隨著被檢測(cè)區(qū)域的輻射強(qiáng)度的增加而增加,閾值也隨著平均輻射的變化而變化,這就難以得到應(yīng)用。圖3顯示了對(duì)粗配準(zhǔn)后的SAR圖像利用這些傳統(tǒng)的光學(xué)圖像梯度算子提取邊緣，可見輪廓不清晰，效果都比較差，這些方法并不適用于SAR圖像。

圖3 傳統(tǒng)光學(xué)圖像邊緣算子提取SAR圖像特征

邊緣提取的質(zhì)量直接影響著同名點(diǎn)匹配的精度，因此對(duì)特征算子的要求是定位精確可靠，對(duì)噪聲不敏感，在對(duì)實(shí)時(shí)性要求高的場(chǎng)合還要有較小的計(jì)算量。通過對(duì)幾種檢測(cè)效果的比較，本文采用基于比值邊緣檢測(cè)算法ROA[11]，它可以較好地克服斑點(diǎn)噪聲，較準(zhǔn)確地反映SAR圖像上的邊緣信息。

計(jì)算過程根據(jù)圖4所示的模板進(jìn)行，模板被分為3個(gè)區(qū)域。對(duì)于4個(gè)常用方向：水平、垂直、左傾、右傾，以感興趣像素為中心，大小為w×w的滑動(dòng)窗口被其分成相鄰且無重疊的兩個(gè)區(qū)域。某個(gè)像素的第k對(duì)區(qū)域邊緣強(qiáng)度定義為

Rk=min(Pk/Qk,Qk/Pk)k=1,2,3,4

(2)

式中,Pk和Qk分別是當(dāng)前像素第k對(duì)邊側(cè)區(qū)域像素的平均灰度值。采用4個(gè)方向的最值作為當(dāng)前像素的梯度強(qiáng)度值，計(jì)算公式為

R=min(R1,R2,R3,R4)

(3)

圖4 ROA計(jì)算原理

(3) 基于Sobel和ROA特征的互相關(guān)匹配

提取同名特征點(diǎn)時(shí)，首先以光學(xué)圖像所提取的特征點(diǎn)位為中心，截取鄰近一定范圍的光學(xué)Sobel梯度圖像(如64×64)定義為匹配窗口模板；然后根據(jù)獲取的粗匹配系數(shù)計(jì)算該特征點(diǎn)在SAR原始圖像上所對(duì)應(yīng)的大概位置，并將該位置鄰近稍大范圍的ROA梯度圖像定義為搜索窗口；最后使用匹配窗口在搜索窗口內(nèi)滑動(dòng)，并逐點(diǎn)計(jì)算歸一化互相關(guān)系數(shù)，相關(guān)系數(shù)最大的位置就被認(rèn)為是該光學(xué)特征點(diǎn)在SAR圖像上對(duì)應(yīng)的同名點(diǎn)。

由于兩幅圖像是不同載荷在不同時(shí)間獲取的，相同地區(qū)的圖像差異仍然較大，導(dǎo)致相關(guān)系數(shù)不高，因此即使存在相關(guān)系數(shù)最大的位置，也并不能確保一定是同名點(diǎn)。在之前的研究工作中，經(jīng)過多次數(shù)據(jù)對(duì)的試驗(yàn)，發(fā)現(xiàn)最大相關(guān)系數(shù)一般優(yōu)于0.2才能確保是同名點(diǎn)，而小于0.2就可能顯示兩幅圖像在此區(qū)域內(nèi)不存在同名點(diǎn)。筆者統(tǒng)計(jì)了最大相關(guān)系數(shù)的直方圖分布情況，可以得到相關(guān)系數(shù)大部分位于0.2～1.0之間，而在0～0.1之間出現(xiàn)高值尖峰是由于地物變化或成像差異而沒有同名點(diǎn)導(dǎo)致。

三、試驗(yàn)結(jié)果與分析

本文試驗(yàn)采用TerraSar衛(wèi)星和SPOT光學(xué)衛(wèi)星圖像數(shù)據(jù),其空間分辨率分別為3 m和1 m，成像時(shí)間分別為2009年11月9日和2009年3月13日。圖幅覆蓋范圍為北京城郊地區(qū)，影像中顯示局部區(qū)域地勢(shì)錯(cuò)落有致，有一定范圍的丘陵分布，包含了大量的道路、建筑、溝渠等地物。兩幅圖像的粗匹配和光學(xué)圖像的特征點(diǎn)提取過程不再贅述，這里重點(diǎn)給出同名點(diǎn)自動(dòng)精匹配過程的試驗(yàn)結(jié)果。

圖5給出了對(duì)同一地區(qū)的光學(xué)和SAR圖像分別用Sobel算子和ROA算子進(jìn)行邊緣提取的結(jié)果。經(jīng)過反復(fù)試驗(yàn)發(fā)現(xiàn)，ROA滑動(dòng)窗口的大小設(shè)置為9像素×9像素時(shí)邊緣提取的效果較好。與圖3其他幾種算子比較可以看出，ROA算子提取的邊界清楚、細(xì)致，且與光學(xué)圖像輪廓非常相似。雖然得到的圖像灰度值較低，同時(shí)損失了一些邊緣細(xì)節(jié)，但這對(duì)特征點(diǎn)的選取及配準(zhǔn)計(jì)算都是不重要的。

圖5 兩種算子提取結(jié)果

本文試驗(yàn)的運(yùn)行環(huán)境為Intel Xeon W3520主頻2.67 GHz。光學(xué)特征點(diǎn)模板的大小設(shè)置為128×128，SAR圖像的候選區(qū)域大小設(shè)置為192×192，共測(cè)試了15個(gè)點(diǎn)的平均匹配時(shí)間為1～2 s左右。實(shí)際上獲取的特征點(diǎn)數(shù)隨試驗(yàn)區(qū)域的大小和地物密集程度的不同而改變，1 m分辨率200 km2范圍單幅光學(xué)圖像的點(diǎn)數(shù)通常在幾千個(gè)左右，因此整個(gè)同名點(diǎn)匹配過程占用的時(shí)間基本在1 h以內(nèi)。

從一維上分別觀察兩幅原始圖像和特征提取圖像的某一行如圖6所示，可以看出，在原始灰度圖像上除了在150附近的兩個(gè)峰值有重合跡象外，其他部分都是混亂的。而在經(jīng)過邊緣提取后，得到相同行的邊緣圖趨勢(shì)基本一致。

圖6 SAR和光學(xué)圖像一維剖面灰度對(duì)比和梯度對(duì)比

為了檢驗(yàn)同名點(diǎn)提取的精度，筆者拼合光學(xué)和SAR圖像包含特征點(diǎn)的相鄰子塊得到一幅新圖像(如圖7所示)，進(jìn)行定性的目視檢驗(yàn)。圖中深淺程度不同的方塊分別是SAR和光學(xué)圖像子塊，相鄰圖像塊接縫處的地物自然程度和平滑程度就顯示了兩幅圖像的匹配精度?？梢钥闯?，各種地物(河流、道路、建筑和田塊)在子塊交叉處連接十分自然平滑，沒有明顯的拼接錯(cuò)位現(xiàn)象，表明匹配精度是很高的。

圖7 分別提取SAR和光學(xué)子塊交錯(cuò)拼合的圖像

此外筆者進(jìn)行了匹配精度的定量檢查，共檢查了30個(gè)同名特征點(diǎn)，這些點(diǎn)一般都是位于道路河流的交叉點(diǎn)或地物角點(diǎn)上，量取的像素坐標(biāo)數(shù)值見表1。已知光學(xué)參考圖像上的特征點(diǎn)，表中第1列為精匹配獲得的SAR圖像同名點(diǎn)坐標(biāo)，第2列為真實(shí)對(duì)應(yīng)的同名點(diǎn)坐標(biāo)。根據(jù)這些像素坐標(biāo)計(jì)算出坐標(biāo)差，進(jìn)而可以計(jì)算圖像的匹配精度。經(jīng)過統(tǒng)計(jì)，兩幅圖像的匹配精度達(dá)到了0.863 8個(gè)像素。

四、結(jié)束語

光學(xué)和SAR圖像的灰度、紋理等特征都有很大的不同，而邊緣是這兩種類型圖像中的共有特征，因此基于邊緣的圖像匹配算法可以較好地解決這個(gè)問題。此外由粗到精的分級(jí)搜索策略符合人類視覺系統(tǒng)選擇性搜集信息的機(jī)理，且能節(jié)省系統(tǒng)的計(jì)算開銷，是一種有應(yīng)用前途的匹配方法。因此本文提出的算法不僅解決了大圖幅影像同名點(diǎn)搜索的計(jì)算復(fù)雜度問題，而且獲取的同名點(diǎn)對(duì)分布更趨均勻，精度更高。試驗(yàn)結(jié)果也有力地證明了這一點(diǎn)，在地物特征較為明顯地區(qū)的精度可達(dá)到1個(gè)像素以內(nèi)，每個(gè)同名點(diǎn)對(duì)的平均匹配時(shí)間為1～2 s左右。

致謝：特別感謝中國科學(xué)院電子學(xué)研究所尤紅建研究員在本文研究過程中的幫助和指導(dǎo)。

表1 同名點(diǎn)匹配精度檢查 像素

參考文獻(xiàn)：

[1] ZITOVA B, FLUSSER J. Image Registration Methods: A Survey[J]. Image and Vision Computing, 2003, 21(11): 977-1000.

[2] LEILA M, MANJUNATH B. Registration Techniques for Multisensor Remotely Sensed Imagery [J]. PE&RS, 1996, 62(9): 1049-1052.

[3] TON J. Registering Landsat Images by Point Matching [J]. IEEE Transactions on Geoscience and Remote Sensing, 1989, 27(5): 642-652.

[4] 王小睿, 吳信才. 遙感多圖像的自動(dòng)配準(zhǔn)方法[J]. 中國圖象圖形學(xué)報(bào), 1997, 2(10): 735-739.

[5] DAI X, KHORRAM S. A Feature-based Image Registration Algorithm using Improved Chain Code Representation Combined with Invariant Moments [J]. IEEE Transactions on Geoscience and Remote Sensing, 1999, 37(5): 2351-2362.

[6] 王磊, 張鈞萍, 張曄. 基于特征的SAR圖像與光學(xué)圖像自動(dòng)配準(zhǔn)[J]. 哈爾濱工業(yè)大學(xué)學(xué)報(bào), 2005, 37(1): 22-25.

[7] ZHANG Z, BLUM R. A Hybrid Image Registration Technique for Digital Camera Image Fusion Application [J]. Information Fusion, 2000(2): 135-149.

[8] JOHNSON H, CHRISTENSEN G. Consistent Landmark and Intensity-based Image Registration [J]. IEEE Transactions on Medical Imaging, 2002, 21(5): 250-461.

[9] PAUL D. An Improved Model for Automatic Feature-based Registration of SAR and SPOT Images [J]. ISPRS Journal of Photogrammetry & Remote Sensing, 2001, 56(1): 13-28.

[10] 彭文. 基于特征的醫(yī)學(xué)圖像配準(zhǔn)中若干關(guān)鍵技術(shù)的研究[D]. 杭州： 浙江大學(xué), 2007.

[11] TOUZI R, LOPES A, BOUSQUET P. A Statistical and Geometrical Edge Detector for SAR Images[J]. IEEE Transactions on Geoscience and emote Sensing, 1988, 26(6): 764-773.