1.相机标定基本原理
1.1 简介
在图像测量过程以及机器视觉应用中，为确定空间物体表面某点的三维几何位置与其在图像中对应点之间的相互关系，必须建立摄像机成像的几何模型,这些几何模型参数就是摄像机参数。在大多数条件下这些参数必须通过实验与计算才能得到，这个求解参数的过程就称之为相机标定。简单来说是从世界坐标系换到图像坐标系的过程，也就是求最终的投影矩阵P P P的过程。
无论是在图像测量或者机器视觉应用中，摄像机参数的标定都是非常关键的环节，其标定结果的精度及算法的稳定性直接影响摄像机工作产生结果的准确性。因此，做好摄像机标定是做好后续工作的前提，是提高标定精度是科研工作的重点所在。其标定的目的就是为了相机内参、外参、畸变参数。
1.2 基本的坐标系
世界坐标系：用户定义的三维世界的坐标系，为了描述目标物在真实世界里的位置而被引入。
相机坐标系：在相机上建立的坐标系，为了从相机的角度描述物体位置而定义，作为沟通世界坐标系和图像/像素坐标系的中间一环。
图像坐标系：为了描述成像过程中物体从相机坐标系到图像坐标系的投影透射关系而引入，方便进一步得到像素坐标系下的坐标。
- 一般来说，标定的过程分为两个部分：
第一步是从世界坐标系转换为相机坐标系，这一步是三维点到三维点的转换，包括R R R，t t t（相机外参）等参数；
第二步是从相机坐标系转为图像坐标系，这一步是三维点到二维点的转换，包括 K K K（相机内参）等参数；
- 同步标定内部参数和外部参数，一般包括两种策略s:
光学标定: 利用已知的几何信息(如定长棋盘格)实现参数求解。
自标定: 在静态场景中利用 structure from motion估算参数
; 1. 3 畸变参数
理想的小孔成像模型，物和像满足相似三角形的关系。实际上由于相机光学系统制造工艺的误差，实际成像与理想成像存在几何失真，称为畸变。畸变主要分为径向畸变和切向畸变。
(1)径向畸变（枕形、桶形）：
①透镜质量原因
②光线在远离透镜中心的地方比靠近中心的地方更加弯曲。
径向畸变可以用如下公式修正：
(2)切向畸变（薄透镜畸变和离心畸变）：
切向畸变是由于透镜制造上的缺陷使得透镜本身与图像平面不平行而产生的。
切向畸变可以用如下公式修正：
其中：
x d i s x_{dis}x d i s ,y d i s y_{dis}y d i s 表示有畸变的坐标；
x c o r r x_{corr}x c o r r ,y c o r r y _{corr}y c o r r 表示修复后的坐标；
k 1 k_1 k 1 ,k 2 k_2 k 2 ,k 3 k_3 k 3 表示径向畸变参数；
p 1 p_1 p 1 ,p 1 p_1 p 1 表示切向畸变参数
- 可知畸变有k 1 k_1 k 1 ,k 2 k_2 k 2 ,k 3 k_3 k 3 ,p 1 p_1 p 1 ,p 2 p_2 p 2 五个参数，对于质量比较好的相机来说，切向畸变很小，可忽略，径向畸变系数k 3 k_3 k 3 也可忽略，只计算k 1 k_1 k 1 ,k 2 k_2 k 2 两个参数。张正友标定中就默认为p 1 p_1 p 1 ,p 2 p_2 p 2 为0
1.4 相机标定
通过空间中已知坐标的(特征)点 ( X i X_i X i ,Y i Y_i Y i ,Z i Z_i Z i ) ，以及它们在图像中的对应坐标( u i u_i u i , v i v_i v i ) ，直接估算 11 个待求解的内部和外部参数。
x x x∼K [ R ∣ t ] X K[R|t]X K [R ∣t ]X=M X MX M X
; 1.4.1 线性标定(最小二乘)
线性标定不考虑相机的畸变而只考虑空间坐标转换。
每个坐标点有X,Y两个变量，可列两个方程，相机内参有5个未知数，外参平移和旋转各3个，共有11个变量，因此至少需要6个特征点来求解。
将x x x∼K [ R ∣ t ] X K[R|t]X K [R ∣t ]X=M X MX M X表示为
可以表示为矩阵形式：
便可用最小二乘法求解。
给定超定方程（超定方程组是指方程个数大于未知量个数的方程组。对于方程组R a Ra R a=y y y，R为n × m n×m n ×m矩阵，如果R R R列满秩，且n>m。则方程组没有精确解，此时称方程组为超定方程组）：A X = b A_X=b A X =b
其中x的解为等式两边的误差平方和最小化。
1.4.2 非线性标定
当镜头畸变明显时必须引入畸变模型，将线性标定模型转化为非线性标定模型，
通过非线性优化的方法求解相机参数：
用概率的视角去看最小二乘问题
- 特征点投影方程为
- 给定{(ui,vi)}，标定参数矩阵 M M M 的概率为：
- 给定{(ui,vi)}，标定参数矩阵 M 的似然函数为：
相应求解策略: 牛顿方法、高斯-牛顿方法、Levenberg-Marquardt算法等
; 1.4.3 张正友标定法（只考虑径向畸变，不考虑切向畸变）
张正友标定法利用棋盘格标定板可以利用相应的角点提取算法（如Harris角点）得到每一个角点的像素坐标( u , v ) (u,v)(u ,v )。
张正友标定法将世界坐标系固定于棋盘格上，则棋盘格上任一点的物理坐标 W = 0 W=0 W =0，由于标定板的世界坐标系是人为事先定义好的，标定板上每一个格子的大小是已知的，我们可以计算得到每一个角点在世界坐标系下的物理坐标( U , V , W ) = 0 (U,V,W)=0 (U ,V ,W )=0。
我们将利用这些信息：每一个角点的像素坐标( u , v ) (u,v)(u ,v ) 、每一个角点在世界坐标系下的物理坐标( U , V , W ) = 0 (U,V,W)=0 (U ,V ,W )=0，来进行相机的标定，获得相机的内外参矩阵、畸变参数。
1.4.3.1 基本原理
(1) 基本概念
其中：
- s s s: 世界坐标系到图像坐标系的尺度因子
- A A A: 相机内参矩阵
- (u 0 u_0 u 0 ,v 0 v_0 v 0 ): 像主点坐标
- α, β: 焦距与像素横纵比的融合
- γ: 径向畸变参数
(2)求解Homographic矩阵
不妨设棋盘格位于Z Z Z = 0
H H H矩阵可以根据特征点/棋盘格角点的空间坐标，以及其图像坐标用最小二乘法很容易求解。
Homography 有 8 个自由度，
由r1和r2正交，且r1和r2的模相等，可以得到如下约束：
（3）计算内参数矩阵
（4）计算外参数矩阵
（5）极大似然估计
给定n n n张棋盘格图像，每张图像有m m m个角点
上述非线性优化问题可以利用Levenberg-Marquardt 算法求解
需要初值𝐀,{𝐑𝑖,t𝑖 |𝑖=1…𝑛}
; 1.4.3.2 主要流程
- 打印一张棋盘格A4纸张（黑白间距已知），并贴在一个平板上
- 针对棋盘格拍摄若干张图片（一般10-20张）
- 在图片中检测特征点（Harris特征）
- 利用解析解估算方法计算出5个内部参数，以及6个外部参数
- 根据极大似然估计策略，设计优化目标并实现参数的refinement
2.相机标定具体实现
2.1 实验数据集
棋盘格(建议使用x方向和y方向不等的棋盘格)：
在电脑上绘制了7*5 每小正方格的长宽为1.8cm的黑白棋盘格，然后再使用手机从不同角度和距离对棋盘进行拍摄。
数据集如下：
- 共拍摄15张图片
; 2.2 实验代码
import cv2
import numpy as np
import glob
np.set_printoptions(suppress=True)
criteria = (cv2.TERM_CRITERIA_EPS + cv2.TERM_CRITERIA_MAX_ITER, 30, 0.001)
w = 6
h = 4
objp = np.zeros((w * h, 3), np.float32)
objp[:, :2] = np.mgrid[0:w, 0:h].T.reshape(-1, 2)
objpoints = []
imgpoints = []
images = glob.glob('p2/*.jpg')
i = 0
for fname in images:
img = cv2.imread(fname)
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
ret, corners = cv2.findChessboardCorners(gray, (w, h), None)
if ret == True:
i += 1
cv2.cornerSubPix(gray, corners, (11, 11), (-1, -1), criteria)
objpoints.append(objp)
imgpoints.append(corners)
cv2.drawChessboardCorners(img, (w, h), corners, ret)
cv2.imshow('findCorners', img)
cv2.imwrite('p2/h' + str(i) + '.jpg', img)
cv2.waitKey(10)
cv2.destroyAllWindows()
ret, mtx, dist, rvecs, tvecs = cv2.calibrateCamera(objpoints, imgpoints, gray.shape[::-1], None, None)
print(("ret:"), ret)
print(("mtx:\n"), mtx)
print(("dist:\n"), dist)
print(("rvecs:\n"), rvecs)
print(("tvecs:\n"), tvecs)
img2 = cv2.imread('p2/1.jpg')
h, w = img2.shape[:2]
newcameramtx, roi = cv2.getOptimalNewCameraMatrix(mtx, dist, (w, h), 0, (w, h))
dst = cv2.undistort(img2, mtx, dist, None, newcameramtx)
cv2.imwrite('p2/calibresult.jpg', dst)
total_error = 0
for i in range(len(objpoints)):
imgpoints2, _ = cv2.projectPoints(objpoints[i], rvecs[i], tvecs[i], mtx, dist)
error = cv2.norm(imgpoints[i], imgpoints2, cv2.NORM_L2) / len(imgpoints2)
total_error += error
print(("total error: "), total_error / len(objpoints))
2.3 实验结果及分析
（1）角点检测结果
使用 findChessboardCorners函数提取角点，这里的角点专指的是标定板上的内角点，保存了每一张棋盘格的角点检测结果，每一张均可检测到角点。简单列出上图可以看出共检测到24个角点，其它图的角点检测结果均与上图相同。
（2）重投影误差
重投影误差定义为一个特征点在归一化相机坐标系下的估计值与观测值的差，网上均重投影误差说小于0.5 就算效果良好(还不太清楚是为什么)，可以看出我们计算出来的重投影误差为0.32，效果较好。
（3）内参数矩阵
（4）畸变系数
由于distortion cofficients = (k 1 k_1 k 1 ,k 2 k_2 k 2 ,p 1 p_1 p 1 ,p 2 p_2 p 2 ,k 3 k_3 k 3 )
可以看出我手机的径向畸变为k 1 k_1 k 1 =0.308,k 2 = − 2.496 k_2=-2.496 k 2 =−2 .4 9 6，k 3 = 5.167 k_3=5.167 k 3 =5 .1 6 7 切向畸变为p 1 = 0.001 p_1=0.001 p 1 =0 .0 0 1,p 2 = 0.001 p_2=0.001 p 2 =0 .0 0 1
（5）外部参数：
①外部参数（旋转向量）
rvecs:
[array([[ 0.16629272],
[-0.05250479],
[ 0.01478173]]), array([[ 0.18137582],
[-0.32449346],
[ 0.02120859]]), array([[ 0.36048591],
[-0.02939182],
[ 0.02555599]]), array([[ 0.41331229],
[ 0.10562771],
[-0.15170865]]), array([[ 0.33406759],
[-0.00680696],
[ 0.0758683 ]]), array([[ 0.12600849],
[-0.02437899],
[ 0.00584542]]), array([[ 0.01192604],
[-0.02320687],
[-0.00056763]]), array([[ 0.14358827],
[ 0.21669299],
[-0.018957 ]]), array([[0.39108078],
[0.03463762],
[0.02774009]]), array([[-0.12300712],
[-0.01348449],
[ 0.00011424]]), array([[ 0.10045939],
[-0.3586192 ],
[ 0.01283334]]), array([[ 0.14020908],
[-0.37377349],
[ 0.01964494]]), array([[ 0.01841512],
[ 0.10686787],
[-0.02817669]]), array([[ 0.35815865],
[-0.00262295],
[ 0.02612322]]), array([[-0.19383748],
[-0.04212595],
[-0.0351099 ]])]
②外部参数（平移向量）
tvecs:
[array([[-2.2324222 ],
[-1.02512617],
[ 6.52257936]]), array([[-1.64776615],
[-1.644847 ],
[ 6.70674111]]), array([[-2.33111107],
[-0.69591389],
[ 6.74235331]]), array([[-2.19518285],
[ 0.21609157],
[ 8.44286089]]), array([[-2.278302 ],
[-1.14775155],
[ 6.97844932]]), array([[-2.38326812],
[-1.38000887],
[ 6.31510077]]), array([[-2.33657618],
[-1.42579754],
[ 6.7801348 ]]), array([[-2.74363267],
[-1.06954044],
[ 7.7096763 ]]), array([[-2.39996295],
[-0.72121477],
[ 6.63976727]]), array([[-2.72550172],
[-1.85128945],
[ 7.41543224]]), array([[-2.40023387],
[-1.51597381],
[ 6.42318618]]), array([[-1.43093631],
[-1.21085161],
[ 6.8980046 ]]), array([[-2.29867829],
[-1.62369345],
[ 9.52627184]]), array([[-2.24827193],
[-0.90869927],
[ 7.33784455]]), array([[-2.79625849],
[-1.29122356],
[ 8.01478454]])]
（6）反向投影误差
可以利用反向投影误差对我们找到的参数的准确性进行估计，得到的结果越接近 0 越好。通过之前计算的内参数矩阵、畸变系数、旋转矩阵和平移向量，使用 cv2.projectPoints()计算三维点到二维图像的投影，然后计算反投影得到的点与图像上检测到的点的误差，最后计算一个对于所有标定图像的平均误差。控制台显示我们得到的结果为0.07，与0较接近，效果较好。
（7）畸变矫正
原图:
矫正后：
肉眼观察不太能观察出两张图片有什么太大的差别，相机镜头不存在太大的畸变。
3.实验中遇到的问题
1.畸变参数计算错误，畸变矫正效果离谱。
原图
畸变矫正后
原因：一开始使用的标定板为 5×5 的棋盘格，可以由上图看出畸变矫正的结果非常糟糕。仔细检查后发现矫正算法不存在问题，畸变参数计算的结果是错误的。在观察每一张标定板图像的角点检测结果后可以发现，如下图所示，在角点检测时，由于棋盘格长宽是相同的，算法误判图像旋转了一个角度，x方向和y方向颠倒导致计算出的畸变参数就有一定的误差。我最后重新制作了7×5的标定板进行实验后，问题得到解决。由于一开始我们在代码中给定了棋盘格的模板规格，给定了不同方向上的内角点个数。虽然在拍摄图像需要有一定的角度变化，但不能将x方向和y方向完全颠倒，会导致角点检测错误。
角点检测结果1
角点检测结果2
2.标定图片拍摄需要注意的问题
- 标定图片的数量应该在10~25张之间，图像数量太少，容易导致标定参数不准确。
- 标定过程，摄像机的光圈、焦距不能发生改变，改变需要重新标定。
参考：
最详细、最完整的相机标定讲解
相机标定(Camera calibration)
Original: https://blog.csdn.net/weixin_43890415/article/details/117084459
Author: 机智的小陈今天学习了吗
Title: 计算机视觉-相机标定(Camera Calibration)
原创文章受到原创版权保护。转载请注明出处：https://www.johngo689.com/518693/
转载文章受原作者版权保护。转载请注明原作者出处！