Ｌ－Ｍ算法在變換矩陣計算中的應用

摘 要：基于特征的圖像拼接是最受關注的一類拼接方法，它分為特征提取、特征匹配、變換矩陣的計算和圖像融合4步。如何利用獲得的多對特征點來計算較高精度的變換矩陣是圖像拼接中的重要問題。介紹L-M (Levenberg-Marquardt，利弗博格-馬夸特)算法，并將L-M算法運用于獲得較高精度的變換矩陣的計算中。另外，對于應用該算法的過程中出現的矩陣元素過大而無法求逆的問題，也給出了解決的辦法。實驗結果證明，該方法實用可行，能有效提高拼接質量。

關鍵詞：圖像拼接；L-M算法；變換矩陣；投影變換

中圖分類號：TP391.41文獻標識碼：B

文章編號：1004-373X(2008)24-099-03

Computation of Transformation Matrix Using L-M Algorithm

CAO Hongxing1，2，LIU Jiahang1，RUAN Ping1

(1. Xi′an Institute of Optics and Precision Mechanism，Chinese Academy of Sciences，Xi′an，710119，China;

2.Graduate University，Chinese Academy of Sciences，Beijing，100048，China)

Abstract：Image stitching based on features is one of most attended methods.It divides into the feature extraction，matching feature points，transformation matrix computation and blending image.How to get full use of features point pair to get transformation matrix with higher precision is a crucial question.This article introduces the L-M (Levenberg-Marquardt)algorithm and uses the L-M algorithm to obtain higher precision transformation matrix.In addition，the problem that matrix elements are too large to inverse matrix appears in process of applying the algorithm，the solution to the problem is also given in this article.Experiments show that this method is practical and can improve quality of image stitching effectively.

Keywords：image sitching;L-M;transformation matrix;project transformation

1 引 言

圖像拼接在航天，醫學，工業等諸多領域都有廣泛應用。而基于特征點的圖像拼接方法是研究最多的一類拼接方法。特征點不易受光照、旋轉等變換而發生改變，具有較強的穩定性，而且少量的特征點不僅能反映出圖像的重要信息，還有利于提高圖像拼接速度，因此，有許多人在從事這方面的研究^［1-8］。基于特征點的圖像拼接分為特征提取、特征匹配、變換矩陣的計算、圖像融合4步。其中變換矩陣的計算需要4對匹配點，但事實上，經過特征匹配得到的匹配點遠遠多于4對。如何充分利用得到的多對特征點，獲得更高精度的變換矩陣就成為一個問題。

首先介紹L-M算法，然后詳細介紹如何將其應用到變換矩陣的求解中，并解決了迭代過程矩陣元素過大無法求逆的問題。實驗結果證明，通過L-M算法，可以有效提高變換矩陣的精度，從而提高拼接質量。

2 L-M算法簡介

最小二乘法就是要尋找x*，使得x=x*時，函數F(x)=12∑ni=12取得最小值，屬于最優化問題。當參數與函數之間不滿足線性關系，即為非線性最小二乘法。最優化問題是通過不斷迭代逼近最優解。不同的算法的迭代公式不同。

L-M算法又稱為阻尼最小二乘法，是解決非線性最小二乘法的有力工具。它是馬夸特(Marquardt)在利弗博格(Levenberg)工作的基礎上提出了一種非常精致的迭代策略^［9］，是對高斯-牛頓法的改良。為了更好理解該算法，本文首先簡單介紹高斯-牛頓法的迭代公式的推導過程^［10］，然后再介紹L-M算法。

H高斯-牛頓法的迭代公式推導如下：

F(x)=12∑ni=1(fi(x))2=12f(x)Tf(x)(1)

f(x+h)f(x)+J(x)h=l(h)(2)

式(2)中的J(x)為f對x的雅克比行列式:

F(x+h)L(h)=12l(h)Tl(h)(3)

將式(2)代入式(3)，得：

L(h)=F(x)+hTJ(x)Tf(x)+1/2hTJ(x)TJ(x)h(4)

L′(h)=J(x)Tf(x)+J(x)TJ(x)h(5)

L″(h)=J(x)TJ(x)(6)

F(x+h)取得最小值，則L′(h)=0，即:

hgn=J(x)TJ(x)〗-1J(x)Tf(x)(7)

高斯牛頓法使用公式xk+1=xk+hgn=xk-(J(x)TJ(x))-1J(x)Tf(x)不斷迭代，最終得到最優解。

L-M算法迭代公式為：

hLM=-J(x)TJ(x)+μI〗-1J(x)Tf(x)

(μ>0，I為單位陣)(8)

它比高斯牛頓法的迭代公式多出了μI項。該項的存在，可以保證J(x)TJ(x)+μI為正定，從而可逆;而且當μ較小時，該方法近似為高斯-牛頓法；而μ較大時，hLM的方向為最速下降法的方向。μ的存在使得L-M算法可以更接近全局最優解^［10］ 。

3 應用L-M算法進行變換矩陣的計算

本文采用的是投影變換，即兩幅待拼接圖像之間的匹配點之間滿足式(9)，其中(xi，yi)，(x′i，y′i)分別為兩幅待拼接圖像的對應匹配點。

x′iy′i1〗=m11m12m13

m21m22m23m31m321〗xiyi1〗(9)

由式(9)得：

x′i=m11xi+m12yi+m13m31xi+m32yi+1

y′i=m21xi+m22yi+m23m31xi+m32yi+1(10)

本文的目標是求解變換矩陣的m11~m32八個參數的最優解，使得所有特征點的像點與其匹配點之間的距離和最小，即F(m)取得最小值。

F(m)=∑im11xi+m12yi+m13m31xi+m32yi+1-x′i2+

m21xi+m22yi+m23m31xi+m32yi+1-y′i2〗1/2(11)

令 fi(m)=m11xi+m12yi+m13m31xi+m32yi+1-x′i2+

m21xi+m22yi+m23m31xi+m32yi+1-y′i2〗1/4(12)

則式(11)形式同式(1)，因此可以使用L-M算法進行求解。

本文首先使用SIFT算法，得到了兩幅待拼接圖像的 ni對初始匹配點^［11］ ，再使用RANSAC方法求得變換矩陣的初始解以及較符合變換矩陣的n對匹配點^［12］。然后對n對匹配點，使用L-M算法求得精度較高的變換矩陣。L-M具體計算流程如圖1所示。

圖1 L-M算法流程圖

其中J(m)的形式為：

f₁m11f₁m12…f₁m32

f₂m11f₂m12…f₂m32

…… …

fnm11fnm12…fnm32

令：

error xi=m11xi+m12yi+m13m31xi+m32yi+1-x′i

error yi=m21xi+m22yi+m23m31xi+m32yi+1-y′i

errori=m11xi+m12yi+m13m31xi+m32yi+1-x′i2+

(m21xi+m22yi+m23m31xi+m32yi+1-y′i2〗1/4

deverrori=1/［4*(errori)3］

則：

fim11=2deverrori·error xi*yi/(m31xi+m32yi+1)

fim12=2deverrori·error xi*yi/(m31xi+m32yi+1)

fim13=2deverrori·error xi*1/(m31xi+m32yi+1)

fim21=2deverrori·error yi*xi/(m31xi+m32yi+1)

fim22=2deverrori·error yi*yi/(m31xi+m32yi+1)

fim23=2deverrori·error yi*1/(m31xi+m32yi+1)

fim31=2*deverrori*xi*-error xi*(m11xi+m12yi+m13)-error yi*(m21xi+m22yi+m23)(m31xi+m32yi+1)2

fim32=2*deverrori*yi*-error xi*(m11xi+m12yi+m13)-error yi*(m21xi+m22yi+m23)(m31xi+m32yi+1)2

如果直接采用該算法計算，由于J(m)矩陣中的元素特別大，導致求逆矩陣存在問題。而J(m)矩陣中的最主要是最后2列元素特別大，這是因為：一般情況下，m31和m32較接近0，因此fim31和fim32為xi，yi的平方量級。而1幅圖像的長或寬的像素至少也有上百個像素，導致J(m)中最后兩列值大(104量級)，導致后續矩陣J(m)TJ(m)+μI數據更大，最終無法計算J(m)TJ(m)+μI〗-1。為了減小xi，yi量級，本文在迭代之前，先將圖1上的特征點的橫、縱坐標xi，yi均縮小10 000倍，為了進行此改動時，不影響F(m)的值，將初始單應矩陣中m11，m12，m21，m22，m31和m32放大10 000倍，而m13和m23不變。迭代過程不做任何改動。迭代結束后，將m11，m12，m21，m22，m31和m32縮小10 000倍。此時得到的m11 - m32就是所求的最優解。由于xi，yi的縮小，J(m)矩陣元素值大大減小，使(J(m)TJ(m)+μI)-1計算得以順利進行。

本文使用L-M算法得到的變換矩陣后，將兩幅待拼接圖像進行融合，最終得到拼接圖像，如圖2所示。

圖2 待拼接圖像及拼接結果

4 實驗與結果分析

本試驗的待拼接圖像由kodak Z650相機手持拍攝圖片，為了處理的方便，先將數碼照片轉換為灰度圖。拼接程序的開發平臺為VC++6.0。

本文先將圖2中的(a)和(b)拼接，得到(d)，再將(d)與(c)拼接。 圖2中的(e)為使用L-M算法進行變換矩陣計算后的拼接結果。 (f)為未使用L-M算法進行拼接結果。在(f)中，水池左側小路在接縫處存在一定的錯位，而(e)中幾何線條連貫自然，因此拼接質量好于(f)。圖3為半透明拼接圖(即重疊區，兩幅圖都可以顯示出來，便于比對)，圖中在重疊區域內，除了水中倒影外無明顯重影，道路磚塊線條清晰，這都說明變換矩陣精度較好。

圖3 半透明拼接圖

5 結 語

從圖像拼接效果來看，使用L-M算法進行變換矩陣的求解，可以顯著提高了圖像拼接質量。

參考文獻

［1］陳利軍.圖像角點檢測和匹配算法得研究［D］.西安:西安電子科技大學，2005.

［2］陳靜.圖像配準特征點提取算法研究［D］.南京:南京理工大學，2006.

［3］汪華琴.基于特征點匹配的圖像拼接方法研究［D］.武漢:華中師范大學，2007.

［4］Yao Li，Lizhuang Ma.A Fast and Robust Image Stitching Algorithm［A］.Proceedings of the 6th World Congress on Intelligent Control and Automation.2006:21-23.

［5］胡社教，葛西旺，陳宗海.基于角點特征的KLT跟蹤全景圖像拼接算法［J］.系統仿真學報，2007，19(8):1 742-1 745.

注：本文中所涉及到的圖表、注解、公式等內容請以PDF格式閱讀原文