《增强现实：原理、算法与应用》读书笔记（1）

入坑增强现实，学长让我先把这本书看懂，看了一小半，确实有很多对数学要求挺高的地方，所以打算写个读书笔记，分享一些体会。

第三章：实景的三维结构恢复与重建

双视图几何原理：基础矩阵与本质矩阵

首先，双视图几何原理需要理解，左图中的一个点会在右图中对应一条直线，这是因为从相机光心到图像上该点可以连出一条射线，该点的空间位置就在这条射线上，而这条射线在另一台相机的图像上的投影就是一条直线，称为极线。

两台相机光心的连线称为基线，基线与两个视图平面的交点称为极点eee,e′e'e′，所有极线lll,l′l'l′必过极点（很好理解，因为极线是光心出发的射线的投影，光心的投影就是极点，极线必过极点）。两个光心与空间中的一个三维点构成极平面πππ。

假设左右视图的投影矩阵分别为PPP和P′P'P′，相机内参矩阵分别为KKK和K′K'K′，以左视图为参考矩阵，相机光心位置C=0C=0C=0，则其投影矩阵可以表示为P=K[I∣0]P=K[I|0]P=K[I∣0]

右视图投影矩阵P′=K′[R∣t]P'=K'[R|t]P′=K′[R∣t]

相机光心位置C′=−R⊤tC'=-R^{\top }tC′=−R⊤t

R和t分别是右视图相机的旋转矩阵和平移向量。

这里可以顺带回顾下相机的投影矩阵P和内参矩阵K，其中，考虑线性畸变的内参矩阵K的一般形式如下：
K=[fkxfkytanαcx0fkycy001]K=\begin{bmatrix} \frac{f}{k_{x}} & \frac{f}{k_{y}}tan\alpha & c_{x}\\ 0 & \frac{f}{k_{y}} & c_{y}\\ 0 & 0 & 1 \end{bmatrix}K=⎣⎢⎡kxf00kyftanαkyf0cxcy1⎦⎥⎤

投影矩阵P=K[R∣t]P=K[R|t]P=K[R∣t]，其中t=−RC~t=-R\tilde{C}t=−RC~

透视投影矩阵描述了三维空间点与其在图像空间上的投影像素点之间的关系，即：x^∼PX^\hat{x}\sim P\hat{X}x^∼PX^

假设三维点X在左右视图上的投影点分别为xxx和x′x'x′，齐次坐标表示为x^\hat{x}x^和x^′\hat{x}'x^′，则有x^∼K[I∣0]X^\hat{x}\sim K[I|0]\hat{X}x^∼K[I∣0]X^x^′∼K′[R∣t]X^\hat{x}'\sim K'[R|t]\hat{X}x^′∼K′[R∣t]X^

xxx和x’x’x’对应的极线lll,l′l'l′，根据点在线上建立数学约束，可得：x^⊤l=0\hat{x}^{\top }l=0x^⊤l=0x^′⊤l′=0\hat{x}'^{\top }l'=0x^′⊤l′=0

证明：x^=(x,y,1)⊤\hat{x}=\left ( x,y,1 \right )^{\top }x^=(x,y,1)⊤，l:ax+by+c=0l:ax+by+c=0l:ax+by+c=0，剩下的应该不用写了吧。

极点eee是相机光心C′C'C′在左视图的投影点，极点e′e'e′是相机光心CCC在右视图的投影点，带入上面可得：
e=K[I∣0]⋅−R⊤t=−KR⊤te=K[I|0]\cdot-R^{\top}t=-KR^{\top}te=K[I∣0]⋅−R⊤t=−KR⊤t
e′=K′[R∣t]⋅R⊤t=K′te'=K'[R|t]\cdot R^{\top}t=K'te′=K′[R∣t]⋅R⊤t=K′t

对于图像上的任一点xxx，所有在光心CCC与xxx的射线上的点都会成像于xxx，光心位于投影矩阵的零向量，因此PC=0PC=0PC=0。光心CCC与xxx的射线上假定一个自由度，方向有P∗x^P^{*}\hat{x}P∗x^定义，其中P∗P^{*}P∗是PPP的广义逆，即P∗=P⊤(PP⊤)−1P^{*}=P^{\top }\left ( PP^{\top } \right )^{-1}P∗=P⊤(PP⊤)−1

其对应的三维点齐次坐标可以写成X^=sP∗x^+C\hat{X}=sP^{*}\hat{x}+CX^=sP∗x^+C

两边乘以矩阵PPP：

PX^=sPP∗x^+PC=sPP⊤(PP⊤)−1x^+0=sx^P\hat{X}=sPP^{*}\hat{x}+PC=sPP^{\top}\left ( PP^{\top} \right )^{-1}\hat{x}+0=s\hat{x}PX^=sPP∗x^+PC=sPP⊤(PP⊤)−1x^+0=sx^

与之前的x^∼PX^\hat{x}\sim P\hat{X}x^∼PX^相对应。

由于PPP的特殊性（即P=K[I∣0]P=K[I|0]P=K[I∣0]以及PP∗=IPP^{*}=IPP∗=I），可知
P∗=[K−10]P^{*}=\begin{bmatrix} K^{-1}\\ 0 \end{bmatrix}P∗=[K−10]

左视图的光心CCC在右视图的投影为极点e′e'e′，故此，x^\hat{x}x^与x^′\hat{x}'x^′之间存在几何关系：
x^′∼P′X^=P′(sP∗x^+C)=sP′P∗x^+e′\hat{x}'\sim P'\hat{X}=P'\left (sP^{*}\hat{x}+C \right )=sP'P^{*}\hat{x}+e'x^′∼P′X^=P′(sP∗x^+C)=sP′P∗x^+e′

令H=P′P∗H=P'P^{*}H=P′P∗，既得x^′∼sHx^+e′\hat{x}'\sim sH\hat{x}+e'x^′∼sHx^+e′

另外，由于极线l′l'l′经过e′e'e′和l′l'l′，因此有l′=e′×x^′l'=e'\times \hat{x}'l′=e′×x^′

向量叉乘的规则比较特殊，为了保证算法表述的一致性，尽量用矩阵算子来表达，设有两个三维向量a=(ax,ay,az)⊤a=\left( a_{x},a_{y},a_{z} \right )^{\top }a=(ax,ay,az)⊤与b=(bx,by,bz)⊤b=\left( b_{x},b_{y},b_{z} \right )^{\top }b=(bx,by,bz)⊤，定义3×33\times 33×3的斜对称矩阵[a]×[a]_{\times}[a]×为：
[a]×=[0−azayaz0−ax−ayax0][a]_{\times}=\begin{bmatrix} 0 & -a_{z} & a_{y}\\ a_{z} & 0 & -a_{x}\\ -a_{y} & a_{x} & 0 \end{bmatrix}[a]×=⎣⎡0az−ay−az0axay−ax0⎦⎤

则向量叉乘表示为矩阵乘积：a×b=[a]×ba\times b=[a]_{\times}ba×b=[a]×b

由aaa的形式易得[a]×[a]_{\times}[a]×的行列式为0，[a]×[a]_{\times}[a]×是秩为2的退化矩阵。将x^′∼sHx^+e′\hat{x}'\sim sH\hat{x}+e'x^′∼sHx^+e′与l′=[e′]×x^′l'=[e']_{\times}\hat{x}'l′=[e′]×x^′联合可得：
l′=[e′]×(sHx^+e′)=s[e′]×Hx^l'=[e']_{\times}\left(sH\hat{x}+e'\right)=s[e']_{\times}H\hat{x}l′=[e′]×(sHx^+e′)=s[e′]×Hx^

结合x^′⊤l′=0\hat{x}'^{\top }l'=0x^′⊤l′=0，消去sss可得
x^′⊤[e′]×Hx^=0\hat{x}'^{\top }[e']_{\times}H\hat{x}=0x^′⊤[e′]×Hx^=0

再令F=[e′]×HF=[e']_{\times}HF=[e′]×H，可得x^′⊤[e′]×Hx^=0\hat{x}'^{\top }[e']_{\times}H\hat{x}=0x^′⊤[e′]×Hx^=0

FFF就是著名的基础矩阵，它描述了两张图像上对应点的几何约束，即xxx，x′x'x′与两个相机的光心CCC，C′C'C′四点共面。

容易证明：
（1）矩阵FFF的各元素同时乘以一个非零常数，等式依然成立；
（2）FFF的秩为2，所以FFF是自由度为7、秩为2的3×33\times33×3矩阵。

之前已经推导到e′=K′te'=K'te′=K′t，H=P′P∗=K′[R∣t][K−10]=K′RK−1H=P'P^{*}=K'[R|t]\begin{bmatrix}K^{-1}\\ 0\end{bmatrix}=K'RK^{-1}H=P′P∗=K′[R∣t][K−10]=K′RK−1

因此若内参已知，可进一步得到如下等式：

F=[K′t]×K′RK−1=K′−⊤[t]×RK−1F=[K't]_{\times}K'RK^{-1}=K'^{-\top}[t]_{\times}RK^{-1}F=[K′t]×K′RK−1=K′−⊤[t]×RK−1

式中，E=[t]×RE=[t]_{\times}RE=[t]×R又被称为本质矩阵。

单应性矩阵

单应性矩阵是另一个比较常用的概念，尤其在快速匹配与定位时。

假设场景中所有点都落在一个平面上，那么不同视角拍摄的两张图像之间的映射关系可以用一个平面变换来表示（Malis et al., 2007）。而且，对于某平面π=(n⊤,d)⊤\pi=\left(n^{\top},d\right)^{\top}π=(n⊤,d)⊤上的任意一点XXX，都能满足以下约束等式：
n⊤X+d=0n^{\top}X+d=0n⊤X+d=0
这里的参数是在相机坐标系中给出的，ddd是相机光心至平面的距离，由于可见的平面一定不经过光心（经过光心的平面一定不可见），因此有d≠0d\neq0d=0，整理后可得约束式：
−n⊤Xd=1-\frac{n^{\top}X}{d}=1−dn⊤X=1

假设平面π\piπ上存在一个三维点X1X_{1}X1，被相机在两个不同的位置观测到，X1X_{1}X1分别投影到两个像素点x1=(u1,v1)⊤x_{1}=\left(u_{1},v_{1}\right)^{\top}x1=(u1,v1)⊤和x2=(u2,v2)⊤x_{2}=\left(u_{2},v_{2}\right)^{\top}x2=(u2,v2)⊤

假设图像1为参考帧，其相机坐标系为世界坐标系，图像2的投影矩阵为P=K2[R∣t]P=K_{2}[R|t]P=K2[R∣t]，由双视图几何关系可以得到：
x^2∼K2(RX1+t)\hat{x}_{2}\sim K_{2}\left(RX_{1}+t\right)x^2∼K2(RX1+t)

由于齐次坐标乘以比例因子不会改变物理意义，故此用∼\sim∼表示左右两边的相似关系，又因为X1X_{1}X1在平面上，因此满足约束式，代入可得：
x^2∼K2(RX1+t(−n⊤X1d))\hat{x}_{2}\sim K_{2}\left(RX_{1}+t\left(-\frac{n^{\top}X_{1}}{d}\right)\right)x^2∼K2(RX1+t(−dn⊤X1))

利用X1∼K1−1x^1X_{1}\sim K_{1}^{-1}\hat{x}_{1}X1∼K1−1x^1，稍作整理后可得：
x^2∼K2(R−tn⊤d)K1−1x^1\hat{x}_{2}\sim K_{2}\left(R-\frac{tn^{\top}}{d}\right)K_{1}^{-1}\hat{x}_{1}x^2∼K2(R−dtn⊤)K1−1x^1

令H21∼K2(R−tn⊤d)K1−1H_{21}\sim K_{2}\left(R-\frac{tn^{\top}}{d}\right)K_{1}^{-1}H21∼K2(R−dtn⊤)K1−1，得到一个3×33\times33×3的矩阵，称为单应性矩阵，这样就得到了 x^2∼H21x^1\hat{x}_{2}\sim H_{21}\hat{x}_{1}x^2∼H21x^1

对于同一平面上的任一点，单应性矩阵将一幅图像上的任意点映射到另一副图像上的对应点，也就是说，两张图像上处于同一平面上的所有匹配点之间的坐标变换，都可以用同一个单应性矩阵来表示。