魚眼鏡頭標定及畸變校正
魚眼攝像頭畸變校正的方法：
- 棋盤矯正法
- 經緯度矯正法。
相機為什麼會出現畸變？
當前相機的畸變主要分為徑向畸變和切向畸變兩種。
徑向畸變產生的原因：相機的光學鏡頭厚度不均勻，離鏡頭越遠場景的光線就越彎曲從而產生徑向畸變。
切向畸變產生的原因：鏡頭與影像感測器不完全平行造成的。
相機引數有哪些？
相機內參：主要包括相機矩陣（包括焦距，光學中心，這些都是相機本身屬性）和畸變係數（畸變數學模型的5個引數 D = {K1,K2,K3,P1,P2}）。
相機外參：通過旋轉和平移將實際場景3D對映到相機的2D座標過程中的旋轉和平移就是外參。（他描述的是世界座標轉化成相機座標的過程）
相機的標定流程
相機標定流程就是4個座標系在轉換過程中求出計算機的外參和內參的過程。四個座標系分別是：世界座標系（真實的實際場景），相機座標系（攝像頭鏡頭中心），影像座標系（影像感測器成像中心，圖片中心，影布中心），畫素座標系（影像左上角為原點）。如圖三所示，O1是影像座標系，O0 是畫素座標系，兩者之間的區別只是原點發生了變化。
世界座標系 ->相機座標系 求解外參（旋轉和平移）
相機座標系 ->影像座標系 求解內參（攝像頭矩陣，畸變係數）
影像座標系 ->畫素座標系 求解畫素轉化矩陣（可簡單理解為原點從圖片中心到左上角，單位釐米變行列）
一、基於棋盤的相機畸變校正方法
列印棋盤並採集魚眼攝像頭下的棋盤圖片：
1.棋盤獲取：連結: https://pan.baidu.com/s/14qB3kQ_MbWORay1i0GFm1A 提取碼: ksqw 複製這段內容後開啟百度網盤手機App，操作更方便哦
2.採集圖片，採集圖片如下圖所示，可以多采集一些，標註的會更加準確。
3.使用採集的圖片求出相機的內參和矯正係數（DIM, K, D），然後使用得到的（DIM, K, D）再進行測試，程式碼如下。
import cv2
import numpy as np
import glob
def get_K_and_D(checkerboard, imgsPath):
CHECKERBOARD = checkerboard
subpix_criteria = (cv2.TERM_CRITERIA_EPS+cv2.TERM_CRITERIA_MAX_ITER, 30, 0.1)
calibration_flags = cv2.fisheye.CALIB_RECOMPUTE_EXTRINSIC+cv2.fisheye.CALIB_CHECK_COND+cv2.fisheye.CALIB_FIX_SKEW
objp = np.zeros((1, CHECKERBOARD[0]*CHECKERBOARD[1], 3), np.float32)
objp[0,:,:2] = np.mgrid[0:CHECKERBOARD[0], 0:CHECKERBOARD[1]].T.reshape(-1, 2)
_img_shape = None
objpoints = []
imgpoints = []
images = glob.glob(imgsPath + '/*.png')
for fname in images:
img = cv2.imread(fname)
if _img_shape == None:
_img_shape = img.shape[:2]
else:
assert _img_shape == img.shape[:2], "All images must share the same size."
gray = cv2.cvtColor(img,cv2.COLOR_BGR2GRAY)
ret, corners = cv2.findChessboardCorners(gray, CHECKERBOARD,cv2.CALIB_CB_ADAPTIVE_THRESH+cv2.CALIB_CB_FAST_CHECK+cv2.CALIB_CB_NORMALIZE_IMAGE)
if ret == True:
objpoints.append(objp)
cv2.cornerSubPix(gray,corners,(3,3),(-1,-1),subpix_criteria)
imgpoints.append(corners)
N_OK = len(objpoints)
K = np.zeros((3, 3))
D = np.zeros((4, 1))
rvecs = [np.zeros((1, 1, 3), dtype=np.float64) for i in range(N_OK)]
tvecs = [np.zeros((1, 1, 3), dtype=np.float64) for i in range(N_OK)]
rms, _, _, _, _ = cv2.fisheye.calibrate(
objpoints,
imgpoints,
gray.shape[::-1],
K,
D,
rvecs,
tvecs,
calibration_flags,
(cv2.TERM_CRITERIA_EPS+cv2.TERM_CRITERIA_MAX_ITER, 30, 1e-6)
)
DIM = _img_shape[::-1]
print("Found " + str(N_OK) + " valid images for calibration")
print("DIM=" + str(_img_shape[::-1]))
print("K=np.array(" + str(K.tolist()) + ")")
print("D=np.array(" + str(D.tolist()) + ")")
return DIM, K, D
def undistort(img_path,K,D,DIM,scale=0.6,imshow=False):
img = cv2.imread(img_path)
dim1 = img.shape[:2][::-1] #dim1 is the dimension of input image to un-distort
assert dim1[0]/dim1[1] == DIM[0]/DIM[1], "Image to undistort needs to have same aspect ratio as the ones used in calibration"
if dim1[0]!=DIM[0]:
img = cv2.resize(img,DIM,interpolation=cv2.INTER_AREA)
Knew = K.copy()
if scale:#change fov
Knew[(0,1), (0,1)] = scale * Knew[(0,1), (0,1)]
map1, map2 = cv2.fisheye.initUndistortRectifyMap(K, D, np.eye(3), Knew, DIM, cv2.CV_16SC2)
undistorted_img = cv2.remap(img, map1, map2, interpolation=cv2.INTER_LINEAR, borderMode=cv2.BORDER_CONSTANT)
if imshow:
cv2.imshow("undistorted", undistorted_img)
return undistorted_img
if __name__ == '__main__':
# 開始使用圖片來獲取內參和畸變係數
DIM, K, D = get_K_and_D((6,9), '')
# 得到內參和畸變係數畸變矯正進行測試
'''
DIM=(2560, 1920)
K=np.array([[652.8609862494474, 0.0, 1262.1021584894233], [0.0, 653.1909758659955, 928.0871455436396], [0.0, 0.0, 1.0]])
D=np.array([[-0.024092199861108887], [0.002745976275100771], [0.002545415522352827], [-0.0014366825722748522]])
img = undistort('../imgs/pig.jpg',K,D,DIM)
cv2.imwrite('../imgs/pig_checkerboard.jpg', img)
'''
二、基於經緯度的矯正方法
1.演算法原理：經緯度矯正法， 可以把魚眼圖想象成半個地球， 然後將地球展開成地圖，經緯度矯正法主要是利用幾何原理， 對影像進行展開矯正。（此演算法操作簡單不需要使用棋盤來進行資料採集）
2.程式碼實現：
import cv2
import numpy as np
import time
# 魚眼有效區域擷取
def cut(img):
img_gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
(_, thresh) = cv2.threshold(img_gray, 20, 255, cv2.THRESH_BINARY)
contours, hierarchy = cv2.findContours(thresh, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
cnts = sorted(contours, key=cv2.contourArea, reverse=True)[0]
x,y,w,h = cv2.boundingRect(cnts)
r = max(w/ 2, h/ 2)
# 提取有效區域
img_valid = img[y:y+h, x:x+w]
return img_valid, int(r)
# 魚眼矯正
def undistort(src,r):
# r： 半徑， R: 直徑
R = 2*r
# Pi: 圓周率
Pi = np.pi
# 儲存對映結果
dst = np.zeros((R, R, 3))
src_h, src_w, _ = src.shape
# 圓心
x0, y0 = src_w//2, src_h//2
# 陣列， 迴圈每個點
range_arr = np.array([range(R)])
theta = Pi - (Pi/R)*(range_arr.T)
temp_theta = np.tan(theta)**2
phi = Pi - (Pi/R)*range_arr
temp_phi = np.tan(phi)**2
tempu = r/(temp_phi + 1 + temp_phi/temp_theta)**0.5
tempv = r/(temp_theta + 1 + temp_theta/temp_phi)**0.5
# 用於修正正負號
flag = np.array([-1] * r + [1] * r)
# 加0.5是為了四捨五入求最近點
u = x0 + tempu * flag + 0.5
v = y0 + tempv * np.array([flag]).T + 0.5
# 防止陣列溢位
u[u<0]=0
u[u>(src_w-1)] = src_w-1
v[v<0]=0
v[v>(src_h-1)] = src_h-1
# 插值
dst[:, :, :] = src[v.astype(int),u.astype(int)]
return dst
if __name__ == "__main__":
t = time.perf_counter()
frame = cv2.imread('../imgs/pig.jpg')
cut_img,R = cut(frame)
t = time.perf_counter()
result_img = undistort(cut_img,R)
cv2.imwrite('../imgs/pig_vector_nearest.jpg',result_img)
print(time.perf_counter()-t)
影像演算法中會經常用到攝像機的畸變校正，有必要總結分析OpenCV中畸變校正方法，其中包括普通針孔相機模型和魚眼相機模型fisheye兩種畸變校正方法。
普通相機模型畸變校正函式針對OpenCV中的cv::initUndistortRectifyMap()，魚眼相機模型畸變校正函式對應OpenCV中的cv::fisheye::initUndistortRectifyMap()。兩種方法算出對映Mapx和Mapy後，統一用cv::Remap()函式進行插值得到校正後的影像。
1. FishEye模型的畸變校正。
方便起見，直接貼出OpenCV原始碼，我在裡面加了註釋說明。建議參考OpenCV官方文件看畸變模型原理會更清楚：https://docs.opencv.org/3.0-beta/modules/calib3d/doc/camera_calibration_and_3d_reconstruction.html#fisheye-initundistortrectifymap
簡要流程就是：
1.求內參矩陣的逆，由於攝像機座標系的三維點到二維影像平面，需要乘以旋轉矩陣R和內參矩陣K。那麼反向投影回去則是二維影像座標乘以K*R的逆矩陣。
2.將目標影像中的每一個畫素點座標(j,i)，乘以1中求出的逆矩陣iR，轉換到攝像機座標系（_x,_y,_w）,並歸一化得到z=1平面下的三維座標(x,y,1)；
3.求出平面模型下畫素點對應魚眼半球模型下的極座標(r, theta)。
4.利用魚眼畸變模型求出擁有畸變時畫素點對應的theta_d。
5.利用求出的theta_d值將三維座標點重投影到二維影像平面得到(u,v)，(u,v)即為目標影像對應的畸變影像中畫素點座標
6.使用cv::Remap（）函式，根據mapx,mapy取出對應座標位置的畫素值賦值給目標影像，一般採用雙線性插值法，得到畸變校正後的目標影像。
#include <opencv2\opencv.hpp>
void cv::fisheye::initUndistortRectifyMap( InputArray K, InputArray D, InputArray R, InputArray P,
const cv::Size& size, int m1type, OutputArray map1, OutputArray map2 )
{
CV_Assert( m1type == CV_16SC2 || m1type == CV_32F || m1type <=0 );
map1.create( size, m1type <= 0 ? CV_16SC2 : m1type );
map2.create( size, map1.type() == CV_16SC2 ? CV_16UC1 : CV_32F );
CV_Assert((K.depth() == CV_32F || K.depth() == CV_64F) && (D.depth() == CV_32F || D.depth() == CV_64F));
CV_Assert((P.empty() || P.depth() == CV_32F || P.depth() == CV_64F) && (R.empty() || R.depth() == CV_32F || R.depth() == CV_64F));
CV_Assert(K.size() == Size(3, 3) && (D.empty() || D.total() == 4));
CV_Assert(R.empty() || R.size() == Size(3, 3) || R.total() * R.channels() == 3);
CV_Assert(P.empty() || P.size() == Size(3, 3) || P.size() == Size(4, 3));
//從內參矩陣K中取出歸一化焦距fx,fy; cx,cy
cv::Vec2d f, c;
if (K.depth() == CV_32F)
{
Matx33f camMat = K.getMat();
f = Vec2f(camMat(0, 0), camMat(1, 1));
c = Vec2f(camMat(0, 2), camMat(1, 2));
}
else
{
Matx33d camMat = K.getMat();
f = Vec2d(camMat(0, 0), camMat(1, 1));
c = Vec2d(camMat(0, 2), camMat(1, 2));
}
//從畸變係數矩陣D中取出畸變係數k1,k2,k3,k4
Vec4d k = Vec4d::all(0);
if (!D.empty())
k = D.depth() == CV_32F ? (Vec4d)*D.getMat().ptr<Vec4f>(): *D.getMat().ptr<Vec4d>();
//旋轉矩陣RR轉換資料型別為CV_64F，如果不需要旋轉，則RR為單位陣
cv::Matx33d RR = cv::Matx33d::eye();
if (!R.empty() && R.total() * R.channels() == 3)
{
cv::Vec3d rvec;
R.getMat().convertTo(rvec, CV_64F);
RR = Affine3d(rvec).rotation();
}
else if (!R.empty() && R.size() == Size(3, 3))
R.getMat().convertTo(RR, CV_64F);
//新的內參矩陣PP轉換資料型別為CV_64F
cv::Matx33d PP = cv::Matx33d::eye();
if (!P.empty())
P.getMat().colRange(0, 3).convertTo(PP, CV_64F);
//關鍵一步：新的內參矩陣*旋轉矩陣，然後利用SVD分解求出逆矩陣iR，後面用到
cv::Matx33d iR = (PP * RR).inv(cv::DECOMP_SVD);
//反向對映，遍歷目標影像所有畫素位置，找到畸變影像中對應位置座標(u,v)，並分別儲存座標(u,v)到mapx和mapy中
for( int i = 0; i < size.height; ++i)
{
float* m1f = map1.getMat().ptr<float>(i);
float* m2f = map2.getMat().ptr<float>(i);
short* m1 = (short*)m1f;
ushort* m2 = (ushort*)m2f;
//二維影像平面座標系->攝像機座標系
double _x = i*iR(0, 1) + iR(0, 2),
_y = i*iR(1, 1) + iR(1, 2),
_w = i*iR(2, 1) + iR(2, 2);
for( int j = 0; j < size.width; ++j)
{
//歸一化攝像機座標系，相當於假定在Z=1平面上
double x = _x/_w, y = _y/_w;
//求魚眼半球體截面半徑r
double r = sqrt(x*x + y*y);
//求魚眼半球面上一點與光心的連線和光軸的夾角Theta
double theta = atan(r);
//畸變模型求出theta_d，相當於有畸變的角度值
double theta2 = theta*theta, theta4 = theta2*theta2, theta6 = theta4*theta2, theta8 = theta4*theta4;
double theta_d = theta * (1 + k[0]*theta2 + k[1]*theta4 + k[2]*theta6 + k[3]*theta8);
//利用有畸變的Theta值，將攝像機座標系下的歸一化三維座標，重投影到二維影像平面，得到(j,i)對應畸變影像中的(u,v)
double scale = (r == 0) ? 1.0 : theta_d / r;
double u = f[0]*x*scale + c[0];
double v = f[1]*y*scale + c[1];
//儲存(u,v)座標到mapx,mapy
if( m1type == CV_16SC2 )
{
int iu = cv::saturate_cast<int>(u*cv::INTER_TAB_SIZE);
int iv = cv::saturate_cast<int>(v*cv::INTER_TAB_SIZE);
m1[j*2+0] = (short)(iu >> cv::INTER_BITS);
m1[j*2+1] = (short)(iv >> cv::INTER_BITS);
m2[j] = (ushort)((iv & (cv::INTER_TAB_SIZE-1))*cv::INTER_TAB_SIZE + (iu & (cv::INTER_TAB_SIZE-1)));
}
else if( m1type == CV_32FC1 )
{
m1f[j] = (float)u;
m2f[j] = (float)v;
}
//這三條語句是上面 ”//二維影像平面座標系->攝像機座標系“的一部分，是矩陣iR的第一列，這樣寫能夠簡化計算
_x += iR(0, 0);
_y += iR(1, 0);
_w += iR(2, 0);
}
}
}
opencv-python 實現魚眼矯正 棋盤矯正法
魚眼矯正有很多的方法， 比較常用的有：
1、棋盤標定法
2、經緯度法
opencv自帶魚眼矯正演算法， 也就是第一種， 棋盤矯正法。
第一步：製作棋盤格
用A4紙列印一張棋盤格， 固定到硬紙板上， 然後用魚眼鏡頭對著拍攝。 保留拍到的圖片， 如下圖所示：
可以從不同角度拍攝， 多儲存一些。
第二步： 計算內參和矯正係數
棋盤標定法， 必須要先計算出魚眼的內參和矯正係數， 可直接呼叫以下函式計算
import cv2
assert cv2.__version__[0] == '3'
import numpy as np
import os
import glob
def get_K_and_D(checkerboard, imgsPath):
CHECKERBOARD = checkerboard
subpix_criteria = (cv2.TERM_CRITERIA_EPS+cv2.TERM_CRITERIA_MAX_ITER, 30, 0.1)
calibration_flags = cv2.fisheye.CALIB_RECOMPUTE_EXTRINSIC+cv2.fisheye.CALIB_CHECK_COND+cv2.fisheye.CALIB_FIX_SKEW
objp = np.zeros((1, CHECKERBOARD[0]*CHECKERBOARD[1], 3), np.float32)
objp[0,:,:2] = np.mgrid[0:CHECKERBOARD[0], 0:CHECKERBOARD[1]].T.reshape(-1, 2)
_img_shape = None
objpoints = []
imgpoints = []
images = glob.glob(imgsPath + '/*.png')
for fname in images:
img = cv2.imread(fname)
if _img_shape == None:
_img_shape = img.shape[:2]
else:
assert _img_shape == img.shape[:2], "All images must share the same size."
gray = cv2.cvtColor(img,cv2.COLOR_BGR2GRAY)
ret, corners = cv2.findChessboardCorners(gray, CHECKERBOARD,cv2.CALIB_CB_ADAPTIVE_THRESH+cv2.CALIB_CB_FAST_CHECK+cv2.CALIB_CB_NORMALIZE_IMAGE)
if ret == True:
objpoints.append(objp)
cv2.cornerSubPix(gray,corners,(3,3),(-1,-1),subpix_criteria)
imgpoints.append(corners)
N_OK = len(objpoints)
K = np.zeros((3, 3))
D = np.zeros((4, 1))
rvecs = [np.zeros((1, 1, 3), dtype=np.float64) for i in range(N_OK)]
tvecs = [np.zeros((1, 1, 3), dtype=np.float64) for i in range(N_OK)]
rms, _, _, _, _ = \
cv2.fisheye.calibrate(
objpoints,
imgpoints,
gray.shape[::-1],
K,
D,
rvecs,
tvecs,
calibration_flags,
(cv2.TERM_CRITERIA_EPS+cv2.TERM_CRITERIA_MAX_ITER, 30, 1e-6)
)
DIM = _img_shape[::-1]
print("Found " + str(N_OK) + " valid images for calibration")
print("DIM=" + str(_img_shape[::-1]))
print("K=np.array(" + str(K.tolist()) + ")")
print("D=np.array(" + str(D.tolist()) + ")")
return DIM, K, D
# 計算內參和矯正係數
'''
# checkerboard： 棋盤格的格點數目
# imgsPath: 存放魚眼圖片的路徑
'''
get_K_and_D((6,9), 'fisheyeImage')
第三步： 根據計算的K和D， 矯正魚眼圖
儲存第二步計算的DIM， K， D， 利用一下程式， 矯正魚眼圖
def undistort(img_path):
img = cv2.imread(img_path)
img = cv2.resize(img, DIM)
map1, map2 = cv2.fisheye.initUndistortRectifyMap(K, D, np.eye(3), K, DIM,cv2.CV_16SC2)
undistorted_img = cv2.remap(img, map1, map2, interpolation=cv2.INTER_LINEAR,borderMode=cv2.BORDER_CONSTANT)
cv2.imwrite('unfisheyeImage.png', undistorted_img)
DIM， K， D是固定不變的，因此，map1和map2也是不變的， 當你有大量的資料需要矯正時， 應當避免map1和map2的重複計算， 只需要計算remap即可。
矯正效果
後續
矯正結束， 對於上面的矯正結果， 應該適用於很多場景了， 但是矯正之後， 很多有效區域被截掉了， 能否矯正出更大的有效面積呢？ 答案當然是有的， 具體的細節， 去github上看專案裡面的pdf吧， 這裡不再介紹了。
https://github.com/HLearning/fisheye
相關文章
- 魚眼相機成像模型和畸變模型模型
- 影像畸變校正解決方案，無畸變的高質量視覺體驗視覺
- 手機中的計算攝影：超廣角畸變校正
- 影像工具(1)照片去畸變
- 基於yolov2深度學習網路模型的魚眼鏡頭中人員檢測演算法matlab模擬YOLO深度學習模型演算法Matlab
- 熱像儀長焦鏡頭、廣角鏡頭和微距鏡頭的應用
- 定價讓玩家大跌眼鏡，RX 5500 XT真“智商檢測卡”
- 鏡頭品牌加入，鐳射電視迎來質變？
- extern 用法，全域性變數與標頭檔案（重複定義）變數
- Unity 鏡頭控制 第三人稱自由鏡頭Unity
- 巨集_變數_函式_指標_標頭檔案變數函式指標
- 相機camera針孔模型與畸變矯正模型
- Timeline 鏡頭技巧（一）
- 8.10 標頭檔案剖析：變數宣告變數
- PTGUI全景合成軟體使用教程之鏡頭引數設定GUI
- 鏡頭光暈外掛Optical Flares ae鏡頭光暈耀斑外掛
- pycharm設定python標頭檔案模版PyCharmPython
- 眼鏡 進銷存 ERP設計（skycto JEEditor）
- OpenCV開發筆記（八十一）：透過棋盤格使用魚眼方式標定相機內參矩陣矯正攝像頭影像OpenCV筆記矩陣
- Unity製作特寫鏡頭Unity
- maven中properties標籤定義變數Maven變數
- 為什麼不在標頭檔案做定義
- 艾瑞諮詢：2021年中國眼鏡鏡片行業白皮書（附下載）行業
- 艾瑞諮詢：2020年中國眼鏡鏡片行業白皮書（附下載）行業
- 如何通過鏡頭運用提升遊戲體驗？ 《死亡擱淺》鏡頭藝術賞析遊戲
- 彙總 | AR眼鏡開發資源合集（建議收藏）
- 拓端：資料下的隱形眼鏡初體驗
- 易優CMS模板標籤assign定義變數模板檔案中定義變數，可在其他標籤裡使用該變數變數
- mac CLion cmake 呼叫自己定義的標頭檔案Mac
- OpenCV開發筆記（七十七）：相機標定（二）：透過棋盤標定計算相機內參矩陣矯正畸變攝像頭影像OpenCV筆記矩陣
- ODG釋出兩款AR眼鏡，將在CES 2017上展示
- 這款新筆記本沒有螢幕，只有VR眼鏡筆記VR
- C語言中的標頭檔案中的巨集定義C語言
- win10設定護眼模式的方法_win10如何設定護眼模式Win10模式
- Unity第一人稱鏡頭模板Unity
- AR眼鏡，機器學習領域的殺手級應用？機器學習
- ODG釋出兩款AR眼鏡，將在CES2017上展示
- IDEA-idea設定護眼模式Idea模式