一、相机标定基本理论
1、相机成像系统介绍
图中总共有4个坐标系:
- 图像坐标系：Op 坐标表示方法(u,v) Unit：Dots(个)
- 成像坐标系：Oi 坐标表示方法(x\',y\',z\') Unit：mm(毫米)
- Camera坐标系：Oc 坐标表示方法(x,y,z) Unit：mm(毫米)
- World世界坐标系：Ow 坐标表示方法(X,Y,Z) Unit：mm(毫米)
图中所示的坐标转换关系：
{World世界坐标系：Ow }---Extrinsics [R|t]---{Camera坐标系：Oc}---Intrinsics A---{成像坐标系：Oi }---Liner Convert---{图像坐标系：Op }
这样就完成了世界坐标系中的点到图像坐标系的映射关系！
2、相机内参介绍
1、标准的针孔模型（针孔模型实际上没有焦距f，这里只是为了方便）：
在实际计算过程中，为了计算方便，我们将像平面翻转平移到针孔前，从而得到一种数学上更为简单的等价形式（方便相似三角形的计算）
图1.标定全过程解析
这里的Q点实际上是O点对面的X点的位置，q实际为x的位置
注意：由于世界坐标系的选定不是唯一的，因此为了方便起见，我们在Opencv使用的过程中，Ow世界坐标系选取的位置为标定板平面及其法线组成的三维空间坐标系，因此，世界坐标系中的所有的标定板上的点的Zw=0，所以在下面的python-Opencv程序中，我们角点的objpoints的取值Zw都为0，如旁边的小图所示：红色为Zw、蓝色为Yw、绿色为Xw
以Oc点为原点建立摄像机坐标系。点Q(x,y,z)为摄像机坐标系空间中的一点，该点被光线投影到图像平面Oi上的q(x\',y\',z\')点，由于图像平面的位置再相机坐标系的焦点的位置，所以这里z‘=f，q点的坐标为(x\',y\',f)。
图像平面与光轴z轴垂直，和投影中心距离为f （f是相机的焦距）。按照三角比例关系可以得出：
Formulation-001
注意这里的Zc实际上的值为s，下面的过程我们将会带入s数值
以上将坐标为(x,y,z)的Q点映射到投影平面上坐标为(x\',y\')的q点的过程称作投影变换
上述Q点到q点的变换关系用3*3的矩阵可表示为：q = A1Q ，其中
Formulation-002
s*A1称为摄像机的内参数矩阵，单位均为物理尺寸（Units：mm）
通过上面，可以把相机坐标系转换到像图像坐标系的物理单位，即{Camera坐标系：Oc}---Intrinsics A1---{成像坐标系：Oi }3、成像平面坐标系→像素坐标系
通过下图，可以把像平面坐标系物理单位到像素单位，即{成像坐标系：Oi }---Liner Convert---{图像坐标系：Op }
以图像平面的左上角或左下角为原点建立坐标系Op。假设像平面坐标系原点位于图像左下角，水平向右为u轴，垂直向上为v轴，均以像素为单位
以图像平面与光轴的交点Oi 为原点建立坐标系，水平向右为x轴，垂直向上为y轴。原点Oi 一般位于图像中心处，Oi 在以像素为单位的图像坐标系中的坐标为(u0, v0)
像平面坐标系和像素坐标系虽然在同一个平面上，但是原点并不是同一个。
设每个像素的物理尺寸大小为 dx * dy (mm) ( 由于单个像素点投影在图像平面上是矩形而不是正方形，因此可能dx != dy)，
图像平面上某点在成像平面坐标系中的坐标为(x\', y\')，在像素坐标系中的坐标为(u, v)，则二者满足如下关系：[即(x\', y\')→(u, v)]
Formulation-004
用齐次坐标与矩阵形式表示为：
Formulation-005
将等式两边都乘以点Q(x,y,zc)坐标中的Zc(也就是s)可得：
Formulation-006
将摄像机坐标系中的Formulation-003 式代入上式可得：
Formulation-007
则右边第一个矩阵和第二个矩阵的乘积亦为摄像机的内参数矩阵A（单位为像素），相乘后可得：
Formulation-008
和 Formulation-003 式相比，此内参数矩阵中f/dx, f/dy, u0, v0 的单位均为像素。令内参数矩阵为A，则上式可写成：
Formulation-009
三.相机内参A（与棋盘所在空间的3D几何无关）
在计算机视觉中，摄像机内参数矩阵
Formulation-010
其中 f 为摄像机的焦距，单位一般是mm;dx,dy 为像元尺寸;u0,v0 为图像中心。
fx = f/dx, fy = f/dy,分别称为x轴和y轴上的归一化焦距.
为更好的理解f、dx、dy、fx、fy、u0、v0，举个实例：
现以NiKon D700相机为例进行求解其内参数矩阵：
焦距 f = 35mm 最高分辨率：4256×2832 传感器尺寸：36.0×23.9 mm
根据以上定义可以有：
u0= 4256/2 = 2128 v0= 2832/2 = 1416 dx = 36.0/4256 dy = 23.9/2832
fx = f/dx = 4137.8 fy = f/dy = 4147.3
分辨率可以从显示分辨率与图像分辨率两个方向来分类:
- 显示分辨率（屏幕分辨率）是屏幕图像的精密度，是指显示器所能显示的像素有多少。由于屏幕上的点、线和面都是由像素组成的，显示器可显示的像素越多，画面就越精细，同样的屏幕区域内能显示的信息也越多，所以分辨率是个非常重要的性能指标之一。可以把整个图像想象成是一个大型的棋盘，而分辨率的表示方式就是所有经线和纬线交叉点的数目。
- 显示分辨率一定的情况下，显示屏越小图像越清晰，反之，显示屏大小固定时，显示分辨率越高图像越清晰。图像分辨率则是单位英寸中所包含的像素点数，其定义更趋近于分辨率本身的定义。
3、相机外参介绍
旋转向量R=[Rx,Ry,Rz]（大小为1×3的矢量或旋转矩阵3×3）和平移向量 t=(Tx,Ty,Tz)
旋转向量:旋转向量是旋转矩阵紧凑的变现形式，旋转向量为1×3的行矢量。
Formulation-011
r就是旋转向量，旋转向量的方向是旋转轴 ,旋转向量的模为围绕旋转轴旋转的角度，如图中的e向量，即旋转向量表现形式。
通过上面的公式，我们就可以求解出旋转矩阵R。同样的已知旋转矩阵，我们也可以通过下面的公式求解得到旋转向量：
Fmulation-012
根据前面的 图1.相机标定全过程解析 中的变换方式，标定过程中的外参Extrinsics = [R|t]，如下所示为[R|t]:
Formulation-013
则从世界坐标系变换到Q点的方式如下：
Formulation-014
即简写为如下形式：
Formulation-015
4、畸变参数介绍
采用理想针孔模型，由于通过针孔的光线少，摄像机曝光太慢，在实际使用中均采用透镜，可以使图像生成迅速，但代价是引入了畸变。有两种畸变对投影图像影响较大： 径向畸变和切向畸变。
1、径向畸变
对某些透镜，光线在远离透镜中心的地方比靠近中心的地方更加弯曲，产生“筒形”或“鱼眼”现象，称为径向畸变。
一般来讲，成像仪中心的径向畸变为0，越向边缘移动，畸变越严重。不过径向畸变可以通过下面的泰勒级数展开式来校正：
xcorrected = x(1+k1r2+k2r4+k3r6)
ycorrected = y(1+k1r2+k2r4+k3r6)
这里(x, y)是畸变点在成像仪上的原始位置，r为该点距离成像仪中心的距离，(xcorrected ，ycorrected )是校正后的新位置。
对于一般的摄像机校正，通常使用泰勒级数中的前两项k1和k2就够了；对畸变很大的摄像机，比如鱼眼透镜，可以使用第三径向畸变项k3
无畸变 桶形畸变k1>0 枕形畸变k1<0
2、切向畸变
当成像仪被粘贴在摄像机的时候，会存在一定的误差，使得图像平面和透镜不完全平行，从而产生切向畸变。也就是说，如果一个矩形被投影到成像仪上时，可能会变成一个梯形。切向畸变可以通过如下公式来校正：
xcorrected = x + [ 2p1y + p2 (r2 + 2x2) ]
ycorrected = y + [ 2p2x + p1 (r2 + 2y2) ]
这里（x, y）是畸变点在成像仪上的原始位置，r为该点距离成像仪中心的距离，（xcorrected ，ycorrected ）是校正后的新位置。
3、参数表示
x、y、xcorrected 、ycorrected 、r参数的表示如下图所示：
畸变系数的作用在于计算矫正之后的实际的成像平面位置，将得到的矫正公式带入到上述的计算公式Formulation-015式中即可得到矫正后的图像。
二、相机标定方法实践
1、采样标准pattern模板图像若干
1 import cv2 2 import glob 3 import time 4 import threading 5 import numpy as np 6 7 cap = cv2.VideoCapture(0) 8 while not cap.isOpened(): # 检查摄像头是否打开成功 9 time.sleep(100) 10 print(\'Camera is Initialize...\') 11 12 width = int(cap.get(3)) # 读取摄像头分辨率参数 13 height = int(cap.get(4)) 14 15 frame = np.zeros((width,height,3),dtype=np.uint8) # 创建图像模板 16 ret, frame = cap.read() 17 18 Key_val = 0 # 保存键值 19 20 def Keybo_Moni(): # 按键测试函数 21 count = 0 22 while True: 23 global Key_val, frame, process_flag, cap 24 if Key_val == ord(\'r\'): 25 Key_val= 0 26 cv2.imwrite(\'ResPic\' + str(count) + \'.jpg\', frame) # 保存图像 27 count += 1 28 print(\'Get new pic %d\' % count) 29 30 try: 31 32 Keybo_Moni_Thread = threading.Thread(target=Keybo_Moni, name=\'Keyboard-Thread\') # 创建键盘监控线程 33 Keybo_Moni_Thread.start() # 启动键盘监测线程 34 except: 35 print(\'Error:uqnable to start the thread!\') 36 37 while True: 38 ret, frame = cap.read() # 读取视频帧 39 cv2.imshow(\'Video_Show\', frame) # 显示图像 40 Key_val = cv2.waitKey(1) # 获取键值，不加此句，无法运行程序！(Ref:https://www.cnblogs.com/kissfu/p/3608016.html) 41 if Key_val == ord(\'q\'): 42 cv2.destroyAllWindows() 43 cap.release() 44 print(\'Pic Sample Finished!\') 45 break 46
1)系统实际拍摄场景 2)硬件设备Logitech_Camera
3)捕获素材如下所示(部分图像)
2、使用MATLAB工具箱进行相机标定
1)使用流程：
step1:打开Matlab 2015软件，上方栏目有一个应用程序，在应用程序下有一个工具箱，Camera Calibrator.
step2:打开工具箱，单击左上角的Add Images，选中捕获的所有的ResPic?.jpg，打开之后需要选择方块大小，取d=30mm.
step3:点击工具箱中间的绿色按钮Calibrate，单击开始我们的内参以及外参的计算和仿真,之后控制台输出Camera参数.
step4:我们可以选择到处MATLAB Script来分析Camera的标定流程.
2)仿真结果：
1 cameraParams = 2 cameraParameters (具有属性): 3 Camera Intrinsics 4 IntrinsicMatrix: [3x3 double] part1：内参矩阵 5 FocalLength: [725.6633 724.1423] 焦距F 6 PrincipalPoint: [336.0815 233.9824] 主点P 7 Skew: 0 8 Lens Distortion 镜头畸变 9 RadialDistortion: [0.0391 -0.3128] 径向畸变参数(k1,k2) xcorrected = x(1+k1r2+k2r4+k3r6) ycorrected = y(1+k1r2+k2r4+k3r6) 10 TangentialDistortion: [0 0] 切向畸变参数(p1,p2) xcorrected = x+[2p1y+p2(r2+2x2)] ycorrected = y+[2p2x+p1(r2+2y2)] 11 Camera Extrinsics 相机外参 12 RotationMatrices: [3x3x11 double] part2：世界外参旋转矩阵 13 TranslationVectors: [11x3 double] part3：世界外参平移矩阵 14 Accuracy of Estimation 计算精度 15 MeanReprojectionError: 0.1804 16 ReprojectionErrors: [48x2x11 double] 17 ReprojectedPoints: [48x2x11 double] 18 Calibration Settings 19 NumPatterns: 11 20 WorldPoints: [48x2 double] 21 WorldUnits: \'mm\' 22 EstimateSkew: 0 23 NumRadialDistortionCoefficients: 2 24 EstimateTangentialDistortion: 0
对上述脚本进行分析如下所示：
图中已经用红色标注出来了Lens的焦距f=725mm,主点分别为P=336mm和P`=234mm
part1：内参矩阵
图像参数 实际图像
从上图中的参数可知
a)相机CCD的点阵总个数为640*480(units：Dots)，人工计算U0以及V0如下：(U0,V0)=(width/2,height/2)=(320,240)
b)注意这里的分辨率dpi是指显示器的dpi，和相机的像素尺寸无关，有关ppi和dpi的区别参考线面的链接。
part2：世界外参旋转矩阵
1 val(:,:,1) = -0.0005 -0.9542 0.2991 0.9816 -0.0575 -0.1819 0.1908 0.2935 0.9367 2 val(:,:,2) = -0.0406 -0.9795 0.1971 0.9472 0.0250 0.3195 -0.3179 0.1997 0.9269 3 val(:,:,3) = 0.0217 -0.9162 0.4002 0.9898 -0.0368 -0.1378 0.1410 0.3991 0.9060 4 val(:,:,4) = 0.5801 -0.7895 0.2002 0.7831 0.4730 -0.4039 0.2242 0.3911 0.8926 val(:,:,5) = 0.9954 -0.0478 -0.0836 0.0074 0.9032 -0.4291 0.0960 0.4265 0.8994 5 val(:,:,6) = 0.4233 0.7063 -0.5674 -0.7637 0.6151 0.1959 0.4874 0.3504 0.7998 val(:,:,7) = -0.0189 -0.9922 -0.1236 0.9970 -0.0095 -0.0762 0.0745 -0.1247 0.9894 6 val(:,:,8) = 0.8786 -0.0926 -0.4685 -0.0071 0.9784 -0.2068 0.4775 0.1851 0.8589 7 val(:,:,9) = 0.0320 -0.9789 0.2020 0.7774 -0.1027 -0.6206 0.6282 0.1769 0.7577 8 val(:,:,10) = 0.0167 -0.9229 0.3846 0.9836 -0.0541 -0.1723 0.1799 0.3811 0.9069 9 val(:,:,11) = -0.5271 -0.7808 0.3356 0.8481 -0.4583 0.2658 -0.0537 0.4247 0.9038
外参旋转矩阵的来源在原理中已经介绍，每一个轴都需要旋转x，y，z旋转矩阵相乘之后就得到了如上所示的旋转矩阵参数！3*3的矩阵，最后任然为3*3.
part3：世界外参平移矩阵
外参平移矩阵主要用来对准棋盘图像的中心，T3=[dx,dy,dz]
MATLAB导出的脚本代码：Camera_Calibrator.m
1 % Auto-generated by cameraCalibrator app on 13-Mar-2019 2 %------------------------------------------------------- 3 4 5 % Define images to process 6 imageFileNames = {\'D:\PythonDevelopment\Opencv_Test\data\ResPic0.jpg\',... 7 \'D:\PythonDevelopment\Opencv_Test\data\ResPic1.jpg\',... 8 \'D:\PythonDevelopment\Opencv_Test\data\ResPic2.jpg\',... 9 \'D:\PythonDevelopment\Opencv_Test\data\ResPic3.jpg\',... 10 \'D:\PythonDevelopment\Opencv_Test\data\ResPic4.jpg\',... 11 \'D:\PythonDevelopment\Opencv_Test\data\ResPic5.jpg\',... 12 \'D:\PythonDevelopment\Opencv_Test\data\ResPic6.jpg\',... 13 \'D:\PythonDevelopment\Opencv_Test\data\ResPic7.jpg\',... 14 \'D:\PythonDevelopment\Opencv_Test\data\ResPic8.jpg\',... 15 \'D:\PythonDevelopment\Opencv_Test\data\ResPic9.jpg\',... 16 \'D:\PythonDevelopment\Opencv_Test\data\ResPic10.jpg\',... 17 }; 18 19 % Detect checkerboards in images 20 [imagePoints, boardSize, imagesUsed] = detectCheckerboardPoints(imageFileNames); 21 imageFileNames = imageFileNames(imagesUsed); 22 23 % Generate world coordinates of the corners of the squares 24 squareSize = 30; % in units of \'mm\' 25 worldPoints = generateCheckerboardPoints(boardSize, squareSize); 26 27 % Calibrate the camera 28 [cameraParams, imagesUsed, estimationErrors] = estimateCameraParameters(imagePoints, worldPoints, ... 29 \'EstimateSkew\', false, \'EstimateTangentialDistortion\', false, ... 30 \'NumRadialDistortionCoefficients\', 2, \'WorldUnits\', \'mm\', ... 31 \'InitialIntrinsicMatrix\', [], \'InitialRadialDistortion\', []); 32 33 % View reprojection errors 34 h1=figure; showReprojectionErrors(cameraParams, \'BarGraph\'); 35 36 % Visualize pattern locations 37 h2=figure; showExtrinsics(cameraParams, \'CameraCentric\'); 38 39 % Display parameter estimation errors 40 displayErrors(estimationErrors, cameraParams); 41 42 % For example, you can use the calibration data to remove effects of lens distortion. 43 originalImage = imread(imageFileNames{1}); 44 undistortedImage = undistortImage(originalImage, cameraParams); 45 46 % See additional examples of how to use the calibration data. At the prompt type: 47 % showdemo(\'MeasuringPlanarObjectsExample\') 48 % showdemo(\'StructureFromMotionExample\')
3、采用Python及Opencv进行相机标定
注意标定板的方块大小Opencv是已知的30mm，根据Opencv官网给出的标定模板，打印在A4的纸上。
代码如下：
1 import cv2 2 import glob 3 import time 4 import threading 5 import numpy as np 6 7 print(\'Starting Capture...\') 8 cap = cv2.VideoCapture(0) 9 while not cap.isOpened(): # 检查摄像头是否打开成功 10 time.sleep(100) 11 print(\'Camera is Initialize...\') 12 13 width = int(cap.get(3)) # 读取摄像头分辨率参数 14 height = int(cap.get(4)) 15 16 frame = np.zeros((width,height,3),dtype=np.uint8) # 创建图像模板 17 ret, frame = cap.read() 18 19 Key_val = 0 # 保存键值 20 21 def Keybo_Moni(): # 按键测试函数 22 count = 0 23 while True: 24 global Key_val, frame, process_flag, cap 25 if Key_val == ord(\'r\'): 26 Key_val= 0 27 cv2.imwrite(\'ResPic\' + str(count) + \'.jpg\', frame) # 保存图像 28 count += 1 29 print(\'Get new pic %d\' % count) 30 31 try: 32 Keybo_Moni_Thread = threading.Thread(target=Keybo_Moni, name=\'Keyboard-Thread\') # 创建键盘监控线程 33 Keybo_Moni_Thread.start() # 启动键盘监测线程 34 except: 35 print(\'Error:uqnable to start the thread!\') 36 37 while True: 38 ret, frame = cap.read() # 读取视频帧 39 cv2.imshow(\'Video_Show\', frame) # 显示图像 40 Key_val = cv2.waitKey(1) # 获取键值，不加此句，无法运行程序！(Ref:https://www.cnblogs.com/kissfu/p/3608016.html) 41 if Key_val == ord(\'q\'): 42 cv2.destroyAllWindows() 43 cap.release() 44 print(\'Pic Sample Finished!\') 45 break 46 47 print(\'Starting CalibrateCamera...\') 48 49 # 标定算法开始 50 51 # 设置寻找亚像素角点的参数，采用的停止准则是最大循环次数30和最大误差容限0.001 52 criteria = (cv2.TERM_CRITERIA_EPS | cv2.TERM_CRITERIA_MAX_ITER, 30, 0.001) 53 # 获取标定板角点的位置 54 ojbp = np.zeros((6*7,3),np.float32) 55 # reshape(-1,N)转换矩阵的为N列的矩阵，-1表示可以代表任意行 56 # 比如有10个元素： 57 # reshape(-1,1)==>10*1的矩阵 58 # reshape(-1,2)===>5*2的矩阵 59 # T 表示将矩阵当中的所有的元素全部顺序反转一遍 60 # 将世界坐标系建在标定板上，所有点的Z坐标全部为0，所以只需要赋值x和y 61 ojbp[:,:2] = np.mgrid[0:7,0:6].T.reshape(-1,2) 62 63 objpoints = [] # 存储世界坐标系中的3D点(实际上Zw在标定板上的值为0) 64 imgpoints = [] # 存储图像坐标系中的2D点 65 66 images = glob.glob(\'*.jpg\') 67 for fname in images: 68 img = cv2.imread(fname) 69 gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY) 70 cv2.imshow(\'gray\',gray) 71 cv2.waitKey(300) 72 # Find the chess board corners 73 ret, corners = cv2.findChessboardCorners(gray, (7,6), flags=3) 74 # If found, add object points, image points(after refining them) 75 if ret == True: 76 if fname.find(\'ResPic\') != -1: 77 time_name = str(int(time.time())) 78 cv2.imwrite(time_name+\'.jpg\',img) 79 os.remove(fname) 80 objpoints.append(ojbp) 81 corners2 = cv2.cornerSubPix(gray,corners,(11,11),(-1,-1),criteria) # 亚像素级焦点检测，基于提取的角点，进一步提高精确度 82 imgpoints.append(corners2) 83 # Draw and display the corners 84 img = cv2.drawChessboardCorners(img,(7,6),corners2,ret) 85 cv2.imshow(\'img\',img) 86 cv2.waitKey(500) 87 else: 88 os.remove(fname) 89 90 cv2.destroyAllWindows() 91 92 gray = cv2.cvtColor(cv2.imread(images[0]),cv2.COLOR_RGB2GRAY) 93 ret, mtx, dist, rvecs, tvecs = cv2.calibrateCamera(objpoints, imgpoints, gray.shape[::-1],None,None) 94 95 print(\'retval:\n\',ret) 96 print(\'Camera-Intrinsics:\n\',mtx) 97 print(\'Camera-Distortion:\n\',dist) 98 print(\'Camera-Extrinsics-R:\') 99 N = 1 100 for r in rvecs: 101 print(\'r\'+str(N)+\':\',r[0],r[1],r[2]) 102 N += 1 103 print(\'Camera-Extrinsics-t:\') 104 N = 1 105 for t in tvecs: 106 print(\'t\'+str(N)+\':\',t[0],t[1],t[2]) 107 N += 1 108 # Test The Result 109 110 # The Picture need to Correct 111 img = cv2.imread(\'./Calibsource/test2.jpg\')
请发表评论