相机标定
相机标定
相机标定是什么
相机的操作是将我们所处的三维空间中的某个场景通过某个模型转换成二维图像，并且这个3D→2D的过程是不可逆的。而相机标定的目标是寻找一个合适的数学模型，求出这个模型的参数，这样我们能够近似这个三维到二维的过程，使这个三维到二维的过程的函数找到反函数。整个逼近的过程就是相机标定，用简单的数学模型来表达复杂的成像过程，并且求出成像的反过程。标定之后的相机，可以进行三维场景的重建，即深度的感知。
在图像测量过程以及机器视觉应用中，为确定空间物体表面某点的三维几何位置与其在图像中对应点之间的相互关系，必须建立相机成像的几何模型，这些几何模型参数就是相机参数。在大多数条件下这些参数必须通过实验与计算才能得到，这个求解参数的过程就是相机标定。无论是在图像测量或者机器视觉应用中，相机参数的标定都是非常关键的环节，其标定结果的精度及算法的稳定性直接影响相机工作产生结果的准确性。
坐标系
世界坐标系
世界坐标系（World Coordinate）$(x_w,y_w,z_w)$，也称为测量坐标系，三维直角坐标系，以其为基准可以描述相机和待测物体的空间位置。世界坐标系的位置可以根据实际情况自由确定
相机坐标系
相机坐标系（Camera Coordinate）$(x_c,y_c,z_c)$，三维直角坐标系，原点位于镜头光心处，$x,y$轴分别与相面的两边平行，$z$轴为镜头光轴，与像平面垂直
像素坐标系
像素坐标系$uov$是一个二维直角坐标系，反映了相机CCD/CMOS芯片中像素的排列情况。原点$o$位于图像的左上角，$u,v$轴分别于像面的两边平行。像素坐标系中坐标轴的单位是像素（整数）
图像坐标系
像素坐标系不利于坐标变换，因此需要建立图像坐标系$XOY$，其坐标轴的单位通常为毫米，原点是相机光轴与相面的交点（称为主点），即图像的中心点，$X,Y$轴分别与$u,v$轴平行。故两个坐标系实际是平移关系，即可以通过平移就可得到
世界坐标系转换为相机坐标系
$$
\left[\begin{aligned}&x_c\&y_c\&z_c\&1\end{aligned}\right]=\left[\begin{aligned}&R &t\&0 &1\end{aligned}\right]\left[\begin{aligned}&x_w\&y_w\&z_w\&1\end{aligned}\right]
$$其中$R$为$3\times3$的旋转矩阵，$t$为$3\times1$的平移矢量，$[x_c,y_c,z_c,1]^T$为相机坐标系的齐次坐标，$[x_w,y_w,z_w,1]^T$为世界坐标系的齐次坐标
图像坐标系转换为像素坐标系
$$
\left[\begin{aligned}&u\&v\&1\end{aligned}\right]=\left[\begin{aligned}&\frac{1}{dX} &0& \ &u_0\&0 &\frac{1}{dY}& \ &v_0\&0 &0& \ &1\end{aligned}\right]\left[\begin{aligned}&X\&Y\&1\end{aligned}\right]
$$其中，$dX,dY$分别为像素在$X,Y$轴方向上的物理尺寸，$u_0,v_0$为主点（图像原点）坐标
针孔相机模型
相机成像模型简化为针孔模型，如下图所示。
空间任意一点$P$与其图像点$p$之间的关系，$P$与相机光心$o$的连线为$oP$，$oP$与像面的交点$p$即为空间点$P$在图像平面上的投影。该过程为透视投影，如下矩阵表示
$$
s\left[\begin{aligned}&X\&Y\&1\end{aligned}\right]=\left[\begin{aligned}&f &0& &0& &0\&0 &f& &0& &0 \&0 &0& &1& &0\end{aligned}\right]\left[\begin{aligned}&x\&y\&z\&1\end{aligned}\right]
$$其中$s$为缩放因子，$f$为有效焦距（光心到图像平面的距离），$[x,y,z,1]^T$是空间点$P$在相机坐标系$oxyz$中的齐次，$[X,Y,1]^T$是图像点$p$在图像坐标系$OXY$中的齐次坐标
世界坐标系转换为像素坐标系
$$
s\left[\begin{aligned}&u\&v\&1\end{aligned}\right]=\left[\begin{aligned}&\frac{1}{dX} &0& \ &u_0\&0 &\frac{1}{dY}& \ &v_0\&0 &0& \ &1\end{aligned}\right]\left[\begin{aligned}&f &0& &0& &0\&0 &f& &0& &0 \&0 &0& &1& &0\end{aligned}\right]\left[\begin{aligned}&R &t\&0 &1\end{aligned}\right]\left[\begin{aligned}&x_w\&y_w\&z_w\&1\end{aligned}\right]\=\left[\begin{aligned}&\alpha_x &0& &u_0& &0\&0 &\alpha_y& &v_0& &0 \&0 &0& &1& &0\end{aligned}\right]\left[\begin{aligned}&R &t\&0 &1\end{aligned}\right]\left[\begin{aligned}&x_w\&y_w\&z_w\&1\end{aligned}\right]=M_1M_2X_w=MX_w
$$其中，$\alpha_x=f/dX,\alpha_y=f/dY$，称为$u,v$轴的尺度因子，$M_1$称为相机的内部参数矩阵，$M_2$称为相机的外部参数矩阵，$M$称为投影矩阵
畸变参数
畸变（distortion）是对直线投影（rectilinear projection）的一种偏移。简单来说直线投影是场景内的一条直线投影到图片上也保持为一条直线。畸变简单来说就是一条直线投影到图片上不能保持为一条直线了，这是一种光学畸变（optical aberration）。
畸变一般可以分为径向畸变和切向畸变。径向畸变来自于透镜形状，切向畸变来自于整个摄像机的组装过程。畸变还有其他类型的畸变，但是没有径向畸变、切向畸变显著。
实际摄像机的透镜总是在成像仪的边缘产生显著的畸变，这种现象来源于“筒形”或“鱼眼”的影响。光线在原理透镜中心的地方比靠近中心的地方更加弯曲。对于常用的普通透镜来说，这种现象更加严重。对于径向畸变，成像仪中心（光学中心）的畸变为0，随着向边缘移动，畸变越来越严重。
切向畸变是由于透镜制造上的缺陷使得透镜本身与图像平面不平行而产生的。
张正友标定法
张正友标定法利用如下图所示的棋盘格标定板，在得到一张标定板的图像之后，可以利用相应的图像检测算法得到每一个角点的像素坐标$(u,v)$。
张正友标定法将世界坐标系固定于棋盘格上，则棋盘格上任一点的物理坐标$W=0$，由于标定板的世界坐标系是人为事先定义好的，标定板上每一个格子的大小是已知的，我们可以计算得到每一个角点在世界坐标系下的物理坐标$(U,V,W=0)$。
我们将利用这些信息：每一个角点的像素坐标$(u,v)$ 、每一个角点在世界坐标系下的物理坐标$(U,V,W=0)$，来进行相机的标定，获得相机的内外参矩阵、畸变参数。
标定相机的内外参数的思路
- 求解内参矩阵与外参矩阵的积
- 求解内参矩阵
- 求解外参矩阵
求解内参矩阵与外参矩阵的积
将世界坐标系固定于棋盘格上，则棋盘格上任一点的物理坐标$W=0$，因此，原单点无畸变的成像模型可以化为下式。其中， $R_1,R_2$为旋转矩阵$R$的前两列。为了简便，将内参矩阵记为$A$
$$
s\left[\begin{aligned}&u\&v\&1\end{aligned}\right]=\left[\begin{aligned}&\frac{f}{dX} &-\frac{f\cot\theta}{dX}& \ &u_0\&0 &\frac{f}{dY\sin\theta}& \ &v_0\&0 &0& \ &1\end{aligned}\right]\left[\begin{aligned}&R_1& R_2&&T\end{aligned}\right]\left[\begin{aligned}&U\&V\&1\end{aligned}\right]=A[R_1\ \ R_2\ \ T]\left[\begin{aligned}&U\&V\&1\end{aligned}\right]
$$对于不同的图片，内参矩阵$A$为定值；对于同一张图片，内参矩阵$A$，外参矩阵$[R_1\ \ R_2\ \ T]$为定值；对于同一张图片上的单点，内参矩阵$A$，外参矩阵，尺度因子$s$为定值
我们将$A[R_1\ \ R_2\ \ T]$记为矩阵$H$， $H$即为内参矩阵和外参矩阵的积，记矩阵$H$的三列为$[H_1,H_2,H_3]$，则有
$$
\left[\begin{aligned}&u\&v\&1\end{aligned}\right]=\frac{1}{s}H\left[\begin{aligned}&U\&V\&1\end{aligned}\right]=\frac{1}{s}\left[\begin{aligned}&H_{11} &H_{12}& &H_{13}\&H_{21} &H_{22}& &H_{23}\&H_{31} &H_{32}& &H_{33}\end{aligned}\right]\left[\begin{aligned}&U\&V\&1\end{aligned}\right]
$$利用上式，消去尺度因子$s$，可得
$$
\begin{aligned}u&=\frac{H_{11}U+H_{12}V+H_{13}}{H_{31}U+H_{32}V+H_{33}}\v&=\frac{H_{21}U+H_{22}V+H_{23}}{H_{31}U+H_{32}V+H_{33}}\end{aligned}
$$此时，尺度因子$s$已经被消去，因此上式对于同一张图片上所有的角点均成立。$(u,v)$是像素坐标系下的标定板角点的坐标，$(U,V)$是世界坐标系下的标定板角点的坐标。通过图像识别算法，我们可以得到标定板角点的像素坐标$(u,v)$，又由于标定板的世界坐标系是人为定义好的，标定板上每一个格子的大小是已知的，我们可以计算得到世界坐标系下的$(U,V)$
由这里的$H$是齐次矩阵，有8个独立未知元素。每一个标定板角点可以提供两个约束方程（$v,U,V$的对应关系、$u,U,V$的对应关系提供了两个约束方程），因此，当一张图片上的标定板角点数量等于4时，即可求得该图片对应的矩阵$H$。当一张图片上的标定板角点数量大于4时，利用最小二乘法回归最佳的矩阵$H$
求解内参矩阵
我们已知了矩阵$H
=A[R_1\ \ R_2\ \ T]$，接下来需要求解相机的内参矩阵$A$。我们利用$R_1,R_2$作为旋转矩阵$R$的两列，存在单位正交的关系，即：
$$
R_1^TR_2=0\R_1^TR_1=R_2^TR_2=I
$$则由$H$和$R_1,R_2$的关系，可知
$$
R_1=A^{-1}H_1\R_2=A^{-1}H_2
$$带入得到
$$
H_1^TA^{-T}A^{-1}H_2=0\H_1^TA^{-T}A^{-1}H_1=H_2^TA^{-T}A^{-1}H_2=I
$$另外，我们发现，上述两个约束方程中均存在矩阵$A^{-T}A^{-1}$。因此，我们记$A^{-T}A^{-1}=B$，则$B$为对称阵。我们试图先求解出矩阵 $B$，通过矩阵 $B$再求解相机的内参矩阵 $A$。同时，为了简便，我们记相机内参矩阵$A$为
$$
A=\left[\begin{aligned}&\frac{f}{dX} &-\frac{f\cot\theta}{dX}& \ &u_0& &0\&0 &\frac{f}{dY\sin\theta}& \ &v_0& &0\&0 &0& \ &1& &0\end{aligned}\right]=\left[\begin{aligned}&\alpha &\gamma& &u_0\&0 &\beta& &v_0\&0 &0& &1\end{aligned}\right]
$$则有
$$
A^{-1}=\left[\begin{aligned}&\frac{1}{\alpha} &-\frac{\gamma}{\alpha\beta}& &\frac{\gamma v_0-\beta u_0}{\alpha\beta}\&0 &\frac{1}{\beta}& &-\frac{v_0}{\beta}\&0 &0& &1\end{aligned}\right]
$$而$B$为
$$
B=A^{-T}A^{-1}=\left[\begin{aligned}&B_{11} &B_{12}& &B_{13}\&B_{21} &B_{22}& &B_{23}\&B_{31} &B_{32}& &B_{33}\end{aligned}\right]
$$通过$R_1,R_2$正交的关系，需要计算$H_i^TBH_j$来得到最终的矩阵$B$，令
$$
\begin{aligned}v_{ij}&=[H_{1i}H_{1j}\ \ H_{1i}H_{2j}+H_{2i}H_{1j}\ \ H_{1i}H_{3j}+H_{3i}H_{1j}\ \ H_{2i}H_{3j}+H_{3i}H_{2j}\ \ H_{3i}H_{3j}]^T\ b&=[B_{11}\ \ B_{12}\ \ B_{22}\ \ B_{13}\ \ B_{23}\ \ B_{33}]^T\end{aligned}
$$有$H_i^TBH_j=v_{ij}b$，$R_1,R_2$正交方程简化为
$$
v_{12}^Tb=0\v_{11}^Tb=v_{22}^Tb=I
$$即
$$
\left[\begin{aligned}v_{12}^T&\v_{11}^T-&v_{22}^T\end{aligned}\right]b=0
$$由于矩阵$H$已知，矩阵$v=\left[\begin{aligned}v_{12}^T&\v_{11}^T-&v_{22}^T\end{aligned}\right]$$\left[\begin{aligned}v_{12}^T&\v_{11}^T-&v_{22}^T\end{aligned}\right]$又全部由矩阵$H$的元素构成，因此矩阵$v$已知。
此时，我们只要求解出向量$b$，即可得到矩阵$B$。每张标定板图片可以提供一个$vb=0$的约束关系，该约束关系含有两个约束方程。但是，向量$b$有6个未知元素。因此，单张图片提供的两个约束方程是不足以解出来向量$b$。因此，我们只要取3张标定板照片，得到3个$vb=0$的约束关系，即6个方程，即可求解向量$b$。当标定板图片的个数大于3时（事实上一般需要15到20张标定板图片），可采用最小二乘拟合最佳的向量$b$，并得到矩阵$B$。
根据矩阵$B$的元素和相机内参$\alpha,\beta,\gamma,u_0,v_0$的对应关系，可得到
$$
v_0=\frac{B_{12}B_{13}-B_{11}B_{23}}{B_{11}B_{22}-B_{12}^2}\ \alpha=\sqrt{\frac{1}{B_{11}}}\ \beta=\sqrt{\frac{B_{11}}{B_{11}B_{22}-B_{12}^2}}\ \gamma=-B_{12}\alpha^2\beta \ u_0=\frac{\gamma v_0}{\beta}-B_{13}\alpha^2
$$即可求得内参矩阵$A=\left[\begin{aligned}&\frac{1}{dX} &-\frac{f\cot\theta}{dX}& \ &u_0& &0\&0 &\frac{1}{dY\sin\theta}& \ &v_0& &0\&0 &0& \ &1& &0\end{aligned}\right]=\left[\begin{aligned}&\alpha &\gamma& &u_0\&0 &\beta& &v_0\&0 &0& &1\end{aligned}\right]$
求解外参矩阵
对于同一个相机，相机的内参矩阵取决于相机的内部参数，无论标定板和相机的位置关系是怎么样的，相机的内参矩阵不变。这也正是在求解内参矩阵中，我们可以利用不同的图片（标定板和相机位置关系不同）获取的矩阵$H$，共同求解相机内参矩阵$A$的原因。
但是，外参矩阵反映的是标定板和相机的位置关系。对于不同的图片，标定板和相机的位置关系已经改变，此时每一张图片对应的外参矩阵都是不同的。
在关系$H
=A[R_1\ \ R_2\ \ T]$中，我们已经求解得到了矩阵$H$（对于同一张图片相同，对于不同的图片不同）、矩阵$A$（对于不同的图片都相同）。通过公式$[R_1\ \ R_2\ \ T]=A^{-1}H$，即可求得每一张图片对应的外参矩阵$[R_1\ \ R_2\ \ T]$。注意，这里值得指出，完整的外参矩阵为$\left[\begin{aligned}&R &T\&0 &1\end{aligned}\right]$。但是，由于张正友标定板将世界坐标系的原点选取在棋盘格上，则棋盘格上任一点的物理坐标$W=0$，将旋转矩阵的$R$的第三列$R_3$消掉，因此，$R_3$在坐标转化中并没有作用。但是$R_3$要使得$R$满足旋转矩阵的性质，即列与列之间单位正交，因此可以通过向量$R_1,R_2$的叉乘，即$R_3=R_1\times R_2$计算得到。
此时，相机的内参矩阵和外参矩阵均已得到。
注：以上推导都是假设不存在畸变参数的情况下成立的。但是事实上，相机是存在畸变参数的，因此，张正友标定法还需要通过L-M算法对于参数进行迭代优化。
标定步骤
- 准备一个张正友标定法的棋盘格，棋盘格大小已知，用相机对其进行不同角度的拍摄，得到一组图像
- 对图像中的特征点如标定板角点进行检测，得到标定板角点的像素坐标值，根据已知的棋盘格大小和世界坐标系原点，计算得到标定板角点的物理坐标值
- 求解内参矩阵与外参矩阵。根据物理坐标值和像素坐标值的关系，求出$H$矩阵，进而构造$v$矩阵，求解$B$矩阵，利用$B$矩阵求解相机内参矩阵$A$，最后求解每张图片对应的相机外参矩阵$\left[\begin{aligned}&R& t\&0&1\end{aligned}\right]$
- 求解畸变参数。利用$\hat{u},u,\hat{v},v$构造$D$矩阵，计算径向畸变参数
- 利用L-M（Levenberg-Marquardt）算法对上述参数进行优化