本文深入探讨双目摄像机定标和校正过程，详细解释了如何通过标定获得每个摄像机的内部参数以及两台摄像机之间的相对位置，包括旋转矩阵和平移向量。同时，介绍了如何利用极线约束减少匹配搜索范围，并阐述了OpenCV中cvStereoCalibrate函数的使用方法。此外，文章还讨论了如何计算Essential Matrix和Fundamental Matrix，以及双目校正的目的和实现方式，最后给出了双目校正后的显示效果改进方法。
http://blog.youkuaiyun.com/chenyusiyuan/article/details/5963256
三、雙目定標和雙目校正
雙目攝像頭定標不僅要得出每個攝像頭的內部參數，還需要通過標定來測量兩個攝像頭之間的相對位置（即右攝像頭相對於左攝像頭的三維平移 t 和旋轉 R 參數）。
圖6
要計算目標點在左右兩個視圖上形成的視差，首先要把該點在左右視圖上兩個對應的像點匹配起來。然而，在二維空間上匹配對應點是非常耗時的，為了減少匹配搜索范圍，我們可以利用極線約束使得對應點的匹配由二維搜索降為一維搜索。
圖7
而雙目校正的作用就是要把消除畸變後的兩幅圖像嚴格地行對應，使得兩幅圖像的對極線恰好在同一水平線上，這樣一幅圖像上任意一點與其在另一幅圖像上的對應點就必然具有相同的行號，只需在該行進行一維搜索即可匹配到對應點。
圖8
1． 關於cvStereoCalibrate的使用
如果按照 Learning OpenCV 的例程，直接通過cvStereoCalibrate來實現雙目定標，很容易產生比較大的圖像畸變，邊角處的變形較厲害。最好先通過cvCalibrateCamera2() 對每個攝像頭單獨進行定標，再利用cvStereoCalibrate進行雙目定標。這樣定標所得參數才比較准確，隨後的校正也不會有明顯的畸變。我使用的程序主要基於Learning OpenCV 的例程ch12_ex12_3.cpp，其中主要部分如下：
////////////////////////////////////////////////////////////////////////// // 是否首先進行單目定標？ cvCalibrateCamera2(&_objectPoints, &_imagePoints1, &_npoints, imageSize, &t_M1, &t_D1, NULL, NULL, CV_CALIB_FIX_K3); cvCalibrateCamera2(&_objectPoints, &_imagePoints2, &_npoints, imageSize, &t_M2, &t_D2, NULL, NULL, CV_CALIB_FIX_K3); ////////////////////////////////////////////////////////////////////////// // 進行雙目定標 cvStereoCalibrate( &_objectPoints, &_imagePoints1, &_imagePoints2, &_npoints, &t_M1, &t_D1, &t_M2, &t_D2, imageSize, &t_R, &t_T, &t_E, &t_F, cvTermCriteria(CV_TERMCRIT_ITER+ CV_TERMCRIT_EPS, 100, 1e-5)); // flags為默認的CV_CALIB_FIX_INTRINSIC
上面的t_M1(2), t_D1(2) 分別是單目定標後獲得的左（右）攝像頭的內參矩陣（3*3）和畸變參數向量（1*5）；t_R, t_T 分別是右攝像頭相對於左攝像頭的旋轉矩陣（3*3）和平移向量（3*1）, t_E是包含了兩個攝像頭相對位置關系的Essential Matrix（3*3），t_F 則是既包含兩個攝像頭相對位置關系、也包含攝像頭各自內參信息的Fundamental Matrix（3*3）。
圖9
2． cvStereoCalibrate 是怎樣計算 Essential Matrix 和 Fundamental Matrix 的？
首先我們以Learning OpenCV第422頁為基礎，討論下 Essential Matrix 和 Fundamental Matrix 的構造過程，再看看 OpenCV 是如何進行計算的。
圖10
注：原文中對pl、pr 和ql、qr 物理意義和計算公式的表述有誤，已修正。（2011-04-12）
（1）Essential Matrix
如上圖所示，給定一個目標點P，以左攝像頭光心Ol為原點。點P相對於光心Ol的觀察位置為Pl，相對於光心Or的觀察位置為Pr。點P在左攝像頭成像平面上的位置為pl，在右攝像頭成像平面上的位置為pr。注意Pl、Pr、pl、pr都處於攝像機坐標系，其量綱是與平移向量T相同的（pl、pr 在圖像坐標系中對應的像素坐標為 ql、qr）。
假設右攝像頭相對於左攝像頭的相對位置關系由旋轉矩陣R和平移向量T表示，則可得：Pr = R(Pl-T)。
現在我們要尋找由點P、Ol和Or確定的對極平面的表達式。注意到平面上任意一點x與點a的連線垂直於平面法向量n，即向量 (x-a) 與向量 n 的點積為0：(x-a)·n = 0。在Ol坐標系中，光心Or的位置為T，則P、Ol和Or確定的對極平面可由下式表示：(Pl-T)T·(Pl×T) = 0。
由Pr = R(Pl-T) 和 RT=R-1 可得：(RTPr)T·(Pl×T) = 0。
另一方面，向量的叉積又可表示為矩陣與向量的乘積，記向量T的矩陣表示為S，得：Pl×T = SPl。
圖11
那麼就得到：(Pr)TRSPl = 0。這樣，我們就得到Essential Matrix：E = RS。
通過矩陣E我們知道Pl和Pr的關系滿足：(Pr)TEPl = 0。進一步地，由 pl = fl*Pl/Zl 和 pr = fr*Pr/Zr 我們可以得到點P在左右兩個攝像機坐標系中的觀察點 pl 和 pr 應滿足的極線約束關系為：(pr)TEpl = 0。
注意到 E 是不滿秩的，它的秩為2，那麼 (pr)TEpl = 0 表示的實際上是一條直線，也就是對極線。
（2）Fundamental Matrix
由於矩陣E並不包含攝像頭內參信息，且E是面向攝像頭坐標系的。實際上我們更感興趣的是在圖像像素坐標系上去研究一個像素點在另一視圖上的對極線，這就需要用到攝像機的內參信息將攝像頭坐標系和圖像像素坐標系聯系起來。在（1）中，pl和pr是物理坐標值，對應的像素坐標值為ql和qr，攝像頭內參矩陣為M，則有：p=M-1q。從而：(pr)TEpl = 0 à qrT(Mr-1)TE Ml-1ql = 0。這裡，我們就得到Fundamental Matrix：F = (Mr-1)TE Ml-1。並有 qrTFql = 0。
（3）OpenCV的相關計算
由上面的分析可見，求取矩陣E和F關鍵在於旋轉矩陣R和平移向量T的計算，而cvStereoCalibrate的代碼中大部分（cvcalibration.cpp的第1886-2180行）也是計算和優化R和T的。在cvcalibration.cpp的第1913-1925行給出了計算R和T初始估計值的基本方法：
/* Compute initial estimate of pose For each image, compute: R(om) is the rotation matrix of om om(R) is the rotation vector of R R_ref = R(om_right) * R(om_left)' T_ref_list = [T_ref_list; T_right - R_ref * T_left] om_ref_list = {om_ref_list; om(R_ref)] om = median(om_ref_list) T = median(T_ref_list) */
具體的計算過程比較繁雜，不好理解，這裡就不討論了，下面是計算矩陣E和F的代碼：
if( matE || matF ) { double* t = T_LR.data.db; double tx[] = { 0, -t[2], t[1], t[2], 0, -t[0], -t[1], t[0], 0 }; CvMat Tx = cvMat(3, 3, CV_64F, tx); double e[9], f[9]; CvMat E = cvMat(3, 3, CV_64F, e); CvMat F = cvMat(3, 3, CV_64F, f); cvMatMul( &Tx, &R_LR, &E ); if( matE ) cvConvert( &E, matE ); if( matF ) { double ik[9]; CvMat iK = cvMat(3, 3, CV_64F, ik); cvInvert(&K[1], &iK); cvGEMM( &iK, &E, 1, 0, 0, &E, CV_GEMM_A_T ); cvInvert(&K[0], &iK); cvMatMul(&E, &iK, &F); cvConvertScale( &F, matF, fabs(f[8]) > 0 ? 1./f[8] : 1 ); } }
3． 通過雙目定標得出的向量T中，Tx符號為什麼是負的？
「@scyscyao：這個其實我也不是很清楚。個人的解釋是，雙目定標得出的T向量指向是從右攝像頭指向左攝像頭(也就是Tx為負)，而在OpenCV坐標系中，坐標的原點是在左攝像頭的。因此，用作校准的時候，要把這個向量的三個分量符號都要換一下，最後求出的距離才會是正的。
圖12
但是這裡還有一個問題，就是Learning OpenCV中Q的表達式，第四行第三列元素是-1/Tx，而在具體實踐中，求出來的實際值是1/Tx。這裡我和maxwellsdemon討論下來的結果是，估計書上Q表達式裡的這個負號就是為了抵消T向量的反方向所設的，但在實際寫OpenCV代碼的過程中，那位朋友卻沒有把這個負號加進去。」
圖13
scyscyao 的分析有一定道理，不過我覺得還有另外一種解釋：如上圖所示，攝像機C1（C2）與世界坐標系相對位置的外部參數為旋轉矩陣R1（R2）和平移向量 t1（t2），如果下標1代表左攝像機，2代表右攝像機，顯然在平移向量的水平分量上有 t1x > t2x；若以左攝像機C1為坐標原點，則可得到如上圖所示的旋轉矩陣R和平移向量t，由於t1x > t2x，則有 tx < 0。
為了抵消Tx為負，在矩陣Q中元素(4,3)應該加上負號，但在cvStereoRectify代碼中並沒有添加上，這就使得我們通過 cvReprojectImageTo3D 計算得到的三維數據都與實際值反號了。
if( matQ ) { double q[] = { 1, 0, 0, -cc_new[0].x, 0, 1, 0, -cc_new[0].y, 0, 0, 0, fc_new, 0, 0, 1./_t[idx], (idx == 0 ? cc_new[0].x - cc_new[1].x : cc_new[0].y - cc_new[1].y)/_t[idx] }; CvMat Q = cvMat(4, 4, CV_64F, q); cvConvert( &Q, matQ ); }
為了避免上述反號的情況，可以在計算得到Q矩陣後，添加以下代碼更改 Q[3][2] 的值。
// Q 是 Mat 類型矩陣，OpenCV2.1 C++ 模式 Q.at<double>(3, 2) = -Q.at<double>(3, 2); // Q 是 double 數組定義時 double Q[4][4]; CvMat t_Q = cvMat(4, 4, CV_64F, Q ); cvStereoRectify(…); Q[3][2]=-Q[3][2];
4． 雙目校正原理及cvStereoRectify 的應用。
圖14
如圖14所示，雙目校正是根據攝像頭定標後獲得的單目內參數據（焦距、成像原點、畸變系數）和雙目相對位置關系（旋轉矩陣和平移向量），分別對左右視圖進行消除畸變和行對准，使得左右視圖的成像原點坐標一致（CV_CALIB_ZERO_DISPARITY 標志位設置時發生作用）、兩攝像頭光軸平行、左右成像平面共面、對極線行對齊。在OpenCV2.1版之前，cvStereoRectify 的主要工作就是完成上述操作，校正後的顯示效果如圖14(c) 所示。可以看到校正後左右視圖的邊角區域是不規則的，而且對後續的雙目匹配求取視差會產生影響，因為這些邊角區域也參與到匹配操作中，其對應的視差值是無用的、而且一般數值比較大，在三維重建和機器人避障導航等應用中會產生不利影響。
因此，OpenCV2.1 版中cvStereoRectify新增了4個參數用於調整雙目校正後圖像的顯示效果，分別是 double alpha, CvSize newImgSize, CvRect* roi1, CvRect* roi2。下面結合圖15-17簡要介紹這4個參數的作用：
（1）newImgSize：校正後remap圖像的分辨率。如果輸入為(0,0)，則是與原圖像大小一致。對於圖像畸變系數比較大的，可以把newImgSize 設得大一些，以保留圖像細節。
（2）alpha：圖像剪裁系數，取值范圍是-1、0~1。當取值為 0 時，OpenCV會對校正後的圖像進行縮放和平移，使得remap圖像只顯示有效像素（即去除不規則的邊角區域），如圖17所示，適用於機器人避障導航等應用；當alpha取值為1時，remap圖像將顯示所有原圖像中包含的像素，該取值適用於畸變系數極少的高端攝像頭；alpha取值在0-1之間時，OpenCV按對應比例保留原圖像的邊角區域像素。Alpha取值為-1時，OpenCV自動進行縮放和平移，其顯示效果如圖16所示。
（3）roi1, roi2：用於標記remap圖像中包含有效像素的矩形區域。對應代碼如下：
02433 if(roi1) 02434 { 02435 *roi1 = cv::Rect(cvCeil((inner1.x - cx1_0)*s + cx1), 02436 cvCeil((inner1.y - cy1_0)*s + cy1), 02437 cvFloor(inner1.width*s), cvFloor(inner1.height*s)) 02438 & cv::Rect(0, 0, newImgSize.width, newImgSize.height); 02439 } 02440 02441 if(roi2) 02442 { 02443 *roi2 = cv::Rect(cvCeil((inner2.x - cx2_0)*s + cx2), 02444 cvCeil((inner2.y - cy2_0)*s + cy2), 02445 cvFloor(inner2.width*s), cvFloor(inner2.height*s)) 02446 & cv::Rect(0, 0, newImgSize.width, newImgSize.height); 02447 }
圖15
圖16
圖17
在cvStereoRectify 之後，一般緊接著使用 cvInitUndistortRectifyMap 來產生校正圖像所需的變換參數（mapx, mapy）。
////////////////////////////////////////////////////////////////////////// // 執行雙目校正 // 利用BOUGUET方法或HARTLEY方法來校正圖像 mx1 = cvCreateMat( imageSize.height, imageSize.width, CV_32F ); my1 = cvCreateMat( imageSize.height, imageSize.width, CV_32F ); mx2 = cvCreateMat( imageSize.height, imageSize.width, CV_32F ); my2 = cvCreateMat( imageSize.height, imageSize.width, CV_32F ); double R1[3][3], R2[3][3], P1[3][4], P2[3][4], Q[4][4]; CvMat t_R1 = cvMat(3, 3, CV_64F, R1); CvMat t_R2 = cvMat(3, 3, CV_64F, R2); CvMat t_Q = cvMat(4, 4, CV_64F, Q ); CvRect roi1, roi2; // IF BY CALIBRATED (BOUGUET'S METHOD) CvMat t_P1 = cvMat(3, 4, CV_64F, P1); CvMat t_P2 = cvMat(3, 4, CV_64F, P2); cvStereoRectify( &t_M1, &t_M2, &t_D1, &t_D2, imageSize, &t_R, &t_T, &t_R1, &t_R2, &t_P1, &t_P2, &t_Q, CV_CALIB_ZERO_DISPARITY, 0, imageSize, &roi1, &roi2); // Precompute maps for cvRemap() cvInitUndistortRectifyMap(&t_M1,&t_D1,&t_R1,&t_P1, mx1, my1); cvInitUndistortRectifyMap(&t_M2,&t_D2,&t_R2,&t_P2, mx2, my2);
5． 為什麼cvStereoRectify求出的Q矩陣cx, cy, f都與原來的不同？
「@scyscyao：在實際測量中，由於攝像頭擺放的關系，左右攝像頭的f, cx, cy都是不相同的。而為了使左右視圖達到完全平行對准的理想形式從而達到數學上運算的方便，立體校准所做的工作事實上就是在左右像重合區域最大的情況下，讓兩個攝像頭光軸的前向平行，並且讓左右攝像頭的f, cx, cy相同。因此，Q矩陣中的值與兩個instrinsic矩陣的值不一樣就可以理解了。」
注：校正後得到的變換矩陣Q，Q[0][3]、Q[1][3]存儲的是校正後左攝像頭的原點坐標（principal point）cx和cy，Q[2][3]是焦距f。
（待續……）
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
