一文搞懂各种形式的单应性矩阵表示
单应性矩阵 H H H表征一些落在同一平面上的一些点在两平面间的映射关系,英文名字叫Homography. 单应性变换在计算机视觉中多视图几何中特别常见,一般应用为(1)已知相机内外参数的情况下,计算两视图之间像素点像素坐标的变换关系 H H H,然后计算坐标变换;(2)或已知两视图之间部分像素点的匹配关系(如通过SIFT、ORB做好了特征点匹配),通过构建多个 x ′ = H x x'=Hx x′=Hx的线性方程解出单应矩阵 H H H,从而恢复出相机的姿态。
最近阅读论文一般为(1)方向的应用,即在多视图稠密立体匹配上用得较多。现在见到以下几种计算 H H H的表达式,分别为:
- H = K 2 ( R − t n T d ) K 1 − 1 H = K_2(R-\frac{tn^{T}}{d})K_1^{-1} H=K2(R−dtnT)K1−1
- H = K 2 R 2 [ R 1 T − ( t 1 − t 2 ) n T d ] K 1 − 1 H = K_2R_2[R_1^T-\frac{(t_1-t_2)n^T}{d}]K_1^{-1} H=K2R2[R1T−d(t1−t2)nT]K1−1
- H = K 2 [ R 2 R 1 T + R 2 ( t 1 − t 2 ) n T n T P 1 ] K 1 − 1 H = K_2[R_2R_1^T+\frac{R_2(t_1-t_2)n^T}{n^TP_1}]K_1^{-1} H=K2[R2R1T+nTP1R2(t1−t2)nT]K1−1
- H = K 2 R 2 [ I − ( t 1 − t 2 ) n T d ] R 1 T K 1 − 1 H = K_2R_2[I-\frac{(t_1-t_2)n^T}{d}]R_1^TK_1^{-1} H=K2R2[I−d(t1−t2)nT]R1TK1−1
- H = d K 2 T 2 T 1 − 1 K 1 − 1 H = dK_2T_2T_1^{-1}K_1^{-1} H=dK2T2T1−1K1−1
- …
这么多的表达式足够让人眼花撩乱,眼冒金星了。不着急,一个一个解释他们为啥长得不一样,虽然他们都是Homography Matrix这个东西。
1. H = K 2 ( R − t n T d ) K 1 − 1 H = K_2(R-\frac{tn^{T}}{d})K_1^{-1} H=K2(R−dtnT)K1−1
变量解释:H表示 p 2 = H p 1 p_2=Hp_1 p2=Hp1的变换, K 1 , K 2 K_1,K_2 K1,K2分别为相机1,2的内参, R , t R,t R,t表示相机坐标系{Camera1}到{Camera2}的旋转和平移.
我们已知相机1和相机2有一对匹配点记作匹配对(
p
1
,
p
2
p_1,p_2
p1,p2),它们在真实场景中落在同一平面
π
(
n
,
P
1
)
\pi(n, P_1)
π(n,P1)上,其中
n
n
n为此平面的单位法向量(对于fronto-Parallel Plane的情况,它便是相机主光轴),
P
1
P_1
P1是该平面上的一点(该点为相机1成像平面上点
p
1
p_1
p1对应的相机坐标),且该平面是定义在相机坐标系{Camera1}下的。该平面的平面方程写作
n
T
P
1
=
−
d
n^TP_1=-d
nTP1=−d
d为
P
1
P_1
P1到光心的距离在z轴投影。已知
p
2
=
K
1
(
R
P
1
+
t
)
p_2 = K_1(RP_1+t)
p2=K1(RP1+t)
该等号认为左右尺度意义下等价,或者认为
p
2
=
[
d
2
u
,
d
2
v
,
d
2
]
T
p_2=[d_2u,d_2v, d_2]^T
p2=[d2u,d2v,d2]T.
接着往下推,我们有
p
2
=
K
2
(
R
P
1
+
t
(
−
n
T
d
)
P
1
)
=
K
2
(
R
−
t
n
T
d
)
P
1
=
K
2
(
R
−
t
n
T
d
)
K
1
−
1
p
1
\begin{aligned} p_2 &= K_2(RP_1+t(-\frac{n^T}{d})P_1) \\ &= K_2(R-\frac{tn^T}{d})P_1 \\ &= K_2(R-\frac{tn^T}{d})K_1^{-1}p_1 \end{aligned}
p2=K2(RP1+t(−dnT)P1)=K2(R−dtnT)P1=K2(R−dtnT)K1−1p1
其中认为
p
1
p_1
p1中也包含了深度
d
d
d.
2. H = K 2 R 2 [ R 1 T − ( t 1 − t 2 ) n T d ] K 1 − 1 H = K_2R_2[R_1^T-\frac{(t_1-t_2)n^T}{d}]K_1^{-1} H=K2R2[R1T−d(t1−t2)nT]K1−1
变量解释:H表示 p 2 = H p 1 p_2=Hp_1 p2=Hp1的变换, K 1 , K 2 K_1,K_2 K1,K2分别为相机1,2的内参, R i , t i R_i, t_i Ri,ti( i = 1 , 2 i=1,2 i=1,2)分别表示{World}到{Camera_i}的旋转,以及相机i光心在{World}下的坐标
由1中推导的公式:
H
=
K
2
(
R
−
t
n
T
d
)
K
1
−
1
H = K_2(R-\frac{tn^T}{d})K_1^{-1}
H=K2(R−dtnT)K1−1,我们可以知道其中的旋转和平移按照2中的声明可以表示为
R
=
R
2
R
1
T
R=R_2R_1^T
R=R2R1T
t
=
R
2
(
t
1
−
t
2
)
t=R_2(t_1-t_2)
t=R2(t1−t2)
上述
R
R
R表示相机坐标系{Camera1}到{Camera2}的旋转,同时也用{Camera1}->{World},{World}->{Camera2}两步来表示,与之前不同的是坐标变换使用了中间的一个坐标系{World}.除此之外,平移
t
t
t比较难理解,其实不难。1中的平移t是从齐次变换矩阵
T
T
T中取出来的,它实质上表示{Camera1}到{Camera2}的平移变换,或者表示{Camera1}的原点在{Camera2}中的投影坐标(往哪里投影就相当于往哪个坐标变换)。这样一来,
t
t
t就相当于一个由{Camera2}原点指向{Camera1}原点的向量在坐标系{Camera2}下的表示。既然
t
t
t可以这么解释那就好办了!{Camera2}原点指向{Camera1}原点的向量可以用
t
1
−
t
2
t_1-t_2
t1−t2表示,但是这个表示是在{World}下的,为了统一到{Camera2}下的表示,需要再左乘一个旋转
R
2
R_2
R2即可。
将上述两个式子代入到1中的单应矩阵表示:
H
=
K
2
(
R
−
t
n
T
d
)
K
1
−
1
=
K
2
[
R
2
R
1
T
−
R
2
(
t
1
−
t
2
)
n
T
d
]
K
1
−
1
=
K
2
R
2
(
R
1
T
−
(
t
1
−
t
2
)
n
T
d
)
K
1
−
1
\begin{aligned} H &= K_2(R-\frac{tn^T}{d})K_1^{-1} \\ &= K_2[R_2R_1^T-\frac{R_2(t_1-t_2)n^T}{d}]K_1^{-1} \\ &= K_2R_2(R_1^T-\frac{(t_1-t_2)n^T}{d})K_1^{-1} \end{aligned}
H=K2(R−dtnT)K1−1=K2[R2R1T−dR2(t1−t2)nT]K1−1=K2R2(R1T−d(t1−t2)nT)K1−1
3. H = K 2 [ R 2 R 1 T + R 2 ( t 1 − t 2 ) n T n T P 1 ] K 1 − 1 H = K_2[R_2R_1^T+\frac{R_2(t_1-t_2)n^T}{n^TP_1}]K_1^{-1} H=K2[R2R1T+nTP1R2(t1−t2)nT]K1−1
式子3和式子2基本是统一的,即把式子2左侧的
R
2
R_2
R2乘到括号里头,然后把深度
d
d
d用平面公式替代即可:
H
=
K
2
(
R
−
t
n
T
d
)
K
1
−
1
=
K
2
[
R
2
R
1
T
−
R
2
(
t
1
−
t
2
)
n
T
d
]
K
1
−
1
=
K
2
[
R
2
R
1
T
+
R
2
(
t
1
−
t
2
)
n
T
n
T
P
1
]
K
1
−
1
\begin{aligned} H &= K_2(R-\frac{tn^T}{d})K_1^{-1} \\ &= K_2[R_2R_1^T-\frac{R_2(t_1-t_2)n^T}{d}]K_1^{-1} \\ &= K_2[R_2R_1^T+\frac{R_2(t_1-t_2)n^T}{n^TP_1}]K_1^{-1} \end{aligned}
H=K2(R−dtnT)K1−1=K2[R2R1T−dR2(t1−t2)nT]K1−1=K2[R2R1T+nTP1R2(t1−t2)nT]K1−1
4. H = K 2 R 2 [ I − ( t 1 − t 2 ) n T d ] R 1 T K 1 − 1 H = K_2R_2[I-\frac{(t_1-t_2)n^T}{d}]R_1^TK_1^{-1} H=K2R2[I−d(t1−t2)nT]R1TK1−1
这个式子和式子2、3如出一辙。唯一不同的是变量的表示。式子1,2,3中的平面
π
\pi
π的单位法向量是表示在相机坐标系{Camera1}下的,而式子4的单位法向量是表示在世界坐标系{World}下的。故平面方程应当变换为
(
R
1
n
)
T
P
1
=
−
d
(R_1n)^TP_1=-d
(R1n)TP1=−d
所以
H
=
K
2
(
R
−
t
n
T
R
1
T
d
)
K
1
−
1
=
K
2
[
R
2
R
1
T
−
R
2
(
t
1
−
t
2
)
R
1
T
n
T
d
]
K
1
−
1
=
K
2
[
R
2
R
1
T
−
R
2
(
t
1
−
t
2
)
R
1
T
n
T
d
]
K
1
−
1
=
K
2
R
2
[
I
−
(
t
1
−
t
2
)
n
T
d
]
R
1
T
K
1
−
1
\begin{aligned} H &= K_2(R-\frac{tn^TR_1^T}{d})K_1^{-1} \\ &= K_2[R_2R_1^T-\frac{R_2(t_1-t_2)R_1^Tn^T}{d}]K_1^{-1} \\ &= K_2[R_2R_1^T-\frac{R_2(t_1-t_2)R_1^Tn^T}{d}]K_1^{-1} \\ &= K_2R_2[I-\frac{(t_1-t_2)n^T}{d}]R_1^TK_1^{-1} \end{aligned}
H=K2(R−dtnTR1T)K1−1=K2[R2R1T−dR2(t1−t2)R1TnT]K1−1=K2[R2R1T−dR2(t1−t2)R1TnT]K1−1=K2R2[I−d(t1−t2)nT]R1TK1−1
5. H = d K 2 T 2 T 1 − 1 K 1 − 1 H = dK_2T_2T_1^{-1}K_1^{-1} H=dK2T2T1−1K1−1
变量解释:H表示 p 2 = H p 1 p_2=Hp_1 p2=Hp1的变换, K 1 , K 2 , T 1 , T 2 K_1,K_2,T_1,T_2 K1,K2,T1,T2分别为相机1,2的内参和外参({Wolrd}->{Camera})
因此,
p
2
=
d
K
2
T
2
T
1
−
1
K
1
−
1
p
1
=
K
2
T
2
(
T
1
−
1
K
1
−
1
d
p
1
)
=
K
2
T
2
P
1
,
W
\begin{aligned} p_2 &= dK_2T_2T_1^{-1}K_1^{-1}p_1 \\ &=K_2T_2(T_1^{-1}K_1^{-1}dp_1) \\ &=K_2T_2P_{1,W} \end{aligned}
p2=dK2T2T1−1K1−1p1=K2T2(T1−1K1−1dp1)=K2T2P1,W
显然左边等于右边,唯一值得注意的是,坐标
p
2
p_2
p2为三维向量,隐含了像素坐标及其在对应相机下的深度。而
p
1
p_1
p1虽然也为三维,但第三个维度为1,仅仅表示像素坐标的其次表示形式,并且显示地表示出了在相机1下的深度值
d
d
d.另外,该单应矩阵的计算式中不显示包含点对(
p
1
,
p
2
p_1,p_2
p1,p2)所在的空间平面的信息。此空间平面实际上为平行于相机frustum的平面,距离光心为
d
d
d,法向量为Priciple Axis. 因此,该单应性变换表达式仅能表达对应两视图匹配点对在深度值已知的情况下的匹配关系,或者表达匹配点在各自视图下的具有相同深度或在同一个Fronto-Parallel Plane的领域像素点的匹配关系。那就意味着,如果邻域像素点与匹配点不在Fronto-Parallel平面上则无法通过此单应变换式计算和表示。