利用单个球及平行圆切线性质标定拋物折反射摄像机的制作方法

本发明属于计算机视觉领域，涉及一种利用空间中一个球及平行圆切线的性质求解抛物折反射摄像机内参数的方法。

背景技术：

计算机视觉的中心任务就是对图像进行理解，而它的最终目标是使计算机具有通过二维图像认知三维环境信息的能力。这种能力将不仅使机器能感知包括形状、姿态、运动等在内的三维环境中物体的几何信息，而且能对它们进行描述、存储、识别与理解。摄像机标定就是确定从三维空间点到它的二维图像点之间的映射关系，它是许多计算机视觉应用必不可少的步骤。为了确定这一映射过程，需要建立摄像机的几何成像模型，几何模型的参数称为摄像机参数，摄像机参数可分为内参数和外参数两类。内参数描述成像系统的成像几何特性，外参数描述成像系统关于世界坐标系的方向和位置。摄像机标定可分为传统标定、自标定和基于几何实体的标定。无论哪种标定方法，都旨在建立二维图像与摄像机内参数之间的约束关系，特别是线性约束关系，这是目前摄像机标定所追求的目标，也是目前计算机视觉领域研究的热点之一。

抛物折反射摄像机由一个抛物镜面和一个正交摄像机组成，它的成像视野大，是全景视觉领域研究的热点之一。文献“catadioptricself-calibration”，(kangs.b.,proceedingsofieeeconferenceoncomputervisionandpatternrecognition,vol.1,pp.201-207,2000.)提出了一种折反射摄像机自标定方法，这类方法的优点是不需要使用标定块，缺点是必须获得图像之间的对应点。而在计算机视觉中，实现一个十分有效的寻找对应点的方法是很困难的。文献“geometricpropertiesofcentralcatadioptriclineimagesandtheirapplicationincalibration”,(barretoj.p.,araujoh.,ieeetransactionsonpatternanalysisandmachineintelligence,vol.27,no.8,pp.1327-1333,2005)研究了中心折反射摄像机下直线的像的几何性质，并将这些性质应用于中心折反射摄像机的标定。文献“anewlinearalgorithmforcalibratingcentralcatadioptriccameras”,(wuf.,duanf.,huz.etal.,patternrecognition,vol.41,no.10,pp.3166-3172,2008)介绍了对拓点和对拓像点的定义，导出了空间中的一个点在视球上的投影和它的折反射图像点之间的关系，使用这个关系建立了中心折反射摄像机内参数的线性约束，通过此线性约束即可获得中心折反射摄像机内参数。文献“calibrationofcentralcatadioptriccamerasusingadlt-likeapproach”,(puigl.,bastanlary.,sturmp.,etal.internationaljournalofcomputervision,vol.93,no.1,pp.101-114,2011)提出了一种基于三维控制点的标定方法，通过使用veronese映射对三维点和其图像点的坐标进行了扩展，在扩展坐标的基础上基于dlt(直接线性变换)——相似方法实现了中心折反射摄像机的标定，但是这类方法需要已知三维点的位置，并且容易从图像中提取其图像点。

球作为一种常见的几何体，其最重要的优点在于无自身遮挡，从任何一个方向看空间中一个球的封闭轮廓线总是一个圆，并且它的投影轮廓线可全部提取。由于球具有丰富的视觉几何特性，因此利用球进行摄像机标定已成为近年来的一个热点。文献“catadioptriccameracalibrationusinggeometricinvariants”，(yingx.,huz.,ieeetransactionsonpatternanalysisandmachineintelligence,vol.26,no.10,pp.1260-1271,2004)首次提出了利用球标定中心折反射摄像机。他们证明了球在中心折反射摄像机的单位球投影模型下的像为椭圆，并且在非退化情况下一个球的投影二次曲线提供两个不变量。为了降低求解的复杂度，他们提出了一种分步标定方法，该方法至少需要4个球的投影才能完成摄像机的标定。但是该文献提出的标定方法是非线性的，计算的复杂度较高，并且该标定方法只能标定抛物折反射摄像机的部分内参数。文献“identicalprojectivegeometricpropertiesofcentralcatadioptriclineimagesandsphereimageswithapplicationstocalibration”，(yingx.,zhah.,internationaljournalofcomputervision,vol.78,no.1,pp.89-105,2008)介绍了修正绝对二次曲线的像(miac)在中心折反射摄像机标定中的作用。他们通过研究球在中心折反射摄像机下的像与miac的几何与代数关系提出了两种线性标定算法。它们得出的结论对于对偶形式也是成立的。但是这篇文献中的理论和标定方法对于抛物折反射摄像机的情况是退化的。文献“acalibrationmethodforparacatadioptriccamerafromsphereimages”，(duanh.,wuy.,patternrecognitionletters,vol.33,no.6,pp.677-684,2012)基于圆环点理论提出了一种利用对拓球像标定抛物折反射摄像机的线性方法。但是这篇文献中关于圆环点的像的选取比较复杂。

技术实现要素：

本发明提供了一种制作简单，适用广泛，稳定性好的利用靶标求解抛物折反射摄像机内参数的方法，该靶标由空间中一个球构成。在求解抛物折反射摄像机内参数的过程中，需使用抛物折反射摄像机从不同的位置拍摄靶标的3幅图像线性求解出抛物折反射摄像机的5个内参数。

本发明采用如下技术方案：

用抛物折反射摄像机从不同的位置拍摄3幅含有一个球的图像。本发明是利用空间中一个球作为靶标用于求解抛物折反射摄像机内参数的方法，其特征在于仅利用球元素。首先，分别从3幅图像中提取镜面轮廓投影的边缘点和靶标图像边缘点，使用最小二乘法拟合获得镜面轮廓投影和球像的投影。其次，根据像点和其对拓像点的关系获得对拓像点，从而拟合出球像的对拓球像。在球像上取一个点，并求出对拓像点。由对拓像点的定义和平行圆的性质，上述一组对拓像点提供一个消失点。在抛物折反射系统下，球第一次投影到单位视球上形成两个平行小圆。异面的两个平行圆有两对共轭的虚交点，其中与圆所在平面的无穷远点共线的两个共轭虚交点为圆环点。在像平面上根据圆环点的像与消失点的共线性,获得圆环点的像，三幅图像提供三组圆环点的像。最后，利用圆环点的像对绝对二次曲线像的约束求解摄像机内参数。具体的步骤包括：拟合出镜面轮廓投影方程及靶标投影方程，估计球像的对拓球像，确定圆环点的像，求解抛物折反射摄像机内参数。

1.拟合镜面轮廓投影方程及靶标投影方程

利用matlab程序中的edge函数提取镜面轮廓投影边缘点和靶标图像边缘点的像素坐标，并用最小二乘法拟合获得镜面轮廓投影方程和球像的方程。

2.估计球像的对拓球像

空间中的球q，在拋物折反射摄像机的单位球模型下的投影分为两步。第一步，将球q投影为以o为中心的单位视球上的平行小圆sn+和sn-(n＝1,2,3表示拍摄的第n幅图像)，称sn-为sn+的对拓圆，并且sn-以单位视球中心o对称于sn+，单位视球直径的两个端点互为对拓。第二步，以单位视球表面上的一点oc为投影中心，这里oc可看作一个摄像机的光心，将平行小圆sn+和sn-分别投影为抛物折反射图像平面上的二次曲线cn+,cn-。其中称可见的二次曲线cn+为球q的像，不可见的二次曲线cn-为球像cn+的对拓球像。对应于单位视球上的sn-称为对拓圆，抛物折反射图像平面与直线oco垂直。对应于单位视球直径的一个端点，该端点在抛物折反射图像平面不可见称为对拓像点，则单位视球直径的这个端点称为对拓点。令以oc为光心的摄像机的内参数矩阵为其中rc是纵横比，fc是有效焦距，s是倾斜因子,[u0v01]^t是摄像机主点p的齐次坐标矩阵形式，rc，fc，s，u0，v0为抛物折反射摄像机的5个内参数。利用matlab中的函数提取第一幅图像中的镜面轮廓投影边缘点和3幅图像靶标图像边缘点的像素坐标，通过最小二乘法拟合得到相应的二次曲线方程。这里用c0表示第一幅图像镜面轮廓投影曲线的系数矩阵，cn+分别表示第n幅图像中的球像的系数矩阵。本文为了简化表述，用相同字母表示曲线和它的系数矩阵。通过c0可获得摄像机内参数矩阵kc的一个初始矩阵值从而得到绝对二次曲线的像ω的初始值这里:其中是纵横比的初始值，是有效焦距的初始值，是倾斜因子的初始值，是摄像机主点的初始齐次坐标矩阵表示，记取cn+上的一组点则与它相对应的一组对拓像点可由关系确定，用齐次坐标矩阵表示。根据对拓像点的定义，点在球像cn+的对拓球像cn-上，因此可用最小二乘法拟合得到对拓球像cn-的方程。

3.确定圆环点的像

在球q在单位视球上的投影小圆s1+上取一个点a1+，用a1-表示a1+关于单位视球中心o的对称点，即对拓点，则点a1-在s1+的关于单位视球中心o对称的圆，即对拓圆s1-上。l1+表示a1+关于圆s1+的切线，l1+表示a1-关于圆s1-的切线。根据对拓点的定义及平行圆的性质，l1+//l1-，于是l1+和l1-具有相同的无穷远点，这里用v1∞表示l1+,l1-上的无穷远点。小圆s1+和s1-为两个平行圆，它们相交于两对共轭的虚交点。因为平面上的每个圆都经过圆环点，所以圆与无穷远直线的交点为两个圆环点，圆环点是两个共轭虚点i(1,i,0)和j(1,-i,0)，即复点。因为平行圆所在的平行平面相交于相同的无穷远直线，因此平行的圆都交于圆环点，所以圆环点i,j为平行圆s1+和s1-共轭交点。空间圆与其所在平面的无穷远直线均交于圆环点，即圆环点与该平面无穷远点是共线的，故无穷远点v1∞与圆环点i,j是共线的。在像平面上，用c1+,c1-分别表示s1+,s1-的像；用a1+,a1-,分别表示a1+,a1-的像，则a1+,a1-为一对对拓像点。记过点a1+关于二次曲线c1+的切线为l1+，过点a1-关于二次曲线c1-的切线为l1-。则根据射影变换的性质，直线l1±为直线l1±的像，于是可通过直线l1+,l1-确定s1+所在平面上的一个消失点v1，即v1∞的像。二次曲线c1+,c1-有两组共轭的虚交点，即s1+和s1-的四个虚交点的像，其中一组必为圆环点的像m1i,m1j。透视变换保持点线的结合性，故在像平面上圆环点的像m1i,m1j与消失点v1也是共线的。根据共线性的判断，在四个虚交点中可获得圆环点的像。对于小圆s2+,s3+所在平面上的圆环点的像m2i,m2j},{m3i,m3j}可获得。下标i,j表示对应于圆环点i,j。

4.求解抛物折反射摄像机内参数

由圆环点的像mni,mnj(n＝1,2,3)对绝对二次曲线的像ω的线性约束获得ω，即：其中re,im分别表示复数的实部和虚部。最后，对ω进行cholesky分解，再求逆可获得内参数矩阵kc，即获得射线机的内参数。

本发明优点：

(1)该靶标制作简单，只需一个球。

(2)对该靶标的物理尺度没有要求，无需知道球心在世界坐标系下的坐标。

(3)该靶标的图像边界点几乎可以全部提取，这样可以提高曲线拟合的精确度，从而提高标定精度。

附图说明

图1是用于求解抛物折反射摄像机内参数的靶标在单位视球上的示意图。

图2是靶标在抛物折反射图像平面上的投影。

具体实施方式

本发明提供了一种利用靶标求解抛物折反射摄像机内参数的方法，靶标是由空间中的一个球构成，如图1。用此靶标完成抛物折反射摄像机内参数的求解需要经过以下步骤：从图像中提取镜面轮廓投影边缘点和靶标图像边缘点，使用最小二乘法拟合获得镜面轮廓投影和球的像。根据像点和其对拓像点的关系获得对拓像点，从而拟合出球像的对拓球像。在球像上取一个点，并求得对拓像点，由对拓像点的定义，上述一组拓像点提供一个消失点。球像与对拓球像有两组共轭的虚交点，其中一组是圆环点的像，根据圆环点的像与消失点的共线性获得圆环点的像。从三个不同的方位对球拍摄，得到三组圆环点的像，利用圆环点的像对绝对二次曲线的像的约束求解摄像机内参数。利用本发明中的方法对用于实验的拋物折反射摄像机进行标定，具体步骤如下：

1.拟合镜面轮廓投影方程及靶标投影方程

利用matlab程序中的edge函数提取镜面轮廓投影边缘点和靶标图像边缘点的像素坐标，并用最小二乘法拟合获得镜面轮廓投影方程和球像的方程。

2.估计球像的对拓球像

空间中的球q(如图1)，在拋物折反射摄像机的单位球模型下的投影分为两步。第一步，将球q投影为以o中心的单位视球上的平行小圆sn+,sn-，称sn-为sn+的对拓圆(n＝1,2,3)，(如图1所示，n＝1为例)。第二步，通过摄像机的光心oc将平行小圆sn+,sn-分别投影为抛物折反射图像平面上的二次曲线cn+,cn-，这里称可见的二次曲线cn+为球q的像，不可见的二次曲线cn-为球像cn+的对拓球像(如图1所示，n＝1为例)。利用matlab中的edge函数分别提取3幅图像中的靶标图像边缘点和第1幅图像的镜面轮廓投影的像素坐标，通过最小二乘法拟合得到相应的二次曲线方程，这里用c0表示第1幅图像镜面轮廓投影曲线的系数矩阵，cn+表示第n幅图像中的球像的系数矩阵。通过c0可获得摄像机内参数矩阵kc的一个初始矩阵值具体如公式(1)：

这里，c0(p,q)(p＝1,2；q＝1,2,3)表示矩阵c0的第p行第q列元素，φ为摄像机视场角的一半，ρ为抛物面镜轮廓投影椭圆c0的长半轴长。在获得的基础上可得到ω的初始值如公式(2)：

取cn+上的一组点则与它相对应的一组对拓像点可由关系式(3)确定：

根据对拓像点的定义，点在球像cn+的对拓球像cn-上，因此可用最小二乘法拟合得到对拓球像cn-的方程。

3.确定圆环点的像

在球q的投影小圆s1+上取一个点a1+，用a1-表示a1+的对拓点(如图1所示)，则点a1-在s1+的对拓圆s1-上。l1+为a1+关于圆s1+的切线，l1-为a1-关于圆s1-的切线，根据对拓点的定义及平行圆的性质，l1+//l1-，于是l1+和l1-具有相同的无穷远点，这里用v1∞表示l1+,l1-上的无穷远点。小圆s1+和s1-为两个平行圆，有两对共轭的虚交点。因为平行圆都交无穷远直线于相同的圆环点，所以圆环点i,j为平行圆s1+和s1-共轭交点。空间圆与其所在平面的无穷远直线均交于圆环点，即圆环点与该平面无穷远点是共线的，故无穷远点v1∞与圆环点i,j是共线的。如图2所示，在像平面上，用c1+,c1-分别表示s1+,s1-的像，用a1+,a1-分别表a1+,a1-的像，则a1+,a1-为一组对拓像点。记过点a1+关于二次曲线c1+的切线为l1+，过点a1-关于二次曲线c1-的切线为l1-。设a1±的齐次坐标矩阵为[ua1±va1±1]^t，直线l1±的齐次线坐标矩阵为[u1±v1±1]^t,则：

λ1+[u1+v1+1]^t＝c1+[ua1+va1+1]^t，(4)

λ1-[u1-v1-1]^t＝c1-[ua1-va1-1]^t，(5)

其中λ1±是非零常数因子，下标±表示+和-的简写，+表示可见，-表示不可见，即对拓含义。则根据射影变换的性质，直线l1±为直线l1±的像，于是可通过直线{l1+,l1-}确定小圆s1+所在平面上的一个消失点v1，通过联立l1+和l1-的方程可得：

其中[uv1]^t为抛物折反射图像平面上的像素齐次坐标矩阵。

设v1的齐次坐标矩阵为[uv1vv11]^t，则v1＝[uv1vv11]^t是方程组(6)的解。c1+,c1-有两对共轭的虚交点，用matlab中的solve函数获得虚交点，其中一组必为圆环点的像m1i,m1j。透视变换保持点线的结合性，故在像平面上圆环点的像m1i,m1j与消失点v1也是共线的。三点共线的充要条件是点的齐次坐标构成矩阵的行列式为零。故在两组共轭虚交点中满足等式(7)的点组为圆环点的像：

其中[xiyi1]^t(i＝1,2,3)表示空间点的齐次坐标矩阵。分别把两组共轭虚交点同v1代入(7)进行检验，满足条件的点即为圆环点的像。对于小圆s2+,s3+所在平面上的圆环点的像m2i,m2j和m3i,m3j用同样的方法获得。

4.求解抛物折反射摄像机内参数

由圆环点的像对绝对二次曲线的像的线性约束有：

其中re,im分别表示复数的实部和虚部。可用最小二乘法优化求解(8)获得ω，最后对进行cholesky分解得再求逆可获得内参数矩阵kc，即获得摄像机5个内参数。

实施例

本发明提出了一种利用空间的一个球作为靶标线性确定抛物折反射摄像机内参数的方法。本发明采用的实验模板如示意图1中所示的球。下面以一实例对本发明的实施方案做出更为详细的描述。

基于空间中球的抛物折反射摄像机标定采用的实验模板是空间中的一个球，如图1所示，球为q。利用本发明中的方法对用于实验的抛物折反射摄像机进行标定，具体步骤如下：

1.拟合图像边界及靶标曲线方程

本发明采用的图像大小为690×1198。用抛物折反射摄像机拍摄靶标的3幅实验图像，读入图像，利用matlab中的edge函数提取第1幅图像镜面轮廓投影边缘点和3幅图像靶标图像边缘点的像素坐标，并用最小二乘法拟合获得镜面轮廓投影方程和球像的方程。第1幅图像镜面轮廓投影方程的系数矩阵为c0，3幅球像的方程的系数矩阵分别为cn+(n＝1,2,3)，结果如下：

2.估计球像的对拓球像

将(9)代入(1)和(2)可得矩阵结果如下：

先分别在球像c1+、c2+和c3+上至少取5个互异的点，再将所取的点和(13)分别代入(3)得对拓像点的坐标，通过最小二乘法拟合获得对拓球像c1-、c2-和c3-的估计，系数矩阵结果如下：

3.确定圆环点的像

在球像c1+上任取一个点a1+，它的齐次坐标矩阵为：

a1+＝[2306.81595265189-1103.757884164641]^t,(17)

根据对拓像点的性质，可获得与点a1+对应的对拓像点a1-，结果如下：

a1-＝[12.170928237332754.5795250408911]^t，(18)

将(17)和(10)代入(4)可得过点a1+关于二次曲线c1+的切线l1+的齐次线坐标矩阵，结果如下：

l1+＝[-0.000208460481840.000470320567981]^t，(19)

将(18)和(14)代入(5)可得过点a1-关于二次曲线c1-的切线l1-的齐次线坐标矩阵，结果如下：

l1-＝[-0.00425493484596-0.001503983560751]^t；(20)

将(19)和(20)代入(6)可得消失点v1，结果如下：

v1＝[0.85294002718423-1.748160227140451]^t；(21)

用matlab中solve函数求解c1+和c1-的交点，把求得的两对共轭的虚交点和(21)代入(7)验证求得s1+所在平面上圆环点的像m1i,m1j，结果的齐次坐标如下：

m1i＝[-134.6434007080-638.3999342633i133.6434007080+638.3999342633i1]^t，(22)

m1j＝[-134.6434007080+638.3999342633i133.6434007080-638.3999342633i1]^t；(23)

对于小圆s2+所在平面上的圆环点的像m2i,m2j用类似方法代入(7)验证获得，结果的齐次坐标如下：

m2i＝[874.018799097-675.567729827i1174.610328832+0666.004896463i1]^t，(24)

m2j＝[874.018799097+675.567729827i1174.610328832-0666.004896463i1]^t；(25)

对于小圆s3+所在平面上的圆环点的像m3i,m3j可用类似方法代入(7)验证获得，结果的齐次坐标如下：

m3i＝[48.5579120270-358.0765980870i679.0263343998-528.7511368867i1]^t,(26)

m3j＝[48.5579120270+358.0765980870i679.0263343998+5.7511368867i1]^t。(27)

4.求解抛物折反射摄像机内参数

将(22)、(24)和(26)代入(8)得到ω中元素的线性方程组，使用最小二乘法求解该线性方程组得到ω的系数矩阵。结果如下：

最后，对(28)中的ω进行cholesky分解,再求逆可获得kc，有

其中纵横比rc＝kc(1,1)/kc(2,2)(kc(1,1)表示矩阵kc的第1行第1列的元素，kc(2,2)表示矩阵kc的第2行第2列的元素)，故抛物折反射摄像机的5个内参数分别为：rc＝1.00455475191844，fc＝573.261169054269，s＝0.519283357335，u0＝343.330665425737，v0＝625.176990789005。