楊 毅，李大成，于 杰，韓啟金

(1.太原理工大學(xué) 光電工程學(xué)院，山西 太原 030024; 2.太原理工大學(xué) 礦業(yè)工程學(xué)院，山西 太原 030024; 3.中國資源衛(wèi)星應(yīng)用中心，北京 100094)

0 引言

圖像配準(zhǔn)是通過一定的方法獲得基準(zhǔn)圖像與待配準(zhǔn)圖像之間的映射關(guān)系實(shí)現(xiàn)圖像的幾何校準(zhǔn)，是諸多遙感應(yīng)用(如圖像融合[1-2]和變化檢測[3-4]等)的關(guān)鍵前置環(huán)節(jié)，對系統(tǒng)性能有著舉足輕重的影響。圖像配準(zhǔn)方法大致可以分為2大類：基于區(qū)域的配準(zhǔn)方法和基于特征的配準(zhǔn)方法。其中，基于區(qū)域的配準(zhǔn)方法[5-8]一般采用模板匹配策略，利用圖像灰度信息尋找基準(zhǔn)圖像和待配準(zhǔn)圖像中的同名點(diǎn)實(shí)現(xiàn)圖像配準(zhǔn)。這類方法存在計(jì)算量大，對圖像旋轉(zhuǎn)、縮放等仿射變換魯棒性較差的問題?；谔卣鞯呐錅?zhǔn)方法首先在基準(zhǔn)圖像和待配準(zhǔn)圖像上提取局部特征點(diǎn)，之后使用圖像局部特征描述方法對局部特征進(jìn)行描述，接著根據(jù)一定的評價(jià)準(zhǔn)則如歐拉距離尋找同名點(diǎn)，最后根據(jù)同名點(diǎn)計(jì)算基準(zhǔn)圖像與待配準(zhǔn)圖像之間的映射關(guān)系實(shí)現(xiàn)圖像配準(zhǔn)。與基于區(qū)域的配準(zhǔn)方法相比，這類方法的計(jì)算量更小，對圖像旋轉(zhuǎn)、縮放等仿射變換有更好的魯棒性。例如，Lowe[9]于1999年提出并于2004年完善的尺度不變特征變換(Scale-Invariant Feature Transform，SIFT)算法，該算法通過構(gòu)建高斯差分金字塔提取多尺度圖像局部特征點(diǎn)，建立局部特征點(diǎn)支持區(qū)域的三維梯度直方圖完成對圖像局部特征的描述，使用歐拉距離作為準(zhǔn)則尋找同名點(diǎn)。由于SIFT在提取圖像局部特征時(shí)考慮了多尺度問題，因此在自然圖像配準(zhǔn)場景中(如圖像縮放、旋轉(zhuǎn))表現(xiàn)出了良好的性能和魯棒性。之后出現(xiàn)了很多SIFT的變種[10-13]，如Mikolajczyk等[10]使用對數(shù)極坐標(biāo)系替換SIFT算法中使用的笛卡爾坐標(biāo)系，進(jìn)一步提升了圖像局部特征描述方法的魯棒性。Bay等[11]使用積分圖像和箱式濾波器提高了算法的運(yùn)行速度和配準(zhǔn)精度。

由于SIFT算法在自然圖像配準(zhǔn)中展示出優(yōu)良的性能，開始被引入遙感圖像配準(zhǔn)領(lǐng)域，衍生出多個(gè)適用于遙感圖像配準(zhǔn)場景的變種[14-18]。如Dellinger等[15]為降低合成孔徑雷達(dá)(Synthetic Aperture Radar，SAR)圖像斑點(diǎn)噪聲的影響，在SIFT算法基礎(chǔ)上改變梯度計(jì)算方法提出了SAR-SIFT用于SAR圖像配準(zhǔn)。Sedaghat等[16]于2011年在SIFT的基礎(chǔ)上提出了UR-SIFT算法，該算法通過改變SIFT提取特征點(diǎn)的策略，使特征點(diǎn)在圖像中的分布更加均勻，提高了遙感圖像配準(zhǔn)的精度。

遙感圖像中全色圖像通常是覆蓋整個(gè)可見光波段(380～780 nm)的混合圖像，多光譜圖像是對多個(gè)單波段單獨(dú)成像，獲取相應(yīng)波段的光譜信息。因此，全色圖像與多光譜圖像之間輻射差異并不是線性關(guān)系，造成相似的圖像局部特征點(diǎn)在不同波段圖像中占據(jù)比例和可重復(fù)性降低，使得圖像局部特征的誤匹配率上升，最終導(dǎo)致圖像配準(zhǔn)質(zhì)量和精度下降。通過SIFT算法進(jìn)行全色圖像與多光譜圖像配準(zhǔn)時(shí)，由于不同波段遙感圖像之間存在復(fù)雜輻射差異，使用SIFT算法進(jìn)行圖像配準(zhǔn)時(shí)無法保證配準(zhǔn)質(zhì)量。因此，本文提出一種基于迭代更新策略的SIFT算法變種，該算法通過遍歷基準(zhǔn)圖像與待配準(zhǔn)圖像之間的同名點(diǎn)對，迭代求解基準(zhǔn)圖像與待配準(zhǔn)圖像之間的映射關(guān)系(即單應(yīng)性矩陣)，得到最優(yōu)單應(yīng)性矩陣，對待配準(zhǔn)圖像進(jìn)行重采樣完成圖像配準(zhǔn)。

1 本文方法

基于特征的遙感圖像配準(zhǔn)方法實(shí)施圖像配準(zhǔn)的一般流程是：① 在基準(zhǔn)圖像和待配準(zhǔn)圖像上檢測并提取圖像局部特征；② 通過圖像局部特征描述方法對圖像局部特征進(jìn)行描述，使用一定準(zhǔn)則尋找同名點(diǎn)；③ 根據(jù)同名點(diǎn)計(jì)算得出基準(zhǔn)圖像與待配準(zhǔn)圖像之間的映射關(guān)系，對待配準(zhǔn)圖像進(jìn)行重采樣完成圖像配準(zhǔn)。其中，通過基準(zhǔn)圖像與待配準(zhǔn)圖像間的同名點(diǎn)對計(jì)算單應(yīng)性矩陣時(shí)，同名點(diǎn)對會直接影響單應(yīng)性矩陣求解的質(zhì)量，進(jìn)而影響圖像配準(zhǔn)精度。本文提出一種迭代更新策略的SIFT算法，該算法的核心是通過遍歷基準(zhǔn)圖像與待配準(zhǔn)圖像之間的同名點(diǎn)對，迭代求解單應(yīng)性矩陣，得到最優(yōu)單應(yīng)性矩陣，根據(jù)單應(yīng)性矩陣對待配準(zhǔn)圖像進(jìn)行重采樣完成圖像配準(zhǔn)，算法流程如圖1所示。

圖1 本文算法流程Fig.1 Flowchart of the proposed method

2 單應(yīng)性矩陣質(zhì)量評價(jià)

根據(jù)基準(zhǔn)圖像和待配準(zhǔn)圖像之間同名點(diǎn)對計(jì)算獲得的單應(yīng)性矩陣的質(zhì)量直接影響圖像配準(zhǔn)的精度和效果，本文采用同名點(diǎn)對的RMSE作為單應(yīng)性矩陣質(zhì)量評價(jià)指標(biāo)。為計(jì)算同名點(diǎn)的RMSE，首先使用單應(yīng)性矩陣H計(jì)算待配準(zhǔn)圖像中的局部特征點(diǎn)xw和yw在基準(zhǔn)圖像中的坐標(biāo)xb′=H·xw和yb′=H·yw，接著根據(jù)式(1)和式(2)分別計(jì)算其與SIFT算法確定的同名點(diǎn)坐標(biāo)在xb和yb方向的殘差Vx和Vy。之后根據(jù)式(3)和式(4)計(jì)算xb和yb方向的均方根誤差mx和my，其中n代表特征點(diǎn)的數(shù)量。最后根據(jù)式(5)計(jì)算得到同名點(diǎn)之間的RMSE(x,y)。同名點(diǎn)集中的每對同名點(diǎn)均通過以上步驟獲得其RMSE，RMSE值越小單應(yīng)性矩陣質(zhì)量越好。

(1)

(2)

(3)

(4)

(5)

3 實(shí)驗(yàn)及結(jié)果分析

為驗(yàn)證所提算法的有效性和魯棒性，在網(wǎng)通一號衛(wèi)星和高分二號某地區(qū)的遙感影像數(shù)據(jù)上進(jìn)行了圖像配準(zhǔn)實(shí)驗(yàn)。用于圖像配準(zhǔn)實(shí)驗(yàn)的網(wǎng)通遙感圖像如圖2所示。

圖2 網(wǎng)通衛(wèi)星圖像Fig.2 Images of WT satellite

這一景圖像共有8個(gè)波段，其中波段0的影像為12 000 pixel×12 000 pixel的全色圖像，其余波段影像為3 000 pixel×3 000 pixel的多光譜圖像。

高分二號某地區(qū)遙感圖像如圖3所示，其中全色圖像尺寸為29 200 pixel×27 620 pixel，多光譜圖像尺寸為7 300 pixel×6 908 pixel，限于篇幅本文取全色圖像及其多光譜圖像的波段1～3用于圖像配準(zhǔn)實(shí)驗(yàn)。從圖2和圖3所示的遙感圖像可以看出，全色圖像和多光譜圖像之間不僅僅是分辨率不同，還存在復(fù)雜的輻射差異，造成不同波段的圖像存在明顯且非線性的灰度變化。

圖3 高分二號衛(wèi)星圖像Fig.3 Images of GF-2 satellite

本文在Windows 10系統(tǒng)上基于QT64位開發(fā)平臺，結(jié)合GDAL(Geospatial Data Abstraction Library)和OpenCV開源庫搭建了遙感圖像配準(zhǔn)實(shí)驗(yàn)環(huán)境。由于遙感圖像的位深是16 bit，在進(jìn)行圖像配準(zhǔn)實(shí)驗(yàn)前采用最大值最小值拉伸算法將遙感圖像位深轉(zhuǎn)換為8 bit以方便之后實(shí)驗(yàn)。

首先驗(yàn)證基于迭代更新策略的SIFT算法獲取基準(zhǔn)圖像與待配準(zhǔn)圖像之間最優(yōu)映射關(guān)系的有效性。通過SIFT算法分別從基準(zhǔn)圖像和待配準(zhǔn)圖像提取1 000個(gè)特征點(diǎn)，接著使用歐拉距離作為評價(jià)準(zhǔn)則，確定基準(zhǔn)圖像和待配準(zhǔn)圖像之間的同名點(diǎn)。由于SIFT算法完成粗匹配后可能會存在誤匹配點(diǎn)，在計(jì)算基準(zhǔn)圖像與待配準(zhǔn)圖像間的單應(yīng)性矩陣前使用RANSAC算法對同名點(diǎn)進(jìn)行初步提純。首先以網(wǎng)通衛(wèi)星影像波段0作為基準(zhǔn)圖像，其他波段1～7作為待配準(zhǔn)圖像進(jìn)行實(shí)驗(yàn)。限于篇幅，僅呈現(xiàn)波段0與波段1、波段0與波段2、波段0與波段4、波段0與波段7的實(shí)驗(yàn)結(jié)果，如圖4～圖7所示。

(a)整體實(shí)驗(yàn)結(jié)果

(b)迭代次數(shù)10～74的結(jié)果圖4 波段0與波段1的RSME之和與迭代次數(shù)的實(shí)驗(yàn)結(jié)果Fig.4 Experimental results of the sum of RMSE and number of iterations for bands 0 and 1

(a)整體實(shí)驗(yàn)結(jié)果

(b)迭代次數(shù)10～216的結(jié)果圖5 波段0與波段2的RSME之和與迭代次數(shù)的實(shí)驗(yàn)結(jié)果Fig.5 Experimental results of the sum of RMSE and number of iterations for the bands 0 and 2

(a)整體實(shí)驗(yàn)結(jié)果

(b)迭代次數(shù)10～177的結(jié)果圖6 波段0與波段4的RSME之和與迭代次數(shù)的實(shí)驗(yàn)結(jié)果Fig.6 Experimental results of the sum of RMSE and number of iterations for the bands 0 and 4

(a)整體實(shí)驗(yàn)結(jié)果

(b)迭代次數(shù)10～141的結(jié)果圖7 波段0與波段7的RSME之和與迭代次數(shù)的實(shí)驗(yàn)結(jié)果Fig.7 Experimental results of the sum of RMSE and number of iterations for bands 0 and 7

接著，本文以高分二號全色圖像作為基準(zhǔn)圖像，多光譜波段1～3作為待配準(zhǔn)圖像進(jìn)行了圖像配準(zhǔn)實(shí)驗(yàn)，實(shí)驗(yàn)結(jié)果如圖8～圖10所示。

(a)整體實(shí)驗(yàn)結(jié)果

(b)迭代次數(shù)10～66的結(jié)果圖8 全色圖像與多光譜波段1的RSME之和與迭代次數(shù)的實(shí)驗(yàn)結(jié)果Fig.8 Experimental results of the sum of RMSE and number of iterations for panchromatic image and multi-spectral band 1

(a)整體實(shí)驗(yàn)結(jié)果

(b)迭代次數(shù)10～114的結(jié)果圖9 全色圖像與多光譜波段2的RSME之和與迭代次數(shù)的實(shí)驗(yàn)結(jié)果Fig.9 Experimental results of the sum of RMSE and number of iterations for panchromatic image and multi-spectral band 2

(a)整體實(shí)驗(yàn)結(jié)果

(b)迭代次數(shù)10～115的結(jié)果圖10 全色圖像與多光譜波段3的RSME之和與迭代次數(shù)的實(shí)驗(yàn)結(jié)果Fig.10 Experimental results of the sum of RMSE and number of iterations for panchromatic image and multi-spectral band 3

從圖4～圖10可以看出，迭代策略剛開始階段RMSE之和會有非常明顯的下降，說明同名點(diǎn)對的數(shù)量對求取基準(zhǔn)圖像與待配準(zhǔn)圖像之間的映射關(guān)系有非常大的影響。換句話說，雖然基準(zhǔn)圖像和待配準(zhǔn)圖像間的單應(yīng)性矩陣僅有8個(gè)自由變量，理論上4對同名點(diǎn)便可將8個(gè)自由變量全部求解出來，但此時(shí)獲得的單應(yīng)性矩陣質(zhì)量往往不太好。隨著迭代次數(shù)進(jìn)一步增加，RMSE之和下降的速度明顯變緩，但總體仍是下降，說明單應(yīng)性矩陣的質(zhì)量在進(jìn)一步提升，當(dāng)遍歷完所有同名點(diǎn)對后便可獲得RMSE之和最小時(shí)基準(zhǔn)圖像與待配準(zhǔn)圖像之間的最優(yōu)映射關(guān)系，即單應(yīng)性矩陣。

獲得基準(zhǔn)圖像與待配準(zhǔn)圖像之間的最優(yōu)單應(yīng)性矩陣后，利用單應(yīng)性矩陣對待配準(zhǔn)圖像進(jìn)行重采樣和插值運(yùn)算，生成分辨率與基準(zhǔn)圖像一致的配準(zhǔn)圖像。本部分仍將網(wǎng)通衛(wèi)星圖像的波段0作為基準(zhǔn)圖像，波段1～7作為待配準(zhǔn)圖像，通過之前部分獲得的基準(zhǔn)圖像與待配準(zhǔn)圖像間最優(yōu)映射關(guān)系對波段1～7的圖像進(jìn)行重采樣插值運(yùn)算，生成與基準(zhǔn)圖像波段0同樣尺寸的配準(zhǔn)結(jié)果圖。限于篇幅，圖11～圖14分別展示波段0-1、波段0-2、波段0-4和波段0-7的圖像配準(zhǔn)結(jié)果圖與基準(zhǔn)圖像制成的棋盤圖。

圖11 波段0與波段1的圖像配準(zhǔn)結(jié)果Fig.11 Image registration results of band 0 and band 1

圖12 波段0與波段2的圖像配準(zhǔn)結(jié)果Fig.12 Image registration results of band 0 and band 2

圖13 波段0與波段4的圖像配準(zhǔn)結(jié)果Fig.13 Image registration results of band 0 and band 4

圖14 波段0與波段7的圖像配準(zhǔn)結(jié)果Fig.14 Image registration results of band 0 and band 7

同樣，本部分仍將之前高分二號的全色圖像作為基準(zhǔn)圖像，使用得到的最優(yōu)單應(yīng)性矩陣對多光譜的波段1～波段3圖像進(jìn)行重采樣和插值計(jì)算，生成與全色圖像尺寸一致的配準(zhǔn)結(jié)果圖，配準(zhǔn)結(jié)果棋盤圖如圖15～圖17所示。從配準(zhǔn)結(jié)果棋盤圖可以看出，雖然待配準(zhǔn)圖像與基準(zhǔn)圖像之間的輻射差異不盡相同，但通過迭代更新策略的SIFT算法均能實(shí)現(xiàn)良好的圖像配準(zhǔn)，進(jìn)一步提升了SIFT算法在不同遙感圖像配準(zhǔn)場景下的魯棒性。

圖15 全色圖像與多光譜波段1的圖像配準(zhǔn)結(jié)果Fig.15 Image registration results of panchromatic image and multi-spectral band 1

綜上所述，實(shí)驗(yàn)結(jié)果表明，基于迭代更新策略的SIFT算法不僅可以得到基準(zhǔn)圖像與待配準(zhǔn)圖像之間的最優(yōu)映射關(guān)系，實(shí)現(xiàn)高精度的圖像配準(zhǔn)，還進(jìn)一步提升了SIFT算法應(yīng)對不同遙感圖像配準(zhǔn)場景的魯棒性。

4 結(jié)束語

本文提出一種基于迭代更新策略的SIFT算法用于遙感圖像配準(zhǔn)，該算法通過迭代更新求解基準(zhǔn)圖像與待配準(zhǔn)圖像之間的映射關(guān)系，即求解單應(yīng)性矩陣，最后得到RMSE之和最小的單應(yīng)性矩陣即基準(zhǔn)圖像與待配準(zhǔn)圖像間的最優(yōu)映射關(guān)系矩陣。實(shí)驗(yàn)結(jié)果表明，基于迭代更新策略的SIFT算法不僅可以得到基準(zhǔn)圖像與待配準(zhǔn)圖像之間的最優(yōu)映射關(guān)系矩陣，同時(shí)也進(jìn)一步提升了SIFT算法在不同遙感圖像配準(zhǔn)場景下的魯棒性。由于圖像特征點(diǎn)提取的均勻性也會對圖像配準(zhǔn)的精度產(chǎn)生影響，以后將對本文算法進(jìn)行改進(jìn)和優(yōu)化，提升特征點(diǎn)提取的均勻度，進(jìn)一步提升圖像配準(zhǔn)的精度。