多视图射影重构.doc

第三章 多视图射影重构
多视图（三幅以上）分层重构的基本思路：
（1）射影重构
　　　　①估计射影深度；
　　　　②构建一个关于射影深度和图像坐标的测量矩阵；
　　　　③通过矩阵分解的方法，求解射影空间下的摄像机运动参数和空间物体三维几何形状。
　　（2）欧氏重构
　　　　①摄像机自标定（求解摄像机内参数）；
　　　　②求解射影空间到欧氏空间的变换；
③恢复欧氏空间下的摄像机运动参数和空间物体三维几何形状。
由此，我们知道在实现分层重构的过程中，射影重构是非常重要的第一步，它对欧氏重构结果的精度起着决定性的作用，因此引起许多研究者极大的兴趣，涌现出了许多方法和技术[1-24]，其中矩阵分解方法经理论和实验证明是一种非常有效的方法，是近年来研究的焦点之一[1]。
本章从摄像机模型出发，详细推导了透视模型与仿射模型的关系，讨论了射影深度的性质，为迭代估计射影深度时合理选择初始值，提供了理论依据。介绍了基于奇异值分解的射影重构算法的一般框架，分析并实现了基于基本矩阵和极点的射影深度估计算法。以测量矩阵的秩为4作为约束，以仿射投影逼近透视投影，提出以下迭代估计射影深度的方法：①基于共轭梯度法的估计射影深度算法，②基于遗传算法的估计射影深度算法。在获得正确的射影深度后，通过奇异值分解将测量矩阵分解为射影空间下的摄像机运动和物体三维几何形状（射影重构）。实验证明：相对于基于基本矩阵和极点方法来计算射影重构，我们所提出的方法对噪声具有更好的鲁棒性。
3.1多视图摄像机模型
根据第二章所述，我们可以写出多视图摄像机模型。
假设摄像机从F个视点对P个三维空间物体点拍摄F幅图像，如图3-1所示。其中（j=1,2,…,P）为第j个三维空间物体点，它在第i（i=1,2,…,F）幅图像上的像点为。这样，我们可以给出摄像机的多视图投影方程为：
（3.1）
其中是一个3维向量，是一个3×3旋转矩阵，分别表示第i个摄像机的位置和方向，是一个3×3上三角矩阵，表示摄像机的内参数，表示空间点到第i个视点的射影深度。令，即为一个2×3矩阵,表示的第一、二行，为一个3维向量，表示的第三行，则
（3.2）
从（j=1,2,…,P）中任取一点作为参考点，令它在第i幅图像上的像点为。将世界坐标系原点平移到，同样，将每一个图像坐标系原点平移到，记：

令，则式（3.2）可改写为：

即

（3.3）
由于是在第i幅图像上的像点，所以

进而得

代入式（3.3）可得
（3.4）
式中 为参考点的深度、为深度的变化。式(3.4)称为透视摄像机模型。
如果相对于充分小，则由式(3.4)可近似得到:
　　（35）
式(3.5)称为仿射摄像机模型。

xj
C1
C2
CF
m1j
m2j
mFj
Ow
······
Q1
Q2
QF

图3-1多视图投影示意图
3.2基于奇异值分解的射影重构
3.2.1多视图测量矩阵
引入齐次坐标向量和，分别表示空间物体点和其像点的齐次坐标，同时，令

这样，根据式（3.4）,我们可以得到以下投影方程：
　 （3.6）
其中表示第i个摄像机的投影矩阵，是一个3×4秩3的齐次矩阵。
对于F幅图像和P个空间物体点，则有：
　 （3.7）
其中矩阵称为测量矩阵，式（3.7）表明可以分解为一个表示摄像机运动的矩阵和一个表示空间物体形状的矩阵。在一般情况下，由于未知，从而也未知，其分解式（3.7）更不得而知。下面我们先来讨论一下射影深度的性质。
3.2.2射影深度
射影深度是一个比例因子，它完全取决于、和所采用的归一化约束。为了便于说明，将式（3.1）改写为：
（3.8）
当图像点坐标、空间点坐标以及投影矩阵被仿射归一化时，即对应的齐次坐标的最后一个元素取1。由式（3.8）可得：
（3.9）
其中是的第三行，表示一个通过摄像机光心且平行于像平面的平面，称为摄像机光平面。如果是一个无穷远平面，则任何有限远点的射影深度均变为：

此时摄像机为仿射摄像机。因此，也就不需要恢复射影深度了，通过对测量矩阵进行SVD分解，即可获得射影意义下的摄像机的运动和空间物体的形状。
另外，如果是一个有限远平面，并且在欧氏意义下将的前三个元素归一化为一个单位向量，即取

其中表示向量的欧氏范数，在这种情况下

即，射影深度就变成了传统的光学深度，也就是到光平面的正交距离。
我们还可以引入一个3D参考平面和一个参考点，对和进行归一化处理，使

这种归一化是可以实现的。先选择和，使之满足，由于和都是齐次化的，在相差一个非零常数因子的情况下，其非齐次元素是相等的。所以，可以对选择一个适当的常数因子，使，对选择一个适当的常数因子，使，因此可以实现上述归一化处理。 经过这种归一化处理之后，我们可得到：
(3.10)
其中是一个射影不变量，称为交比，可用来表示。式（3.10）表明：经过这种归一化处理之后，射影深度可以利用射影不变量—交比来确定。
如果或者，则。当时，我们有

即相当于仿射归一化的情况。当时，摄像机即变为仿射摄像机。
结合3.1节的分析，我们可以得出这样的结论：在迭代估计射影深度时，我们可以选择的初值为1，即可以用仿射投影来逼近透视投影。
3.2.3基于奇异值分解的射影重构
如果由某种方法给出了正确的射影深度，则将作奇异值分解（SVD）得：
　　　　　　 (3.11)
式中为的奇异值，，、分别为、的次正交矩阵。如果图像坐标无噪声（即Ｆ幅图像精确匹配），则的秩为４，即，从而
(3.12)
式(3.12)为的一个满秩分解，其中和分别为和的前4列。结合式（3.7）可得：
　　(3.13)
但此分解不是唯一，因为式（3.7）在相差任意非零常数因子的情况下可以写为：
　　　 　 　　(3.14)
其中为任意非奇异4×4矩阵，因此均为分解。这就表明可以在一个未知射影变换下恢复摄像机运动和空间物体形状，并可以任意选择一个3D射影坐标系来描述这种运动和形状，这种恢复称为射影重构。
由于构成测量矩阵的图像数据之间的不平衡，优化过程有时是非常不稳定的，为了避免这种不稳定，对图像坐标进行归一化处理是非常有效的。本章归一化处理的方法是：每一幅图像中的乘上一公共因子(i=1,…F)，使它们的平均范数为。
由于式（3.7）在相差一个非零常数因子的情况下仍然成立，所以，对于式（3.7）来说，在射影重构未完成之前，我们可以让图像点坐标、空间点坐标以及投影矩阵的归一化保持浮动，即可以任意地以代替、代替，而相应地代之以。也就是说，射影深度矩阵[]的每一行或每一列可以任意地引入非零常数因子，而不会破坏测量矩阵W的秩为4的性质。这样，在求解射影深度的过程中，我们可以保持射影深度矩阵[]的某一行和某一列的元素为1，于是射影深度矩阵[]的自由度变为：

剩余的的元素仍可给出射影重构所需的信息。在具体实验中，我们保持射影深度矩阵[]的第一行和第一列的元素为1。
由以上分析可知，在透视分解算法中，最关键的问题是如何估计射影深度。下面先分析基于基本矩阵和极点的射影深度估计算法，然后再介绍我们所提出的基于共轭梯度法的射影深度估计和基于遗传算法的射影深度估计算法。
3.3射影深度估计及射影重构
3.3.1基于基本矩阵和极点的射影深度估计及射影重构
对于摄自不同视点的未标定图像，基本矩阵是两幅图像间极几何的代数刻划，反映了两幅图像之间的对应关系。我们假设从不同视点摄取三幅图像u、v、w，是图像v和图像u之间的基本矩阵，是摄像机v在图像u中的极点，（j=1,2,…,P）为第j个三维空间物体点，它在第i（i=u,v,w）幅图像上的像点为。
因为图像u中的像点和图像v中的像点是一对对应点，且在图像u中对应的极线与经过其对应点和极点的直线是同一条直线，所以我们有
　　　 　　　　（3.15）
其中，、 是、的齐次坐标表示。由式（3.15）可解得：
　　 　　　　　（3.16）
其中表示向量的欧氏范数。同样，我们可以获得：
　　　　　　　（3.17）
另外，由于摄像机u在图像v中的极点与其在图像w中的极点是一对对应点，因此在相差一个非零常数因子的情况下，可得下列方程：
　　　 　　　　（3.18）
这是３幅图像u、v、w之间的极点约束。
　　根据式（3.16）、（3.17）、（3.18），我们可以利用基本矩阵和极点来计算射影深度，计算过程中我们取(j=1,…,P)为初始值。
基于基本矩阵和极点的射影重构算法如下:
①在不改变摄像机光学系统（如焦距、光圈等）的情况下，从F个不同的视点拍摄F幅图像；
②从F幅图像中选择P个特征对应点并提取它们的图像坐标 (i=1,…,F, j=1,…,P)，在其中任取一点作为图像坐标原点，计算；
③对图像进行归一化，即每一幅图像中的乘上一公共因子(i=1,…F)，使它们的平均范数为。可以推出；
④取初始值（j=1,…,P）；
⑤利用式（3.16）、（3.17）、（3.18）计算射影深度(i=2,…,F, j=1,…,P)；
⑥根据式（3.7）构建测量矩阵；
⑦按式(3.12)、(3.13)将分解为射影空间下的摄像机运动和空间物体形状。构建一个对角矩阵，则就是关于未曾归一化的原始图像坐标的摄像机运动参数。
3.3.2基于共轭梯度法的射影深度估计及射影重构
由式（3.7），我们不难看出，在图像坐标无噪声且射影深度为正确值的理想情况下，测量矩阵的秩，当有噪声但不大时，。因此可以定义一个的测量指标为：
　 　　　　 　　（3.19）
是依赖于的非负实数，是的充要条件，于是求使取最小值的{}将给出所需要的射影深度。由于不是的线性函数，故最小化是一个非线性优化问题。下面我们采用共轭梯度法来解决这一问题。
对于非线性优化问题，在逼近时尽可能使初始估计接近最优解是非常关键的。根据3.1和3.2.2节的分析,我们可以用仿射投影来逼近透视投影，即取作为初始值。
使用共轭梯度法的前提是计算关于的梯度。由式(3.11)知，，其中、分别为次正交矩阵和的第ｎ列。所以

又由可推得
，
则
。
由此得关于的梯度为：
（3.20）
其中、，和均为三维向量。
基于共轭梯度法的射影重构算法如下:
①在不改变摄像机光学系统（如焦距、光圈等）的情况下，从F个不同的视点拍摄F幅图像；
②从F幅图像中选择P个特征对应点并提取它们的图像坐标(i=1,…,F, j=1,…,P)，在其中任取一点作为图像坐标原点，计算；
③对图像进行归一化，即每一幅图像中的乘上一公共因子(i=1,…F)，使它们的平均范数为。可以推出；
④取（i=1,…,F, j=1,…,P）；
⑤根据式（3.7）构建测量矩阵W并对其进行奇异值分解(SVD)；
⑥按式(3.19)计算，若足够小（如小于10-8），则转至⑦；否则按共轭梯度法计算关于的梯度共轭方向的最优迭代增量，将，返回⑤继续迭代；
⑦按式(3.12)、(3.13)将分解为射影空间下的摄像机运动和空间物体形状。构建一个对角矩阵，则就是关于未曾归一化的原始图象坐标的摄像机运动参数。
在迭代估计射影深度时，我们以测量指标作为结束循环迭代的判据，比一般情况下以迭代增量作为判据更为合理和有效，因为在噪声存在的条件下，我们发现会出现某一步迭代增量变得很小而却相对很大的情况。
3.3.3 基于遗传算法的射影深度估计及射影重构
1.遗传算法原理
遗传算法是模拟生物进化过程的计算模型。遗传算法作为一种新的全局优化搜索算法，以其简单通用、鲁棒性强、适于并行处理以及应用范围广等显著特点奠定了它作为21世纪关键智能计算之一的地位。
在问题解决的过程中，遗传算法是一种迭代算法，它从一组变量中搜寻一个最优解。对于一个具体的问题，用遗传算法搜寻最优解的过程有4个主要步骤：
（1）编码：在候选解的范围内，随机选择一个候选解作为初始解。每一个候选解都可用一个序列表示，它们可以用二进制或者十进制编码，例如，8被编码成00001000，被编码后的候选解被称为个体。
（2）选择：根据每个个体的适应度函数的值，基于一个确定的概率选择一些个体作为父代，一个个体的适应度函数值越大，被选择的概率就越大。
（3）交叉：从当前交配池中，基于一个确定的概率随机地选择一些个体作为下一代的父代，每两个父代个体按一定的概率通过交叉产生两个新的个体，这两个新的个体保留了它们父代的部分基因（它们的一部分编码未变），也引进了一些新的基因。
（4）变异：基于一个确定的概率，随机地选择一些个体的部分基因并且把它们变成新的基因。例如，对于二进制编码，一个初始位元为1的位变成0，初始位元为0的位变成1；对于十进制编码，即在设定的范围内，随机地选择一个新的数值替代原来的基因。变异在遗传算法中非常重要，虽然变异的概率很小，但却能够在解决问题的过程中防止早熟收敛所引起的局部最优解。
与传统优化算法相比较，遗传算法具有如下特征：
（1）最优解搜寻的过程不影响变量，但是变量的编码对寻优过程会产生影响；
（2）搜寻最优解的过程是从一组解到另外一组解的过程，而不是从一个解到另外一个解的过程，这个过程能够防止局部收敛，因此它的全局收敛能力较强；
（3）对搜寻空间没有特别的要求，仅仅利用适应值而不需要辅助的信息，因此它的应用范围很广泛。
2.基于遗传算法的射影深度估计及射影重构
实现射影重构的关键是获得正确的射影深度。下面利用遗传算法并结合奇异值分解技术来估计射影深度并实现射影重构。
(Ⅰ) 解的编码与适应度函数
对于参数（i=1,…,F, j=1,…,P）采用十进制编码。令是射影深度的一个解，以作为一个个体进行编码。由于透视投影模型可以利用仿射投影模型来迭代逼近，而仿射投影模型下的仿射深度是一个常量（通常取1），因此可以在（0,2]区间内随机产生。
在图像坐标无噪声且射影深度为正确值的理想情况下，测量矩阵的秩，当有噪声但不大时，。因此我们可以定义遗传算法的适应度函数为：
　 　 　　　　（3.21）
是依赖于的非负实数，为了使必有，于是求使取最大值的将给出所需要的正确射影深度。由于不是线性函数，故最大化是一个非线性优化问题。
（Ⅱ）基于遗传算法的射影深度估计
①在不改变摄像机光学系统（如焦距、光圈等）的情况下，从F个不同的视点拍摄F幅图像；
②从F幅图象中选择P个特征对应点并提取它们的图象坐标(i=1,…,F, j=1,…,P)，在其中任取一点作为图像坐标原点，计算；
③对图像进行归一化，即每一幅图像中的乘上一公共因子 (i=1,…F)，使它们的平均范数为；
④在（0,2]区间内随机选取个体，取交配池的容量为20；
⑤对每一个个体所对应的测量矩阵()进行奇异值分解(SVD)，并由适应度函数公式计算每个个体的适应度,；
⑥每次进行遗传操作,以概率复制；
⑦以概率用赌轮法随机选取两个个体，对其进行交叉；以概率进行变异操作；
⑧把经过遗传操作后得到的个体都放在交配池中。当个体的个数超过交配池容量时，就将适应度小的个体从交配池中删去；
⑨进行上述遗传操作至第代后(是预先给定的常数),以(也是预先给定的常数)对交配池中适应值小的个体进行删除，代之以新的个体；
⑩进行上述遗传操作至第代后(是预先给定的常数)，在第代的交配池中取适应度最大的个体，即为所求的个体。
（III）按式(3.12)、(3.13)将分解为射影空间下的摄像机运动和空间物体形状。构建一个对角矩阵，则就是关于未曾归一化的原始图像坐标的摄像机运动参数。
3.4仿真实验
在仿真实验中，我们采取以下实验步骤：
①随机模拟P个三维空间物体点；
②将P个三维空间物体点通过F个模拟的投影矩阵（i=1,2,…,F）投射到F幅图像上；
③在图像坐标上分别加上0.0、0.2、0.4、0.6、0.8、1.0、1.2、1.4、1.6、1.8、2.0（像素）的高斯噪声；

④利用上述方法，完成P个三维空间物体点的射影重构；
⑤将P个射影重构点再反投影到图像上，计算图像点的误差，计算公式如下：

其中表示向量的欧氏范数。
我们模拟了三个视点下拍摄的三幅图像，第一、二、三个视点的摄像机内参数矩阵的理论值均为：

第一、二、三个视点的摄像机坐标系关于世界坐标系运动的理论值分别为：
，
，
，
我们模拟了14个空间点，利用以上参数将其分别投影到三幅图像上，仿真的三幅图像如图3-2所示。实验结果如图3-3所示，其中(a)为基于基本矩阵和极点估计射影深度，利用SVD实现射影重构；(b)为基于共轭梯度法估计射影深度，利用SVD实现射影重构；(c)为基于遗传算法估计射影深度，利用SVD实现射影重构。
实验中，在每个噪声量级上，重复实验100次，然后求平均。从结果上可以看出，后两种算法的鲁棒性明显好于第一种算法，这是因为后两种算法不涉及对噪声非常敏感的基本矩阵和极点。

图3-2 仿真图像

(a)

(b)

(c)
图3-3 2D误差曲线
参考文献
[1]C.Rothwell, O.Faugeras, A comparison of projective reconstruction methods for pairs of views. Computer vision and image understanding,Vol.68,No.1 pp.37-58,1997.
[2]A.Shashua, Projective depth: A geometric invariant for 3D reconstruction from two perspective/ orthographic views and for visual reconstruction, In Proc. ECCV’92, 1992.
[3]R.Koch, M.Pollefeys, Multi viewpoint stereo from uncalibrated video sequences, Proc. ECCV’98, Vol.I, pp.55-71, 1998.
[4]S.M.Seitz, Image－based transformation of viewpoint and Scene appearance, Ph.D. Thesis, Wisconsin University, 1997
[5]C.Tomasi, T.Kanade, Shape and motion from image streams under orthography: a factorization method. International Journal of Computer Vision, Vol.9 No.2, pp.137-154, 1992.
[6]C.J.Poelman,T.Kanade, A paraperspective factorization method for shape and motion recovery. Proc.ECCV’3, Vol.2, pp.97-108, 1994.
[7]C.Christy, Horaud.R. Euclidean reconstruction:from paraperspective to perspective. In Proc. ECCV’4, Vol.2, pp.129-140, 1996.
[8]B.Triggs, Plane+parallax, tensor and factorization. Proc. ECCV’2000, 6/7, pp.523- 538. 2000.
[9]B.Triggs, Linear projective reconstruction from matching tensors. Image &Vision Computing, Vol15, No.8,pp.617-626,1997.
[10]B.Triggs, Factorization methods for projective structure and motion. In Proc. CVPR’96, pp.845- 851, 1996.
[11]A.Shashua, The rank 4 constraint in multiple (³3) views geometry. In Proc. ECCV’96, pp.196- 206, 1996.
[12]B.Triggs, The geometry of projective reconstruction I:Matching constraints and the joint image. Int.Conf. Computer Vision, pp.338-343,1995.
[13]T.Morita, T.Kanade, A sequential factorization method for recovering shape and motion from image streams. In ARPA Image Understanding Workshop, pp.1177-1188, 1994.
[14]R.Hartley, Projective Reconstruction and Invariants from Multiple Images, IEEE-T PAMI 16:10, pp.1036-1041, 1994
[15]Long QUAN and Takeo Kanade, Affine structure from line correspondences with uncalibrated affine cameras, IEEE-Trans. PAMI, Vol.19, No.8, pp.834-845, 1997.
[16]M.Berry, Large scale singular value computations. Int.J. Supercomputer Applica- tions, Vol.6, No. 1,pp.13-49, 1992.
[17]U.Toshio,T.Fumiaki, A factorization method for projective and Euclidean recons- truction from multiple perspective views via iterative depth estimation. ECCV’5, Freiburg, Vol 1, pp 296-310,1998.
[18]R.Hartley, A.Zisserman, Multiple View Geometry in Computer Vision. Cambridge University Press 2000.
[19]A.Shashua, Projective Structure from Uncalibrated Images: Structure From Motion and Recognition, IEEE-T PAMI, Vol.16,No.8, pp.778-790, 1994.
[20]A.Shashua, Algebraic functions for recognition. IEEE Transaction on Pattern Analysis and Machine Intelligence, Vol17,No.8,pp.779-789,1995.
[21]A.Heyden, Reconstruction from image sequences by means of relative depths.IEEE Int. Conf.Computer Vision, pp.1058-1063, 1995.
[22]O.Faugeras, B.Mourrain, On the geometry and algebra of the point and line corres- pondences between N images. IEEE Int. Conf. Computer Vision, pp.951-956, 1995.
[23]R.Mohr, B.boufama, Accurate projective reconstruction. In 2nd Europe-U.S.Work- shop on Invariance, pp.257-263, 1993
[24]Q.T.Luong,T.Vieville,Canonical representations for the geometries of multiple pro-
jective views. Technical Report, INRIA, France, 1993.
[25]刘勇等，非数值并行算法----遗传算法，科学出版社，1995。
17