目录
1 导论
计算机视觉系列教程1-1：透视空间与透视变换中提到，透视空间所有变换都是投影变换的特例，本节进一步研究投影变换矩阵(单应性矩阵)的估计。
透视变换的核心是单应性矩阵 H H H或单应性向量 h h h。
H = [ h 11 h 12 h 13 h 21 h 22 h 23 h 31 h 32 h 33 ] ⇔ h = [ h 11 h 12 h 13 h 21 h 22 h 23 h 31 h 32 h 33 ] T H=\left[
h11h21h31h12h22h32h13h23h33
h11h12h13h21h22h23h31h32h33
\right] \Leftrightarrow h=\left[
h11h12h13h21h22h23h31h32h33
h11h12h13h21h22h23h31h32h33
\right] ^T
H=
⎣
⎡
h
11
h
21
h
31
h
12
h
22
h
32
h
13
h
23
h
33
⎦
⎤
⇔h=[
h
11
h
12
h
13
h
21
h
22
h
23
h
31
h
32
h
33
]
T
设 p s r c = [ x y 1 ] T p_{src}=\left[
xy1
xy1
\right] ^T p
src
=[
x
y
1
]
T
与 p d s t = [ x ′ y ′ 1 ] T p_{dst}=\left[
x′y′1
x′y′1
\right] ^T p
dst
=[
x
′
y
′
1
]
T
是二维透视空间 P 2 \mathbb{P}^2 P
2
中，一次透视变换前后的对应点，因此其满足
p d s t = H p s r c ⟺ [ x ′ y ′ 1 ] = [ h 11 h 12 h 13 h 21 h 22 h 23 h 31 h 32 h 33 ] [ x y 1 ] p_{dst}=Hp_{src}\Longleftrightarrow \left[
x′y′1
x′y′1
\right] =\left[
h11h21h31h12h22h32h13h23h33
h11h12h13h21h22h23h31h32h33
\right] \left[
xy1
xy1
\right]
p
dst
=Hp
src
⟺
⎣
⎡
x
′
y
′
1
⎦
⎤
=
⎣
⎡
h
11
h
21
h
31
h
12
h
22
h
32
h
13
h
23
h
33
⎦
⎤
⎣
⎡
x
y
1
⎦
⎤
若将单应性矩阵进行尺度缩放后作用于 p s r c p_{src} p
src
，则
k H p s r c = k [ h 11 h 12 h 13 h 21 h 22 h 23 h 31 h 32 h 33 ] [ x y 1 ] = k p d s t kHp_{src}=k\left[
h11h21h31h12h22h32h13h23h33
h11h12h13h21h22h23h31h32h33
\right] \left[
xy1
xy1
\right] =kp_{dst}
kHp
src
=k
⎣
⎡
h
11
h
21
h
31
h
12
h
22
h
32
h
13
h
23
h
33
⎦
⎤
⎣
⎡
x
y
1
⎦
⎤
=kp
dst
而透视空间中， k p d s t kp_{dst} kp
dst
与 p d s t p_{dst} p
dst
实际上对应同一点，因此 k H kH kH与 H H H相当于同一次透视变换，故单应性矩阵 H H H仅有8个自由度，通常通过设置 h 33 = 1 h_{33}=1 h
33
=1或 ∥ h ∥ 2 2 = 1 \lVert h \rVert _{2}^{2}=1 ∥h∥
2
2
=1来约束冗余的参数。
下面详细阐述单应性矩阵估计方法。
2 基本直接线性变换(Basic DLT)
将上式改写为齐次形式
[ 0 0 0 − x − y − 1 y ′ x y ′ y y ′ x y 1 0 0 0 − x ′ x − x ′ y − x ′ − y ′ x − y ′ y − y ′ x ′ x x ′ y x ′ 0 0 0 ] [ h 11 h 12 h 13 h 21 h 22 h 23 h 31 h 32 h 33 ] = [ 0 0 0 ] \left[
0x−y′x0y−y′y01−y′−x0x′x−y0x′y−10x′y′x−x′x0y′y−x′y0y′−x′0
000−x−y−1y′xy′yy′xy1000−x′x−x′y−x′−y′x−y′y−y′x′xx′yx′000
\right] \left[
h11h12h13h21h22h23h31h32h33
h11h12h13h21h22h23h31h32h33
\right] =\left[
000
000
\right]
⎣
⎡
0
x
−y
′
x
0
y
−y
′
y
0
1
−y
′
−x
0
x
′
x
−y
0
x
′
y
−1
0
x
′
y
′
x
−x
′
x
0
y
′
y
−x
′
y
0
y
′
−x
′
0
⎦
⎤
⎣
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎡
h
11
h
12
h
13
h
21
h
22
h
23
h
31
h
32
h
33
⎦
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎤
=
⎣
⎡
0
0
0
⎦
⎤
其中系数矩阵的秩为2，因此一对变换点仅能确定2个自由度。因此需要无三点共线的四对变换点才能确定单应性矩阵 H H H。
[ 0 0 0 − x 1 − y 1 − 1 y 1 ′ x 1 y 1 ′ y 1 y 1 ′ x 1 y 1 1 0 0 0 − x 1 ′ x 1 − x 1 ′ y 1 − x 1 ′ 0 0 0 − x 2 − y 2 − 1 y 2 ′ x 2 y 2 ′ y 2 y 2 ′ x 2 y 2 1 0 0 0 − x 2 ′ x 2 − x 2 ′ y 2 − x 2 ′ ⋮ ] [ h 11 h 12 h 13 h 21 h 22 h 23 h 31 h 32 h 33 ] = [ 0 0 0 ] ⇔ A h = 0 \left[
0x10x20y10y20101−x10−x20−y10−y20⋮−10−10y′1x1−x′1x1y′2x2−x′2x2y′1y1−x′1y1y′2y2−x′2y2y′1−x′1y′2−x′2
000−x1−y1−1y1′x1y1′y1y1′x1y11000−x1′x1−x1′y1−x1′000−x2−y2−1y2′x2y2′y2y2′x2y21000−x2′x2−x2′y2−x2′⋮
\right] \left[
h11h12h13h21h22h23h31h32h33
h11h12h13h21h22h23h31h32h33
\right] =\left[
000
000
\right] \Leftrightarrow {Ah=0}
⎣
⎢
⎢
⎢
⎢
⎢
⎡
0
x
1
0
x
2
0
y
1
0
y
2
0
1
0
1
−x
1
0
−x
2
0
−y
1
0
−y
2
0
⋮
−1
0
−1
0
y
1
′
x
1
−x
1
′
x
1
y
2
′
x
2
−x
2
′
x
2
y
1
′
y
1
−x
1
′
y
1
y
2
′
y
2
−x
2
′
y
2
y
1
′
−x
1
′
y
2
′
−x
2
′
⎦
⎥
⎥
⎥
⎥
⎥
⎤
⎣
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎡
h
11
h
12
h
13
h
21
h
22
h
23
h
31
h
32
h
33
⎦
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎤
=
⎣
⎡
0
0
0
⎦
⎤
⇔Ah=0
对于形如 A x = b Ax=b Ax=b的非齐次线性方程可以通过伪逆形式计算 x x x；对于形如 A x = 0 Ax=0 Ax=0的齐次线性方程的解则对应矩阵 A A A右奇异向量 v p v_p v
p
，其中 v p v_p v
p
对应的奇异值 σ p ≈ 0 \sigma _p\approx 0 σ
p
≈0或不对应奇异值。
3 归一化直接线性变换(Normalized DLT)
基本DLT估计算法的缺陷是：
(1) 单应性估计不具有相似不变性
假设在第一次估计下有 x d s t = H x s r c x_{dst}=Hx_{src} x
dst
=Hx
src
。现对两张图像分别进行相似变换并重新进行单应性估计，得到 ( T d s t x d s t ) = H ′ ( T s r c x s r c ) \left( T_{dst}x_{dst} \right) =H'\left( T_{src}x_{src} \right) (T
dst
x
dst
)=H
′
(T
src
x
src
)，改写为 x d s t = ( T d s t − 1 H ′ T s r c ) x s r c x_{dst}=\left( T_{dst}^{-1}H'T_{src} \right) x_{src} x
dst
=(T
dst
−1
H
′
T
src
)x
src
，大部分情况下 H ≠ T d s t − 1 H ′ T s r c H\ne T_{dst}^{-1}H'T_{src} H
=T
dst
−1
H
′
T
src
，这表明基本DLT算法无法抵抗相似变换的干扰。
(2) 估计的单应性矩阵容易产生病态条件，鲁棒性差
由于默认透视空间的尺度变换因子 w = 1 w=1 w=1，所以齐次坐标下很可能产生分量幅度差异大的情况，例如某特征点 X = [ 100 101 1 ] T X=\left[
1001011
1001011
\right] ^T X=[
100
101
1
]
T
。在这种情况下估计出的单应性矩阵，各个元素数值数量级可能会相差 1 0 4 10^4 10
4
以上，导致病态条件——特征点的轻微变化都会造成单应性矩阵的剧变。
基于(1)(2)两种缺陷，需要将基本DLT算法进行优化,优化的核心就是特征点坐标的归一化。设原图像特征点集合为 ，目标图像特征点集合为 ，则具体的算法为：
① 将特征点集合 X s r c X_{src} X
src
、 X d s t X_{dst} X
dst
归一化
使用相似变换矩阵 T s r c T_{src} T
src
、 T d s t T_{dst} T
dst
将特征点集合中心移至原点，且与原点平均距离为 2 \sqrt{2}
2
。由于默认尺度因子为 w = 1 w=1 w=1，所以归一化到 2 \sqrt{2}
2
可以保持齐次坐标的三个分量有相同的幅度，例如 X = [ 100 100 1 ] T ⇒ X n o r m a l = [ 1 1 1 ] T X=\left[
1001001
1001001
\right] ^T\Rightarrow X^{normal}=\left[
111
111
\right] ^T X=[
100
100
1
]
T
⇒X
normal
=[
1
1
1
]
T
。
这里给出一种相似变换矩阵的计算方式，设
X n o r m a l = T X ⇔ [ x ~ i y ~ i 1 ] = [ s t x s t y 1 ] [ x i y i 1 ] X^{normal}=TX\Leftrightarrow \left[
x~iy~i1
x~iy~i1
\right] =\left[
sstxty1
stxsty1
\right] \left[
xiyi1
xiyi1
\right]
X
normal
=TX⇔
⎣
⎡
x
~
i
y
~
i
1
⎦
⎤
=
⎣
⎡
s
s
t
x
t
y
1
⎦
⎤
⎣
⎡
x
i
y
i
1
⎦
⎤
令
{ 1 N ∑ i = 1 N x ~ i = 1 N ∑ i = 1 N ( s x i + t x ) = 0 1 N ∑ i = 1 N y ~ i = 1 N ∑ i = 1 N ( s y i + t y ) = 0 1 N ∑ i = 1 N x ~ i 2 + y ~ i 2 = 2
⎧⎩⎨⎪⎪⎪⎪⎪⎪1N∑Ni=1x~i=1N∑Ni=1(sxi+tx)=01N∑Ni=1y~i=1N∑Ni=1(syi+ty)=01N∑Ni=1x~2i+y~2i−−−−−−√=2–√
{1N∑i=1Nx~i=1N∑i=1N(sxi+tx)=01N∑i=1Ny~i=1N∑i=1N(syi+ty)=01N∑i=1Nx~i2+y~i2=2
⎩
⎪
⎨
⎪
⎧
N
1
∑
i=1
N
x
~
i
=
N
1
∑
i=1
N
(sx
i
+t
x
)=0
N
1
∑
i=1
N
y
~
i
=
N
1
∑
i=1
N
(sy
i
+t
y
)=0
N
1
∑
i=1
N
x
~
i
2
+
y
~
i
2
=
2
解得
{ t x = − s 1 N ∑ i = 1 N x i = − s x ˉ t y = − s 1 N ∑ i = 1 N y i = − s y ˉ s = 2 1 N ∑ i = 1 N x ~ i 2 + y ~ i 2 = 2 1 N ∑ i = 1 N ( x i − x ˉ ) 2 + ( y i − y ˉ ) 2
⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪tx=−s1N∑Ni=1xi=−sx¯ty=−s1N∑Ni=1yi=−sy¯s=2√1N∑Ni=1x~2i+y~2i√=2√1N∑Ni=1(xi−x¯)2+(yi−y¯)2√
{tx=−s1N∑i=1Nxi=−sx¯ty=−s1N∑i=1Nyi=−sy¯s=21N∑i=1Nx~i2+y~i2=21N∑i=1N(xi−x¯)2+(yi−y¯)2
⎩
⎪
⎪
⎨
⎪
⎪
⎧
t
x
=−s
N
1
∑
i=1
N
x
i
=−s
x
ˉ
t
y
=−s
N
1
∑
i=1
N
y
i
=−s
y
ˉ
s=
N
1
∑
i=1
N
x
~
i
2
+
y
~
i
2
2
=
N
1
∑
i=1
N
(x
i
−
x
ˉ
)
2
+(y
i
−
y
ˉ
)
2
2
② 运用基本DLT算法由 与 估计单应性矩阵
③ 解归一化，映射回实际图像
4 鲁棒单应性估计(Robust Homography Estimation)
结合基本DLT算法、归一化DLT算法及计算机视觉系列教程6-2：图像匹配算法概述.的RANSAC算法进行单应性矩阵估计，具体流程如下：
(1) 设置迭代次数 K = ∞ K=\infty K=∞，内点集 S i n = ⊘ S_{in}=\oslash S
in
=⊘，模型参数 H = H 0 H=H_0 H=H
0
；
(2) 随机从样本数据集 S S S中选取4对特征点，并通过基本DLT算法确定测试模型 H t e s t H_{test} H
test
；
(3) 用 H t e s t H_{test} H
test
遍历样本数据集 S S S，估计误差 ε \varepsilon ε在距离阈值 t t t内的点加入内点集 S i n S_{in} S
in
。其中阈值 t = 5.99 σ t=\sqrt{5.99}\sigma t=
5.99
σ， σ \sigma σ为估计不确定度，估计误差 ε \varepsilon ε主要有两种度量方式：
① 代数误差 ε i = ∥ A i h t e s t ∥ \varepsilon _i=\left\| A_ih_{test} \right\| ε
i
=∥A
i
h
test
∥，其中
A i = [ 0 0 0 − x − y − 1 y ′ x y ′ y y ′ − y ′ x − y ′ y − y ′ x ′ x x ′ y x ′ 0 0 0 ] A_i=\left[
0−y′x0−y′y0−y′−xx′x−yx′y−1x′y′x0y′y0y′0
000−x−y−1y′xy′yy′−y′x−y′y−y′x′xx′yx′000
\right]
A
i
=[
0
−y
′
x
0
−y
′
y
0
−y
′
−x
x
′
x
−y
x
′
y
−1
x
′
y
′
x
0
y
′
y
0
y
′
0
]
② 几何误差(二次投影误差) ε i = ∥ H X s r c , i , X d s t , i ∥ 2 2 + ∥ X s r c , i , H − 1 X d s t , i ∥ 2 2 \varepsilon _i=\left\| HX_{src,i}, X_{dst, i} \right\| _{2}^{2}+\left\| X_{src,i}, H^{-1}X_{dst, i} \right\| _{2}^{2} ε
i
=∥HX
src,i
,X
dst,i
∥
2
2
+
∥
∥
X
src,i
,H
−1
X
dst,i
∥
∥
2
2
，可以视作交叉检验。
(4) 若 S i n S_{in} S
in
的大小小于阈值 T T T，则放弃该模型，重复(2)；若 S i n S_{in} S
in
的大小大于阈值 t t t，则通过归一化DLT算法或Levenberg Marquardt等迭代优化算法，利用 S i n S_{in} S
in
中的所有点重新估计模型 H t e s t ∗ H_{test}^{*} H
test
∗
；
(5) 计算当前模型 H t e s t ∗ H_{test}^{*} H
test
∗
下的内点率 ω = ∣ S i n ∣ ∣ S ∣ \omega =\frac{|S_{in}|}{|S|} ω=
∣S∣
∣S
in
∣
，根据 K = ln ( 1 − p ) ln ( 1 − ω n ) K=\frac{\ln \left( 1-p \right)}{\ln \left( 1-\omega ^n \right)} K=
ln(1−ω
n
)
ln(1−p)
更新迭代次数；
(6) 至此完成一次迭代，若 H t e s t ∗ H_{test}^{*} H
test
∗
下内点率为最大，则更新 H = H t e s t ∗ H=H_{test}^{*} H=H
test
∗
，重复(2) ~ (5)直至迭代次数满足要求。
🚀 计算机视觉基础教程说明
🔥 更多精彩专栏：
