龐戰1,2，滕奇志1，何海波2

　?。?. 四川大學 電子信息學院圖像信息研究所，四川 成都 610064；2. 成都西圖科技有限公司，四川 成都 610064）

       摘要：為了解決傳統的基于尺度不變特征變換（SIFT）算法的圖像拼接方法無法使用SIFT特征匹配算法對巖石薄片圖像中的空白區域圖像進行特征提取，進而影響整個巖石薄片圖像拼接的問題，在尺度不變特征變換匹配算法的基礎上，提出了一種將巖石薄片圖像的行列中心位置作為拼接基準點的圖像拼接方法。該方法改進了特征點的提取方式、變換矩陣的計算順序以及矩陣優化的順序，并利用相鄰圖像之間的位置關系估計出空白區域圖像的變換矩陣，實現了整個巖石薄片圖像的拼接。實驗結果表明，該方法可以較好地完成巖石薄片圖像中空白圖像的拼接，能夠較為完整地保留整個巖石薄片的紋理信息，具有一定的實際應用價值。

　　關鍵詞：尺度不變特征變換；空白區域；巖石薄片；估計；圖像拼接

　　中圖分類號：TP391.4文獻標識碼：ADOI： 10.19358/j.issn.1674-7720.2017.06.015

　　引用格式：龐戰，滕奇志，何海波. 基于SIFT的巖石薄片圖像拼接［J］.微型機與應用，2017,36（6）：46-50.

0引言

　　*基金項目：國家自然科學基金項目（61372174）圖像拼接是指將在同一環境或者條件下拍攝的一組相互之間具有一定重合量的圖像，通過提取特征點、圖像配準等處理后，得到一幅包含巖石薄片中各個圖像紋理信息的全景圖。目前，圖像拼接技術在很多方面都有廣泛的應用價值［1］。

　　在石油地質分析中，由于顯微鏡視域的限制，無法一次性采集得到寬視域巖石薄片的全景圖［2］。為了從整體上分析巖石薄片，只能先對巖石薄片分區域采集得到整個薄片的圖像，然后將圖像拼接成一幅大視域圖像。目前，圖像拼接主要采用SIFT算法，通過圖像預處理、圖像配準以及融合完成圖像拼接［35］ 。SIFT算法主要應用在相鄰的具有明顯特征點和一定重疊量的圖片拼接中。然而，由于巖石薄片的不規則性，巖石薄片通常包括空白區域(如圖1所示)，由于無法進行特征提取，導致整個巖石薄片圖像的拼接失敗。

　　本文在SIFT算法的基礎上，通過改變SIFT特征提取的順序與方式、配準的基準點位置和變換矩陣的計算以及優化方式，實現整個巖石薄片圖像的拼接。實驗結果表明，該方法能夠得到比較完整并且包含了各個圖像信息的巖石薄片全景圖。

1相關理論

　　1.1SIFT算法的基本原理

　　SIFT算法是由哥倫比亞大學的LOWE D G教授于1999年提出，并在2004年完善的［6］。該算法的具體流程如下：

　?。?）構建尺度空間，檢測極值點

　　尺度空間是模擬圖像數據的多尺度特征，檢測圖像中具有尺度不變性的位置。1994年，LINDEBERG T發現唯一有可能的尺度空間內核是高斯函數［7］。假設I(x,y)表示一幅圖像，尺度可變的高斯函數用G(x,y,σ)表示，則其尺度空間的定義如式（1）所示。

　　L(x，y,σ)=G(x,y,σ)*I(x,y)（1）

　　其中，σ代表尺度因子，其大小表示平滑程度，*為卷積運算符，高斯函數的定義用式（2）表示：

　　G(x,y,σ)=12πσ2exp(－(x2+y2)/2σ2)（2）

　　LOWE D G使用高斯差分算子(DoG)主要有兩個原因：一是能夠更方便地檢測出穩定的關鍵點；二是考慮到計算的復雜程度。

　　D(x,y,σ)=(G(x,y,kσ)－G(x,y,σ))*I(x,y)

　　=L(x,y,kσ)－L(x,y,σ)（3）

　　其中，k為閾值。LOWE D G在實驗中，通過對原始圖像進行降采樣構造若干階尺度空間，如圖2（a）所示，然后在每一階尺度空間內利用高斯函數卷積得到若干層高斯圖像，再差分生成高斯差分尺度空間，如圖2（b）所示。在同一個尺度內，將每一個像素點同時與該像素點同尺度的8個相鄰點以及上下兩個相鄰尺度對應的18個像素點的值進行比較，只有當該像素點的值是最大或者最小值時，才稱該像素點是圖像在該尺度下的一個特征點，也叫關鍵點。最后，對其余的每一階尺度空間都重復以上操作，直到檢測完所有尺度空間中的關鍵點為止。

　?。?）定位關鍵點，去除不穩定點

　　LOWE D G使用擬合三維二次函數的方法準確定位檢測出來的關鍵點的位置和尺度。該方法可以有效去除對比度低以及穩定性較差的關鍵點。使用的尺度空間函數的泰勒公式展開式如式（4）所示：

　　其中，X=(x,y,σ)T，在采樣點處，計算D和D的導數，當其導數為零時，得到極值點的X1，如式（5）所示：

　　利用尺度空間函數在極值點處的函數值去除對比度較低的響應點。將式（5）代入式（4），計算得到極值點處函數值，如式（6）。

　　通常情況下，若|D(X1)|≥0.03，該特征點就保留，否則就去除。

　　因為高斯差分算子在邊緣處會產生比較強烈的邊緣響應，在該算子中，差的極值點在邊緣處具有較大的曲率，在垂直方向上具有較小的曲率。因此，可以利用該性質去除邊緣響應。主曲率可以通過一個2×2的Hessian矩陣H求出，如式（7）所示。

　　D的主曲率和H的特征值成正比，若假設α為最大特征值，β為最小特征值，則矩陣H主對角線的值和行列式的值用式（8）和式（9）表示。

　　如果曲率大于式（10）結果的比值，則舍棄。LOWE D G在實驗中通過對比發現，當r=10時，能夠最大程度去除不穩定的邊緣響應點［6］。

　　（3）確定關鍵點的大小和方向

　　在初步確定了關鍵點的位置之后，可以利用圖像的局部特性為每一個關鍵點分配一個基準點，使該關鍵點的特征描述符具有旋轉不變性。對于高斯金字塔檢測出來的關鍵點，使用其高斯金字塔鄰域窗口內像素的梯度特征和方向特征對描述符進行描述。(x,y)處的梯度值m(x,y)和方向θ(x,y)的表達式如式（11）和式（12）所示。

　　其中，L所用的尺度為每個關鍵點各自所在的尺度。為了進一步準確給出關鍵點的方向，一般采用梯度直方圖統計的方法，以關鍵點為中心的鄰域窗口內的像素的梯度和方向。直方圖峰值所代表的方向即為關鍵點的方向。通過上述各個步驟之后，檢測出的包含位置、所在尺度以及方向的關鍵點稱為圖像的特征點。

　　（4）生成特征描述子

　　為了確保關鍵點不受光線等因素的影響，可以為每一個關鍵點建立一個描述子，使其具有獨特性，進而可以提高特征點的匹配率。為了保證描述子具有旋轉不變性，首先將方向旋轉到關鍵點的方向。LOWE D G的實驗中，在關鍵點所在的尺度空間中，以關鍵點為中心，選取一個8×8的窗口，然后在每個4×4的窗口內計算8個方向的梯度方向直方圖，繪制每個梯度方向的累加可形成一個種子點，最終獲得一個4×4×8共128維的向量特征。

　　（5）特征點的匹配

　　SIFT特征向量生成后，一般使用最近鄰域算法對特征點進行匹配。利用特征點的最短歐氏距離和次最短歐氏距離的比值與所設定的閾值的大小關系作為判斷兩個特征點是否為一匹配對的依據。如果結果小于閾值，則認為這一對特征點是匹配的。LOWE D G在實驗中發現，將閾值設定為0.8比較合適，在計算特征點的歐氏距離之后，可以先利用BBF（BestBinFirst）算法對特征向量進行快速處理，再使用RANSAC算法去除誤匹配的特征點［7］。

　　1.2變換矩陣計算的基本原理

　　描述圖像之間位置關系的矩陣稱為變換矩陣。圖像的變換矩陣可以通過特征匹配對計算得到［89］。在同一個空間下，兩幅圖之間的變換關系用式（13）表示，其中，(x1,y1)為參考圖像上的點，(x2,y2)為目標圖像上與(x1,y1)對應的點，M稱為變換矩陣:

　　在式（13）中，線性方程組中有8個未知變量，則至少需要4對特征點匹配對才能求解方程組的解。整理并將結果進行向量化用式（14）表示:

　　A1=MA2（14）

　　若矩陣A2AT2可逆，則可以求解得到變換矩陣M，如式（15）所示：

　　M=(A1AT2)(A2AT2)-1（15）

　　只要有足夠多可用的特征點匹配對就可以通過式（15）計算得到圖像之間的變換矩陣，從而確定圖像之間的相對位置關系。

2空白區域圖像拼接方法

　　巖石薄片中空白區域的圖像（見圖1）雖然可以使用SIFT算法對相鄰圖片進行特征匹配，但由于空白區域圖像的表面特征比較單一，紋理信息不豐富，檢測到的特征點匹配對大多都是錯誤匹配，在使用RANSAC算法［10］去除錯誤匹配之后，可用的特征匹配數目很少。如圖3所示，從實驗對比結果可以看出，巖石薄片中空白區域的相鄰圖像的特征匹配數目在使用RANSAC算法去除錯誤匹配之后，沒有可用的特征點匹配對，無法利用式（15）求出圖像的變換矩陣，從而無法確定兩張圖片的相對位置，影響整個巖石薄片圖像的拼接。

　　圖4拼接流程圖在實際應用中，為了解決該問題，在SIFT算法的基礎上，提出了一種選取巖石薄片圖像的中心位置為基準點的拼接方法。并利用空白區域圖像與其相鄰圖像之間的位置關系估計出空白區域圖像的變換矩陣，實現了整個巖石薄片的拼接。流程圖如圖4所示。

　　2.1特征點提取

　　為了解決SIFT算法無法提取空白區域圖像的特征點的問題，本文在SIFT算法的基礎上，改變了特征點提取的方式，并在RANSAC算法篩選時設定閾值，如篩選結果小于閾值則視為不能正常匹配，并標記。本實驗設定的閾值為30，該閾值也可以根據圖片特征點的具體情況進行設定。如圖5所示，該示意圖表示3行7列圖片。其中，橫向箭頭表示每一行圖片的特征點提取順序，豎向箭頭表示中間一列圖片的特征點提取順序。具體步驟如下：

　?。?）獲取待拼接巖石薄片圖像的行和列的數目。

　?。?）分別對巖石薄片圖像中的每一行和中間列進行特征提取，并對不能正常匹配的圖片進行標記。

　　2.2計算變換矩陣

　　由于無法得到空白區域圖片的特征點匹配對，從而無法確定空白區域圖片的相對位置。在前文特征提取順序的基礎上，對變換矩陣的計算順序作了相應的調整。如圖6所示，其中每一格代表一張圖片，圓形代表本文選取的中心位置，五角星代表空白圖片，箭頭表示變換矩陣的計算方向，每一行圖片變換矩陣的計算方向都與示意圖中箭頭表示的方向一致。文中，“左（右）側”是表示分別對左側和右側都進行相同的操作。具體步驟如下：

　?。?）選取巖石薄片圖像的中心位置（圖6中圓形所在位置），并初始化其變換矩陣為單位矩陣，作為變換矩陣計算的基準點。

　　（2）從中心位置所在行開始，到圖像的起始行，分別從每一行的中間位置開始，對左（右）側的圖片計算變換矩陣，如果遇到標記圖像（圖6中圓形左（右）側五角星代表的圖片）則停止計算當前圖片的變換矩陣，并對標記圖片左（右）側的所有圖片的變換矩陣賦值為空。

　　（3）從中心位置所在行開始到圖像的最后一行，分別從每一行的中間位置開始重復步驟（2）。

　?。?）從圖像中間列的中間位置開始，對中間列圖像的上下兩部分的圖片重復步驟（1）和（2）。

　　2.3優化變換矩陣

　　計算完所有圖片的變換矩陣之后，本文在上一節提出的變換矩陣的計算順序的基礎上，對變換矩陣的優化方式作了相應的調整。具體優化方式如下：

　　（1）從中心位置所在行開始到圖像的起始行，分別從每一行的中間位置左（右）側的圖片進行變換矩陣優化，如果遇到標記圖像則停止優化。直至每一行圖像都重復該操作為止。

　　（2）從中心位置所在行開始到圖像的最后一行，分別從每一行的中間位置開始重復步驟（1）。

　　（3）從薄片圖像中間列的中間位置開始，對中間列圖像的上下兩部分的圖片重復步驟（1）和（2）。

　　2.4估計空白區域圖片的變換矩陣

　　變換矩陣優化完成之后，為了獲得完整的拼接結果，需要估計空白區域圖片的變換矩陣。具體方法如下。

　?。?）從中心位置所在行開始到圖像的起始行，分別從每一行的中間位置開始，對其左（右）側的圖片進行標志位掃描，遇到標記圖片，則利用該標記圖片右（左）側相鄰兩張圖片的變換矩陣與標記圖像之間的相對位置關系估計出標記圖片的變換矩陣，同時位于標記圖片左（右）側的所有圖片的變換矩陣均通過該方法依次得到。直至每一行圖像都重復該操作為止。

　?。?）從中心位置所在行開始到圖像的最后一行，分別從每一行的中間位置開始重復步驟（1）。

　?。?）從圖像中間列的中間位置開始，對中間列圖像的上下兩部分的圖片則重復步驟（1）和（2）。

　　完成上述處理之后，通過相應的融合算法對結果圖進行融合，消除拼接縫［11］，完成整個巖石薄片圖像的拼接。

3實驗對比結果及分析

　　本文使用6行6列共36張圖片進行實驗，圖片尺寸都是5 184×3 456，實驗所用的圖像包含了實驗所需要的薄片信息。實驗結果如圖7~9所示，其中，圖7是實驗所用的圖像，圖8是基于SIFT算法的拼接結果圖，圖9是本文提出的拼接方法的結果圖。

　　從實驗的對比結果中可以明顯看出，當采集得到的巖石薄片圖像中包含空白區域圖片時，如圖7實驗使用的圖像中“C0001.jpg”和“C0033.jpg”等圖片，因為無法使用SIFT算法對巖石薄片中空白區域圖片完成配準，因此，只能在拼接前剔除這些圖片所在的行或者列的所有圖片，對比可以看到基于SIFT算法的拼接結果圖中，第一張圖片所在的列和下方最后一行都存在嚴重缺失（結果圖中的黑色缺失部分）。而本文提出的圖像拼接方法則較為完整地完成了整個巖石薄片圖像的拼接，較好地反映出整個巖石薄片中顆粒的分布情況，保證了巖石薄片的完整性。因此，本文提出的圖像拼接方法能較好地完成巖石薄片空白區域圖像的拼接，能較為完整地展示巖石薄片的全貌，具有一定的實際應用價值。

4結論

　　本文首先介紹了SIFT算法和變換矩陣計算的基本原理，同時對巖石薄片中空白區域圖像不能使用SIFT算法進行特征點提取的原因作了簡要分析，并在此基礎上分別對本文提出的圖像拼接方法中的特征點的提取方式、變換矩陣的計算方式、變換矩陣的優化方式以及估計空白區域圖片的變換矩陣的具體步驟進行了描述，在不改變巖石薄片整體性基礎上，完成了包含空白區域的巖石薄片圖像的拼接，保證了巖石薄片紋理信息的完整性。本文的不足之處在于，若在高倍顯微鏡下巖石薄片的中心區域存在大量空白區域，則無法具體確定拼接的基準位置，對該方法的使用具有一定的限制。

　　參考文獻

　?。?］ 徐鑫, 孫韶媛, 沙鈺杰，等. 一種基于RANSAC的紅外圖像拼接方法［J］. 激光與光電子學進展，2014，51(11):129-134.

　?。?］ 岳永娟, 苗立剛, 彭思龍. 大規模顯微圖像拼接算法［J］. 計算機應用，2006, 5(26):1012-1014.

　?。?］ 鄒承明, 侯小碧, 馬靜. 基于幾何學圖像配準的SIFT圖像拼接算法［J］. 華中科技大學學報（自然科學版），2016，44(4)：32-36.

　　［4］ 汪松, 王俊平, 萬國挺，等. 基于SIFT算法的圖像匹配方法［J］. 吉林大學學報（工學版），2013，43(S1):279-282.

　?。?］ 白廷柱, 侯喜報. 基于SIFT算子的圖像匹配算法研究［J］. 北京理工大學學報，2013，33(6):622-627.

　?。?］ LOWE D G. Distinctive image features from scaleinvariant keypoints［J］. International Journal of Computer Vision，2004，60(2): 91110.

　?。?］ LINDEBERG T. Scalespace theory: a basic tool for analyzing structures at different scales ［J］. Journal of Applied Statistics, 1994, 21(2):225-270.

　　［8］ 郭紅玉, 王鑒. 一種基于RANSAC基本矩陣估計的圖像匹配方法［J］. 紅外，2008，29（2）: 5-8.

　?。?］ 曲天偉, 安波, 陳桂蘭. 改進的RANSAC算法在圖像配準中的應用［J］. 計算機應用，2010，30（7）：1849-1851.

　?。?0］ FISCHER M A, BOLLES R C. Random sample consensus: a paradigm for model fitting with applications to analysis and automated cartography ［J］. Communication of the ACM,1981, 24 (6):381-395．

　　［11］ 苗啟廣, 王寶樹. 基于改進的拉普拉斯金字塔變換的圖像融合方法［J］.光學學報，2007，9(27): 1605-1610.