MVG（second）学习笔记- 三角形法

 

三角化，是已经知道相机矩阵和匹配点x-x',求出空间点。

线性三角形法

推导：

       线性三角形法，通常，被估计的点不能准确地满足儿何关系，因而不是一个最优估计.
       线性三角形法是直接模仿介绍的 DLT 方法.在每一幅图像上分别有x = PX ， x '= P ' X ，且这些方程可以组合成 AX = 0 的形式，它是关于 X 的线性方程.
首先通过叉乘消去齐次纯盘因子，使每一个图像点给出 三个方程，其中两个是线性独立的 ， 例如，对第一幅图像有 x X (PX) = 0 并把它展开得

其中 piT 是 P 的行 ，这些方程关于 X 的分量是线性 的.
可以组成形如 AX = 0 的方程，其中

       这里从每一像点取两个方程，总共绘出 4 个齐次朱知量的 4 个方程.这是一个冗余方程组，因为解在相差一个尺度因子的意义下被确定. 讨论过的求解形如AX = 0的方程组的两种方法.

（1）齐次方法

（2）非齐次方法

   几何法(见)（大佬博客）

（1）非线性优化法

（2）最优解法
 例子：

Triangulate in ORB-SLAM2

Initializer::Triangulate() 使用2d-2d点对 和 变换矩阵 R，t 三角变换恢复2d点对应的3d点

rianularization: 已知匹配特征点对{p1  p2} 和 
各自相机投影矩阵{P1 P2}, P1 = K*[I 0], P2 = K*[R t], 尺寸为[3,4]
估计三维点 X3D
     p1 = P1 * X3D
     p2 = P2 * X3D
采用直接线性变换DLT的方法(将式子变换成A*X=0的形式后使用SVD奇异值分解求解线性方程组)：
对于 p1 = P1 * X3D: 方程两边 左边叉乘 p1，可使得式子为0
p1叉乘P1*X3D = 0
其叉乘矩阵为:
叉乘矩阵 =  
	  |0  -1  y|
	  |1   0 -x| 
	  |-y  x  0| 
上述等式可写: 
  |0  -1  y|  |P1.row(0)|  
  |1   0 -x| *|P1.row(1)|* X3D = 0
  |-y  x  0|  |P1.row(2)|  

对于第一行 |0  -1  y| 会与P的三行分别相乘 得到四个值 与齐次3d点坐标相乘得到 0
有 (y * P1.row(2) - P1.row(1) ) * X3D = 0
对于第二行 |1   0 -x|有：
有 (x * P1.row(2) - P1.row(0) ) * X3D = 0
得到两个约束，另外一个点 p2 = P2 * X3D,也可以得到两个式子：
(y‘ * P2.row(2) - P2.row(1) ) * X3D = 0
(x’ * P2.row(2) - P2.row(0) ) * X3D = 0
写成 A*X = 0的形式有：
A =(维度4*4)
|y * P1.row(2) - P1.row(1) |
|x * P1.row(2) - P1.row(0) |
|y‘ * P2.row(2) - P2.row(1)|
|x’ * P2.row(2) - P2.row(0)|
对A进行奇异值分解求解X
cv::SVD::compute(A,w,u,vt,cv::SVD::MODIFY_A| cv::SVD::FULL_UV);
x3D = vt.row(3).t();// vt的最后一列,为X的解
x3D = x3D.rowRange(0,3)/x3D.at<float>(3);//  转换成非齐次坐标  归一化

源代码:

void Initializer::Triangulate(const cv::KeyPoint &kp1, const cv::KeyPoint &kp2, const cv::Mat &P1, const cv::Mat &P2, cv::Mat &x3D)
{
    cv::Mat A(4,4,CV_32F);

    A.row(0) = kp1.pt.x*P1.row(2)-P1.row(0);
    A.row(1) = kp1.pt.y*P1.row(2)-P1.row(1);
    A.row(2) = kp2.pt.x*P2.row(2)-P2.row(0);
    A.row(3) = kp2.pt.y*P2.row(2)-P2.row(1);

    cv::Mat u,w,vt;
    cv::SVD::compute(A,w,u,vt,cv::SVD::MODIFY_A| cv::SVD::FULL_UV);
    x3D = vt.row(3).t();
    x3D = x3D.rowRange(0,3)/x3D.at<float>(3);
}

谭平老师的ppt：

Mid‐point Algorithm

    

    参考：

1. https://blog.csdn.net/weixin_44580210/article/details/90679847
2. https://blog.csdn.net/kokerf/article/details/72844455
3. https://blog.csdn.net/u011178262/article/details/86729887
4. Multiple View Geometry in Computer Vision,Second Edition, 12.2 P312 & A5.2.1 P591
5. https://blog.csdn.net/xiaoxiaowenqiang/article/details/80871464
6. 谭平老师ppt