计算机视觉课后作业5——相机标定

                      相机标定

1.实验要求：

2.实验环境:

pyCharm

3.实验简介：

相机标定——顾名思义，是针对照相机的相关实验。那么何为标定？我们需要了解到，相机由于生产过程中光学仪器的不精确性等原因，其参数并不像厂商所标注的那么精确，对于成像的准确度会产生误差。相机标定要做的就是建立图像像素位置与场景点位置之间的关系，根据相机成像模型，由特征点在图像中坐标与世界坐标的对应关系，求解相机模型的内部参数和外部参数。

4.实验原理：

相机的成像过程涉及到四个坐标系——世界坐标系、相机坐标系、图像坐标系、像素坐标系及其转换。

世界坐标系：用户定义的三维世界的坐标系，描述目标物在真实世界里的位置。

相机坐标系：从相机的角度描述物体位置，作为沟通世界坐标系和图像/像素坐标系的中间一环。

图像坐标系：描述成像过程中物体从相机坐标系到图像坐标系的投影透射关系，方便进一步得到像素坐标系下的坐标。

像素坐标系：描述物体成像后的像点在数字图像上（相片）的坐标，是我们真正从相机内读取到的信息所在的坐标系。

针孔照相机模型如图所示：

 

小孔成像基本模型：

在此模型下，物体的世界坐标和图像坐标之间是线性的关系，因而对相机参数的求解就归结到求解线性方程组上。四个坐标系的关系图如下图所示，其中 M 为三维空间点，m 为 M 在图像平面投影成的像点。

首先，根据相机针孔原理，可以推导出空间3D点M的坐标和它的像素坐标之间的关系如下所示

                                                                                                            （1）

式中：，表示像素坐标系下的坐标；，表示世界坐标系下的坐标；s为尺度因子；R为世界坐标系与图像坐标系的旋转矩阵，t为世界坐标系与图像坐标系的平移向量，R和t中的参数均被称为相机的外部参数，简称外参；矩阵A被称为内参矩阵，如下所示：

 

 

矩阵A中，为光轴与像平面焦点在像素坐标系下的坐标，将它成为主点；α和β是图像u轴和v轴上的尺度因子，γ为像素两轴非正交项。

其次，是之前学过的单应性变换矩阵，目的是将一张图中的点映射到另一张图中对应的点，进行特征点匹配，求解Homography 矩阵。

。

世界坐标系（Z=0）平面下的点M及其在像平面中的像坐标关系可以用如下的单应性矩阵H联系起来：

                                          其中，                                （2）

最后对相关相机矩阵求解内部参数和外部参数，涉及到的方法属于线性微分方程部分的内容。

上图所示模型是理想模型，正如前文所述，实际上由于生产过程中光学仪器的不精确性等原因，造成相机图像平面上实际所成的像与理想成像之间会存在畸变现象。我们这边主要讨论径向畸变，径向畸变分为两种，分别是枕形畸变和桶型畸变。枕形畸变是在长焦距下发生的，而桶形畸变是在中短焦距下发生的。如下图所示，前阵是枕形畸变，后者为桶形畸变：

从图像可以看出，径向畸变以某一个中心往外延伸，且越往外，畸变越大；显然畸变与距离成一种非线性的变换关系，可以用多项式来近似。

                  

x,y是归一化的图像坐标，即坐标原点已经移动到主点，并且像素坐标除以焦距。K1,k2,k3是径向畸变系数，r2=x2+y2（2表示2次方）。

如前文所述，我们求解相机参数是在物体的世界坐标和图像坐标之间求解，那我们就必须使用现实世界的三维图像。实验要求中提到的棋盘格则是出自张正友教授提出的棋盘格标定法——通过棋盘格的多视角拍摄多张图片，用二维图像代替三维图像进行相机标定，大大减少了标定的计算量。并且，棋盘格角点突出，有助于我们找到有效的特征点，提高标定准确度。

棋盘格是打印下来的一张由黑白方块间隔组成的A4纸，并把它夹在板上固定，作为相机标定实验的标定物。如前文所示，二维平面比三维立体图形更容易标定，但与此同时，二维平面相对于三维立体图形会缺少一部分信息，于是我们通过拍摄不同角度的棋盘格，来获得更丰富的坐标信息。

5.实验数据——棋盘图（11张）：

手机型号：红米note9Pro

6.实验代码：

import cv2
import numpy as np
import glob

# 找棋盘格角点
# 阈值
criteria = (cv2.TERM_CRITERIA_EPS + cv2.TERM_CRITERIA_MAX_ITER, 30, 0.001)
#棋盘格模板规格
w = 6   #内角点个数，内角点是和其他格子连着的点
h = 4

# 世界坐标系中的棋盘格点,例如(0,0,0), (1,0,0), (2,0,0) ....,(8,5,0)，去掉Z坐标，记为二维矩阵
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('C:/Users/oulia/Pictures/camera/*.jpg')
for fname in images:
    img = cv2.imread(fname)
    gray = cv2.cvtColor(img,cv2.COLOR_BGR2GRAY)
    # 找到棋盘格角点
    # 棋盘图像(8位灰度或彩色图像)  棋盘尺寸  存放角点的位置
    ret, corners = cv2.findChessboardCorners(gray, (w,h),None)
    # 如果找到足够点对，将其存储起来
    if ret == True:
        # 角点精确检测
        # 输入图像 角点初始坐标 搜索窗口为2*winsize+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.waitKey(1000)
cv2.destroyAllWindows()
#标定、去畸变
# 输入：世界坐标系里的位置 像素坐标 图像的像素尺寸大小 3*3矩阵，相机内参数矩阵 畸变矩阵
# 输出：标定结果 相机的内参数矩阵 畸变系数 旋转矩阵 平移向量
ret, mtx, dist, rvecs, tvecs = cv2.calibrateCamera(objpoints, imgpoints, gray.shape[::-1], None, None)
# mtx：内参数矩阵
# dist：畸变系数
# rvecs：旋转向量 （外参数）
# tvecs ：平移向量 （外参数）
print (("ret:"),ret)
print (("mtx:\n"),mtx)        # 内参数矩阵
print (("dist:\n"),dist)      # 畸变系数   distortion cofficients = (k_1,k_2,p_1,p_2,k_3)
print (("rvecs:\n"),rvecs)    # 旋转向量  # 外参数
print (("tvecs:\n"),tvecs)    # 平移向量  # 外参数
# 去畸变
img2 = cv2.imread('C:/Users/oulia/Pictures/camera/1.jpg')
h,w = img2.shape[:2]
# 我们已经得到了相机内参和畸变系数，在将图像去畸变之前，
# 我们还可以使用cv.getOptimalNewCameraMatrix()优化内参数和畸变系数，
# 通过设定自由自由比例因子alpha。当alpha设为0的时候，
# 将会返回一个剪裁过的将去畸变后不想要的像素去掉的内参数和畸变系数；
# 当alpha设为1的时候，将会返回一个包含额外黑色像素点的内参数和畸变系数，并返回一个ROI用于将其剪裁掉
newcameramtx, roi=cv2.getOptimalNewCameraMatrix(mtx,dist,(w,h),0,(w,h)) # 自由比例参数

dst = cv2.undistort(img2, mtx, dist, None, newcameramtx)
# 根据前面ROI区域裁剪图片
x,y,w,h = roi
dst = dst[y:y+h, x:x+w]
cv2.imwrite('calibresult.jpg',dst)

# 反投影误差
# 通过反投影误差，我们可以来评估结果的好坏。越接近0，说明结果越理想。
# 通过之前计算的内参数矩阵、畸变系数、旋转矩阵和平移向量，使用cv2.projectPoints()计算三维点到二维图像的投影，
# 然后计算反投影得到的点与图像上检测到的点的误差，最后计算一个对于所有标定图像的平均误差，这个值就是反投影误差。
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))

7.实验结果截图：

特征点匹配过程：

内参数矩阵：

畸变参数：

旋转向量（第一张图）：

平移向量（第一张图）：

8.实验小结：

使用cv.getOptimalNewCameraMatrix()优化内参数和畸变系数，修改自由自由比例因子alpha=1后得到的实验结果和=0也是一样的，应该说只是多进行了一步ROI的裁剪。在结果部分无法直接体现出优化的结果。但是是提供了一种可以对内部参数合畸变参数进行不同取值，不同操作优化的方法。