多视图几何(11)_结构计算

结构计算

本章介绍 3 D 3D 3D点在两幅视图上的像和这些视图的摄像机矩阵给定时,如何计算该点的位置。因为被测量点 x , x ′ x,x^\prime x,x都有误差,所以从点反向投影额射线不共面。这意味着不存在准确的满足 x = P X , x ′ = P ′ X x=PX,x^\prime = P^\prime X x=PX,x=PX的点 X X X,并且图像点不满足对极几何约束 x ′ T F x = 0 x^{\prime T}F x = 0 xTFx=0

三角测量法的一个理想特性是其在重构的某种适当的变换下保持不变。三角测量法的关键思想是估计一个精确满足给定的摄像机几何的 3 D 3D 3D X ^ \hat X X^,因此它投影为 x ^ = P X ^ , x ^ ′ = P ′ X ^ \hat x = P\hat X,\hat x^\prime = P^\prime \hat X x^=PX^,x^=PX^且目标是由图像测量 x , x ′ x,x^\prime x,x来估计 X X X

  • 线性三角测量法
    由于 x = P X , x ′ = P ′ X x=PX,x^\prime = P^\prime X x=PX,x=PX,构建 x × ( P X ) = 0 x \times (PX)=0 x×(PX)=0构建 A X = 0 AX=0 AX=0的方程求解。但是这种方法不是射影不变的
  • 几何误差代价函数
    寻求最小化函数
    C ( x , x ′ ) = d ( x , x ^ ) 2 + d ( x ′ , x ^ ′ ) 2 , x ^ ′ F x ^ = 0 \mathcal C(x,x^\prime) = d(x,\hat x)^2 + d(x^\prime,\hat x ^\prime)^2,\hat x^\prime F \hat x = 0 C(x,x)=d(x,x^)2+d(x,x^)2,x^Fx^=0
  • Sampson近似
    在这里我们关心的是对测量点进行矫正计算。测量点 X = ( x , y , x ′ , y ′ ) X=(x,y,x^\prime,y^\prime) X=(x,y,x,y) S a m p s o n Sampson Sampson矫正 δ X \delta_X δX
    δ X = − J T ( J J T ) − 1 ε \delta_X = -J^T(JJ^T)^{-1} \pmb \varepsilon δX=JT(JJT)1ε
    矫正的点为
    X ^ = X + δ X = X − J T ( J J T ) − 1 ε \hat X = X+\delta_X = X-J^T(JJ^T)^{-1} \pmb \varepsilon X^=X+δX=XJT(JJT)1ε
    前面已经证明,对 x ′ T F x = 0 x^{\prime T}Fx=0 xTFx=0所定义的族,该误差是 ε = x ′ T F x \pmb \varepsilon =x^{\prime T}Fx ε=xTFx,其 J a c o b i a n Jacobian Jacobian
    J = ∂ ε / ∂ x = [ ( F T x ′ ) 1 , ( F T x ′ ) 2 , ( F x ) 1 , ( F x ) 2 ] J = \partial \pmb \varepsilon / \partial x = [(F^Tx^\prime)_1,(F^Tx^\prime)_2,(Fx)_1,(Fx)_2] J=ε/x=[(FTx)1,(FTx)2,(Fx)1,(Fx)2]
    矫正点的一阶近似为
    [ x ^ y ^ x ^ ′ y ^ ′ ] = [ x y x ′ y ′ ] − x ′ T F x ( F T x ′ ) 1 2 2 + ( F T x ′ ) 2 2 2 + ( F x ) 1 2 2 + ( F x ) 2 2 [ ( F T x ′ ) 1 ( F T x ′ ) 2 ( F x ) 1 ( F x ) 2 ] \begin{bmatrix} \hat x \\ \hat y \\ \hat x^\prime \\ \hat y^\prime \end{bmatrix} = \begin{bmatrix} x \\ y \\ x^\prime \\ y^\prime \end{bmatrix} - \frac{x^{\prime T} Fx}{(F^Tx^\prime)_12^2+(F^Tx^\prime)_22^2+(Fx)_12^2+(Fx)_2^2} \begin{bmatrix} (F^Tx^\prime)_1 \\ (F^Tx^\prime)_2 \\ (Fx)_1 \\ (Fx)_2 \end{bmatrix} x^y^x^y^ = xyxy (FTx)122+(FTx)222+(Fx)122+(Fx)22xTFx (FTx)1(FTx)2(Fx)1(Fx)2
    针对其中最优的几何代价函数,可以通过转换为解一个 6 6 6次多项式的非迭代方法得到

最优解

最小化问题的重新公式化

给定一组测量得到的对应 x ↔ x ′ x \leftrightarrow x^\prime xx,我们要寻找 x ^ , x ^ ′ \hat x,\hat x^\prime x^,x^,它们最小化几何距离平方和并满足对极几何约束 x ^ ′ T F x ^ = 0 \hat x^{\prime T}F\hat x=0 x^TFx^=0。任意一对满足对极约束的点必须在两幅图像中对应的对极线上, x ⊥ , x ⊥ ′ x_\perp,x_\perp^\prime x,x分别为在 I , I ′ I,I^\prime I,I上最接近 x , x ′ x,x^\prime x,x的点,因此,几何代价函数变为
d ( x , I ) 2 + d ( x ′ , I ′ ) 2 d(x,I)^2 + d(x^\prime,I^\prime)^2 d(x,I)2+d(x,I)2
最小化的策略如下:

  1. 用参数 t t t参数化第一幅图像中的对极约束。因此第一幅图像中的对极线可以被记为 I ( t ) I(t) I(t)
  2. 利用基本矩阵 F F F,计算第二幅图像中对应的对极线 I ′ ( t ) I^\prime(t) I(t)
  3. 把距离函数 d ( x , I ( t ) ) 2 + d ( x ′ , I ′ ( t ) ) 2 d(x,I(t))^2+d(x^\prime,I^\prime(t))^2 d(x,I(t))2+d(x,I(t))2显示的表示为 t t t的函数
  4. 求这个函数最小化

外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传

最小化细节

假设两个图像点 x , x ′ x,x^\prime x,x都不是对极点。对每一幅图像都做刚性变换使得两个点 x , x ′ x,x^\prime x,x都是齐次坐标系的原点 ( 0 , 0 , 1 ) T (0,0,1)^T (0,0,1)T。此外,将对极点分别置于 x x x轴上的点 ( 1 , 0 , f ) T (1,0,f)^T (1,0,f)T ( 1 , 0 , f ′ ) T (1,0,f^\prime)^T (1,0,f)T。由于 F ( 1 , 0 , f ) T = ( 1 , 0 , f ′ ) T F = 0 F(1,0,f)^T = (1,0,f^\prime)^T F = 0 F(1,0,f)T=(1,0,f)TF=0,该基本矩阵有以下特殊形式
F = [ f f ′ d − f ′ c − f ′ d − f b a b   − f d c d ] F = \begin{bmatrix} ff^\prime d & -f^\prime c & -f^\prime d \\ -fb & a & b \\\ -fd & c & d \end{bmatrix} F= ffdfb fdfcacfdbd
考虑第一幅图像中过点 ( 0 , t , 1 ) T (0,t,1)^T (0,t,1)T和对极点 ( 1 , 0 , f ) T (1,0,f)^T (1,0,f)T的对极线,记这条极线为 I ( t ) I(t) I(t),这个对极线的向量表示为 ( 0 , t , 1 ) T × ( 1 , 0 , f ) T = ( t f , 1 , − t ) T (0,t,1)^T \times (1,0,f)^T = (tf,1,-t)^T (0,t,1)T×(1,0,f)T=(tf,1,t)T给出,所以由原点到该直线的距离的平方是
d ( x , I ( t ) ) 2 = t 2 1 + ( t f ) 2 d(x,I(t))^2 = \frac{t^2}{1+(tf)^2} d(x,I(t))2=1+(tf)2t2
利用基本矩阵计算另一幅图像的对极线,得到
I ′ ( t ) = F ( 0 , t , 1 ) T = ( − f ′ ( c t + d ) , a t + b , c t + d ) T I^\prime(t) = F(0,t,1)^T = (-f^\prime(ct+d),at+b,ct+d)^T I(t)=F(0,t,1)T=(f(ct+d),at+b,ct+d)T
从原点到该直线的距离的平方为
d ( x ′ , I ′ ( t ) ) 2 = ( c t + d ) 2 ( a t + b ) 2 + f ′ 2 ( c t + d ) 2 d(x^\prime,I^\prime(t))^2 = \frac{(ct+d)^2}{(at+b)^2 +f^{\prime 2}(ct+d)^2} d(x,I(t))2=(at+b)2+f′2(ct+d)2(ct+d)2
因此总的距离的平方是
s ( t ) = t 2 1 + ( t f ) 2 + ( c t + d ) 2 ( a t + b ) 2 + f ′ 2 ( c t + d ) 2 s(t) = \frac{t^2}{1+(tf)^2} + \frac{(ct+d)^2}{(at+b)^2 +f^{\prime 2}(ct+d)^2} s(t)=1+(tf)2t2+(at+b)2+f′2(ct+d)2(ct+d)2
我们的任务就是求这个式子的最小值,求其导数, s ( t ) s(t) s(t)的最大值和最小值将在 s ′ ( t ) = 0 s^\prime(t)=0 s(t)=0处出现。把 s ′ ( t ) s^\prime(t) s(t)中的两项在一个公分母上合并并且使分子等于 0 0 0,给出条件
g ( t ) = t ( ( a t + b ) 2 + f ′ 2 ( c t + d ) 2 ) 2 − ( a d − b c ) ( 1 + f 2 t 2 ) 2 ( a t + b ) ( c t + d ) = 0 g(t) = t((at+b)^2 + f^{\prime 2}(ct+d)^2)^2 - (ad-bc)(1+f^2t^2)^2 (at+b)(ct+d) = 0 g(t)=t((at+b)2+f′2(ct+d)2)2(adbc)(1+f2t2)2(at+b)(ct+d)=0
s ( t ) s(t) s(t)的最大值和最小值将在这个多项式的根处出现,这是一个 6 6 6次多项式,它最多有 6 6 6个实根,对应于 s ( t ) s(t) s(t) 3 3 3个极大值和 3 3 3个最小值

算法流程

目标

已知一组测量得到的点对应 x ↔ x ′ x \leftrightarrow x^\prime xx和基本矩阵 F F F,计算在对极几何约束 x ^ ′ F x ^ = 0 \hat x^\prime F \hat x = 0 x^Fx^=0下最小化几何误差的矫正对应 x ^ ↔ x ^ ′ \hat x \leftrightarrow \hat x^\prime x^x^

算法
  1. 定义变换矩阵
    T = [ 1 − x 1 − y 1 ]  and  T ′ = [ 1 − x ′ 1 − y ′ 1 ] T = \begin{bmatrix} 1 & & -x \\ & 1 & -y \\ & &1 \end{bmatrix} \text{ and } T^\prime = \begin{bmatrix} 1 & & -x^\prime \\ & 1 & -y^\prime \\ & &1 \end{bmatrix} T= 11xy1  and T= 11xy1
    这些变换是将 x = ( x , y , 1 ) T \pmb x=(x,y,1)^T x=(x,y,1)T x ′ = ( x ′ , y ′ , 1 ) T \pmb x^\prime=(x^\prime,y^\prime,1)^T x=(x,y,1)T变到坐标原点的平移

  2. T ′ − 1 F T − 1 T^{\prime -1} F T^{-1} T1FT1代替 F F F。新的 F F F对应于平移后的坐标

  3. 计算满足 e ′ T F = 0 , F e = 0 e^{\prime T}F=0,Fe=0 eTF=0,Fe=0的右、左对极点 e = ( e 1 , e 2 , e 3 ) T e=(e_1,e_2,e_3)^T e=(e1,e2,e3)T e ′ = ( e 1 ′ , e 2 ′ , e 3 ′ ) T e^\prime=(e_1^\prime,e_2^\prime,e_3^\prime)^T e=(e1,e2,e3)T。归一化 e e e使得 e 1 2 + e 2 2 + e 3 2 = 1 e_1^2+e_2^2+e_3^2=1 e12+e22+e32=1,并对 e ′ e^\prime e做同样的处理

  4. 构造矩阵
    R = [ e 1 e 2 − e 2 e 1 1 ]  and  R ′ = [ e 1 ′ e 2 ′ − e 2 ′ e 1 ′ 1 ] R=\begin{bmatrix} e_1 & e_2 &\\ -e_2 & e_1 & \\ && 1 \end{bmatrix} \text{ and } R^\prime = \begin{bmatrix} e_1^\prime & e_2^\prime &\\ -e_2^\prime & e_1^\prime & \\ && 1 \end{bmatrix} R= e1e2e2e11  and R= e1e2e2e11
    注意 R , R ′ R,R^\prime R,R是旋转矩阵,并满足 R e = ( 1 , 0 , e 3 ) T , R ′ e ′ = ( 1 , 0 , e 3 ′ ) T Re=(1,0,e_3)^T,R^\prime e^\prime = (1,0,e_3^\prime)^T Re=(1,0,e3)T,Re=(1,0,e3)T

  5. R ′ F R − 1 R^\prime FR^{-1} RFR1代替 F F F

  6. f = e 3 , f ′ = e 3 ′ , a = F 22 , b = F 23 , c = F 32 , d = F 33 f=e_3,f^\prime =e_3^\prime,a=F_{22},b=F_{23},c=F_{32},d=F_{33} f=e3,f=e3,a=F22,b=F23,c=F32,d=F33

  7. 根据上述分析,构造关于 t t t的多项式 g ( t ) g(t) g(t),得到 t t t 6 6 6个根

  8. 计算代价函数在 s ( t ) s(t) s(t) g ( t ) g(t) g(t)的每个根的实部的值。同时求当 t = ∞ t= \infty t=时几何代价函数的渐近线,选择使代价函数取最小值的 t m i n t_{min} tmin

  9. 计算两条直线 I = ( t f , 1 , − t ) T I=(tf,1,-t)^T I=(tf,1,t)T I ′ ( t ) I^\prime(t) I(t) t m i n t_{min} tmin处的值,并求出这些直线上最接近于原点的点 x ^ , x ^ ′ \hat x,\hat x^\prime x^,x^,对于一般的直线 ( λ , μ , v ) T (\lambda , \mu,v)^T (λ,μ,v)T,直线上最接近原点的点是 ( − λ v , − μ v , λ 2 + μ 2 ) T (-\lambda v,-\mu v,\lambda^2+\mu^2)^T (λv,μv,λ2+μ2)T

  10. T − 1 R T x ^ T^{-1}R^T\hat x T1RTx^替代 x ^ \hat x x^并用 T ′ − 1 R ′ T x ^ ′ T^{\prime -1}R^{\prime T}\hat x^\prime T1RTx^替代 x ^ ′ \hat x^\prime x^变回原来的坐标

  11. 3 3 3维空间点 X ^ \hat X X^由齐次方法线性三角测量法得到

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

愤怒的卤蛋

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值