机器视觉测量技术4

① 场景坐标：OXYZ

绝对坐标或世界坐标 通常情况下x′o′y′为oxyz的一部分，即x′o′y′与xoy重合。

在上述几何模型的基础上，进一步建立摄像机坐标系O1x1y1z1，其中O1y1轴线平行于CCD像面的水平像素行，O1x1平行于CCD像面的竖直像素列，O1z1垂直于CCD像面。如图所示，摄像机坐标系O1x1y1z1可以认为是由世界坐标系Oxyz通过绕y轴线顺时针旋转θ角度，绕x轴线逆时针旋转α角度，绕z轴线逆时针旋转ω角度，再加上一定的平移变换得到[200]。

世界坐标系与摄像机坐标系

在从世界坐标系Oxyz到摄像机坐标系O1x1y1z1的变换中，我们设R′为旋转变换矩阵，M′为平移变换矩阵。可以得到世界坐标系Oxyz与摄像机坐标系O1x1y1z1的变换关系，如下式：

x1xy1RyM z1z在从世界坐标系Oxyz到摄像机坐标系O1x1y1z1的变换中，我们第一步将坐标系Oxyz绕x轴逆时针旋转α角度, 使得Oz1 与Oz2重合，得到过渡坐标Ox2y2z2 ，有下列关系式成立：

x21y020z20cossinxsiny

cosz0再将坐标系Ox2y2z2 绕y轴顺时针旋转θ角度, 使得Ox1 与Ox3重合，得到过渡坐标Ox3y3z3 ，有下列关系式成立：

x3cosy03sinz3010sinx20y

2cosz2 最后，将坐标系Ox3y3z3 绕z轴逆时针旋转ω角度, 使得Oy1 与Oy4重合，得到过渡坐标Ox4y4z4 ，有下列关系式成立：

x4cosysin40z4sincos00x30y3 1z3如果世界坐标系的原点O在摄像机坐标系O1x1y1z1 中的坐标为 M（Mx′, My′,Mz′）,则式(2-5)中的平移变换矩阵MMMx1x4MxyyM 14yz1z4Mzy，可以得到： zx结合式(2-5)～式(2-8)，可以解得旋转变换矩阵R′为：

cos()cos()Rsin()cos()sin()cos()sin()sin()cos()sin()cos()cos()sin()sin()sin()cos()sin()sin()sin()cos()sin()cos()sin()cos()cos()sin()sin()

cos()cos() 我们设世界坐标系Oxyz的坐标原点O在摄像机坐标系O1x1y1z1中的坐标为：（Mx＇,My＇,Mz＇），则可以得到平移矩阵M＇如下：

MxMMy

Mz为了实际计算的方便，将式(2-5)变换为：

2

xx1yRy1M zz1 依照上述，如果摄像机坐标系O1x1y1z1的坐标原点O1在世界坐标系Oxyz中的坐标为：（Mx，My，Mz），则根据上述的公式推导变换，可以得到旋转矩阵R和平移矩阵M分别为：

cos()cos()Rcos()sin()sin()cos()sin()sin()sin()cos()sin()cos()sin()cos()cos()cos()sin()sin()sin()sin()cos()cos()sin()sin()cos()sin()cos()cos()sin()

r1为方便计算，记作 Rr4r7r2r5r8r3r6 r9xMMMMy z 理想的三维空间物体到视平面的投影关系即成像模型是光学中的小孔成像模型[95]。实际的CCD摄像机成像系统通常都是透镜成像,图2-4为透镜成像的原理图。

图透镜成像原理

其中，f微透镜的焦距，v为像距，u为物距，且有如下关系：

1f1u1v

一般地由于u>>f，于是v≈f，这是可以将透镜成像模型近似地用小孔模型代替。在上述模型中，为了方便起见，取坐标系为成正实像的投影变换坐标系，即将视平面的位置与光学中心的位置对调。形成的图像坐标系O＇x′y′与摄像机坐标系O1x1y1z1的成像模型如图所示，其中坐标轴x1、y1分别与坐标轴x′、y′平行，投影中心O1到图像坐标平面O＇x′y′的距离O＇O1即为摄像机的焦距参数f。

3

摄像机成像模型 根据投影变换关系，应有下式成立：

xfx1z1 yfy1z1

为求得x和z，我们联立上列式解得：

xr1x1r2y1r3z1Mx yr4x1r5y1r6z1My

zr7x1r8y1r9z1Mz

结合式(2-17)，(2-18)可以得到： z1xMxfr1xr2yr3f

结合式(2-17)，(2-19)，(2-21)解得：

yr4xr5yr6fr1xr2yr3fr7xr8yr9fr1xr2yr3f(xMx)My

z(xMx)Mz

这样,我们就建立了被测曲面的世界坐标y、z与图像坐标x′、y′以及系统参数α、θ、ω、f、M等的函数关系式。

yfy(x,y,,,,f,Mx,My)

4

zfz(x,y,,,,f,Mx,Mz)

6.1 传统标定

坐标系①世界坐标系(x,y,z) ②摄像机坐标系(xc,yc,zc)

③图像坐标系（u，v） 图像的物理坐标系（x,y）

(以图像与光轴交点为原心，x，y轴分别与u，v平行) 1>坐标变换关系为： xc yczcxwRywzwTT11rr21r31r1222rr32r13r23r33xwywzw T 式① Tt,ty,tz为t界坐标系原点在摄像机坐标系中的坐标。R含有三个独立变量α

θ，w加之tx,ty,tz共六个参数决定了摄像机光轴在一世界坐标系中的空间位置。这6个参数称为摄像机外部参数。 2>图像坐标系与摄像机坐标系

xfxc/zcyfyc/zc

进而可得：

uu0x/dxSxxvv0Y/dySyY

u0，v0是图像中心O的坐标dx，dy分别为一个像素在x，y方向上的物理尺寸，Sx＝1/dx Sy＝1/dy，分别为x，y方向上的采样频率，即单位长度的像素个数。

综合有P与P的变换关系有

uu0fSxxc/zcfxxc/zcvv0fSyyc/zcfyyc/zc

其中

fxfSxfyfSy 分别定义为x和y方向的等效焦距。fx，fy，u0，v0四个参数只与摄

像机内部结构有关，称为内部参数。

5

3>世界坐标系与图像坐标系（共线方程） 式①代入式⑤可得：

xfyfuu0fxvv0fyr11xwr12ywr13zwtxr31xwr32ywr33zwtzr21xwr22ywr23zwtyr31xwr32ywr33zwtz

上式说明：物点，光心和像点必须在同一直线上。

由于成像系统存在着各种误差因素，使得像点，光点和物点不在同一直线上。实际成像模型并不满足线形关系而是一种非线形关系。 中焦镜头图像边界有1～3个像素的畸变 广角镜头更大

要获得高精度必须采用非线形模型来标定。描述图像点的非线形畸变可用下面的公式

XuXdxYuYdy

Xu，Yu为图像点理想坐标 Xd，Yd为图像点实际坐标

x，y为x和y方向上的畸变值，它和图像点的位置有关。将上式代入共线方程中：

(Xdx)/fYdyfr11xwr12ywr13zwtxr31xwr32ywr33zwtzr21xwr22ywr23zwtyr31xwr32ywr33zwtz 即为线形成像模型

根据物点和像点的已知坐标，代入上列方程求解内外参数主要求解方法： ① 直接非线形求解

建立标定点坐标与图解点坐标的投影关系，使用迭代算法对非线形方程求解。 优点：可覆盖所有像差模型，可达到很高精度。 缺点:迭代初始值提供不当，会导致错误稳定性差。 ② 线形求解法

将非线形方程的变量组合成一组新变量——中继参数 原非线形方程转化为中继参数的线形方程

利用最小二乘法求得中继参数值后，再求得原变量的值

6

优点：不需要迭代运算，速度快

缺点：不能考虑系统误差的修正，精度不高

（直接线形变换法DLT直接用共线方程进行标定） ③ 分步法

6.2 Tsais万能摄像机标定法

先获得足够多的空间点和对应像点，用于摄像机标定。认为镜头畸变起主导作用的是径向轴对称畸变。 取轴对称误差模型

xXd(k1k2……)yYd(k1k2……)22424

其中XdYd是图像点到原点的距离

2k1，k2是轴对称误差的系数

将上式代入非线形成像模型中 Xd(1k1k2……)fYd(1k1k2……)f2424r11xwr12ywr13zwtxr31xwr32ywr33zwtzr21xwr22ywr23zwtyr31xwr32ywr33zwtz

第一步：求解系统的外部参数T和R

忽略镜头误差后，标定方程可简化为如下中间变量的线形方程： 1112 [YdiXwi,Ydiywi,Ydi,XdiXwi,Xdiywi]txXdi

2122上式有5个未知数，利用最小二乘方法可以求解。 从这5个中间变量可求得T和R。

第二步：计算内部参数fx,fy和镜头畸变系数k1,k2，令k1,k2初值为零，fx,fy取不考虑畸变估计值，对成像模型进行迭代求解，求得内部参数。 6.3 Weng’s标定法

Tsai’s只考虑了轴对称畸变、偏心畸变、薄棱镜畸变 综合误差模型为

7

x(g1g2)Xdg4XdYdg1Ydk1Xd(XdYd)y(g2g4)Ydg3XdYdg2Xdk1Yd(XdYd)22222222

求解方法与Tsai’s相同，分两步求解。 6.4 几何映射变换

几何映射变换可改变图像中各物体之间的空间位置关系，这种运算可以被看成是将各个图像元素在图像内移动。其效果可形象化地理解为把图像看成一张橡皮膜，对该膜作任意伸展，并在若干点把变形后的橡皮膜固定下来的结果。实际上，几何映射的结果可以比上述比喻宽得多，一个不受约束的几何映射，可以将输入图像中的一个点映射到输出图像中任意位置。也就是说，几何映射可将原图像变得面目全非。但我们所关心的映射是一种保持变换前后图像局部特征相似性的映射。

几何映射首先需要一个算法来定义空间映射本身，用它描述每个像素如何从其初始位置移动到终止位置，即每个像素的运动。

图像几何映射运算的基本定义为：

g(x,y)f(u,v)f(p(x,y),q(x,y))

在上式中，u=p(x,y)，v=q(x,y)唯一地描述了一组几何映射变换，即将输入图像f(u,v)从u-v坐标系变换到x-y坐标系的输出图像g(x,y)。

简单的几何映射有平移、缩放、旋转等，它们的变换矩阵如下： u1v010010x0xy0y 11u1/cv01001/d0x0xy0y 11ucosvsin10sincos00x0y 11平移变换 缩放变换 旋转变换

在上式平移变换中，点(x0,y0)被平移到原点。在缩放变换中，图像在x轴方向上放大c倍，在y轴方向上放大d倍。在旋转变换中，图像被绕原点顺时针旋转角度。

显然，将平移和旋转、平移和缩放结合起来，就可以使图像围绕任意点产生旋转或缩放。利用齐次坐标变换方法可以很简单地确定复合变换公式。

控制点空间映射法是通过测定若干特定坐标点的位移量来确定坐标变换方程的系数的方法。如图4-2所示，若已知输入图像与输出图像上四对对应点对（即控制点对）的位移量，则可利用这四对已知控制点对，求解下列坐标变换方程中的系数。

uaxbycxydvexfygxy

8

系数a～h确定后就可以得到确定的坐标变换关系，用此算法即可确定所有落入矩形框内的输出点。

输入 输出 控制点的映射变换

上式可被看作是下列多项式变换幂函数方程的一个特例。

NNiijuai0j0NNixy

ijvbi0j0ijxy

ij 对一般的N次幂函数，系数个数以及为确定这些系数所需的控制点对数应根据具体N值来确定。一般来讲，适当提高阶次N可提高精度，但随着阶次N的提高，为确定幂函数系数所需的控制点对数量亦随之增多，并且坐标变换所需的计算量亦增加。实际应用中，一般取N=2～3。

一般情况，已知四个像素的实际坐标，就足够建立一组映射函数，但是对于含有径向畸变的图像要实现高精度的几何映射，就需要采用较多的控制点。平均地分布于整张图像中的控制点数愈多，则映射结果愈精确。

对于较多控制点的映射变换，一般是利用如下的二次多项式关系：

ua00a10xa01ya20xa11xya02yvb00b10xb01yb20xb11xyb02y2222

如果选取M个控制点，则有：

u(u1,u2,......,uM)v(v1,v2,......,vM)TT

Ta(a00,a10,a01,a20,a11,a02)b(b00,b10,b01,b20,b11,b02)T

变换矩阵为：

9

11H1x1x2xMy1y2yMx122x1y1x2y2xMyMx22xM2y12y2 2yM以均方误差最小为目标，求a、b的估计值：

E(uHa)(uHa)(vHb)(vHb)

TT令

EaT0,Eb0

T则得到

a(HH)T11HuTb(HH)Hv

6.5 重采样算法

几何映射还需要一个算法用于灰度级的重采样，这是因为在一般情况下，输入图像的每个像素的坐标值为整数，它们经过空间变换映射到目标图像上后，其坐标位置不一定是整数，而数字图像只记录整数位置上灰度值，这就需要将非整数位置上的灰度值经过一定的算法重新分配到整数坐标位置上，反过来也是如此。

1) 最近邻法

最近邻法是将距(u0,v0)点最近的整数坐标(u,v)点的灰度值取为(u0,v0)的灰度值，如下图4-3所示。在(u0,v0)点个相邻像素间灰度变化较小时，这种方法是一种简单快速的方法，但当(u0,v0)点相邻像素间灰度差较大时，这种灰度估值方法会产生较大的误差。表现为在图像中相邻像素间灰度差很大的边缘处出现明显的锯齿，图4-4为最近邻算法重采样效果，4-4(a)为原图像，4-4(b)为旋转后最近邻法重采样图像， 4-4(c)为重采样图像局部放大图。

(u,v)(u0,v0)(u,v+1)(u+1,v)(u+1,v+1) 图4-3 最近邻法 10

(a) (b)

(c)

图4-4 最近邻法重采样效果 2) 双线性插值法

双线性插值法是对近邻法的一种改进，即用线性内插方法，根据(u0,v0)点的4个相邻点的灰度值，插值计算出f(u0,v0)值。具体计算过程如图4-5所示。 f(u+1,v)f(u,v+1)f(u0,v0)f(u,v)(u+1,v)(u+1,v+1)(u,v)(u,v+1)(u0,v0)f(u+1,v+1) 图4-5双线性插值法 1、先根据f(u,v)及f(u+1,v) 插值求f(u0,v)。

f(u0,v)f(u,v)[f(u1,v)f(u,v)]

(3-10)

2、再根据f(u,v+1)及 f(u+1,v+1)插值求的f(u0,v+1)。

f(u0,v1)f(u,v1)[f(u1,v1)f(u,v1)] (3-11)

3、最后根据f(u0,v)和f(u0,v+1)插值求f(u0,v0)。

(3-12)

上述f(u0,v0)的计算过程实际上是根据f(u,v),f(u+1,v),f(u,v+1)和f(u+1,v+1)四个整数点的

11