单应性（Homography）变换与单应性矩阵的求解

单应性（Homography）变换

单应性变换的严格数学定义请参考：

• 《Multiple View Geometry in Computer Vision -2nd Edition》 by Richard Hartley, Andrew Zisserman
• 第2.3节Projective transformations

1. 概念

• 单应性变换又叫投影变换：应用在平面坐标变换中：

  • 平面投影变换是在三元素向量的齐次坐标下进行的线性变换，他由一个3×3的非奇异变换矩阵$H$表示，具体如下：
  • 
  • $x^{{}^{\prime }}$和$x$都是齐次坐标，一般的有：$x_3^{{}^{\prime }}=x_3=1$

• 单应性矩阵：

  • 单应矩阵描述两个平面上的对应点之间的变换关系；
  • 同一个平面在任意坐标系之间都可以建立单应性变换关系；
  • 如（a）：plannar surface上的X点可以通过单应性矩阵$H_1$和$H_2$变换到image1和image2，（b）和（c）同理。
  • 

2. 在CV方面的应用

• 图像校正、图像拼接、相机位姿估计、视觉SLAM等。

• 图像校正

  • 

• 视角变化

  • 

• 图像拼接

  • 
  • 

• AR

  • 

3. 求解单应性矩阵

3.1 假设

• 首先，我们假设两张图像中的对应点对齐次坐标为$P_{xy}=(x,y,1{)}^T$和$P_{uv}=(u,v,1{)}^T$，单应矩阵H定义为：
  • $H=\left[\begin{array}{ccc}{h}_{11}& h_{12}& h_{13}\\ h_{21}& h_{22}& h_{23}\\ h_{31}& h_{32}& h_{33}\end{array}\right]$
  • $\left[\begin{array}{c}u\\ v\\ 1\end{array}\right]=s\left[\begin{array}{ccc}{h}_{11}& h_{12}& h_{13}\\ h_{21}& h_{22}& h_{23}\\ h_{31}& h_{32}& h_{33}\end{array}\right]\left[\begin{array}{c}x\\ y\\ 1\end{array}\right](1)$

3.2 性质

• 这里使用的是齐次坐标系，也就是说可以进行任意尺度的缩放（s为尺度因子），也就是说把H乘以任意一个非零常数k并不改变上式结果，无非就是尺度因子s有所改变。
  • 比如我们把H乘以$\frac{1}{h_{33}}$可以得到：
    • $H^{\prime }=\frac{1}{h_{33}}H=\frac{1}{h_{33}}\left[\begin{array}{ccc}{h}_{11}& h_{12}& h_{13}\\ h_{21}& h_{22}& h_{23}\\ h_{31}& h_{32}& h_{33}\end{array}\right]=\left[\begin{array}{ccc}{h}_{11}^{\prime }& h_{12}^{\prime }& h_{13}^{\prime }\\ h_{21}^{\prime }& h_{22}^{\prime }& h_{23}^{\prime }\\ h_{31}^{\prime }& h_{32}^{\prime }& 1\end{array}\right]$
    • 而上述等式依然成立
  • 同理H乘以$\frac{1}{\sqrt{{\sum }_{i=1}^3{\sum }_{j=1}^3{h}_{ij}^2}}$ 可以得到约束：
    • $h_{11}^{\prime 2}+h_{12}^{\prime 2}+h_{13}^{\prime 2}+h_{21}^{\prime 2}+h_{22}^{\prime 2}+h_{23}^{\prime 2}+h_{31}^{\prime 2}+h_{32}^{\prime 2}+h_{33}^{\prime 2}=1$
    • 即${\sum }_{i=1}^3{\sum }_{j=1}^3{h}_{ij}^{\prime 2}=1$
  • 由此我们可以得出单应性矩阵有8个自由度，并非9个自由度。

3.3 求解

• 由3.1中假设公式（1）可得:
  • $\{\begin{array}{l}u=\frac{h_{11}x+h_{12}y+h_{13}}{h_{31}x+h_{32}y+h_{33}}\\ v=\frac{h_{21}x+h_{22}y+h_{23}}{h_{31}x+h_{32}y+h_{33}}\end{array}$，进一步变换得：
  • $\{\begin{array}{l}u(h_{31}x+h_{32}y+h_{33})=h_{11}x+h_{12}y+h_{13}\\ v(h_{31}x+h_{32}y+h_{33})=h_{21}x+h_{22}y+h_{23}\end{array}$，进一步得到：
  • $\{\begin{array}{l}x{h}_{11}+y{h}_{12}+h_{13}-ux{h}_{31}-uy{h}_{32}-u{h}_{33}=0\\ x{h}_{21}+y{h}_{22}+h_{23}-vx{h}_{31}-vy{h}_{32}-v{h}_{33}=0\end{array}$，化成矩阵形式有：
  • $\left[\begin{array}{ccccccccc}x& y& 1& 0& 0& 0& -ux& -uy& -u\\ 0& 0& 0& x& y& 1& -vx& -vy& -v\\ ...& & & & & & & & \\ ...& & & & & & & & \\ ...& & & & & & & & \\ ...& & & & & & & & \\ ...& & & & & & & & \\ ...& & & & & & & & \\ ...& & & & & & & & \end{array}\right]\left[\begin{array}{c}{h}_{11}\\ h_{12}\\ h_{113}\\ h_{21}\\ h_{22}\\ h_{23}\\ h_{31}\\ h_{32}\\ h_{33}\end{array}\right]=0$，更进一步抽象：
  • $\left[\begin{array}{ccc}{P}_{xy}^T& {\vec{0}}^T& -u{P}_{xy}^T\\ {\vec{0}}^T& P_{xy}^T& -v{P}_{xy}^T\end{array}\right]h=0$
    • 其中$h={\left[\begin{array}{ccccccccc}{h}_{11}& h_{12}& h_{113}& h_{21}& h_{22}& h_{23}& h_{31}& h_{32}& h_{33}\end{array}\right]}^T$
    • ${\vec{0}}^T=\left[\begin{array}{ccc}0& 0& 0\end{array}\right]$
  • 我们发现一对点可以提供两个方程
  • 由于单应矩阵H包含了$||H|{|}_F={\sum }_{i=1}^3{\sum }_{j=1}^3{h}_{ij}^{\prime 2}=1$或者$h_{33}=1$约束，因此根据上式的线性方程组8自由度的H我们至少需要4对点才能计算出单应矩阵。

4. 优化

但是，以上只是理论推导，在真实的应用场景中，我们计算的点对中都会包含噪声。比如点的位置偏差几个像素，甚至出现特征点对误匹配的现象，如果只使用4个点对来计算单应矩阵，那会出现很大的误差。因此，为了使得计算更精确，一般都会使用远大于4个点对来计算单应矩阵。另外上述方程组采用直接线性解法通常很难得到最优解，所以实际使用中一般会用其他优化方法如奇异值分解、Levenberg-Marquarat(LM)算法等进行求解，具体的优化过程可以参考张正友相机标定中得优化方法。

5. 推荐阅读

（超详细）张氏标定法原理及公式推导