相机标定，是图像测量和机器视觉应用时，绕不过去的关键步骤。通过标定，可以获得相机成像几何模型的参数，也就是三维空间中点与二维图像中点的对应关系。

1 相机成像模型

相机成像，实际上是一个光学成像过程。将相机的镜头组看作一个凸透镜，光线通过透镜在感光元件（CCD/CMOS）上成像，感光元件将光电信号转换为数字信号，再经数字信息处理（DSP）成数字图像，存储到存储介质当中。

透镜成像原理如下图，凸透镜的中心为光心，光线平行于主光轴（虚线）穿过透镜时，会汇聚到焦点，然后折射成像。其中，a为实物，u为物距；a1为成像，v为相距；f为焦距，表示焦点到光心的距离。

透镜成像原理

当相机感光元件位于凸透镜焦点附近，焦距与光心到感光元件距离无限接近时，即f≈vf\approx v，相机成像模型就成了我们熟悉的「小孔成像」。

小孔成像原理

2 坐标系转换

相机成像系统中，共包含四个坐标系：世界坐标系、相机坐标系、图像坐标系、像素坐标系。相机标定的主要工作就是求解三维世界坐标系空间的点到二维像平面的转换关系，总的转换公式如下：

zc[uv1]=[fdX0u000fdYv000010][RT01][xcyczc1]z_{c}\left[\begin{array}{l}{u} \\ {v} \\ {1}\end{array}\right]=\left[\begin{array}{cccc}{\frac{f}{d X}} & {0} & {u_{0}} & {0} \\ {0} & {\frac{f}{d Y }} & {v_{0}} & {0} \\ {0} & {0} & {1} & {0}\end{array}\right]\left[\begin{array}{ll}{R} & {T} \\ {0} & {1}\end{array}\right]\left[\begin{array}{l}{x_{c}} \\ {y_{c}} \\ {z_{c}} \\ {1}\end{array}\right] \\

相机成像模型

  1. 世界坐标系到相机坐标系

    世界坐标系到相机坐标系为同一个点从一个坐标系到另一个坐标系的转换（刚体变换），主要涉及到两个坐标系之间的旋转平移矩阵。

    这里的R就是3*3的旋转矩阵，T是3*1的平移矩阵

    把1*3的坐标变成1*4的原因是使用齐次坐标。原因主要是三维笛卡尔坐标与矩阵的乘法只能实现三维坐标的缩放和旋转，无法实现坐标平移。添加一个额外坐标，就可以实现坐标平移了，而且保持了三维向量与矩阵乘法具有的缩放和旋转操作。

    [xcyczc1]=[RT01][xwywzw1]\left[\begin{array}{c} x_c \\ y_c \\ z_c \\ 1 \end{array}\right]=\left[\begin{array}{ll} R & T \\ 0 & 1 \end{array}\right]\left[\begin{array}{c} x_w \\ y_w \\ z_w \\ 1 \end{array}\right]

trans1

  1. 相机坐标系到图像坐标系

    相机坐标系到图像坐标系是三维坐标到二维坐标的透视投影。

    x=f⋅xczc,y=f⋅yczczc[xy1]=[f0000f000010][xcyczc1]\begin{gathered} x=f \cdot \frac{x_c}{z_c}, y=f \cdot \frac{y_c}{z_c} \\ z_c\left[\begin{array}{l} x \\ y \\ 1 \end{array}\right]=\left[\begin{array}{llll} f & 0 & 0 & 0 \\ 0 & f & 0 & 0 \\ 0 & 0 & 1 & 0 \end{array}\right]\left[\begin{array}{c} x_c \\ y_c \\ z_c \\ 1 \end{array}\right] \end{gathered}

trans2

  1. 图像坐标系到像素坐标系

    最后图像坐标系到像素坐标系是两个二维坐标系的平移（一个在像平面中心一个在左上角）和缩放（用到每个像素在x和y方向的物理尺寸）。

    {u=u0+xdxv=v0+ydy[uv1]=[1dx0u001dyv0001][xy1]\begin{gathered} \left\{\begin{array}{l} u=u_0+\frac{x}{d x} \\ v=v_0+\frac{y}{d y} \end{array}\right. \\ {\left[\begin{array}{l} u \\ v \\ 1 \end{array}\right]=\left[\begin{array}{lll} \frac{1}{d x} & 0 & u_0 \\ 0 & \frac{1}{d y} & v_0 \\ 0 & 0 & 1 \end{array}\right]\left[\begin{array}{l} x \\ y \\ 1 \end{array}\right]} \end{gathered}

trans3

  1. 相机畸变

    真实的镜头在成像时，会有畸变存在。由透镜形状引起的畸变是径向畸变，而由透镜安装与成像平面不平行引起的畸变是切向畸变，图像上的点离中心越远畸变也会越严重。

    相机畸变由$(k_1,k_2,k_3,p_1,p_2)5个参数描述，而这些参数都可以使用相机标定来求解。质量较好的相机，切向畸变小到可忽略，径向畸变系数5个参数描述，而这些参数都可以使用相机标定来求解。质量较好的相机，切向畸变小到可忽略，径向畸变系数k_3也可忽略，只计算也可忽略，只计算k_1,k_2$两个参数就可以了。

x^=x(1+k1r2+k2r4+k3r6+2p1y+p2(r2+2x2))y^=y(1+k1r2+k2r4+k3r6+p1(r2+2y2)+2p2x)\hat{x} = x(1+k_1r^2+k_2r^4+k_3r^6+2p_1y+p_2(r^2+2x^2)) \\ \hat{y} = y(1+k_1r^2+k_2r^4+k_3r^6+p_1(r^2+2y^2)+2p_2x)

3 张正友标定法

相机标定过程最常用的就是的张正友标定法（张氏标定法），操作简单、精度较高，可以满足大部分场合。

张正友标定法，是利用棋盘格标定板进行标定，将世界坐标系固定在棋盘格上。棋盘格上的每个格子大小都是已知的，即棋盘格每一个角点在世界坐标系下的坐标都可以计算得到。

张正友标定法标定相机的内外参数的思路如下：1）求解内参矩阵与外参矩阵的积；2）求解内外参矩阵；3）求解畸变参数；4）参数优化。

3.1 求解单应性矩阵

利用平面坐标和空间坐标的对应关系，可以计算单应性矩阵。假设我们的棋盘平面标定模板落在世界坐标系的Zw=0Z_{w}=0的平面上（式子中A为内参矩阵，H为单应性矩阵，这里的单应性矩阵就是内参矩阵与外参矩阵的乘积）。

单应性矩阵求解的是一个平面到另外一个平面的投影映射。

s[uv1]=A[RT][xwyw01]=H[xwyw1]=[H11H12H13H21H22H23H31H32H33][xwyw1]s\left[\begin{array}{c} u \\ v \\ 1 \end{array}\right]=A\begin{bmatrix}R&T\end{bmatrix}\left[\begin{array}{c} x_w \\ y_w \\ 0 \\ 1 \end{array}\right]=H \left[\begin{array}{c} x_w \\ y_w \\ 1 \end{array}\right]=\left[\begin{array}{lll} 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_w \\ y_w \\ 1 \end{array}\right]

消除尺度因子s（尺度因子对于同一张图片的不同角点是不同的）得到下面的式子，这里的式子对同一张图片上的所有角点都是成立的，且(u,v)(u,v)可以通过角点检测算法得到，(xw,yw)(x_w,y_w)是自己定义的棋盘格坐标。在单应性矩阵中有8个独立未知元素，因此至少需要8个方程，每个角点提供2个方程，因此当一张图片上的标定板角点数量等于4时，即可求得该图片对应的矩阵H。当一张图片上的标定板角点数量大于4时，利用最小二乘法回归最佳的矩阵。

u=H11xw+H12yw+H13H31xw+H32yw+H33v=H21xw+H22yw+H23H31xw+H32yw+H33u = \frac{H_{11}x_w+H_{12}y_w+H_{13}}{H_{31}x_w+H_{32}y_w+H_{33}}\\ v = \frac{H_{21}x_w+H_{22}y_w+H_{23}}{H_{31}x_w+H_{32}y_w+H_{33}} \\

3.2 求解内外参矩阵

旋转矩阵R描述了坐标系之间的旋转变换，由3 × 3矩阵描述，是一个正交矩阵。因此有

R1TR2=0R1TR1=R2TR2=1R 1^{T} R 2=0 \\ R 1^{T} R 1=R 2^{T} R 2=1

根据上面单应性矩阵的描述：H=A[R1R2T]H=A\begin{bmatrix}R1&R2&T\end{bmatrix}，因此有R1=A−1H1R1= A^{-1} H1，R2=A−1H2R 2= A^{-1} H 2，代入可得：

H1TA−TA−1H2=0H1TA−TA−1H1=H2TA−TA−1H2=1H 1^{T} A^{-T} A^{-1} H 2=0 \\ H 1^{T} A^{-T} A^{-1} H 1 = H 2^{T} A^{-T} A^{-1} H 2=1

令B=A−TA−1B=A^{-T} A^{-1}，则B为对称阵，有六个未知元素，单张图片只能提供两个约束方程，因此，我们只要取3张标定板照片，即6个方程，即可求解 B。当标定板图片的个数大于3时（事实上一般需要15到20张标定板图片），可采用最小二乘拟合最佳矩阵B。然后根据内参矩阵与矩阵B的对应关系就可以求出内参矩阵A。

对于外参，其反映的是标定板和相机的位置关系。对于不同的图片，标定板和相机的位置关系已经改变，此时每一张图片对应的外参矩阵都是不同的。在H=A[R1R2T]H=A\begin{bmatrix}R1&R2&T\end{bmatrix}中我们已经求解得到了矩阵HH（对于同一张图片相同，对于不同的图片不同）、矩阵AA（对于不同的图片都相同），因此可以方便求出外参[R1R2T]\begin{bmatrix}R1&R2&T\end{bmatrix}。

这里的外参没有R3是因为棋盘格上的点Z坐标为0使得其无效，即R3在坐标转化中并没有作用。但是R3要使得R满足旋转矩阵的性质，即列与列之间单位正交，因此可以通过向量 R1,R2 的叉乘，即 R3=R1×R2 ，计算得到R3。

3.3 求解畸变系数

张正友标定法仅仅考虑了畸变模型中影响较大的径向畸变。根据上面相机畸变的公式，将图像坐标替换为像素坐标后有：

u^=u+(u−u0)(k1r2+k2r4)v^=v+(v−v0)(k1r2+k2r4)\begin{array}{l}{\hat{u}=u+\left(u-u_{0}\right)\left(k_{1} r^{2}+k_{2} r^{4}\right)} \\ {\hat{v}=v+\left(v-v_{0}\right)\left(k_{1} r^{2}+k_{2} r^{4}\right)}\end{array} \\

[(u−u0)r2(u−u0)r4(v−v0)r2(v−v0)r4][k1k2]=[u^−uv^−v]\left[\begin{array}{cc}{\left(u-u_{0}\right) r^{2}} & {\left(u-u_{0}\right) r^{4}} \\ {\left(v-v_{0}\right) r^{2}} & {\left(v-v_{0}\right) r^{4}}\end{array}\right]\left[\begin{array}{l}{k_{1}} \\ {k_{2}}\end{array}\right]=\left[\begin{array}{c}{\hat{u}-u} \\ {\hat{v}-v}\end{array}\right] \\

上式中的u^,u,v^,v\hat{u},u,\hat{v},v可以通过识别标定板的角点获得，每一个角点可以构造两个上述等式。有m幅图像，每幅图像上有n个标定板角点，则将得到的所有等式组合起来，可以得到个mn未知数为的k=[k1,k2]Tk=[ k_1, k_2]^T约束方程将约束方程系数矩阵记为DD，等式右端非齐次项记为dd，可将其记着矩阵形式：

Dk=dD k=d

使用最小二乘法可得：

k=[k1k2]=(DTD)−1DTdk=\left[\begin{array}{l}{k_{1}} \\ {k_{2}}\end{array}\right]=\left(D^{T} D\right)^{-1} D^{T} d \\

4 OpenCV相机标定

OpenCV中给出了基于棋盘格的相机标定方法，其主要步骤为角点检测、亚像素细化、可视化展示（非必须）、构建世界坐标系三维坐标和标定，使用流程主要如下。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
//*******************************************************************
//Fuction Name:cameraCalibration
//Summary:单个相机标定，输出内外参等
//Params:  -vector<String> firImgPathList 输入参数，标定板图像文件名数组
//		   -int rowNums, colNums, blockSize 输入参数，棋盘格标定板每行每列的角点数以及方格尺寸
//		   -vector<String> errImageList 输出参数，角点检测出错文件名
//		   -vector<Mat>& cornersMat 输出参数，标记角点后的图像数组
//		   -Mat& cameraK 输出参数，3*3的相机内参矩阵
//		   -Mat& distCoeffs 输出参数，1*5的相机畸变系数
//		   -vector<Mat>& rotationMat 输出参数，外参中旋转矩阵数组（世界坐标系转换到相机坐标系），每幅图像对应一个
//		   -vector<Mat>& translationMat 输出参数，外参中平移向量数组，每幅图像对应一个
//		   -vector<vector<Point2f>>& allPts 输出参数，文件夹中每幅图像的角点坐标
//		   -
//retrun:None
//********************************************************************
void cameraCalibration(std::vector<cv::String> firImgPathList, int rowNums, int colNums, int blockSize, std::vector<cv::String> errImageList, std::vector<cv::Mat>& cornersMat, cv::Mat& cameraK, cv::Mat& distCoeffs, std::vector<cv::Mat>& rotationMat, std::vector<cv::Mat>& translationMat, std::vector<std::vector<cv::Point2f>>& allPts)
{
	Size cornerSize = Size(rowNums, colNums); //标定板每行每列角点个数
	Size imgSize;//图像尺寸
	vector<cv::Point2f> imgPts;//缓存每幅图检测到的角点

	for (int imgCount = 0;imgCount < firImgPathList.size();imgCount++)
	{
		//读入图片
		cv::Mat imgRaw = cv::imread(firImgPathList[imgCount], 0);

		//获取图像尺寸
		if (!imgCount)
		{
			imgSize.width = imgRaw.cols;  //图像的宽对应着列数
			imgSize.height = imgRaw.rows; //图像的高对应着行数
		}

		//角点检测
		//step1 提取角点
		bool success = cv::findChessboardCorners(imgRaw, cornerSize, imgPts, CV_CALIB_CB_ADAPTIVE_THRESH | CV_CALIB_CB_FAST_CHECK | CV_CALIB_CB_NORMALIZE_IMAGE);
		if (!success)
		{
			errImageList.push_back(firImgPathList[imgCount]);
			continue;
		}
		else
		{
			//step2 亚像素角点
			//cv::find4QuadCornerSubpix(image_gray, points_per_image, cv::Size(5, 5));
			cornerSubPix(imgRaw, imgPts, cv::Size(11, 11), cv::Size(-1, -1), cv::TermCriteria(CV_TERMCRIT_EPS | CV_TERMCRIT_ITER, 30, 0.1));
			allPts.push_back(imgPts); //保存亚像素角点
			//在图中画出角点位置
			//step3 角点可视化
			cv::Mat imgCorners = imgRaw.clone();
			cvtColor(imgCorners, imgCorners, CV_GRAY2BGR);
			cv::drawChessboardCorners(imgCorners, cornerSize, imgPts, success); //将角点连线
			cornersMat.push_back(imgCorners);
		}
	}

	//开始相机标定
	//初始化角点三维坐标,从左到右,从上到下
	std::vector<cv::Point3f>imgPts3D;
	for (int i = 0; i < cornerSize.height; i++)
		for (int j = 0; j < cornerSize.width; j++)
		{
			imgPts3D.push_back(cv::Point3f(blockSize * j, blockSize * i, 0));
		}

	int imgNums = firImgPathList.size();
	vector<vector<cv::Point3f>> allPts3D(imgNums, imgPts3D);        //保存所有图像角点的三维坐标, z=0

	//step4 标定								 					  
	cv::calibrateCamera(allPts3D, allPts, imgSize, cameraK, distCoeffs, rotationMat, translationMat, 0);
}
1
2
3
4
5
6
//1. 棋盘格角点检测函数
bool cv::findChessboardCorners	( InputArray 	image,
Size 	patternSize,
OutputArray 	corners,
int 	flags = CALIB_CB_ADAPTIVE_THRESH+CALIB_CB_NORMALIZE_IMAGE 
)
参数 说明
InputArray image 输入棋盘格图像，必须是8位灰度或彩色图
Size patternSize 棋盘每行和每列的内角数
OutputArray corners 检测到的角点的输出数组
int flags 各种操作标志，可以是0或者下面值的组合：CV_CALIB_CB_ADAPTIVE_THRESH：该函数的默认方式是根据图像的平均亮度值进行图像二值化，设立此标志位的含义是采用变化的阈值进行自适应二值化； CV_CALIB_CB_NORMALIZE_IMAGE：在二值化之前，调用EqualizeHist()函数进行图像归一化处理； CV_CALIB_CB_FILTER_QUADS：二值化完成后，函数开始定位图像中的四边形（这里不应该称之为正方形，因为存在畸变），这个标志设立后，函数开始使用面积、周长等参数来筛选方块，从而使得角点检测更准确更严格。 
1
2
3
4
5
6
7
//2.亚像素细化函数（有博客说cornerSubPix比find4QuadCornerSubpix精度高）
void cv::cornerSubPix	(   InputArray 	image,
							InputOutputArray 	corners,
							Size 	winSize,
							Size 	zeroZone,
							TermCriteria 	criteria 
						 )
参数 说明
InputArray image 输入单通道、8 位或浮点图像。
InputOutputArray corners 检测到的角点，即是输入也是输出
Size winSize 计算亚像素角点时考虑的区域的大小，大小为N * N; N=(winSize*2+1)
Size zeroZone 通常忽略（即Size(-1, -1)）
TermCriteria criteria 表示计算亚像素时停止迭代的标准，可选的值有cv::TermCriteria::MAX_ITER 、cv::TermCriteria::EPS（可以是两者其一，或两者均选），前者表示迭代次数达到了最大次数时停止，后者表示角点位置变化的最小值已经达到最小时停止迭代。二者均使用cv::TermCriteria()构造函数进行指定。 
1
2
3
4
5
6
//3.角点可视化
void cv::drawChessboardCorners	(	InputOutputArray 	image,
									Size 	patternSize,
									InputArray 	corners,
									bool 	patternWasFound 
								 )
参数 说明
InputOutputArray image 用于绘制角点可视化的图像，必须为8为彩色图像
Size patternSize 棋盘每行和每列的内角数
InputArray corners 检测到的角点坐标数组，可以输入findChessboardCorners或者cornerSubPix的结果
bool patternWasFound 该参数指示是否找到了完整棋盘。应在此处传递 findChessboardCorners 的返回值。 
1
2
3
4
5
6
7
8
9
double cv::calibrateCamera	(	InputArrayOfArrays 	objectPoints,
                                InputArrayOfArrays 	imagePoints,
                                Size 	imageSize,
                                InputOutputArray 	cameraMatrix,
                                InputOutputArray 	distCoeffs,
                                OutputArrayOfArrays 	rvecs,
                                OutputArrayOfArrays 	tvecs,
                                int 	flags = 0,
)
参数 说明
InputArrayOfArrays objectPoints 输入的世界坐标系中的点。在使用时，应该输入一个三维点vector的vector，即vector<vector< Point3f>> objectPoints。第一层vector表示每一个视角，第二层vector表示每一个点。需要注意的是点的数量应该和角点数量对应，且三维点的Z坐标值应该为0。
InputArrayOfArrays imagePoints 其对应的图像坐标系中的点。和objectPoints一样，应该输入std::vector< std::vector< cv::Point2f>> imagePoints型的变量。来源应该是findChessboardCorners或者cornerSubPix的结果。
Size imageSize 图像的大小，在计算相机的内参数和畸变矩阵需要用到
InputOutputArray cameraMatrix 3*3的相机内参矩阵，输入一个cv::Mat cameraMatrix即可
InputOutputArray distCoeffs 畸变矩阵。输入一个cv::Mat distCoeffs即可。具体尺寸取决于falgs。
OutputArrayOfArrays rvecs 旋转向量；应该输入一个cv::Mat的vector，即vector< cv::Mat> rvecs因为每个vector< Point3f>会得到一个rvecs。每个vec为3*1，可以用Rodrigues函数转换为3*3的旋转矩阵。
OutputArrayOfArrays tvecs 位移向量；和rvecs一样，也应该为vector< cv::Mat> tvecs。每个vec为3*1。
int flags CV_CALIB_USE_INTRINSIC_GUESS：使用该参数时，在cameraMatrix矩阵中应该有fx,fy,cx,cy的估计值。否则的话，将初始化(cx,cy）图像的中心点，使用最小二乘估算出fx，fy。如果内参数矩阵和畸变居中已知的时候，应该标定模块中的solvePnP()函数计算外参数矩阵。CV_CALIB_FIX_PRINCIPAL_POINT：在进行优化时会固定光轴点。当CV_CALIB_USE_INTRINSIC_GUESS参数被设置，光轴点将保持在中心或者某个输入的值。CV_CALIB_FIX_ASPECT_RATIO：固定fx/fy的比值，只将fy作为可变量，进行优化计算。当CV_CALIB_USE_INTRINSIC_GUESS没有被设置，fx和fy将会被忽略。只有fx/fy的比值在计算中会被用到。CV_CALIB_ZERO_TANGENT_DIST：设定切向畸变参数（p1,p2）为零。CV_CALIB_FIX_K1,…,CV_CALIB_FIX_K6：对应的径向畸变在优化中保持不变。如果设置了CV_CALIB_USE_INTRINSIC_GUESS参数，CV_CALIB_RATIONAL_MODEL：计算k4，k5，k6三个畸变参数。如果没有设置，则只计算其它5个畸变参数。

5 参考

深度科普：一文搞懂相机标定

相机标定数学原理

单应性矩阵

旋转矩阵

亚像素角点