Stitching模塊中focalsFromHomography初步研究

在Stitching模塊中，通過“光束法平差”的時候，有一個步驟為“通過單應矩陣估算攝像頭焦距”，調用的地方為：

 

void focalsFromHomography( const Mat & H,  double  &f0,  double  &f1,  bool  &f0_ok,  bool  &f1_ok)
{
    CV_Assert(H.type()  == CV_64F  && H.size()  == Size( 3,  3));
     const  double * h  = H.ptr < double >();
     double d1, d2;  // Denominators
     double v1, v2;  // Focal squares value candidates
    f1_ok  =  true;
    d1  = h[ 6]  * h[ 7];
    d2  = (h[ 7]  - h[ 6])  * (h[ 7]  + h[ 6]);
    v1  =  -(h[ 0]  * h[ 1]  + h[ 3]  * h[ 4])  / d1;
    v2  = (h[ 0]  * h[ 0]  + h[ 3]  * h[ 3]  - h[ 1]  * h[ 1]  - h[ 4]  * h[ 4])  / d2;
     if (v1  < v2) std : :swap(v1, v2);
     if (v1  >  0  && v2  >  0) f1  = std : :sqrt(std : :abs(d1)  > std : :abs(d2)  ? v1  : v2);
     else  if (v1  >  0) f1  = std : :sqrt(v1);
     else f1_ok  =  false;
    f0_ok  =  true;
    d1  = h[ 0]  * h[ 3]  + h[ 1]  * h[ 4];
    d2  = h[ 0]  * h[ 0]  + h[ 1]  * h[ 1]  - h[ 3]  * h[ 3]  - h[ 4]  * h[ 4];
    v1  =  -h[ 2]  * h[ 5]  / d1;
    v2  = (h[ 5]  * h[ 5]  - h[ 2]  * h[ 2])  / d2;
     if (v1  < v2) std : :swap(v1, v2);
     if (v1  >  0  && v2  >  0) f0  = std : :sqrt(std : :abs(d1)  > std : :abs(d2)  ? v1  : v2);
     else  if (v1  >  0) f0  = std : :sqrt(v1);
     else f0_ok  =  false;
}

本文具體分析 focalsFromHomography，函數的參數定義：

     Tries to estimate focal lengths from the given homography  
     under the assumption that the camera undergoes rotations around its centre only.    
     Parameters  
     H – Homography.  
     f0 – Estimated focal length along X axis.  
     f1 – Estimated focal length along Y axis.  
     f0_ok – True,  if f0 was estimated successfully,  false otherwise.  
     f1_ok – True,  if f1 was estimated successfully,  false otherwise.  

 

可以看到它通過輸入的單應矩陣，最后得到了相機焦距的估計值，計算的過程也比較復雜。那這樣做的理由是什么了？具體計算的時候又是如何實現的了？

 

論文也就是算法的依據為《Construction of Panoramic Image Mosaics with Global and Local Alignment 》, Heung-Yeung Shum ([e-mail]) and Richard Szeliski ([e-mail])  page 17.method "focals from homgraphy matrix" 

我將具體的內容截出來：

 

 

原論文中40-44的推導，分為兩個部分。一個部分是從“8參數”的變換，得出和x軸，y軸兩個方向焦距的關系；一個部分是通過行列式的數學性質，計算出兩個方向的焦距。這兩個部分我目前都沒有掌握足夠的資料來進行證明，如果有能夠證明的同學麻煩聯系我一下。

然后來看算法實現。如果認為論文的表述是正確的，那么依據數學函數來對比c++的實現

代碼中的h0-h8直接對應論文中的m0-m8,僅以f0來觀察，那么f0^2可能有兩種取值（這里x^2 是 x * x的一種簡單表示方法，代表階乘）

f0^2 = - m2*m5/(m0*m3+m1*m4)或 f0^2 = m5^2 - m2^2 /(m0^2 + m1^2 - m3^2 - m4 ^2) 

看代碼

 

d1  = h[ 0]  * h[ 3]  + h[ 1]  * h[ 4];
那么
 v1  =  -h[ 2]  * h[ 5]  / d1  =   -h[ 2]  * h[ 5]  /(h[ 0]  * h[ 3]  + h[ 1]  * h[ 4])

 

 d2  = h[ 0]  * h[ 0]  + h[ 1]  * h[ 1]  - h[ 3]  * h[ 3]  - h[ 4]  * h[ 4];
則
 v2  = (h[ 5]  * h[ 5]  - h[ 2]  * h[ 2])  / d2  = (h[ 5]  * h[ 5]  - h[ 2]  * h[ 2])  / （h[ 0]  * h[ 0]  + h[ 1]  * h[ 1]  - h[ 3]  * h[ 3]  - h[ 4]  * h[ 4]）

 

前后是一一對應的，計算f0^2的兩個值是沒有問題的。但是這里f0有兩個計算結果，最后選擇哪個了？這一點在論文中沒有說，在代碼中采用的方法是首先判斷v1,v2的符號，如果都是負數，那么肯定是計算錯誤了，因為它們所代表的f0^2肯定是非負數；然后判斷v1,v2的大小，取其中比較大的那個來進行計算。

但是在

if (v1  >  0  && v2  >  0) f0  = std : :sqrt(std : :abs(d1)  > std : :abs(d2)  ? v1  : v2);

我認為這樣寫是沒有用的，我也在嘗試聯系一下相關對這個問題比較熟悉的人共同討論。f1的計算方法是同樣的。到這里已經得到f0和f1，分別對應x軸和y軸，為了得到最后的結果，那么會取

f = sqrt(f0 * f1) 

則得到這個當應矩陣的對於焦距的估計值。那么focalsFromHomography的一次運算也就結束了。