實(shí)戰(zhàn)：如何在真實(shí)場景實(shí)現(xiàn)雙目立體匹配獲取深度圖？
點(diǎn)擊下方卡片，關(guān)注“新機(jī)器視覺”公眾號
視覺/圖像重磅干貨，第一時(shí)間送達(dá)
來源：博客園
作者：一度逍遙
文章僅用于學(xué)術(shù)分享。如有侵權(quán)，請聯(lián)系刪除。
雙目立體匹配一直是雙目視覺的研究熱點(diǎn)，雙目相機(jī)拍攝同一場景的左、右兩幅視點(diǎn)圖像，運(yùn)用立體匹配匹配算法獲取視差圖，進(jìn)而獲取深度圖。而深度圖的應(yīng)用范圍非常廣泛，由于其能夠記錄場景中物體距離攝像機(jī)的距離，可以用以測量、三維重建、以及虛擬視點(diǎn)的合成等。
目前大多數(shù)立體匹配算法使用的都是標(biāo)準(zhǔn)測試平臺提供的標(biāo)準(zhǔn)圖像對，比如著名的有如下兩個(gè)：
MiddleBury: http://vision.middlebury.edu/stereo/；
KITTI：http://www.cvlibs.net/datasets/kitti/eval_scene_flow.php?benchmark=stereo
但是對于想自己嘗試拍攝雙目圖片進(jìn)行立體匹配獲取深度圖，進(jìn)行三維重建等操作的童鞋來講，要做的工作是比使用校正好的標(biāo)準(zhǔn)測試圖像對要多的。因此博主覺得有必要從用雙目相機(jī)拍攝圖像開始，捋一捋這整個(gè)流程。
主要分四個(gè)部分講解：
攝像機(jī)標(biāo)定（包括內(nèi)參和外參）
雙目圖像的校正（包括畸變校正和立體校正）
立體匹配算法獲取視差圖，以及深度圖
利用視差圖，或者深度圖進(jìn)行虛擬視點(diǎn)的合成
注：如果沒有雙目相機(jī)，可以使用單個(gè)相機(jī)平行移動(dòng)拍攝，外參可以通過攝像機(jī)自標(biāo)定算出。我用自己的手機(jī)拍攝，拍攝移動(dòng)時(shí)盡量保證平行移動(dòng)。
一、攝像機(jī)標(biāo)定
1.內(nèi)參標(biāo)定
攝像機(jī)內(nèi)參反映的是攝像機(jī)坐標(biāo)系到圖像坐標(biāo)系之間的投影關(guān)系。攝像機(jī)內(nèi)參的標(biāo)定使用張正友標(biāo)定法，簡單易操作，具體原理請拜讀張正友的大作《A Flexible New Technique for Camera Calibration》。當(dāng)然網(wǎng)上也會(huì)有很多資料可供查閱，MATLAB 有專門的攝像機(jī)標(biāo)定工具包，OpenCV封裝好的攝像機(jī)標(biāo)定API等。
攝像機(jī)的內(nèi)參包括，fx, fy, cx, cy,以及畸變系數(shù)[k1,k2,p1,p2,k3]，詳細(xì)就不贅述。我用手機(jī)對著電腦拍攝各個(gè)角度的棋盤格圖像，棋盤格圖像如圖所示：
使用OpenCV3.4+VS2015對手機(jī)進(jìn)行內(nèi)參標(biāo)定。標(biāo)定結(jié)果如下，手機(jī)鏡頭不是魚眼鏡頭，因此使用普通相機(jī)模型標(biāo)定即可：
圖像分辨率為：3968 x 2976。上面標(biāo)定結(jié)果順序依次為fx, fy, cx, cy, k1, k2, p1, p2, k3， 保存到文件中供后續(xù)使用。
2.外參標(biāo)定
攝像機(jī)外參反映的是攝像機(jī)坐標(biāo)系和世界坐標(biāo)系之間的旋轉(zhuǎn)R和平移T關(guān)系。如果兩個(gè)相機(jī)的內(nèi)參均已知，并且知道各自與世界坐標(biāo)系之間的R1、T1和R2，T2，就可以算出這兩個(gè)相機(jī)之間的Rotation和Translation，也就找到了從一個(gè)相機(jī)坐標(biāo)系到另一個(gè)相機(jī)坐標(biāo)系之間的位置轉(zhuǎn)換關(guān)系。攝像機(jī)外參標(biāo)定也可以使用標(biāo)定板，只是保證左、右兩個(gè)相機(jī)同時(shí)拍攝同一個(gè)標(biāo)定板的圖像。外參一旦標(biāo)定好，兩個(gè)相機(jī)的結(jié)構(gòu)就要保持固定，否則外參就會(huì)發(fā)生變化，需要重新進(jìn)行外參標(biāo)定。
那么手機(jī)怎么保證拍攝同一個(gè)標(biāo)定板圖像并能夠保持相對位置不變，這個(gè)是很難做到的，因?yàn)楹罄m(xù)用來拍攝實(shí)際測試圖像時(shí)，手機(jī)的位置肯定會(huì)發(fā)生變化。因此我使用外參自標(biāo)定的方法，在拍攝實(shí)際場景的兩張圖像時(shí)，進(jìn)行攝像機(jī)的外參自標(biāo)定，從而獲取當(dāng)時(shí)兩個(gè)攝像機(jī)位置之間的Rotation和Translation。
比如：我拍攝這樣兩幅圖像，以后用來進(jìn)行立體匹配和虛擬視點(diǎn)合成的實(shí)驗(yàn)。
cv::Mat E = cv::findEssentialMat(tmpPts1, tmpPts2, camK, CV_RANSAC);cv::Mat R1, R2;cv::decomposeEssentialMat(E, R1, R2, t);R = R1.clone();t = -t.clone();
cv::stereoRectify(camK, D, camK, D, imgL.size(), R, -R*t, R1, R2, P1, P2, Q);
cv::initUndistortRectifyMap(P1(cv::Rect(0, 0, 3, 3)), D, R1, P1(cv::Rect(0, 0, 3, 3)), imgL.size(), CV_32FC1, mapx, mapy);cv::remap(imgL, recImgL, mapx, mapy, CV_INTER_LINEAR);cv::imwrite("data/recConyL.png", recImgL);cv::initUndistortRectifyMap(P2(cv::Rect(0, 0, 3, 3)), D, R2, P2(cv::Rect(0, 0, 3, 3)), imgL.size(), CV_32FC1, mapx, mapy);cv::remap(imgR, recImgR, mapx, mapy, CV_INTER_LINEAR);cv::imwrite("data/recConyR.png", recImgR);
int numberOfDisparities = ((imgSize.width / 8) + 15) & -16;cv::Ptr<cv::StereoSGBM> sgbm = cv::StereoSGBM::create(0, 16, 3);sgbm->setPreFilterCap(32);int SADWindowSize = 9;int sgbmWinSize = SADWindowSize > 0 ? SADWindowSize : 3;sgbm->setBlockSize(sgbmWinSize);int cn = imgL.channels();sgbm->setP1(8 * cn*sgbmWinSize*sgbmWinSize);sgbm->setP2(32 * cn*sgbmWinSize*sgbmWinSize);sgbm->setMinDisparity(0);sgbm->setNumDisparities(numberOfDisparities);sgbm->setUniquenessRatio(10);sgbm->setSpeckleWindowSize(100);sgbm->setSpeckleRange(32);sgbm->setDisp12MaxDiff(1);int alg = STEREO_SGBM;if (alg == STEREO_HH)sgbm->setMode(cv::StereoSGBM::MODE_HH);else if (alg == STEREO_SGBM)sgbm->setMode(cv::StereoSGBM::MODE_SGBM);else if (alg == STEREO_3WAY)sgbm->setMode(cv::StereoSGBM::MODE_SGBM_3WAY);sgbm->compute(imgL, imgR, disp);
sgbm->setMinDisparity(-numberOfDisparities);sgbm->setNumDisparities(numberOfDisparities);sgbm->compute(imgR, imgL, disp);disp = abs(disp);
depth = ( f * baseline) / disp
/*函數(shù)作用：視差圖轉(zhuǎn)深度圖輸入： dispMap ----視差圖，8位單通道，CV_8UC1 K ----內(nèi)參矩陣，float類型輸出： depthMap ----深度圖，16位無符號單通道，CV_16UC1*/void disp2Depth(cv::Mat dispMap, cv::Mat &depthMap, cv::Mat K){int type = dispMap.type();float fx = K.at<float>(0, 0);float fy = K.at<float>(1, 1);float cx = K.at<float>(0, 2);float cy = K.at<float>(1, 2);float baseline = 65; //基線距離65mmif (type == CV_8U){const float PI = 3.14159265358;int height = dispMap.rows;int width = dispMap.cols;uchar* dispData = (uchar*)dispMap.data;ushort* depthData = (ushort*)depthMap.data;for (int i = 0; i < height; i++){for (int j = 0; j < width; j++){int id = i*width + j;if (!dispData[id]) continue; //防止0除depthData[id] = ushort( (float)fx *baseline / ((float)dispData[id]) );}}}else{cout << "please confirm dispImg's type!" << endl;cv::waitKey(0);}}
void insertDepth32f(cv::Mat& depth){const int width = depth.cols;const int height = depth.rows;float* data = (float*)depth.data;cv::Mat integralMap = cv::Mat::zeros(height, width, CV_64F);cv::Mat ptsMap = cv::Mat::zeros(height, width, CV_32S);double* integral = (double*)integralMap.data;int* ptsIntegral = (int*)ptsMap.data;memset(integral, 0, sizeof(double) * width * height);memset(ptsIntegral, 0, sizeof(int) * width * height);for (int i = 0; i < height; ++i){int id1 = i * width;for (int j = 0; j < width; ++j){int id2 = id1 + j;if (data[id2] > 1e-3){integral[id2] = data[id2];ptsIntegral[id2] = 1;}}}// 積分區(qū)間for (int i = 0; i < height; ++i){int id1 = i * width;for (int j = 1; j < width; ++j){int id2 = id1 + j;integral[id2] += integral[id2 - 1];ptsIntegral[id2] += ptsIntegral[id2 - 1];}}for (int i = 1; i < height; ++i){int id1 = i * width;for (int j = 0; j < width; ++j){int id2 = id1 + j;integral[id2] += integral[id2 - width];ptsIntegral[id2] += ptsIntegral[id2 - width];}}int wnd;double dWnd = 2;while (dWnd > 1){wnd = int(dWnd);dWnd /= 2;for (int i = 0; i < height; ++i){int id1 = i * width;for (int j = 0; j < width; ++j){int id2 = id1 + j;int left = j - wnd - 1;int right = j + wnd;int top = i - wnd - 1;int bot = i + wnd;left = max(0, left);right = min(right, width - 1);top = max(0, top);bot = min(bot, height - 1);int dx = right - left;int dy = (bot - top) * width;int idLeftTop = top * width + left;int idRightTop = idLeftTop + dx;int idLeftBot = idLeftTop + dy;int idRightBot = idLeftBot + dx;int ptsCnt = ptsIntegral[idRightBot] + ptsIntegral[idLeftTop] - (ptsIntegral[idLeftBot] + ptsIntegral[idRightTop]);double sumGray = integral[idRightBot] + integral[idLeftTop] - (integral[idLeftBot] + integral[idRightTop]);if (ptsCnt <= 0){continue;}data[id2] = float(sumGray / ptsCnt);}}int s = wnd / 2 * 2 + 1;if (s > 201){s = 201;}cv::GaussianBlur(depth, depth, cv::Size(s, s), s, s);}}
—版權(quán)聲明—
僅用于學(xué)術(shù)分享，版權(quán)屬于原作者。
若有侵權(quán)，請聯(lián)系微信號:yiyang-sy 刪除或修改！