北理工贾云德《计算机视觉》chapter15三维运动估计

211第十五章三维运动估计三维运动估计是指从二维图象序列来估计物体三维运动参数以及三维结构。具体地说，假定三维物体上一点M相对于摄象机坐标系从时刻kt的位置),,(kkkzyx运动到时刻1kt的位置),,(111kkkzyx，它在二维图象平面上的投影从),(kkyx运动到),(11kkyx，然后，通过分析二维运动来恢复物体的三维运动及物体上感兴趣点的深度值。这一点类似于立体视觉的深度恢复，不过立体视觉是从立体图象对来恢复深度值，而三维运动分析是从图象序列中恢复参数。三维运动估计有着广泛的应用，如机器人视觉，自主导航，目标跟踪，图象监控，智能车辆高速公路系统，基于物体的图象压缩等。三维运动估计仍然是一个不适定问题，必须增加适当的约束才能得到有效解。三维运动估计和分析也可以是基于场景的深度图像序列，其方法与基于二维图像序列完全不同。基于深度图像序列的三维运动估计是一个适定问题，求解方法要比基于二维图像序列要简单一些，主要问题是数据量十分巨大，因此许多研究人员一开始就以实现大规模集成电路（VLSI）作为三维场景估计的基础。我们将不讨论深度图像序列运动估计问题，感兴趣的读者可以查阅有关文献[Wheeler1996,Jiar1996].15．1基于成象模型的对应点估计首先回顾一下第十二章引入的三维刚体运动方程。假定三维场景中有一个刚性物体，其上一点M从时刻kt的位置),,(kkkzyx经过旋转和平移，运动到时刻1kt的位置),,(111kkkzyx。设旋转矩阵和平移向量分别是kR和kT,则三维刚体运动模型重新表示为kkkkkzyxkkkzzzyzxyzyyyxxzxyxxkkkzyxtttzyxrrrrrrrrrzyxTR111(15.1)用欧拉角的形式表示上述旋转矩阵(见式(12.2),（12.3）)，并假定旋转角较小，则旋转矩阵可以表示为111kR(15.2)其中,,分别表示绕zyx,,轴逆时针旋转小角位移。15．1．1正交投影模型当物体深度变化范围不大时，正交投影是透视投影的一个很好的逼近。其它逼近方法还有弱透视投影，超透视投影，正交透视投影等[Dementhon1992]。设空间点),,(zyx在图象平面上的投影为),(yx。如果成象模型为正交投影（参见图1.8），则有xxyy(15.3)所以，式(15.1)可表示为)()(1111ykyzkyykyxkkxkxzkxykxxkktzryrxryytzryrxrxx(15.4)上述方程包含有6个参数，即xxr，xyr，yyr，yxr，)(xkxztzr和)(ykyztzr，表示第k帧212图象像素),(kkyx到第1k帧图象像素),(11kkyx的仿射映射关系。显然，正交投影模型是无法确定物体点到成象平面的距离，因为垂直于图像平面的一条直线，其上的所有点都将投影到该图象平面上一点（见图1.8）。但是，如果在物体上选择一个参考点，并设该点的深度值为refz,则通过上式有可能估计出物体上其它点相对于这个参考点在垂直图像平面方向上的距离relz。实际上，我们无法得到真实的相对深度值，只能得到关于一个尺度系数的相对深度值，即relrefzzz(15.5)从方程15.4中，我们看到xzr和yzr缩小，kz放大，方程仍然成立，因此产生多义性。文献[Ull79]证明三帧或三帧以上图象上的四点就可以克服这个问题。15．1．2基于正交投影的运动估计将(15.2)小角度旋转矩阵代入(15.1)，得到如下的正交投影模型：ykkkkxkkkktzyxytzyxx11(15.6)在该方程中，对于每一个给定的对应点),(kkyx和),(11kkyx，共有6个未知参数，其中5个是全局运动参数，，，xt，yt，一个是深度参数kz。另外，这个方程是一个双线性方程，因为kz与未知参数和是乘积关系。理论上，给定三点，就可以根据(15.6)列出6个方程，从而解出六个未知参数。但实际上，由于数值计算误差，常常需要多个点，这样有可能得到较好的结果。文献[Aizawa1989]基于上述正交投影模型提出了基于两帧图象的两步迭代法：首先，根据上一次迭代得到的深度估计值，确定运动参数，然后再使用新的运动参数更新深度估计值。具体实现见算法15.1算法15.1基于两帧图像的运动估计两步迭代算法①给定n个对应点坐标对)},{(,,kikiyx)},{(1,1,kikiyx和深度估计值}{,kiz，ni,,2,1，且3n，这样方程(15.6)可重新写为yxkkkkkkkkttxzyzyyxx10001011(15.7)n个对应点对应着n2个方程，而未知参数仅有5个，因此，可以通过最小二乘法来求解这5个运动参数。深度参数的初始估计值可以根据场景的先验模型来设置，深度估计值应在预先设定的范围内选定，这主要是为了避免解的不唯一性。②根据①得到的运动参数估计值，再对深度值}{,kiz进行估计。将式(15.7)重新写为kykkkxkkkztxyytyxx11(15.8)由于每一个深度值对应两个方程，即方程(15.8)是一个超定方程，因此，可以用最小二乘法来求解。③重复上述两步，直到两次迭代值之差小于给定的某一个阈值。请注意，在上述算法中，运动估计误差和深度估计误差有着密切的关系。由方程(15.7)和(15.8)可知，深度估计的随机误差会重复反馈到运动估计上。因此，当深度估计不够准确或深度的初始值设置不当时，都可能导致迭代算法的错误收敛或收敛在一个局部最小值。为了避免这种错误的收敛，[Bozdagi1994]提出了改进的算法，该算法的基本思想是在每一次修正后，在深度估计值上加一个随机扰动。改进的算法如算法15.1所示。213算法15.2基于两帧图像运动估计扰动迭代算法：①初始化深度值}{,kiz，ni,,2,1，置迭代计数器0m。②在给定深度值下根据式(15.6)估计运动参数③根据当前的运动估计和深度参数，由公式(15.6)计算对应点的坐标),()(1,)(1,mkimkiyx④计算预估误差：niimenE11(15.9)其中2)(1,1,2)(1,1,)()(mkikimkikiiyyxxe,1,kix和1,kiy是已知的对应点真实坐标。⑤如果mE小于预定的误差阈值TE，即TmEE，则终止迭代，否则，置1mm。⑥给深度参数赋一个扰动值)()1(,)(,miimkimkizezz(15.10)其中和是常系数，),0()(2)(miimiN是零均值高斯分布函数，其方差imie)(2。⑦回到第②步实验证明，这种改进的迭代算法在初始深度值有50%的误差的情况下，也能很好地收敛到正确的运动参数值。15．1．3透视投影模型设空间点),,(zyx在图象平面上的投影为),(yx。如果成象模型为透视投影，则zxFxzyFy(15.11)根据（15.1）式有zkzzkzykzxykyzkyykyxkkkzkzzkzykzxxkxzkxykxxkkktzryrxrtzryrxrFzyFytzryrxrtzryrxrFzxFx111111(15.12)由于成象系统的焦距是一个常数，因此，不乏一般性，取1F,即规范化透视投影。上式右边分子分母同除以kz，得到图象平面坐标表示式：kzzzkzykzxkyyzkyykyxkkzzzkzykzxkxxzkxykxxkztryrxrztryrxryztryrxrztryrxrx////11(15.13)按照图像平面坐标，透视投影模型（15.13）是一个非线性方程。因为每一点对应的深度值z是一个自由参数，因此，这个模型适合于任意表面形状三维物体的运动估计。15．1．4外极线方程和基本矩阵由方程(15.13)可见，在透视投影情况下，运动和结构参数之间的关系是非线性的。早期人们使用迭代方法求解运动参数，但迭代过程往往不收敛，比如两步迭代法。[Huang1986]证明，使用八个对应点或更多对应点来求取外极线约束方程并估计运动参数是可以改进两步迭代法的收敛性能。本节主要讨论外极线方程及其性质，关于估计运动参数，将在下一节介214绍。将方程15.1缩写为kkkkTMRM1(15.14)或0)(1kkkTkMRTM(15.15)令kkkRTE(15.16)其中“”表示矢量积。为了书写简单，在不引起混淆的情况下，可以把下标k去掉，即：RTE（15.17)方程(15.15)可以重新写为01kTkEMM(15.18)下面引进一个反对称矩阵：000][xyxzyzttttttT因此式(15.17)可表示为RTE][(15.19)E是一个33矩阵，称为基本矩阵(essentialmatrix)，矩阵元素称为基本参数[Huang1986]。下面首先看一下矩阵E的一些性质[Faugeras1993]：1）0detE2）0TET3）TTTTITTEET)(4）222TE,其中22ijeE，称为Frobenius模。在性质1）中，E是反对称矩阵，因此该性质是成立的；对于性质2），0][)]([TTRTTRTRTTETTxTTxT；对于性质3），TTTxTxTxxTTITTTRRTRTRTEET)(])[]([)])([]([；性质4）显然是成立的。用1kkzz除以等式(15.18)的两边得01)1(1111kkkkkkkkzyzxzyzxE(15.20)根据规范化透视投影方程(15.11)，得到外极线方程：01)1(11kkkkyxyxE(15.21)Tkkkyx),(m和Tkkkyx),(111m是物体上一点M在第k帧和第1k帧图像上的投影点坐标，其齐次坐标分别为Tkkkyx)1,,(~m和Tkkkyx)1,,(~111m，则外极线方程(15.21)可以表示为0~~1kTkmEm(15.22)215该方程的几何意义十分明显，从图15.1可见，三条线kkcc1，11kkmc和kkmc是共面的，将这三个矢量变换到第二个坐标系中，其表示式分别为kT，1~km和kkmR~。第k帧图像上一点km对应矢量kkmc表示的射线上所有点的投影，射线上任何一点可以表示为km~，),0[。根据针孔模型，射线上任何一点在第1k帧图像平面上的投影是TmRk~，而整个射线在第1k帧图像平面上的投影是一条直线，表示成为第k帧图像上一点km在第1k帧图像平面上的外极线。这条直线可以用两点来表示，一点是Te1k，0，对应第k帧图像光学中心kc在第1k帧图像平面上的投影，常称为极点（epipole）；另一点是无穷远点，kkmRm~,1。则外极线可以表示为kkkmRTme~,11(15.23)第k帧图像上一点km在第1k帧图像平面上的对应点1km应该位于该直线上。图15.1基于对应点的外极线约束恢复运动参数示意图方程（15.22）是包含9个未知参数的齐次线性方程，齐次线性方程无解或有无穷解。因此，可以令基本矩阵E的一个系数为1，这样待估计的参数有八个。另外，矩阵E是一个斜对称矩阵和一个旋转矩