雙攝像頭立體成像（二）-攝像頭標定

寫在題前：這篇文章磨磨蹭蹭了好久，曾經兩次接近完稿而丟失。我想任何事情在起步時都會有類似的囧境，還好我還有恆心繼續下去。

攝像頭標定的目的有兩個。第一，要還原攝像頭成像的物體在真實世界的位置就需要知道世界中的物體到計算機圖像平面是如何變換的，攝像頭標定的目的之一就是為了搞清楚這種變換關系，求解內外參數矩陣。第二，針孔攝像頭的發明使得攝像頭變成了親民物品，大行於世，但是針孔攝像頭有個很大的問題——畸變。攝像頭標定的另一個目的就是求解畸變系數，然后用於計算求解正確的成像。

• 數學原理

  • 映射矩陣-內部參數(intrinsic parameters)和外部參數(extrinsic parameters)

前面所說的真實世界到成像平面的變換過程牽扯到四個坐標系，變換矩陣可以大致分為兩組參數——內部參數和外部參數。下面依次根據四種坐標系的關系來推導內外參數的形式。

計算機坐標系和成像平面坐標系

計算機坐標系是指數字圖像在計算機中的保存形式-二維數組的坐標形式。成像平面坐標系是攝像機鏡頭的成像平面上建立的坐標系，一般是位於感光器件上（如CCD）所建立的以鏡頭光軸與成像平面交點為原點的二維坐標系。由於二維數組以離散的形式保存，所以計算機坐標系的長度單位為1，一般感光器件大小約為指甲蓋大小，上面又密集的分別很多感光單元，每個感光單元采集的圖像都會計算成計算機坐標系中的像素值，所以一般長度單位為微米。具體二者關系見下圖，

圖中坐標系O-uv為計算機坐標系，其坐標軸方向從右上到左下遞增。坐標系Q-xy為成像平面坐標系，其坐標原點Q在計算機坐標系的坐標為(u0,v0)。假設有一點P，其在計算機坐標系中的坐標為(u,v)，在成像平面坐標系中的坐標為(x,y)，並設像素單元在x,y方向的尺寸分別為a,b。則有如下等式成立，

整理成齊次變換矩陣的形式如下，

成像平面坐標系和攝像機坐標系

攝像機坐標系是一個以攝像機鏡頭的光心為原點而建立的空間三維坐標系。成像平面坐標系可以看成攝像機坐標系在其Zc軸上投影而成的（其x和y軸與成像平面坐標系方向一致）。其與成像平面坐標系的關系如下圖，

這里需要詳細解釋一下上圖怎樣理解。圖中Oc-XcYcZc為攝像機坐標系，M點為物點，而Q-xy為成像平面坐標系。實際情況下成像平面與物體應該分居攝像頭兩側，但是在這里我們將成像平面以攝像頭坐標系原點為中心對稱點對稱過來。這樣的好處是本來倒立的像變成正立的，且大小不變。如上一篇博客寫到的那樣，中心對稱的后的像點和物點連線必定過光心。我們要記住小孔成像的性質是，成像平面到鏡頭的距離始終等於焦距f。由相似三角形，有如下等式成立，

整理成齊次變換矩陣的形式如下，

攝像機坐標系和世界坐標系

世界坐標系是區別於攝像機坐標系的一個空間坐標系，其依附於我們標定時要使用的物點而存在。既我們觀察的物點相對於世界坐標系的位置是固定的，如張正友標定法中就是選的世界坐標系為以chessboard平面為XOY平面的相對於chessboard不動的空間坐標系。下圖是我費了老大勁繪制的攝像機坐標系和世界坐標系之間的位置關系，

我們知道空間中的兩個坐標系可以通過平移旋轉變換而重合，那么現在我們來簡單粗暴解釋一下物點M在兩個坐標系下坐標的代數關系是什么樣的。假設攝像機坐標系可以通過旋轉平移變換與世界坐標系重合，其齊次變換矩陣形式如下，

也就是說攝像機坐標系上的每個點(在世界坐標系中的坐標)通過上述旋轉平移變換后與世界坐標系中的相應點重合。我們可以理解的是，對攝像機坐標系上的每個點和物點M(在世界坐標系下的坐標)做相同的旋轉平移變換后其相對位置不會改變。也就是說物點M在攝像機坐標系中的坐標不會改變。我們考察一下變換完后是什么情況。這個時候世界坐標系和攝像機坐標系重合，且物點M相對於攝像機坐標系的坐標沒有改變，那么這個時候物點M在世界坐標系中的坐標應該變成了變換前其在攝像機坐標系下的坐標。所以有如下關系成立，

綜合上面的三組關系，我們可以推導出世界坐標系和計算機坐標系的直接代數關系，

進一步化簡，

我們記，

A矩陣式刻畫攝像機的內部參數，包括焦距f、成像中心的位置、及成像單元尺寸，因此稱為內參矩陣。M描述的攝像頭的運動關系，所以稱為外參矩陣。其中R有三個相關參數，分別是繞x，y，z軸的旋轉角度。T也有三個相關參數，分別是在三個坐標軸上的平移量Tx、Ty、Tz。

• 攝像頭的畸變參數(distortion parameters)

攝像頭的畸變是由於成像模型的不精確造成的。人們為了提高光通量用透鏡代替小孔來成像，由於這種代替不能完全符合小孔成像的性質，因此畸變就產生了。另外這里再插句額外的話，現在大量使用的透鏡為球面鏡，原因是其廉價易得。但是真正的完全符合理想光學系統的透鏡實際是個四次曲面（很好證明，依據光程不變），制造成本很大哦。

畸變可以分為兩大類，徑向畸變和切向畸變。詳細的畸變介紹可以參考工程光學的相關課程，下面簡單介紹相關畸變及其修正。

徑向畸變（radial distortion）

徑向畸變的效應有兩種，一種是枕形效應，另一種是桶形效應，具體見下圖(圖片來自互聯網)，

徑向畸變可用下面公式修正，

切向畸變（tangential distortion）

徑向畸變是由於透鏡與成像平面不嚴格的平行，其可以用如下公式修正，

這樣又引入了五個畸變參數，

• 小結

我們記fx=f/a,fy=f/b。通過上面的介紹，我們了解到，攝像機標定共有4個內參，6個外參和五個畸變參數要求。下面就介紹怎么基於OpenCV函數庫標定求得這三組參數。其求解原理放在以后的博客中敘述。

• 基於OpenCV的標定程序

OpenCV中有標定實例哦，寫的很好，功能很完善。一個是基於命令行標定參數讀入的標定程序，另一個是基於xml文件參數讀入的標定程序。它們的位置分別為，

...\opencv\sources\samples\cpp\calibration.cpp
...\opencv\sources\samples\cpp\tutorial_code\calib3d\camera_calibration\ camera_calibration.cpp

但是為了更深入的理解OpenCV的標定方法庫，我自己寫了一個簡單粗暴易讀的標定程序。下面簡單介紹一下標定的過程。

• 標定過程簡介

標定過程如下，

1. 圖像獲取
2. 角點檢測，如果沒有檢測到角點重復第一步
3. 亞像素檢測以提高角點檢測精度
4. 標記檢測出的角點
5. 如果成功檢測到角點圖片小於預設的數目，重復上面四個步驟
6. 標定
7. 標定結果保存

我寫的標定程序如下，

/*
    Writer: Wang Xianshun
    Email:    [e-mail]
*/
#include <iostream>
#include <stdio.h>
#include <time.h>
#include <string.h>

#include <cv.hpp>
#include <highgui\highgui.hpp>
#include <calib3d\calib3d.hpp>
#include <imgproc\imgproc.hpp>
#include <core\core.hpp>

using namespace std;
using namespace cv;

static void calcChessboardCorners(Size boardSize, float squareSize, vector<Point3f>& corners)
{
    corners.resize(0);
    for (int i = 0; i < boardSize.height; i++)        //height和width位置不能顛倒
    for (int j = 0; j < boardSize.width; j++)
    {
        corners.push_back(Point3f(j*squareSize, i*squareSize, 0));
    }
}

int main(int argc, char** argv)
{
    int success = 0;
    int cameraId = 0;
    int nFrames = 10;
    int w = 6;
    int h = 9;
    clock_t prevTimestamp = 0;
    int delay = 1000;

    //相關參數初始化
    Size boardSize, imageSize;
    boardSize.width = w;
    boardSize.height = h;
    vector<vector<Point2f>> imagePoints;
    float squareSize = 1.f;
    Mat intrMatrix, distCoeffs;
    vector<Mat> rvecs, tvecs;

    //標定參數讀取
    if (argc < 5)
    {
        cout << "參數不足" << endl;
        return 0;
    }

    for (int i = 1; i < argc; i++)
    {
        if (!strcmp(argv[i], "-w"))
        {
            if (!sscanf(argv[++i], "%u", &boardSize.width))
            {
                return fprintf(stderr, "無效的標定角點寬度\n"), -1;
            }
        }
        else if (!strcmp(argv[i], "-h"))
        {
            if (!sscanf(argv[++i], "%u", &boardSize.height))
            {
                return fprintf(stderr, "無效的標定角點高度\n"), -1;
            }
        }
        else if (!strcmp(argv[i], "-s"))
        {
            if (!sscanf(argv[++i], "%f", &squareSize) != 1 || squareSize <= 0)
            {
                return fprintf(stderr, "無效的方格尺寸\n"), -1;
            }
        }
        else
            return fprintf(stderr, "未知參數\n"), -1;
    }

    //圖像采集
    VideoCapture capture;
    capture.open(cameraId);
    namedWindow("Image View", 1);

    if (!capture.isOpened())
    {
        cout << "無法打開攝像頭，(づ￣3￣)づ╭❤～……" << endl;
        return -1;
    }

    for (int i = 0; success < nFrames; i++)
    {
        string msg = "CAPTURING";
        Mat viewGray, view;
        capture >> view;
        imageSize = view.size();
        vector<Point2f> pointBuf;
        cvtColor(view, viewGray, COLOR_BGR2GRAY);

        //尋找角點
        bool found = findChessboardCorners(view, boardSize, pointBuf,
            CV_CALIB_CB_ADAPTIVE_THRESH | CV_CALIB_CB_FAST_CHECK | CV_CALIB_CB_NORMALIZE_IMAGE);

        if (found)
        {
            //亞像素檢測以提高精度
            cornerSubPix(viewGray, pointBuf, Size(11, 11),
                Size(-1, -1), TermCriteria(CV_TERMCRIT_EPS + CV_TERMCRIT_ITER, 30, 0.1));
            //標記出檢測到的角點
            drawChessboardCorners(view, boardSize, Mat(pointBuf), found);
        }

        //等待用戶改變姿態
        if (found && clock() - prevTimestamp > delay*1e-3*CLOCKS_PER_SEC)
        {
            imagePoints.push_back(pointBuf);
            prevTimestamp = clock();
            success = success + 1;
            bitwise_not(view, view);
        }

        imshow("Image View", view);
        waitKey(20);
    }

    cout << "圖像采集完成，開始標定……" << endl;

    //標定
    vector<vector<Point3f>> ObjectPoints(1);
    calcChessboardCorners(boardSize, squareSize, ObjectPoints[0]);
    ObjectPoints.resize(imagePoints.size(), ObjectPoints[0]);
    calibrateCamera(ObjectPoints, imagePoints, imageSize, intrMatrix,
        distCoeffs, rvecs, tvecs);
    bool ok = checkRange(intrMatrix) && checkRange(distCoeffs);

    if (!ok)
    {
        cout << "標定失敗，再來一次" << endl;
        return -3;
    }

    //標定結果保存
    FileStorage fs("caliResult.xml", FileStorage::WRITE);
    fs << "cameraId" << cameraId;
    fs << "intrinsic_parameters" << intrMatrix;
    fs << "distortion_parametes" << distCoeffs;
    fs.release();

    return 0;
}

 

該程序有三個參數輸入，基於命令行讀入參數

-w    #標定板一個方向上的角點數

-h    #標定板另一個方向上的角點數

-s    #標定板上正方形的邊長，默認為1

 

另，發表下-s參數設置的觀點。之所以該參數在很多標定實例程序中設置為默認1，是因為該參數的改變確實是會影響到標定結果，但是不會影響到攝像頭的矯正。因為標定和矯正類似一個逆運算過程，單位定義對其沒有影響。

• 實驗結果

標定過程，

標定結果，

<?xml version="1.0"?>
<opencv_storage>
<cameraId>0</cameraId>
<intrinsic_parameters type_id="opencv-matrix">
  <rows>3</rows>
  <cols>3</cols>
  <dt>d</dt>
  <data>
    7.7881772950073355e+002 0. 3.1562441595543476e+002 0.
    7.8624564811643825e+002 2.5630331974129393e+002 0. 0. 1.</data></intrinsic_parameters>
<distortion_parametes type_id="opencv-matrix">
  <rows>1</rows>
  <cols>5</cols>
  <dt>d</dt>
  <data>
    -7.2660835182078581e-002 2.0765291395491934e+000
    5.9477659924542790e-004 -8.2981148319346263e-004
    -7.0307616798578119e+000</data></distortion_parametes>
</opencv_storage>