[发明专利]一种基于正交消隐点的新型Tsai’s摄像机标定改进方法

摘要

本发明提出一种基于正交消隐点的新型Tsai’s摄像机标定改进方法，首先拍摄若干幅棋盘格标定板平面图像，并对所有图像进行预处理，对每一幅图像的角点进行亚像素坐标提取，然后将角点坐标划分成四组，应用最小二乘指标函数拟合每个方向上的平行直线束，计算每个方向上的最佳消隐点坐标，接着建立摄像机内部参数的约束方程组，求解内部参数f_x,f_y,u₀,v₀，将这些计算值作为初始值，代入外部参数的方程组中，求解所有的外部参数；最后引入一阶、二阶径向畸变系数，在实际图像平面坐标系上，建立图像残差当量最小化寻优函数，对所有的参数进行优化。本发明大大减少了计算量，速度快、精度高，可广泛应用于机器视觉研究、工业三维测量等多个领域的摄像机标定。

主权项

一种基于正交消隐点的新型Tsai’s摄像机标定改进方法，其特征在于：采用以下步骤：Step1：搭建光学测量系统，建立CCD摄像机模型；Step2：使摄像机光轴与棋盘格标定板平面所成夹角在0°～80°之间，拍摄若干幅不同位姿的棋盘格标定板平面图像：并对拍摄的图像进行中值滤波、锐化预处理；棋盘格标定板的图案规格为m×n，其中m＝n或m＝n‑1；Step3：应用Harris角点提取算法对预处理后的图像提取角点亚像素坐标，角点亚像素坐标为s_i,j＝[u_i,j,v_i,j]^T,i＝1,2,...m,j＝1,2,...n；Step4：将Step3中提取的角点划分成水平、垂直、主对角线、副对角线四组S＝{^gs|g＝hor,ver,mad,vid}；hor对应水平，ver对应垂直，mad对应主对角线，vid对应副对角线：4.1若m＝n，则执行步骤4.2→4.6；若m＝n‑1，则执行步骤4.7→4.11；4.2定义棋盘格标定板的左下角第一个角点为s_1,1，水平向右依次为s_1,1,s_1,2,...,s_1,n，垂直向上依次为s_1,1,s_2,1,...,s_n,1；4.3第一组：水平方向的角点坐标集合为^hors＝{^hors¹,^hors²,...,^horsⁱ,...,^horsⁿ},i＝1,2,...,n，其中子集合^horsⁱ＝{s_i,1,s_i,2,...s_i,j,...s_i,n},j＝1,2,...n；4.4第二组：垂直方向的角点坐标集合为^vers＝{^vers¹,^vers²,...,^vers^j,...,^versⁿ},j＝1,2,...,n，其中子集合^vers^j＝{s_1,j,s_2,j,...s_i,j,...s_n,j},i＝1,2,...n；4.5第三组：主对角线方向的角点坐标集合为^mads＝{^mads¹,^mads²,...,^mads^k,...,^mads²ⁿ‑³},k＝1,2,...,2n‑3；其中：1)当k＝1时，子集合^mads¹＝{s_2,1,s_1,2}，表示第2行第1列的角点和第1行第2列的角点组成的集合；2)当k＝2时，子集合^mads²＝{s_3,1,s_2,2,s_1,3}，表示第3行第1列、第2行第2列和第1行第3列的角点组成的集合；3)当k＝n‑1时，子集合^mads^n‑1＝{s_n,1,s_n‑1,2,...,s_n‑i,j,...,s_1,n},i＝0,1,...,n‑1；j＝1,2,...,n；4)当k＝n时，子集合^madsⁿ＝{s_n,2,s_n‑1,3,...,s_n‑i,j,...,s_2,n},i＝0,1,...,n‑2；j＝2,3,...,n；5)当k＝2n‑4时，子集合^mads^2n‑4＝{s_n,n‑2,s_n‑1,n‑1,s_n‑2,n}；6)当k＝2n‑3时，子集合^mads^2n‑3＝{s_n,n‑1,s_n‑1,n}；4.6第四组：副对角方向的角点坐标集合为^vids＝{^vids¹,^vids²,...,^vids^q,...,^vids^2n‑3},q＝1,2,...,2n‑3；其中：1)当q＝1时，子集合^vids¹＝{s_n‑1,1,s_n,2}，表示第n‑1行第1列的角点和第n行第2列的角点组成的集合；2)当q＝2时，子集合^vids²＝{s_n‑2,1,s_n‑1,2,s_n,3}，表示第n‑2行第1列、第n‑1行第2列和第n行第3列的角点组成的集合；3)当q＝n‑1时，子集合^vids^n‑1＝{s_1,1,s_2,2,...,s_i,j},i,j＝1,2,...,n；4)当q＝n时，子集合^vidsⁿ＝{s_1,2,s_2,3,...,s_i,j},i＝1,2,...,n‑1；j＝2,3,...,n；5)当q＝2n‑4时，子集合^vids^2n‑4＝{s_1,n‑2,s_2,n‑1,s_3,n}；6)当q＝2n‑3时，子集合^vids^2n‑3＝{s_1,n‑1,s_2,n}；4.7定义棋盘格标定板的左下角第一个角点为s_1,1，水平向右依次为s_1,1,s_1,2,...,s_1,n，垂直向上依次为s_1,1,s_2,1,...,s_m,1；4.8第一组：水平方向的角点坐标集合为^hors＝{^hors¹,^hors²,...,^horsⁱ,...,^hors^m},i＝1,2,...,m，其中子集合^horsⁱ＝{s_i,1,s_i,2,...s_i,j,...s_i,n},j＝1,2,...n；4.9第二组：垂直方向的角点坐标集合为^vers＝{^vers¹,^vers²,...,^vers^j,...,^versⁿ},j＝1,2,...,n，其中子集合^vers^j＝{s_1,j,s_2,j,...s_i,j,...s_m,j},i＝1,2,...m；4.10第三组：主对角线方向的角点坐标集合为^mads＝{^mads¹,^mads²,...,^mads^k,...,^mads^m+n‑3},k＝1,2,...,m+n‑3；其中：1)当k＝1时，子集合^mads¹＝{s_2,1,s_1,2}，表示第2行第1列的角点和第1行第2列的角点组成的集合；2)当k＝2时，子集合^mads²＝{s_3,1,s_2,2,s_1,3}，表示第3行第1列、第2行第2列和第1行第3列的角点组成的集合；3)当k＝m‑1时，子集合^mads^m‑1＝{s_m,1,s_m‑1,2,...,s_m‑i,j},i＝0,1,...,m‑1；j＝1,2,...,n‑1；4)当k＝m时，子集合^mads^m＝{s_m,2,s_m‑1,3,...,s_m‑i,j},i＝0,1,...,m‑1；j＝2,3,...,n；5)当k＝m+n‑4时，子集合^mads^m+n‑4＝{s_m,n‑2,s_m‑1,n‑1,s_m‑2,n}；6)当k＝m+n‑3时，子集合^mads^m+n‑3＝{s_m,n‑1,s_m‑1,n}；4.11第四组：副对角方向的角点坐标集合为^vids＝{^vids¹,^vids²,...,^vids^q,...,^vids^m+n‑3},q＝1,2,...,m+n‑3；其中：1)当q＝1时，子集合^vids¹＝{s_m‑1,1,s_m,2}，表示第m‑1行第1列的角点和第m行第2列的角点组成的集合；2)当q＝2时，子集合^vids²＝{s_m‑2,1,s_m‑1,2,s_m,3}，表示第m‑2行第1列、第m‑1行第2列和第m行第3列的角点组成的集合；3)当q＝n‑1时，子集合^vids^n‑1＝{s_1,1,s_2,2,...,s_i,j},i,j＝1,2,...,n；4)当q＝n时，子集合^vidsⁿ＝{s_1,2,s_2,3,...,s_i,j},i＝1,2,...,m；j＝2,3,...,n；5)当q＝m+n‑4时，子集合^vids^m+n‑4＝{s_1,n‑2,s_2,n‑1,s_3,n}；6)当q＝m+n‑3时，子集合^vids^m+n‑3＝{s_1,n‑1,s_2,n}；Step5：应用Step4中角点的划分结果，采用最小二乘法计算hor、ver、mad、vid方向上的平行直线束参数：若m＝n，则平行直线束参数为：^horlⁱ[^horaⁱ,‑1,^horbⁱ],i＝1,2,...,n^verl^j[^vera^j,‑1,^verb^j],j＝1,2,...,n^madl^k[^mada^k,‑1,^madb^k],k＝1,2,...,2n‑3^vidl^q[^vida^q,‑1,^vidb^q],q＝1,2,...,2n‑3若m＝n‑1，则平行直线束参数为：^horlⁱ(^horaⁱ,‑1,^horbⁱ),i＝1,2,...,m^verl^j(^vera^j,‑1,^verb^j),j＝1,2,...,n^madl^k(^mada^k,‑1,^madb^k),k＝1,2,...,m+n‑3^vidl^q(^vida^q,‑1,^vidb^q),q＝1,2,...,m+n‑3Step6：建立hor、ver、mad、vid方向的最佳消隐点目标函数：hor=argminlihorΣi=1mdis[V(uvahor,vvahor),l[aihor,-1,bihor]]=argminlihorΣi=1m|aihoruvahor-vvahor+bihor|(aihor)2+(-1)2;]]>ver=argminljverΣj=1ndis[V(uvaver,vvaver),l[ajver,-1,bjver]]=argminljverΣj=1n|ajveruvaver-vvaver+bjver|(ajver)2+(-1)2;]]>mad=argminlkmadΣk=1m+n-3dis[V(uvamad,vvamad),l[akmad,-1,bkmad]]=argminlkhorΣk=1m+n-3|akmaduvamad-vvamad+bkmad|(akmad)2+(-1)2;]]>vid=argminlqvidΣq=1m+n-3dis[V(uvavid,vvavid),l[aqvid,-1,bqvid]]=argminlqvidΣq=1m+n-3|aqviduvavid-vvavid+bqvid|(aqvid)2+(-1)2;]]>并应用Levenberg‑Marquard算法对最佳消隐点目标函数寻优计算出最佳的消隐点坐标^gV＝[^gu_va,^gv_va]^T；Step7：利用消隐线，验证消隐点计算的是否合理：7.1对Step6中4个消隐点进行随机组合，生成种组合，并拟合出条直线，计算这些直线两两间的夹角，选择出最大的夹角θ_max；7.2判断若θ≤θ_max，则这4个消隐点共线，形成了消隐线；反之，转到Step3，重新提取角点亚像素坐标，其中θ为阈值；Step8：建立摄像机内部参数几何约束方程组：horωV^ver=0V^madωV^vid=0,V^g=[uvag,vvag,1]T]]>其中：$\mathrm{\ω}=K^{-T}{K}^{-1}=\left[\begin{array}{ccc}\frac{1}{f_x^2}& 0& -\frac{u_0}{f_x^2}\\ 0& \frac{1}{f_y^2}& -\frac{v_0}{f_y^2}\\ -\frac{u_0}{f_x^2}& -\frac{v_0}{f_y^2}& \frac{u_0^2}{f_x^2}+\frac{v_0^2}{f_y^2}+1\end{array}\right]=\left[\begin{array}{ccc}{W}_{11}& 0& W_{13}\\ 0& W_{22}& W_{23}\\ W_{13}& W_{23}& W_{33}\end{array}\right]$f_x,f_y为CCD摄像机模型u轴和v轴上的归一化像素焦距，(u₀,v₀)是CCD摄像机模型成像平面的中心坐标；根据摄像机内部参数几何约束方程组，采用奇异值(SVD)分解和Zhang’s的方法计算该内部参数f_x,f_y,u₀,v₀；Step9：计算CCD摄像机模型中摄像机外部参数平移向量T的Z_C轴方向的分量t_z为：9.1建立摄像机外部参数的方程组：W+r2Yi,jW+txr7Xi,jW+r8Yi,jW+tz+u0vi,j=fyr4Xi,jW+r5Yi,jW+tyr7Xi,jW+r8Yi,jW+tz+v0,i=1,2,...,m;j=1,2,...,n]]>其中：摄像机外部参数的平移向量$T=\left[\begin{array}{c}{t}_x\\ t_y\\ t_z\end{array}\right]$摄像机外部参数的旋转矩阵$R=\left[\begin{array}{ccc}{r}_1& r_2& r_3\\ r_4& r_5& r_6\\ r_7& r_8& r_9\end{array}\right]$^WX_i,j和^WY_i,j为棋盘格标定板中角点的空间坐标；9.2将摄像机外部参数的方程组变形为矩阵形式为：WfxYi,jWfx000-UxXi,jW-UxYi,jW000fyXi,jWfyYi,jWfy-VxXi,jW-VxYi,jWb=UxVx]]>式中：向量b＝[b₁ b₂ b₃ b₄ b₅ b₆ b₇ b₈]^T；U_x＝u_i,j‑u₀；V_x＝v_i,j‑v₀；9.3应用最小二乘法计算9.2中矩阵形式，得到向量b＝[b₁ b₂ b₃ b₄ b₅ b₆ b₇ b₈]^T；并根据公式t_z＝(t_z1+t_z2+t_z3)/3$t_{z1}=\sqrt{({\mathrm{Sr}}_1-\sqrt{{{\mathrm{Sr}}_1}^2-{4\mathrm{Se}}_1})/{2\mathrm{Se}}_1},{\mathrm{Sr}}_1=b_1^2+b_2^2+b_4^2+b_5^2,{\mathrm{Se}}_1={(b_1{b}_5-b_2{b}_4)}^2$$t_{z2}=\sqrt{({\mathrm{Sr}}_2-\sqrt{{{\mathrm{Sr}}_2}^2-{4\mathrm{Se}}_2})/{2\mathrm{Se}}_2},{\mathrm{Sr}}_2=b_4^2+b_5^2+b_7^2+b_8^2,{\mathrm{Se}}_2={(b_5{b}_7-b_4{b}_8)}^2$$t_{z3}=\sqrt{({\mathrm{Sr}}_3-\sqrt{{{\mathrm{Sr}}_3}^2-{4\mathrm{Se}}_3})/{2\mathrm{Se}}_3},{\mathrm{Sr}}_3=b_1^2+b_2^2+b_7^2+b_8^2,{\mathrm{Se}}_3={(b_1{b}_8-b_2{b}_7)}^2$得到t_z；Step10：计算旋转矩阵R中的分量r₁,r₂,r₃,r₄,r₅,r₆,r₇,r₈,r₉和平移向量中的分量t_x,t_y；Step11：引入一、二阶径向畸变系数k₁,k₂，对内部参数f_x,f_y,u₀,v₀和外部参数R,T进行非线性全局优化：在摄像机的实际图像平面坐标系UV内建立图像残差当量最小化寻优函数：式中：：3D空间点^WS从世界坐标系X_WY_WZ_W投影到物理图像平面坐标系UV的畸变2D点实际图像平面坐标；：3D空间点^WS从像素图像平面坐标系uv逆向投影到物理图像平面坐标系UV的畸变2D点实际图像坐标。