孫登第 ，羅 斌 ，卜令斌

（1.安徽大學 計算機科學與技術(shù)學院，安徽 合肥 230039；2.安徽省工業(yè)圖像處理與分析重點實驗室，安徽 合肥 230088）

圖像配準是對同一場景在不同條件下（如不同的時間、拍攝環(huán)境、視角和傳感器等）得到的兩幅或多幅圖像尋求某種空間上的變換，使一幅圖像能夠和另一幅圖像上的對應點達到空間上的一致[1]。圖像配準技術(shù)是圖像處理與分析中的基本任務，已經(jīng)在計算機視覺、圖像融合、全景圖像拼接、醫(yī)學診斷與輔助治療等眾多領(lǐng)域得到廣泛的應用[2-3]。

目前，圖像配準方法大體分為基于灰度和基于特征兩大類?；诨叶鹊姆椒ńD像間像素灰度值的目標函數(shù)，如互信息測度[4]，通過對目標函數(shù)的優(yōu)化實現(xiàn)配準。該方法沒有考慮像素的空間信息，在不同成像條件下的圖像配準中，其精度較低、計算量大且配準時間長。

基于特征的配準方法先提取各個圖像中的特征，再完成特征之間的匹配，通過匹配的特征建立圖像間的映射變換，最后求出配準后的圖像。特征點是該類方法最常使用的圖像特征，其中BAY H等人提出的SURF（Speeded-Up Robust Features）算法是一種尺度空間的特征點描述方法[5]，對圖像間的分辨率、旋轉(zhuǎn)、平移和光照變化等保持不變，且時間復雜度低、速度較快?；谔卣鼽c的配準方法減輕了圖像灰度差異和噪聲的影響，縮短了配準時間。然而，特征點匹配問題本身就是一個尚未得到較好解決的難題，特征點的誤匹配直接影響了圖像的最終配準結(jié)果。

為解決上述問題，RANGARAJAN A等提出了一種結(jié)合互信息與特征點的配準方法[6]，定義了特征點集的互信息函數(shù)，通過對該函數(shù)最大化實現(xiàn)圖像配準。這種方法減輕了灰度差異與特征點誤配對配準的影響，但函數(shù)形式復雜，配準時間較長。

最新的研究結(jié)果表明，通過對隨機抽樣構(gòu)建廣義近鄰圖可以估計隨機變量的熵[7]，這已在統(tǒng)計學與信息論的研究中受到廣泛關(guān)注。本文將該理論引入圖像配準中，將圖像配準中的特征點與互信息結(jié)合起來估計特征點的Rényi互信息，提出了一種結(jié)合SURF描述子和廣義近鄰圖醫(yī)學圖像配準方法。該方法了融合圖像空間信息且無需計算概率直方圖。通過與幾種傳統(tǒng)配準算法相比較，結(jié)果表明，該算法在魯棒性、配準時間和配準精度方面提供了更好的綜合性能。

1 SURF檢測及描述

SURF可以在圖像尺度空間中提取特征點，并對每個特征點賦予特征，即SURF描述符。該算法提取的特征點對尺度、旋轉(zhuǎn)、光照、仿射和透視變換等均具有較強魯棒性，并在計算速度上明顯快于以往同類方法。

1.1 特征點的檢測

SURF是一種基于尺度空間的特征點檢測算法。在尺度空間的每一層圖像上，SURF使用快速Hessian矩陣來檢測圖像的極值點。對于尺度為σ的空間中任意一點（x，y），Hessian 矩陣的定義為：

其中，Lxx、Lxy與 Lyy是點（x，y）分別與高斯函數(shù) g（σ）二階偏導卷積的結(jié)果。

為加快SURF檢測速度，BAY H等人在尺度圖像金字塔中用不同尺寸的框狀濾波 Dxx、Dxy與 Dyy代替二階高斯濾波，用積分圖像來加速卷積[8]，進一步求解得到Hessian矩陣的行列式作為點（x，y）在尺度 σ空間中的響應值：

其中，ω為權(quán)重系數(shù)，約等于0.9。

對每個點（x，y），當 Hessian矩陣行列式大于設(shè)定的閾值時，就作為待判定點。對每個待判定點上下層與同層的3×3×3的立體領(lǐng)域中進行非極大值抑制，當響應值為局部極大或極小值時，該點被確定為特征點，并經(jīng)過插值確定位置[9]。

1.2 主方向確定

為保證特征點描述符的旋轉(zhuǎn)不變性，SURF賦予每個特征點主梯度方向。在以特征點為中心，半徑為6σ（σ為特征點對應的尺度）的鄰域內(nèi)，用邊長為4σ的Haar小波模板計算該點在x、y方向的Haar小波響應，并以該點的高斯函數(shù)對這些響應值加權(quán)，然后將60°范圍內(nèi)的響應相加形成新的向量，遍歷整個圓形區(qū)域，選擇最長向量的方向為該特征點的主方向。這樣，對特征點逐個進行計算，就可以得到每一個特征點的主方向。

1.3 描述子生成

為了構(gòu)建描述子向量，首先以特征點為中心，將坐標軸旋轉(zhuǎn)到主方向，在新坐標軸中選取邊長為20σ的正方形區(qū)域，再將該區(qū)域劃分成4×4個子區(qū)域。計算每個子區(qū)域內(nèi)水平方向與垂直方向的Haar小波響應dx與dy，并用特征點為中心的高斯函數(shù)對 dx、dy加權(quán)，以增加對幾何變換的魯棒性。

在每個子區(qū)域中對Haar小波響應以及響應的絕對值求和，形成這樣，得到一個四維向量因此，對每一特征點，把4×4個子區(qū)域的向量連起來就形成64維的向量，即該特征點的SURF描述子。

2 結(jié)合SURF和廣義近鄰圖的圖像配準

2.1 算法流程

本文算法流程圖如圖1所示。

圖1 算法流程圖

2.2 廣義近鄰圖估計 Rényi熵

Rényi熵 Rα是 Shannon熵 H的廣義形式，具有比Shannon熵更為平滑的函數(shù)形式，因此在圖像配準中受到越來越多的關(guān)注。

對于概率分布為 pi的隨機變量 X，Rényi熵為：

式（3）通過樣本概率來計算熵，但在許多實際問題中，樣本概率并不容易獲得。為克服這一缺陷，許多研究者從樣本的隨機分布情況出發(fā)，利用樣本構(gòu)成的圖結(jié)構(gòu)來描述樣本的隨機分布。最新的研究結(jié)果表明，通過對隨機抽樣構(gòu)建廣義近鄰圖可以較精確地估計隨機變量的熵，收斂速度快且魯棒性較好，已在統(tǒng)計學與信息論的研究中受到廣泛關(guān)注。

廣義近鄰圖是一般的K近鄰圖的推廣，一般記作NNS（V）。 對于頂點集 V，|V|=n，S為某一非空的有限正整數(shù)集，k為S中的最大元素值。對每個i∈S和每個頂點x∈V，x與其第i鄰近點構(gòu)成一個邊，這些邊構(gòu)成了NNs（V）的邊集。

考慮V去除x之后余下頂點的集合 V{x}，不妨記為{y1，y2，…，yn-1}。 對 V{x}到頂點 x 的歐式距離排序：則yi就是頂點x的第i鄰近點。那么對S中的每個正整數(shù) i，在 NNS（V）就有一個 x 到 yi的邊。

PAL D等人最先提出用廣義近鄰圖來描述樣本點的空間位置分布，進而估計隨機變量的Rényi熵。對歐式空間 Rd上的隨機抽樣點集 V，用 Lp（V）表示廣義近鄰圖邊的歐氏距離p次冪的和：

其中，E（NNS（V））表示 NNS（V）邊的集合，p≥0 為參數(shù)。在這里，S是固定的有限非空正整數(shù)集合，如S={1，3，4}表示取點的 1鄰近、3鄰近和4鄰近。Lp（V）是度量樣本分布離散程度的測度。PAL D進一步證明了這一度量與Rényi熵中的樣本概率冪的和成正比。通過合理設(shè)置比例，得到隨機樣本集 V 的 Rényi熵 Rα（V）在 α∈（0，1）時的估計：

其中，γ 為常數(shù)，p=d（1-α）。

2.3 目標函數(shù)的構(gòu)造

將待配準的兩幅圖像定義為浮動圖像A和參考圖像B，分別提取特征點集 V1和 V2，用 64維的 SURF向量描述圖像中特征點。每個SURF描述子可看作全像素上的密度函數(shù)在64維空間中的一個隨機抽樣，因此可以用式（5）估計 A 與 B 的 Rényi熵 Rα（A）、Rα（B），并利用點集 V=V1∪V2估計 A 與 B 的聯(lián)合 Rényi熵 Rα（A，B），然后用 Rényi熵和聯(lián)合 Rényi熵計算它們之間 Rényi互信息：

其中，Rα（A）與 Rα（B）是由圖像自身特征點估計的 Rényi熵。對于單幅圖像，事先按照主方向排序的SURF向量經(jīng)過平移和旋轉(zhuǎn)，其元素按照某一順序進行統(tǒng)一重排，而元素以同樣的順序進行重排并不影響向量之間的歐式距離計算。因此，在插值影響不大的情況下，可以近似地認為對每幅單獨的圖像A或B，其特征點集的廣義近鄰圖距 離 Lp（V1）、Lp（V2）在不同的 變換 T 下是不變的。進而，Rényi互信息簡化為 MIα（A，B）=-Rα（A，B），配準問題變?yōu)榍笞儞Q T使 Rα（A，B）最小，互信息最大。

通過式（5）用特征點集V估計 Rα（A，B），在本文中d=64，由 p=d（1-α），式（5）轉(zhuǎn)化為：

最終本文使用式（7）來估計 Rényi熵 Rα（A，B）。 求解最優(yōu)變換矩陣的過程變轉(zhuǎn)化為：

當A和B完全配準時，其對應的特征點的SURF向量重合或接近，此時構(gòu)造的廣義近鄰圖會包含大量的較短的邊，Lp（V）達到最??；而當 A和 B的對齊度差時，特征點將會呈現(xiàn)更加散布的狀態(tài)，廣義近鄰圖中會增加許多很長的邊，Lp（V）的值就增大。

3 實驗結(jié)果與分析

為了驗證本文提出的結(jié)合SURF描述子與廣義近鄰圖的配準算法（GNN-SURF）的有效性，對多幅真實遙感圖像進行配準實驗。將本方法與傳統(tǒng)的基于Shannon互信息 （NMI）的配準算法和形狀特征點互信息配準算法（Point-MI）進行比較，在配準準確度、時間與魯棒性等多項準則上進行實驗。實驗代碼用Matlab編寫，實驗機器配置為 Inter （R）Dual-Core （TM）2 Quad， E5500 2.80 GHz CPU，4 GB內(nèi)存。

3.1 配準準確度

取一對待配準的圖像A與B分別作為模板圖像和浮動圖像。用3種比較方法對圖像A與B進行配準實驗， 得到配準變換參數(shù) （△x′，△y′，△θ′）。 計算配準參數(shù)（△x′，△y′，△θ′）與真實參數(shù)（△x，△y，△θ）之間誤差的絕對值作為配準準確度，結(jié)果如表1所示。

表1 配準準確性的比較

從表 1可知，采用特征點的Point-MI與GNN-SURF都要比NMI方法準確，這主要是因為實驗中采用的遙感圖像成像條件不同，圖像內(nèi)容、灰度差異較大，因此基于像素灰度的NMI方法配準效果最差。此外，Point-MI僅依賴特征點坐標，而GNN-SURF方法融入了圖像的空間信息，配準更加準確。

圖2是兩幅遙感圖像配準例子。這兩幅圖像拍攝位置不同，圖像內(nèi)容僅有部分重合，特征點提取的數(shù)目與位置都不同（圖 2（a）、（b））。 然而在兩幅圖重合的部分中，圖 2（a）中的特征點與圖 2（b）中特征點是完全對應的，因此，配準時（圖 2（d））的廣義近鄰圖在中間重合部分比未配準之前（圖 2（c））的結(jié)構(gòu)簡單、清晰得多，對應點之間的連線使得SURF向量距離和最小，從而達到整體最優(yōu)配準。

圖2 配準前后廣義近鄰圖比較

3.2 配準魯棒性與時間

為了比較配準魯棒性，先對浮動圖像進行隨機平移和旋轉(zhuǎn)變換，平移參數(shù)（△x，△y）分別在[-20 mm，20 mm]中隨機選擇，旋轉(zhuǎn)參數(shù)（△θ）在[-20°，20°]中隨機選擇，共同構(gòu)成初始誤配參數(shù)（△x，△y，△θ），總共選擇 50 個初始誤配參進行試驗。對每次實驗結(jié)果計算平移和旋轉(zhuǎn)的誤差。根據(jù)常用評估標準[10]，當平移誤差小于1個像素、旋轉(zhuǎn)誤差小于1°時，認為配準達到了亞像素級，配準成功。

表2給出了3種比較方法的配準成功率和配準時間。從表2可知，GNN-SURF方法的成功率最高，而NMI方法效果最差。在配準時間上，NMI方法需要計算兩幅圖像的所有像素，耗時最久；Point-MI雖然只需要計算圖像特征點的坐標信息，但其定義的特征點集Rényi互信息函數(shù)形式過于復雜，增加了配準時間。而本文提出的GNN-SURF方法融合了快速計算魯棒特征點空間信息的SURF描述子，同時采用了廣義近鄰圖估計Rényi熵，在配準魯棒性與配準時間上都明顯優(yōu)于其他兩種方法。本文提出了一種結(jié)合SURF描述子和廣義近鄰圖的圖像配準算法。采用SURF快速、魯棒地提取圖像特征點，并形成特征點描述子，再結(jié)合特征點的廣義近鄰圖估計Rényi互信息進行配準。在真實遙感圖像中進行實驗，結(jié)果表明該算法在配準魯棒性、配準準確性和配準時間3個方面都優(yōu)于另外兩種傳統(tǒng)圖像配準方法。

表2 配準魯棒性與配準時間的比較

[1]ZITOVA B， FLUSSER J.Image registration methods： a survey[J].Image and Vision Computing， 2003， 21：977-1000.

[2]趙仕俊，孫林港.基于紋理特征的圖像自動配準方法研究[J].微型機與應用，2011，30（9）：36-38.

[3]李靖宇，穆偉斌，沈煥泉.醫(yī)學圖像配準的優(yōu)化算法改進研究[J].微型機與應用，2010（8）：47-49.

[4]COLLIGNON A， MAES F， DELAERE D， et al.Automated muhimodality medical image registration using information theory[C].Proceedings of Information Processing in Medical Imaging， 1995：263-274.

[5]BAY H，TUVTELLARS T，GOOL L VAN.SURF：speeded up robust features[C].Proceedings of the European Conference on Computer Vision， 2006：404-417.

[6]RANGARAJAN A， Chui Haili， DUNCAN J S， et al.Rigid point feature registration using mutual information[J].Medical Image Analysis，1999，3（4）：425-440.

[7]PAL D， POCZOS B， SZEPESVARI C.Estimation of Rényi entropy and mutual information based on generalized nearestneighbor graphs[C].NIPS，2010.

[8]VIOLA P，JONES M J.Rapid object detection using a boosted cascade of simple features[C].Proceedings of Computer Vision and Pattern Recognition， 2001：511-518.

[9]BROWN M，LOWED.Invariant features from interest point groups[C].BMVC， 2002：1-10.

[10]LUAN H X，QI F H，XUE Z，et al.Multimodality image registration by maximization ofquantitative-qualitative measure of mutual information[J].Pattern Recognition， 2008，41：285-298.