《A Flexible New Technique for Camera Calibration》

论文简述

总结

输入：

标定板物理坐标系中的角点坐标（p×n×3=姿态数×角点数×3 => z=0）
图像像素坐标系中的角点坐标（p×n×2=姿态数×角点数×2）

输出：

相机内参（3×3）
畸变系数
相机外参（p×4×4=姿态数×4×4）

总体步骤：

估算单应性矩阵H（封闭解）
计算相机内参
计算径向畸变k1,k2,相机外参
最大似然估计

概述

动机：从二维图像中提取几何信息
（我的理解：①像素与物理之间的尺度变化②纠正畸变）

要求：相机从不同方位（至少两个）来观测标定板（贴到平面上），不必知道相机或者标定板的位置

张氏标定法介于摄影测量标定法和自标定法之间。

标定方法	利用信息	特点
张氏标定法	利用二维几何信息	灵活、鲁棒性好
摄影测量标定法 Photogrammetric calibration	三维信息	有效，但标定装置昂贵
自标定法 Self-calibration	纯粹的特征信息	灵活，但鲁棒性不好

基础方程

符号

$s\ \widetilde{m} = A\ [R|t]\ \widetilde{M} \tag{1}$ $实数×(3*1)\ =\ (3*3)×(3*4)×(4*1)$

含义：3D点 $M$ 和它的图像投影点 $m$ 之间的关系
$m=[u,v]^T$ ：2D点
$\widetilde{m}=[u,v,1]^T$ ：2D点的齐次形式 —— 图像像素坐标系
$M=[X,Y,Z]^T$ ：3D点
$\widetilde{M}=[X,Y,Z,1]^T$ ：3D点的齐次形式 —— 标定板物理坐标系
$s$ ：任意的比例因子
$[R|t]$ ：相机外参（联系世界坐标系和相机坐标系）
：相机内参
- $(u_0,v_0)$ ：主点的坐标
- $α,β$ ：图像在u轴和v轴的比例因子
- $γ$ ：描述两个坐标轴倾斜角的参数

标定板平面与其图像之间的单应性（Homograph）

$s\ \begin{bmatrix} u \\ v \\ 1 \\ \end{bmatrix} = A\ [r_1\ r_2\ r_3\ t]\ \begin{bmatrix} X \\ Y \\ 0 \\ 1 \\ \end{bmatrix} \\= A\ [r_1\ r_2\ t]\ \begin{bmatrix} X \\ Y \\ 1 \\ \end{bmatrix}$ $实数×(3*1)\ =\ (3*3)×(3*4)×(4*1) =\ (3*3)×(3*3)×(3*1)$

假设标定板平面在世界坐标系中Z=0
旋转矩阵 $R$ 的第i列记为 $r_i$

$简写为： s\ \widetilde{m} = H\ \widetilde{M} \\ H=A\ [r_1\ r_2\ t] \tag{2}$ $实数×(3*1)\ =\ (3*3)×(3*1)$

内参的约束条件

数学基础

旋转矩阵的性质：旋转矩阵是正交矩阵

行列式为1
各行是单位向量且两两正交
各列是单位向量且两两正交[1]
逆 = 转置

由“各列是单位向量且两两正交”可得（ $r_1$ 、 $r_2$ 为旋转矩阵的第1列、第2列）：

$r_1^T\ r_2 = 0 \\ r_1^T\ r_1 = r_2^T\ r_2 = 1$

矩阵的性质（转置、逆）：

$(A\ B)^T = B^T\ A^T \\ (A\ B)^{-1} = B^{-1}\ A^{-1}$

公式推导

把 $H$ 记为 $H=[h_1\ h_2\ h_3]$ ，由公式(2)可得（λ是任意的标量）：

$[h_1\ h_2\ h_3] = λ A [r_1\ r_2\ T]$

可得：

$\begin{cases} h_1 = λ\ A\ r_1 \\ h_2 = λ\ A\ r_2 \\ h_3 = λ\ A\ t \\ \end{cases}$ $\begin{cases} r_1 = \frac{1}{λ}\ A^{-1}\ h_1 \\ r_2 = \frac{1}{λ}\ A^{-1}\ h_2 \\ t = \frac{1}{λ}\ A^{-1}\ h_3 \\ \end{cases}$

代入“单位向量+两两正交”的公式：

$\begin{cases} r_1^T\ r_2 = (\frac{1}{λ}\ A^{-1}\ h_1)^T\ (\frac{1}{λ}\ A^{-1}\ h_2) = 0 \\ r_1^T\ r_1 = (\frac{1}{λ}\ A^{-1}\ h_1)^T\ (\frac{1}{λ}\ A^{-1}\ h_1) = r_2^T\ r_2 = (\frac{1}{λ}\ A^{-1}\ h_2)^T\ (\frac{1}{λ}\ A^{-1}\ h_2) = 1 \\ \end{cases}$

整理可得：

$h_1^T\ A^{-T}\ A^{-1}\ h_2 = 0 \tag{3}$ $h_1^T\ A^{-T}\ A^{-1}\ h_1 = h_2^T\ A^{-T}\ A^{-1}\ h_2 \tag{4}$

几何解释是绝对二次曲线（absolute conic），见论文P4 2.4 Geometric Interpretation。（此部分没看懂）

实现功能

实用代码

matlab的符号计算

% 推导《A Flexible New Technique for Camera Calibration》 P5 3.1 Colsed-form solution 中的公式(5) B = A^(-T) * A^(-1)
syms a b c u v;
A = [a, c, u; 0, b, v; 0, 0, 1]; % 相机内参
inv(A); % 逆
A.'; % 转置
B = inv(A).' * inv(A);

% 输出结果
A % 相机内参
B % 对称矩阵

% 结果：
% A =
%     [ a, c, u]
%     [ 0, b, v]
%     [ 0, 0, 1]
% B =
%     [                1/a^2,                        -c/(a^2*b),                  -(b*u - c*v)/(a^2*b)]
%     [           -c/(a^2*b),             1/b^2 + c^2/(a^2*b^2),     (c*(b*u - c*v))/(a^2*b^2) - v/b^2]
%     [ -(b*u - c*v)/(a^2*b), (c*(b*u - c*v))/(a^2*b^2) - v/b^2, v^2/b^2 + (b*u - c*v)^2/(a^2*b^2) + 1]

% 推导《A Flexible New Technique for Camera Calibration》 P5 3.1 Colsed-form solution 中的公式(7)、公式(8)
syms B11 B12 B22 B13 B23 B33;
syms h11 h12 h13 h21 h22 h23 h31 h32 h33;
B = [B11, B12, B13; B12, B22, B23; B13, B23, B33]; % 维度=(3*3)
b = [B11, B12, B22, B13, B23, B33].'; % 维度=(6*1)
h1 = [h11, h12, h13].'; % 单应性矩阵(3*3)的列向量 维度=(3*1)
h2 = [h21, h22, h23].';
v12 = [h11*h21, h11*h22+h12*h21, h12*h22, h13*h21+h11*h23, h13*h22+h12*h23, h13*h23].'; % 维度=(6*1)
v11 = [h11*h11, h11*h12+h12*h11, h12*h12, h13*h11+h11*h13, h13*h12+h12*h13, h13*h13].'; % 维度=(6*1)
v22 = [h21*h21, h21*h22+h22*h21, h22*h22, h23*h21+h21*h23, h23*h22+h22*h23, h23*h23].'; % 维度=(6*1)

% 第1个等式：i=1, j=2
equation1_left  = h1.' * B * h2;
equation1_right = v12.' * b;

% 第2个等式：i=1, j=1; i=2, j=2
equation2_left  = h1.' * B * h1 - h2.' * B * h2;
equation2_right = (v11-v22).' * b;

% 输出结果
equation1_left  % 第1个等式的左边
equation1_right % 第1个等式的右边
equation2_left  % 第2个等式的左边
equation2_right % 第2个等式的右边

参考

题目	备注
[1]正交矩阵和旋转矩阵之间关系和性质总结	总结旋转矩阵的性质，各列是单位向量且两两正交
[2][图像]摄像机标定(2) 张正友标定推导详解	精华文章，包含张氏相机标定的公式详细推导
[3]张正友标定法	精华文章，包含张氏相机标定的公式详细推导