1、1对于数码相机定位问题的分析摘要:本文针对双目定位的数码相机相对位置的标定问题进行了深入研究。对于问题一,建立了世界坐标系和相机成像坐标系之间的转换矩阵,可以根据相机与靶面任意的相对位置关系求得靶标圆心的像坐标;对于问题二,针对所给靶标的特征提出了质心模型、切线模型、解析法模型三种求靶标圆心在像坐标投影点的模型,并就所给的靶标相片求出了圆心投影的像坐标。对于问题三,应用了问题一中转换矩阵,建立了判定问题二中模型优劣的仿真算法,对问题二中提出的三种模型进行了优劣分析。对于问题四,应用 Roger Y. Tsai 的单部相机内部和外部参数的标定算法,不用求解多元非线性方程组,直接用最小二乘法求解超
2、定线性方程组,即可求得相机坐标系的变换矩阵和像距,从而能够确定两相机的相对方位和位置关系。本文得到的主要结论如下:问题一:求解靶标圆心像坐标的算法为世界坐标系和相机成像坐标系之间的坐标变换;问题二:利用三种模型得出问题所给靶标相片中靶标圆心投影位置的像坐标,如下表所示;质心模型 切线模型 解析法模型X Y X Y X YA 322.8948 189.4935 323.2771 189.9167 323.0474 189.6938B 422.9960 196.9423 423.2943 197.3401 423.0209 197.0611C 639.8994 213.1522 640.1433
3、213.3946 640.1127 213.2289D 582.7329 502.9824 583.0062 503.2098 582.8821 503.1650E 284.6657 501.7731 284.9940 502.0426 284.8034 502.0957问题三:在像平面和物平面夹角不是很大的情况下,三种模型的精度相差不大;像平面和物平面夹角比较大时,切线模型的精度大于解析法模型,质心模型的精度最差。问题四:运用相机内外参数的标定算法和问题二得到的靶标圆心和切点坐标,得题图三对应相机的相机坐标系的旋转矩阵与平移向量分别为, 。0.93765.190.386-24-517.42
4、-.3094.关键词:机器视觉、坐标变换、透视投影、相机标定2对于数码相机定位问题的分析一、 问题重述双目视觉定位系统是机器视觉学科的主要研究内容,系统通过处理放在两个不同位置的摄像机捕捉到的图像,对空间中的某点或某物体进行定位。这种技术在交通、医疗、工业、军事等领域都有广泛的应用。要使定位系统进行准确的定位,必须知道两台相机精确的相对位置关系,可以使用如下方法对两台相机进行标定以获得相机的相对位置:在一块平板上画若干个点, 同时用这两部相机照相,分别得到这些点在它们像平面上的像点,利用这两组像点的几何关系就可以得到这两部相机的相对位置。然而,无论在物平面或像平面上我们都无法直接得到没有几何尺
5、寸的“点” 。实际的做法是在物平面上画若干个圆(称为靶标) ,它们的圆心就是几何的点了。而它们的像一般会变形,所以必须从靶标上的这些圆的像中把圆心的像精确地找到,标定就可实现。本课题要求完成下列任务:(1) 建立数学模型和算法以确定靶标上圆的圆心在该相机像平面的像坐标, 这里坐标系原点取在该相机的焦点,x-y 平面平行于像平面;(2) 由靶标及其像,计算靶标上圆的圆心在像平面上的像坐标(3) 设计一种方法检验你们的模型,并对方法的精度和稳定性进行讨论;(4) 建立用此靶标给出两部固定相机相对位置的数学模型和方法。二、 问题分析该问题实际上是求解不同坐标系中的点的对应关系的问题。对于第一个问题,
6、首先需要解决的是空间圆成像的映射问题。由于实际情况中维度的差距,可以考虑将成像的模型设置为小孔成像的模型,从而省略了成像的畸变对定位造成的影响对问题进行简化,对于由此造成的误差在后面的模型分析中进行讨论。空间圆小孔成像中有一些基本的原则和理论需要进行分析,然后以此进行建模,考虑空间中各个参考系之间的坐标变换关系,通过旋转、平移等方式达到点的映射的目的。对于问题二,即考虑上述过程的反向过程:已知像的信息和原象的信息,要求得出原象的圆心在像上的位置。一个朴素的想法就是考虑质心的对应关系,如果能找到像和原象的质心的对应关系,就求出了圆心在像平面上的坐标。但是此种做法缺乏理论依据,只是根据生活经验得出
7、,但不失为一种参考的方法,3可以将其与其他方法进行对比,选出更优的方法。再者,还可以通过切线的垂线的交点来确定圆心。首先对给出的像的信息进行处理,拟合出像中的曲线方程,接着在物的由物平面中的切点求出它们在相平面中的坐标(用问题一的方法) ,之后就可以由切点求出切线的方程从而得出圆心的坐标。另一个想法是解析的方法。通过空间解析几何和平面上线、角的关系计算出圆心相对于空间坐标系的坐标,并由第一问中各坐标系的坐标变换关系就可以得出圆心在像平面上的坐标。以上方法都在理论推导的过程中进行了一些近似从而会导致系统误差,所以在由上述一些方法对问题进行分析过后还需要对它们的精确程度和稳定性进行分析。此即问题三
8、的要求。而对于不同的方法,应采用不同的误差分析方法以得出他们的精度和稳定性。对于质心的方法,由于其缺乏理论根据,所以只能将其结果与其他的方法的结果进行比较以考虑其精度;对于切线的方法,由于计算机拟合的过程很精确(点足够多的情况下) ,且在之后的处理过程中没有进行近似处理,其精度的讨论应该也着重与其他方法的比较;对于解析的方法,由于其通过公式的推导,经过适当的近似得到结果,因此应着重讨论其近似时产生的误差。而对于稳定性的讨论,应朝着考虑靶标的移动、标定形状、噪声对方法精度的影响。最后一问要求通过靶标的世界坐标系坐标和靶标在两部相机的成像坐标系中的像点坐标来确定两部相机的相对位置。具体来说,即已知
9、靶标上若干物点的世界坐标和成像坐标,寻找标定两个相机坐标系相对世界坐标系的旋转矩阵 和平移向量 以及相机的内部参数的算法。由这四个矩阵能够自然地求得两相机的相对位置。三、 模型假设1本题中不考虑镜头畸变;2本题中相机模型为小孔模型;3各空间坐标轴均为右手系。四、 符号与术语说明透射投影:坐标变换中的旋转矩阵。R:坐标变换中的平移向量。T:像距,光心到像平面的距离。f世界坐标系: ,被拍摄物体所在坐标系。wXYZ 相机坐标系: ,以相机光轴为 ,光轴平行于相机像平面。cc cZ4成像坐标系: ,像平面坐标系。uv 五、 模型的建立与求解问题一该问要求靶标上圆的圆心在相机像平面上的坐标,即世界坐标
10、系() 、相机坐标系( )和成像坐标系( )三个坐标系之间wZYX, ccZYX, vu,的点坐标的转换。1 世界坐标系和相机坐标系之间的转换该转换用旋转矩阵 和平移向量 来描述,世界坐标和相机坐标间满足下RT式: ZYXRwcc其中 (具体确定见下文) , ,即世界坐标先旋转987654321rR zyxT矩阵,再将坐标原点平移 向量成为相机坐标。为方便表达,上式亦可表示T为,111wwTcc ZYXMORZYX 其中 ,而旋转矩阵 的确定有多种方法,如绕坐标轴的旋转变换,TO)0,(欧拉角旋转变换等,考虑到直观易用性,本文选取绕坐标轴的旋转变换来确定矩阵 , 分别为坐标系旋转变换时绕固定的
11、世界坐标系坐标轴的旋转R,夹角,此时 为,cossinsicosincosinsincosii 5图 1 坐标系的旋转又应用 z-x-y 旋转顺序,即先绕世界坐标系的 z 轴旋转 ,再绕其 x 轴旋转 ,最后绕其 y 轴旋转 ,但每次旋转都是以固定的世界坐标系坐标轴为转轴的。图 2 相机坐标系和成像坐标系2 相机坐标系和成像坐标系之间的转换对于该问题,用小孔模型代表相机工作原理,相机坐标系以小孔(光心)为坐标原点,以光轴为 轴正方向,而成像平面坐标系 则平行于相机坐标cZuov系的 平面,如图 2 所示。原理上说,成像平面应当在物点 P 和光心 O 连cXOY线的延长线上,完成小孔成像,但由于
12、成倒像,为处理和讨论方便,将成像平面放置在光心和物点之间,距光心 的位置上, 此处为像距。由小孔成像原ff6理和相似三角形知识,得到相机坐标系和成像坐标系之间的关系为:,ccZfYvZfXu 其中, 为 P 点的成像坐标, 为 P 点在相机坐标系中的空间坐标。),(yx ),(c矩阵形式表示为:。011cc cXufYZvZ综上所述,可得世界坐标系和成像坐标系之间的变换关系为:(1)0111wwc TXuf YRZv MOZ 在相机内部参数都确定的情况下,M 为一常矩阵, 为随着世界坐标系中c点坐标变化的量,为列向量 的第三项元素。因而得到了世界坐标1wTXYROZ系和成像坐标系之间的变换关系
13、。对于本小问,输入为世界坐标系中靶标的 5 个圆心位置以及数码相机的 6个位置参数(旋转矩阵 中的 3 个旋转角度和平移向量 的 3 个变量) ,要求输RT出靶标圆心在像平面上的坐标。对于靶标,定义如图 3 的世界坐标系,5 个圆心放置在 xoy 平面上,其坐标如表 1 所示;相机的位置参数则由外部设定 6 个参数 。用 Matlab 编写程序(代码见附录) ,可得到应用上述zyxT,坐标变换式所得圆心对应的像坐标的 和 。由于题中要求坐标系的原点定在光uv心,因而最终像坐标为 。),(fvu7图 3 世界坐标系表 1 世界坐标系中靶标圆心的坐标圆编号 圆心 x 坐标 圆心 y 坐标 圆心 z
14、 坐标A -50 50 0B -20 50 0C 50 50 0D 50 -50 0E -50 -50 0问题二问题二要求根据靶标图像和靶标的像,靶标上圆的圆心在像平面上的像坐标。1 图像的预处理由于题中靶标图像的照片给出的是灰度图像,并且实际摄像头拍摄的照片也是灰度图像,所以提取圆心工作以前,必须对图像进行预处理工作。图像预处理流程如图 4 所示,读入图像利用软件提取图形文件中像素的颜色信息,由于圆心投影只与边缘有关,所以用软件算法实现边缘的提取,本文利用了 Sobel边界探测器,由于原图像灰度及阈值设定的原因,完成边界探测后会有许多与边界不连续的点,我们称之为噪声,这里我们手工去除这些噪声
15、点,然后将边界点坐标存入数组。同时我们将图像做二值化处理,作为求质心的原始数据。图 5 显示了图像预处理的效果,A 为初始读入图像数据,是灰度图像, B为经过 Sobel 边界探测后的结果,有噪点存在,C 为手工去除噪声点后的图像,D 为原始图像做二值化处理后的图像。8读入图像噪声点去除边缘提取 读入数据A B C二值化处理D图 4 图像预处理流程A B C D图 5 图像预处理阶段示意2 圆空间透视投影的分析本小节介绍推导模型必须证实的几个定理定理 1:空间内圆的投射投影为椭圆证明:如图 6 所示以 为圆心的圆经过小孔投影到像平面上,以透镜的光心为坐rO标原点建立空间三位坐标系,Z 轴垂直于
16、像平面向下,O 到像平面的距离为 ,fO 到 圆所在的平面的距离为 。r l9图 6 空间圆的透视投影设圆 以 为轴旋转角度 ,再像平面形成一个封闭曲线,下面证明该封rOy闭曲线为椭圆。圆 的参数方程为:rcosinirrXYZl根据相似三角形比例原则,圆 上一点与象平面中该点的坐标 有如rO(,)xyz下关系 cossincosinzfxyzrrlr整理得: 2222sinco()1()isinrfxyl rfl此方程为一椭圆的标准方程,由此可知圆通过照相机在像平面上成的像为椭圆(在 时为圆) 。010定理 2:透射投影中,物圆的圆心的像不是所成椭圆像的中心。证明:图 7 XOZ 平面 AB
17、 剖面图过图 6 的 的面作圆锥的剖面如图 7 所示: xOz其中 , 为圆 在像平面上所成的像,过 作 的平行线 分ACBr CABm别交 的反向延长线为 ,交 于 。DBOE假设圆心 的像 为椭圆的中心,则有 A由相似三角形边的关系有: C将上式代入下式有 。分别从 、 向 作垂线垂足分别为 、 ,由于 为 的中点,所以有ABmFGCAB,由图易知: 。与上述矛盾,故假设不成立,即:物圆的圆FCGED心的像不是所成椭圆像的中心定理 3:空间圆圆心的投影在投影椭圆的短轴上证明:首先说明圆 的圆心的像、椭圆中心都在 。rOAB由于 是直径,易知圆 的圆心的像 上。ABr11图 8 平行于转轴的
18、剖面图如图 8 所示,由于圆 以 为轴旋转角度 ,而 轴过圆心,所以椭圆的rOyry中心在圆上的原象一定在垂直于 轴的直径上,在图 8 中的表述即为椭圆的中r心在 的中点,其原象在 上,所以 为过椭圆中心的一条轴。ABAB下面证明该轴 是椭圆的最短轴:设 、 的 坐标分别为 、X12X由相似三角形边的关系及 1 中坐标的关系有: 12cosAOrXB由上两式可得, 12 cos()AABrO在 上 为单调递减,又由于 过椭圆中心,所以转角 最大0,9cosB的直径所对应的像 即为椭圆的最短轴。即物圆 的圆心的像在其像椭圆的最短轴上。rO123 模型一 像素质心模型31 模型建立一般情况下,平面
19、上一薄片的质心计算公式如下,为薄片的范围,为薄片在(x,y) 的面密度, 为质心坐标。(,)xy(,)xy(,),(,)SSxyddxyAA本模型因为是基于像素的模型,像素是一个个离散的点,又由于图形做了二值化处理,这里 为二值量。定义:(,)xy1(,)0(,)fxy点 为 黑点 为 白则在像素意义下的中心坐标定义式如下所示: 2 21 1(,)(,), ,y yx xy yx xfdfdy其中(x1,y1), (x2,y2)为包含此椭圆而不包含其他椭圆的一个矩形的两个对角顶点的坐标。经过 2 的分析,虽然空间圆圆心的投射投影一般不在投影图形的几何中心,而且利用像素所作的离散计算存在着一定的
20、误差,但是作为一种近似的方法,也是可行的,后文将对这种模型造成的误差做定量的分析。32 模型求解按上述模型求得的圆心坐标如表 2 所示,表 2 切线模型求点坐标圆 X YA 322.8948 189.4935B 422.9960 196.9423C 639.8994 213.1522D 582.7329 502.9824E 284.6657 501.77314 模型二 切线模型13在上述的质心模型中,仅利用了一个椭圆自己的性质就给出了该椭圆内圆心的投影,而本问题给出了五个相等大小并且相对位置已知的圆,所以如果利用圆之间相互关系求解各圆圆心的投影坐标将更为精确。41 切点位置不变的性质空间的投射
21、投影有如下的性质,如图 9 所示,如果空间存在圆 A,B ,这两个圆有公切线 ,切点 a,b,c,d。如果将这两个圆做某透视投影如图所示,12,l此时在做投影椭圆 A,B的公切线 ,则得到的四个切点 a,b,c,d分别是原12,l来 a,b,c,d 的投影。图 9 切点位置不变的性质示意图投射投影变换属于仿射变换,点在线上的位置关系不会改变,在本问题的靶标图形中,有四个相等大小的圆(A,C,D,E)正好在正方形的四个定点上,两两做这四个圆的外公切钱,如图,可得相对的切点的连线经过圆心,两条连线的交点就是圆的圆心。所以在投射的象的图形中,只需找到这些切线与切点,连接各个切点形成交线就可以求得圆心
22、。对于问题中的圆 B,由于不存在相互垂直的切线,所以只能够确定一条直径的位置,另外根据 2 的叙述,圆心的投影肯定在投影椭圆的短轴上,这样用已知直径的投影和椭圆的短轴相交,就可以求得 B 圆圆心投影的坐标。42 椭圆方程的拟合14本模型要求作椭圆的曲线,由于所给图像分辨率不高,仅为 1024*768 像素,所以如果直接用图像中图形边界做切线,精度将会变得非常低,会造成很大的误差,所以在本模型中,先要利用所给图像中图形的边界(在 1 中提取)拟合出椭圆的方程。平面中,二次曲线的一般方程是: 22(,) 0Faxbycdxeyf其中 T21bcdefxya将空间中任意一点 带入 F, 被称为空间点
23、 到圆锥曲线x1(,)ax1x的“代数距离” 。拟合一条圆锥曲线可以看作是求所有被拟合点代数(,)0Fa距离平方和的最小值。如下式所示: 21,()()NiFiax寻 找 使 Dist取 最 小 值设用 N 个点拟合一个椭圆, T12nx又根据文献2中描述如果限制 可行解的 向量肯定是椭圆的参241acba数向量。由此可以把椭圆拟合的问题转换为一个带二次限制条件的最小二乘问题,如下: min().1EstTDaC其中02000图中五个椭圆方程拟合的结果如表 3 所示。15表 3 椭圆方程参数22(,) 0Faxbycdxeyf圆 a b c d e fA 0.4908368 0.0456273
24、2 0.5103946 -325.6590 -208.1902 71450.98B 0.5005452 0.07111793 0.5019815 -437.4560 -227.8848 114165.5C 0.5252867 0.1323467 0.4842668 -700.4724 -291.0644 254440.3D 0.4789880 0.1802621 0.5388937 -648.9240 -647.1483 351269.2E 0.4521153 0.1041725 0.5589569 -309.6662 -590.5596 191540.0图 10 为实际过程中椭圆拟合的效果图
25、,从图中可以看出,方程对于椭圆做了很好的拟合,图 10 椭圆拟合效果图表 4 显示了拟合过程中代数距离平方和的最小值,可看出大约在 数量级,410由于每个椭圆用大约 200 个点进行拟合,所以 应该在 20 左右,x,y 坐标()Fix平均的误差应该在几个像素以内。表 4 代数距离平方和 21min()NFix圆 A B C D E平方和 52242.57 45703.23 39525.20 30956.11 40892.521643 公切线的求取图 11 靶标的公切线图设两椭圆方程分别为 22111: 0Aaxbycdxeyf222BA 椭圆在 处的切线斜率为:1(,)xy11(,)2xya
26、bddceB 椭圆在 处的切线斜率为:,)222(,)xyaxbyddce如果 , 连线是切线,那么12(,)121(,)(,)xyxydk联立上述五式可求得 , 和 的解,其中,关于 ,1(,)2(,)xyk1(,)xy是四次方程,既可以求得四组解,这和切线的物理意义(两组内切,两2(,)xy组外切)是一致的。实际求解过程当中,由于计算规模过大,不宜求解,所以17这里我们采用了将最后一个条件改为 的优化条件,然后限定21minykx, 的范围,用 LINGO 软件解优化(程序见附录)问题,经试验,效1(,)xy2(,)果良好,不仅求解速度很快,而且每次优化目标都等于 0,说明求得了正确解。4
27、4 模型求解按上述模型求得的圆心坐标如表 5 所示,表 5 切线模型求点坐标圆 X yA 323.2771 189.9167B 423.2943 197.3401C 640.1433 213.3946D 583.0062 503.2098E 284.9940 502.04265 模型三 解析法模型51 模型建立在图 12 所示在照相机 坐标系中圆的半径是 ,设该圆绕 y 轴旋转Oxyzr角,图中坐标系原点 设在焦距中心,圆心 ,距坐标系原点的距离为 , rOl摄像机焦距为 , 圆的参数方程为f cosinirrXYZl根据相似三角形比例原则,圆上一点与图像 中该点的坐标 有如下(,)xyz关系
28、 cossincosinzfxyzrrlr18整理化成椭圆的标准方程(在像平面上)为: 2222sinco()1()isinrfxyl rfl图 12 空间圆的透视投影图过图 12 中 的面作圆锥的剖面如图 7 所示: , (设xOz OAB2r圆的半径为 ) ,易知 AB 在像平面上的像 为所形成椭圆的短轴。而过圆心r AB且与 AB 垂直的半径 CD 在像平面上的像 为一条平行于圆 在像平面上的CD像(椭圆)的长轴的直线,设其长度为 ,设其像椭圆的长轴长为 ,短轴长aa为 。如图 8 所示,设圆心到原点 的距离为 ,圆心到像平面的距离为 。 则b l f由相似三角形的关系有:(2)rlaf
29、又由小孔成像的原理,对应的直线在成像后大小变化的比例是相等的,反映在上左图中即:(3)cosAEbBa由上面的椭圆标准方程可知圆 所成的像的中心并不是圆 的圆心在像OO平面上的对应点,圆心的坐标 为(,)eXYZ1922sinco0eerfXlYZf与上面(*)(*)式联立可以得出如下关系:(4)2()ebaXaf为了确定圆心的位置 ,由于 的纸无法直接得出,所以取椭圆的半(,)eYZ长轴 代替,可以得出圆心的近似位置。a52 参数求解对于式(*) ,像距 f 已知,需要解得椭圆的半长轴 和半短轴 ,对于任ab意椭圆的一般式 ,总可以将其化成如下形式:220axbycdxeyf2()()()0
30、uvcyf根据系数相等可得: 2dabecvfuf解得 22dcbeavfucvb由此可得椭圆的圆心为 ,(,)uv20椭圆的半长轴为 ,半短轴为fafb椭圆旋转角度为 2rctn()/b53 模型求解通过 5.2 求得拟合的各椭圆的长轴和短轴,转角,中心,然后通过 5.1 中的(2)-(4)式求得空间圆圆心在象平面投影的坐标。表 6 空间圆圆心在象平面投影的坐标圆 X YA 323.0474 189.6938B 423.0209 197.0611C 640.1127 213.2289D 582.8821 502.0957E 284.8034 502.0426问题三1模型的检验方法该方法的核心
31、是基于问题一的结论,在问题一中,已经获得靶标上的点和像平面上像坐标的对应关系,因此可以按以下步骤完成基于仿真的检验:步骤 1 固定靶标,确定相机相对于靶标的姿态,包括旋转矩阵和位移矢量。步骤 2 通过问题一结论中的变换公式,求出靶标上所有点映射到相机像平面的坐标,这里必须使用离散化算法。 (离散化算法解释)这一步实际上就仿真求出了在此相机姿态下,靶标在相机像平面所成的像,同时,可以准确的知道靶标上圆的圆心在相机像平面投影的坐标。步骤 3 用问题二中的多种方法对步骤 2 中求出的靶标的像求圆心步骤 4 用所求得圆心和映射在相机像平面的准确的圆心的像做比较,两者相差越小说明方法越准确。212检验结
32、论(a) (b)图 13 成像检验图表 7 图 13(a)检验结果标准圆心 质心模型 切线模型 解析法模型X y X Y x y x y519.9263 309.5831 519.9351 310.3686 519.9333 310.4732 520.1446 309.9762589.1007 311.0271 589.0412 311.8446 589.1136 312.0014 588.9056 311.4932747.6679 314.3308 747.5765 315.1844 747.5763 315.4256 747.7089 315.1042722.0924 518.4924 5
33、05.9950 518.2396 723.0138 518.3788 721.8761 518.3584506.1380 517.4491 722.1016 519.3891 504.1380 517.7691 503.8780 517.4360误差平方和 3.475416 7.913014 6.229838表 8 图 13(b)检验结果 标准圆心 质心模型 切线模型 解析法模型X y x y x y x y484.6105 361.5820 486.8028 361.5050 485.7325 361.5549 485.6051 361.5320616.7665 348.9821 616.7
34、102 349.1967 616.7449 349.0115 616.7233 349.1549914.1542 320.6321 914.3273 320.8787 914.2246 320.6122 914.2781 320.6257813.8254 506.6421 813.8603 507.2887 813.8434 506.6281 813.8398 506.6377466.4812 544.9864 466.3966 545.6107 466.4572 545.1974 466.4355 545.4902误差平方和 5.768324 1.311919 1.294977表 7-8 显
35、示了用上述方法检验问题二中方法的结果,生成的虚拟照片如图13 所示。其中图 13(a)为像平面与物平面近似平行的情况,从结果中可以看出,在这种情况下,利用质心模型求解投影模型的误差是十分小的,在某些情况下,质心模型求解的误差小于用理论上更精确的切线模型和解析法模型。经过分析,我们认为,这是因为由于有限分辨率造成的量化误差其实可以等效为均值在质心的二维高斯噪声,而求质心相当于对噪声作了平均,一次可以抵消一部分噪声的影响,而切线模型和解析法模型是基于二次曲线最小二乘拟合的,这不能消除高斯噪声的影响。但图 13(b)中,像平面与物平面有较大的角度,这时,质22心模型有较大的误差,而切线模型和解析法模
36、型的误差较小。问题四此小问即要求通过靶标的世界坐标系坐标和靶标在两部相机的成像坐标系中的像点坐标来确定两部相机的相对位置。具体来说,即已知靶标上若干物点的世界坐标和成像坐标,标定两个相机坐标系相对世界坐标系的旋转矩阵和平移向量 以及相机的内部参数。由这四个矩阵能够自然地求得两21,R21,T相机的相对位置。由问题一可以看到,本小问是问题一的反向求解,即已知两种坐标系下点的坐标,求转移矩阵 的过程。但直接通过式(1)来求解 是非常困难的。原M M因有二:第一, 是 3*4 的矩阵,有 12 个未知数, 需要解 12 个方程,即使借助计算机软件也很难完成任务,Matlab 的 Solve 函数求解
37、能力不足以解决 12个方程,Lingo 也应变量过多无法在有限时间内给出结果;第二,通过消去 ,cZ每个靶标上的物点及其对应成像点可提供 2 个方程,因而共需 6 个点才能解出未知数,但本题每个靶标上只有 5 个圆,直接提供 5 个圆心,提供信息不足。因此,需要更直接的标定相机位置的方法,在参考了文献14后,我们发现了不用解多元或非线形方程的标定方法。该算法依然基于问题一的式(1),同时充分利用了靶标的自由设定来减少未知量。现描述标定一部相机的算法。 、1标定单部相机内部和外部参数的算法步骤 1:令靶标上 个物点 的世界坐标为 ,物点在成像坐n)5(wiizyx,标上的坐标为 , ,矩阵 和
38、的定义同问题一,则ivu,.21RT121415 1,2.yiwiiiwiwiyxyrvxyvuxunTr,23当 时即可用最小二乘法由上式解得 ;5n121415yyxyTrrT步骤 2:令 12124545/yyrrCT当 矩阵没有任何一行或一列为全 0 时(有一行或一列为全 0 的情况见文献C1) , 221/2154()rrySrT其中, 。由此可得 的大小;22145yT步骤 3:假设 的符号为正,在 个物点中选择一个成像坐标远离坐标原点yn的物点 ,计算t 12451245 ()yyxyxywttxttyrTrTrT 如果 与 同号并且 与 同号,则 的符号为正,否则为负;xtuy
39、tvy步骤 4:由于旋转矩阵 为单位正交阵,因此由 四个量即可确定出R5421,r, 21/3 21/64545r=(-r)sgn+(-r)而 可由 与 叉乘得到, 的符号的具体确定请请参987,r,321r,654r63,阅文献1;24步骤 5:构造矩阵 ,其行为 ,再构造向量A),(54iyiwii vTrxa,然后解线性方程 , ,即得相机像距 和iwii vyrxb)(87bAvzf f转移矩阵中的 。至此成功标定了相机的内部参数 和全部外部参数。算法zTMATLAB 代码见附录。2标定两台照相机相对位置的方法使用此算法可以得到两部相机各自的 和 ,即得到如下方程组:RT2211Xwc
40、消去 后得到:wX 21211 TRRXcc即得两相机坐标系之间的相对位置关系: 21T3针对本题的算法实例使用问题二得到的成像坐标系中五个圆的圆心坐标,再结合问题一中的世界坐标系的假设,使用上述算法,编制 Matlab 程序得到了拍摄题目中图三的相机的旋转矩阵 和平移向量 :1R1T0.9372.6035.96-82.4 -15.3809247.6T再用问题一的算法,使用世界坐标系中的靶标圆心坐标和矩阵 、 ,运R用 Matlab 反推出其在成像坐标系中的像(重塑值) ,及其与问题二中得到的圆心像坐标(理论值)对比如图 14。可见两者基本重合,1 中所述算法可认为有效。25图 14 靶标圆心
41、点的理论值和重塑值的像坐标对比尽管五个圆心可以标定出相机的位置,但上述算法的特点是参与标定的物点越多,得到的变换矩阵越精确。因而,我们将第二问中得到的四个圆的十六个切点的两套坐标也代入算法中,此时得到的 和 为:2RT0.93765.190.3862-24-517.4.7.2R -15.3072894.可以发现此时的变换矩阵和只用 5 个圆心得到的矩阵近似相等,但显然和 更为精确,因此将其作为题中拍摄出图 3 的相机的标定位置。图 15 是T圆心和切点的理论值和重塑值对比图。由图中两者的近似重合亦可验证算法的有效性。图 15 靶标圆心点和切点的理论值和重塑值的像坐标对比26六、 模型评价和模型
42、改进1关于问题一,完全按照坐标变换和几何关系进行推导,若不考虑相机镜头的畸变,则模型是精确的,因而改进方法为尽量研究抗量化噪声的算法。但对于问题一的坐标变换模型,世界坐标系到相机坐标系的变换为先旋转后平移,不如先平移后旋转直观。因为旋转后的平移相对于新坐标轴,这个平移量比基于原始坐标轴,即未旋转的坐标轴,要来得复杂和抽象。因此,先平移后旋转的坐标变换方法能使得用户或研究人员更容易地确定相机坐标系方位。2关于问题二提出的三种计算圆心投影的方法中,都遇到了相片分辨率有限存在量化误差的问题,这会造成圆心投影计算的误差。其中质心法和解析法均已做了近似处理,而切线法的精度仅取决于二次曲线拟合的精确与否。
43、一般情况下,精度切线模型 解析法模型 质心模型,运算复杂性也是切线模型 解析法模型 质心模型,所以在精度要求不高,并且物像平面基本平行的情况下,质心模型也不失为一种好的选择。3对于问题四的相机参数标定算法,参量符号的确定方法不够严谨而比较随意,特别是转移矩阵中的 和 的符号有两种组合方式,要根据后面算出的 的符3r6 f号来确定,不具有优秀算法所必须的确定性。参考文献:1Roger Y. Tsai,A Versatile Camera Calibration Technique for High-Accuracy 3D Machine Vision Metrology Using Off-th
44、e-shelf TV Cameras and Lenses,IEEE JOURNAL OF ROBOTICS AND AUTOMATION, VOL. RA-3, NO. 4, AUGUST 19872 Andrew Fitzgibbon, Maurizio Pilu, and Robert B. Fisher. Direct Least Square Fitting of Ellipses, PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. 21, NO. 5, MAY 19993刘红宁,基于双目立体视觉拟人机器人定位关键技术的研究,河北工业大学
45、硕士学位论文,p7-28,2006.124方超,双目 CCO 成像目标器识别算法研究,南京理工大学硕士学位论文,p4-28,2007.65Zhang, Z,A Flexible New Technique for Camera Calibration,Pattern Analysis and Machine Intelligence,Volume 22, Issue 11,p2-8,Nov. 20006韩庆瑶 万书亭,机器视觉中空间圆透视投影特性研究,华北电力大学学报,第 25 卷第1 期,p79-83,1998.127附录:1.问题一 Matlab 程序代码clcclear%ques1 us
46、ing coordinates trans.f=1577/3.78; %mmtrans1=f 0 0 0;0 f 0 0;0 0 1 0;a=0; %rotate angle inputb=pi/2;c=-pi/2;%rotation matrixR= cos(b)*cos(c)+sin(a)*sin(b)*sin(c), cos(a)*sin(c), -cos(c)*sin(b)+sin(a)*sin(c)*cos(b);-cos(b)*sin(c)+sin(a)*sin(b)*cos(c), cos(a)*cos(c), sin(c)*sin(b)+sin(a)*cos(c)*cos(b)
47、;cos(a)*sin(b), -sin(a), cos(a)*cos(b);t=0 0 -400; %displacementtrans2=R t;0 0 0 1;targ1=0 -50 50 1;targ2=0 -20 50 1;targ3=0 50 50 1;targ4=0 -50 -50 1;targ5=0 50 -50 1;targ=targ1 targ2 targ3 targ4 targ5;img=zeros(5,2); %image pointss coordi.for i=1:5mytemp=trans2*targ(:,i);zc=mytemp(3);w=trans1*trans2*targ(:,i)./zc;img(i,1)=w(1);img(i,2)=w(2);endplot(img(:,1),img(:,2),.);grid on;axis(-520 520 -520 520);% px=-50 -20 50 -50 50;% py=50 50 50 -50 -50;% plot(px,py,.);282LINGO 拟合椭圆方程的代码MODEL:SETS:POINT/R1R184/:X,Y;ENDSETSDATA:X=