孫海迅,羅健欣,潘志松,張艷艷,鄭義桀
(陸軍工程大學 指揮控制工程學院,江蘇 南京 210001)
單應性矩陣反映了空間中的點在兩個投影平面的坐標變換關系,在圖像拼接[1]、三維重建[2]等計算視覺領域有廣泛應用。
在多視圖幾何中,單應性矩陣求解的一種經典方法是基于圖像特征點匹配的隨機抽樣一致性算法(RANSAC)[3]與總體最小二乘(Total Least Squares,TLS)[4]相結合的方法。其中,特征點匹配主要是找到不同視角的圖像上特征點的對應關系,特征點檢測可采用SIFT[5]、SURF[6]等方法;RANSAC算法的作用主要是從特征點對中消除誤匹配,準確找出正確的匹配對(內點),根據內點的對應關系構造線性方程組用TLS求解。
由于經典RANSAC方法存在迭代次數多、需要人為定義內點閾值等問題,研究人員提出了多種魯棒方法[7-14]。例如將RANSAC求解過程與貝葉斯方法結合起來的MAPSAC[8],應用圖分割方法進行局部優化的GC-RANSAC[9],通過加權最小二乘擬合優化求解的MAGSAC[10]及其改進版本MAGSAC++[11],將各種現有RANSAC類方法進行了有機融合的USACv20[12]。另外,還有利用特征點的幾何位置關系進行誤匹配濾除方法,如文獻[7]提出了一種基于改進RANSAC的方法,利用了特征點的距離信息;文獻[13]利用四個特征點構成的幾何形狀,提出了一種基于保序性約束和相似性度量的估計方法;文獻[14]提出了一種單應性矩陣自適應估計方法,克服了RANSAC需要人為設置閾值的缺點。
當前求解單應性矩陣各種方法的研究重點大都在于如何濾除誤匹配點得到內點,在得到內點后采用傳統的TLS方法求解線性方程組,而對于內點噪聲普遍較大的情況下,如何精確求解單應性矩陣研究較少。事實上,當內點坐標比較準確時,用TLS方法能準確求解單應性矩陣,而當內點噪聲較大時,用TLS求解精度會有所下降。本文針對以上問題,對方程組中噪聲矩陣的結構進行分析研究,充分考慮到噪聲列之間的相關性,提出了一種基于約束總體最小二乘(CTLS,Constrained Total Least Squares)[15]的單應性矩陣優化求解方法。
如圖1所示,若空間中一平面π方程為:nTX+d=0,其中,n∈R3×1表示平面法向量,X∈R3×1表示平面上任一點的三維坐標,d∈R3×1為一常數向量,文中上標T均表示矩陣或向量的轉置,上標-1均表示矩陣的逆;圖中O1、O2為兩個相機的光心,R∈R3×3、t∈R3×1為相機的旋轉矩陣和平移向量。

圖1 空間中平面在不同視角投影示意圖
則空間點X對應的兩幅圖像上的投影像素x1,x2∈R3×1有如下關系:
(1)


(2)
當匹配點對大于4對時,令H33=1,并構建方程組Ax=b,即:

(3)
其中,上標(i)表示第i對匹配點。求解以上方程的經典方法是TLS,即認為在系數矩陣A∈R2n×8和向量b∈R2n×1同時含有高斯噪聲,對應的噪聲矩陣和噪聲向量分別為E和e,于是將方程重新構造為:
(A+E)x=b+e
(4)

B=UΣVT
(5)
最小奇異值對應的右奇異向量V8=(v0,v1,…,v8)T∈R9×1,則方程組的解為:
(6)
以上方法雖然速度較快,但單純假設系數矩陣A上存在高斯噪聲并不符合方程的原型,通過觀察發現,A的某些元素(如第3列與第6列)是常數,不存在噪聲。另外,A的最后兩列是兩個坐標相乘的形式,疊加在其上的噪聲并不是高斯噪聲。
約束總體最小二乘法(CTLS)常被作為對經典TLS和普通最小二乘法(Ordinary Least Squares,OLS)的一種改進方法,在國內也有相關的研究[16-20]。文獻[17]提出了一種基于CTLS的到達時差到達頻差無源定位算法,降低了定位誤差;文獻[20]還將CTLS應用于點云拼接,提高了參數求解精度。
對于單應性矩陣的求解,更合理的方法是考慮到A與b上噪聲列的相互關系,將方程求解問題重新構造為CTLS問題,即假定高斯噪聲存在于像素坐標(u2,v2)T和(u1,v1)T上,所有的噪聲形成一個基本的噪聲向量w,而疊加在B上的噪聲矩陣D的每一列都可以由w生成,它們具有一定的線性關系,即可寫成如下形式:
D=[F0w,F1w,…,F8w]
(7)
其中,F0,F1,…,F8為不同的系數矩陣,則根據文獻[15],這一CTLS問題的求解即變成如下的優化問題:

(8)
并且可以通過最小化以下損失函數得到方程解:
(9)
其中,
(10)
具體實現方法如下:
設第i對匹配點的像素坐標為:
i=0,1,…,n-1
(11)
疊加在其上的高斯噪聲為:
(12)
于是,可得基本的噪聲向量:
(13)
設bj為矩陣B的第j行,則bj對si的雅各比矩陣(一階偏導數)為:
Jsibj=
當j=2i,i=0,1,…,n-1
(14)
Jsibj=
當j=2i+1,i=0,1,…,n-1
(15)
于是,疊加在B的第k列的一階噪聲擾動可表示為:
Dk=Fkw∈R2n×1,k=0,1,…,8
(16)
(17)



(18)
(19)
(20)

(21)

于是,便可以利用梯度下降法進行優化求解,損失E對z的梯度為:

Jsibj·z)
(22)

Jsibj·z)
(23)

另外,為增強計算的魯棒性,常采用歸一化圖像坐標進行計算,歸一化的方式為:
(24)
其中,u、v表示歸一化之前的像素坐標,l表示圖像的像素寬度,u0、v0表示圖像中心點的像素坐標,歸一化之后的圖像坐標介于[-0.5,0.5]。

為了對比實驗效果,在此加入其他兩種方法——OLS和數據最小二乘法(Data Least Squares,DLS)。OLS的求解可在通用矩陣教材中找到,DLS的求解由文獻[21]給出。
對于最小二乘問題Ax=b各種求解方法總結如下:
(1)普通最小二乘法(OLS)。
xOLS=(ATA)-1ATb
(25)
(2)數據最小二乘法(DLS)。
(26)
其中,vDLS是矩陣Pb·A最小奇異值對應的右奇異向量,Pb=(I-b(bTb)-1bT)是一投影矩陣,I是單位矩陣。
(3)總體最小二乘法(TLS)。
(27)
其中,vTLS是矩陣B=(A|b)最小奇異值對應的右奇異向量,vlast是vTLS的最后一位元素。
(4)約束總體最小二乘法(CTLS)。
(28)
相關符號說明見公式(10)。
評價標準為:
(1)平均坐標轉換損失。
Transfer Loss=
(29)
其中,n表示匹配點對的數量,上標(i)表示第i對匹配點,合成數據采用n=50。
(2)與真實H0矩陣之差的Frobenius范數平方。
(30)

從圖2中可以看出,當噪聲較小時,四種方法都能求解準確,噪聲標準差在1~5像素時,四種方法的平均坐標轉換損失在4個像素以下,且互相之間的求解結果相差不到1像素;但隨著噪聲標準增大,CTLS的平均坐標轉換損失小于其他幾種方法,且相同條件下δH的Frobenius 范數也較小。
單應性矩陣中不同元素的偏差對坐標轉換的影響不盡相同,所以δH的F范數與相應的坐標轉換損失沒有嚴格的正相關關系,如圖2所示,但在噪聲相同的情況下,δH的F范數仍能從一個側面反映出四種方法的求解精度。
單應性矩陣常用在圖像拼接上,下面選取在不同視角下拍攝的圖像。首先利用SURF提取特征點,用RANSAC去除誤匹配,再分別用以上OLS、TLS、CTLS三種方法求單應性矩陣(從實驗3.1發現DLS與TLS求解結果接近,下面不再單獨展示),用單應性矩陣對圖像拼接,拼接的結果可以從直觀上驗證H矩陣的求解精度。

圖2 合成數據下4種方法的平均坐標轉換損失對比與δH的F范數平方對比
由于SIFT、SURF等方法在清晰圖像上的特征點檢測比較準確,包含在特征點圖像坐標上的噪聲較小,經典的TLS可以將圖像拼接的較為準確。但對于模糊圖像,由于噪聲影響,特征點檢測變得不準確,拼接效果會打折扣,下面選取兩組圖像,第一組中相機位置只進行旋轉,第二組中相機位置變化主要是平移。對每組中第二幅圖采用不同的滑動窗口進行均值模糊(5-10),在特征點匹配的基礎上,分別用求單應性矩陣并拼接,為對比實驗效果,未對圖像接縫處采取平滑處理。
對第一組圖像,實驗發現,OLS對圖像的拼接結果存在錯位(框處),TLS雖能拼接準確,但右側建筑出現了傾斜和形變,CTLS拼接效果較好,如圖3、圖4所示。

圖3 待測圖像(右圖經過均值模糊(滑動窗口10*10))

圖4 不同方法對圖3的拼接結果(上:OLS,中:TLS,下:CTLS)
第二組圖像拼接對比結果如圖5與圖6所示,用OLS和TLS拼接結果都出現錯位,CTLS能對圖像拼接準確。

圖5 待測圖像(右圖經過均值模糊(滑動窗口5*5)(原圖來自ETH3D數據集[22]))

圖6 不同方法對圖5的拼接結果(上:OLS,中:TLS,下:CTLS)
在計算機視覺領域,單應性矩陣常用來表示兩個平面的坐標變換關系。傳統的基于RANSAC和TLS的求解方法在某些情況下的求解精度會下降,而當前關于特征點坐標噪聲較大時的求解方法研究較少。該文在文獻[4]和文獻[15]的基礎上,提出了一種基于CTLS的單應性矩陣求解方法,將像素坐標上存在高斯噪聲作為基本假設,考慮了噪聲矩陣列之間的相關性,并給出了具體求解步驟。分別在合成數據和真實圖像上驗證算法的求解精度,證明了在特征點噪聲誤差較大的情況下,基于CTLS的單應性矩陣求解要優于經典的TLS與OLS、DLS方法。另一方面,由于CTLS問題目前還沒有閉式解,只能通過梯度下降等優化方法求解,速度不如求解TLS等其他方法快。未來對于CTLS求解速度的提升具有一定的研究價值。