论文阅读整理笔记
一、StructVIO : Visual-inertial Odometry with Structural Regularity of Man-made Environments
I. 摘要
- 基于Atlanta World模型对人工环境规则性进行描述。
-
- Atlanta World 定义为包含有几个不同方向的局部Manhattan World世界。
-
- 每个局部Manhattan World可以被实时检测,且它们的方向会在新观测到来时由状态估计器不断进行修正。
介绍
- 检测到的局部Manhattan世界无需是个真实的箱型世界,允许结构线只对准到同一个水平方向上。
- 对准到局部Manhattan世界的水平线可以提供航向约束,保证航向误差不会在局部区域增加。
- 即使没有检测到局部的Manhattan世界,垂直线仍然可以指明重力方向,使横滚角和俯仰角可观。
- 主要贡献
-
- 提出了一个新颖的线特征最小参数化表示法,融合了Atlanta World假设和结构化特征。每个局部Manhattan世界用一个航向方向来表示,并且会被状态估计器是为一个状态变量不断修正。
-
- 结构线(和主方向一致)和局部Manhattan World可以实时在线检测。即使没有检测到Manhattan World,垂直线(vertical line)依旧可以提供位姿信息观测。
-
- 对状态估计器和线跟踪方法进行了改进。提出的信息积累方法能够处理长期特征跟踪的观测丢失,保证更好的特征估计;提出的多采样点的线特征跟踪方法和延迟EKF更新让线特征跟踪更加可靠。
II. 相关工作
- 绝大多数VI-SLAM以点特征为主的原因:1. 点特征广泛存在。2. 点特征易于检测和跟踪。3. 线特征的自由度为4,点特征只有3,因此线特征更难被初始化(尤其是其方向),更难被估计(通常需要6个参数来表示,如普吕克坐标表示法)。在SLAM系统中额外添加线特征的应用,偶尔会带来更糟糕的表现。解决方法有:使用立体相机或多等几个关键帧来进行延迟初始化
- 有许多提取灭点来改善VI- SLAM性能的工作,但这些研究中的线特征仅仅是灭点提取的过渡结果,提取完灭点后先特征即被抛弃了。因此,若是像点特征一样,把线特征加进估计器里,理应能提升性能。
- 大多数现有方法只利用了结构特征的部分信息——要么是没有考虑它们的先验方向,要么是仅使用了先验方向而没有作为额外的观测来改善估计结果。
- 本文所提出的方法不需要在初始化过程中提取任何灭点。而且,新的线参数化方法把每条线特征锚定在了一个局部的manhattan世界上,从而减少了线参数的一个自由度,沿着Manhattan世界方向的线特征方向也可以被修正,作为更可靠的观测被收集。
III. 结构线和Atlanta World
Atlanta World是多个具备不同水平方向
ϕ
i
∈
[
0
,
π
/
2
)
\phi_i \in [0,\pi/2)
ϕi∈[0,π/2)的Manhattan Worlds的组合,如图所示。
我们在里程计起点建立
Z
Z
Z 轴朝上(重力方向相反)的全局世界坐标系
W
{W}
W 。
每条线被锚定在第一次被观测到的图像局部坐标系中,称这个坐标系统为“起始帧” S {S} S,其方向与这条直线所属的局部Manhattan World的 ϕ i \phi_i ϕi 相同,原点设置为该帧相机的位置。该起始帧原点位置 W p s ^Wps Wps, 会在滤波器中继续进行优化估计。
对于锚定在起始帧
S
S
S 上的一条给定直线,我们可以找到一个旋转
L
S
R
^S_LR
LSR 来将这条线从参数空间
L
L
L 转换到
S
S
S 。该条直线被对齐到
L
L
L 的
Z
Z
Z 轴上,如下图所示。
在参数空间
L
L
L 中,该条线可以用其与
X
Y
XY
XY平面的交点来表示,即
L
l
p
=
(
a
,
b
,
0
)
T
^Ll_p=(a,b,0)^T
Llp=(a,b,0)T。这里我们使用逆深度重新表示该条直线,即
(
θ
,
ρ
,
0
)
T
(\theta,\rho,0)^T
(θ,ρ,0)T, 其中
ρ
=
1
/
a
2
+
b
2
\rho=1/\sqrt{a^2+b^2}
ρ=1/a2+b2 和
θ
=
atan2
(
b
,
a
)
\theta=\text{atan2}(b,a)
θ=atan2(b,a)。
该条线段在起始帧中可以通过以下公式计算:
L
S
R
L
l
p
=
a
L
S
R
(
:
,
1
)
+
b
L
S
R
(
:
,
2
)
=
cos
θ
ρ
L
S
R
(
:
,
1
)
+
sin
θ
ρ
L
S
R
(
:
,
2
)
\begin{aligned} ^S_LR ^Ll_p &=a{^S_LR(:,1)}+b^S_LR(:,2) \\ &= {\text{cos}\theta \over \rho}{^S_LR(:,1)}+{\text{sin}\theta \over \rho}{^S_LR(:,2)} \end{aligned}
LSRLlp=aLSR(:,1)+bLSR(:,2)=ρcosθLSR(:,1)+ρsinθLSR(:,2)
对于每个结构线,无论对齐到局部Manhattan World的哪一个轴上,旋转
L
S
R
^S_LR
LSR 一定是下面三个矩阵之一:
[ 0 0 1 0 1 0 − 1 0 0 ] \begin{bmatrix}0 & 0 & 1\\ 0 &1 & 0\\ -1 & 0 & 0 \end{bmatrix} ⎣ ⎡00−1010100⎦ ⎤, [ 1 0 0 0 0 1 0 − 1 0 ] \begin{bmatrix}1 & 0 & 0\\ 0 &0 & 1\\ 0 & -1 & 0 \end{bmatrix} ⎣ ⎡10000−1010⎦ ⎤ 和 [ 1 0 0 0 1 0 0 0 1 ] \begin{bmatrix}1 & 0 & 0\\ 0 &1 & 0\\ 0 & 0 & 1 \end{bmatrix} ⎣ ⎡100010001⎦ ⎤,分别是对齐到 X X X, Y Y Y, Z Z Z轴。
从起始帧
S
S
S 到世界坐标系下的变换由旋转
S
W
R
(
ϕ
i
)
^W_SR(\phi_i)
SWR(ϕi) 和 相机位置
W
P
s
^WPs
WPs决定,其中
S
W
R
(
ϕ
i
)
^W_SR(\phi_i)
SWR(ϕi)是绕重力方向轴的旋转(重力方向可观),即
S
W
R
(
ϕ
i
)
=
[
cos
ϕ
i
sin
ϕ
i
0
-sin
ϕ
i
cos
ϕ
i
0
0
0
1
]
^W_SR(\phi_i)=\begin{bmatrix} \text{cos}\phi_i & \text{sin}\phi_i & 0\\\text{-sin}\phi_i & \text{cos}\phi_i & 0\\ 0 & 0& 1\end{bmatrix}
SWR(ϕi)=⎣
⎡cosϕi-sinϕi0sinϕicosϕi0001⎦
⎤
对于竖直线段,它们的起始帧坐标轴和世界坐标系保持一致,即
S
W
R
=
I
3
×
3
^W_SR=I_{3\times3}
SWR=I3×3。我们使用
S
W
R
(
ϕ
0
)
^W_SR(\phi_0)
SWR(ϕ0),
(
ϕ
0
=
0
)
(\phi_0=0)
(ϕ0=0) 来表示这类起始帧。
以下部分属于本文比较关键的部分,结构线参数化和3d空间线与2d观测的关系构建。
为了获取结构线在图像中的投影,可以如下计算:
W
l
p
=
S
W
R
(
ϕ
i
)
L
S
R
L
l
p
+
W
p
S
^Wl_p= {^W_SR}(\phi_i){^S_L}R^Ll_p+{^Wp_S}
Wlp=SWR(ϕi)LSRLlp+WpSand
C
l
p
=
W
C
R
W
l
p
+
C
p
W
^Cl_p={^C_WR}^Wl_p+^Cp_W
Clp=WCRWlp+CpW
用逆深度表示法
(
θ
,
ρ
,
0
)
T
(\theta,\rho,0)^T
(θ,ρ,0)T 来取代
L
l
p
^Ll_p
Llp,则该 交点(注意这里还是参数空间的平面交点) 在图像平面上的2D投影齐次坐标为:
C
l
p
∼
W
C
R
S
W
R
(
ϕ
i
)
L
S
R
⋅
r
+
(
W
C
R
W
p
S
+
C
p
W
)
⋅
ρ
^Cl_p \sim {^C_WR}{^W_SR(\phi_i)}{^S_LR}\cdot r+({^C_WR^Wp_S+{^Cp_W}})\cdot \rho
Clp∼WCRSWR(ϕi)LSR⋅r+(WCRWpS+CpW)⋅ρ,其中
r
=
[
cos
θ
,
sin
θ
,
0
]
T
r=[\text{cos}\theta, \text{sin}\theta,0]^T
r=[cosθ,sinθ,0]T。经参数空间
Z
Z
Z 轴投影后的消影点齐次坐标计算如下:
C
v
∼
W
C
R
S
W
R
(
ϕ
i
)
L
S
R
(
:
,
3
)
^Cv \sim {^C_WR} {^W_SR(\phi_i)} {^S_L}R(:,3)
Cv∼WCRSWR(ϕi)LSR(:,3)
考虑相机内参,则图像中的直线表示为:
i
m
l
=
(
K
−
T
)
(
C
l
p
×
C
v
)
^{im}l=(K^{-T})(^Cl_p \times {}^Cv)
iml=(K−T)(Clp×Cv)
通过上述定义,给定逆深度表示
l
=
(
θ
,
ρ
)
T
l=(\theta,\rho)^T
l=(θ,ρ)T,局部Manhattan World 方向
ϕ
i
\phi_i
ϕi ,以及该直线对齐的轴(用
L
S
R
^S_LR
LSR 描述)等参数,我们在3D线和2D观测之间建立了联系。
直线投影可以写作这些参数的函数,
i
m
l
=
Π
(
l
,
ϕ
i
,
L
S
R
,
C
I
τ
,
I
W
τ
)
^{im}l=\Pi(l,\phi_i,^S_LR,{^I_C\tau},{^W_I\tau})
iml=Π(l,ϕi,LSR,CIτ,IWτ)
其中
C
I
τ
^I_C\tau
CIτ 为IMU和相机之间的相对位姿变换,可以继续被滤波器优化。
与基于点的EKF方法相比,虚线框为Struct VIO提出的新框架,包括线段检测和分类、Manhattan世界检测和融合、特征边缘化和重新三角化等。
系统概述
状态向量定义为:
-
- 其中 x I k x_{I_k} xIk为当前的IMU状态向量,包括位姿、速度和陀螺仪加速度计的偏差等。
-
- 采用Atlanta世界模型,本方法实时检测每个局部的manhattan世界并将其方向 ϕ i \phi_i ϕi 加入状态向量中并不断优化,共计 N N N个。
-
- I i w τ ^w _{I_i}\tau Iiwτ为历史IMU位姿,共计 M M M个。当超出窗口阈值时则从状态向量中排除。
-
- 所有图像特征,包括点、结构线,都没有包含在状态向量中。会在滤波器外被单独估计且仅用来推导IMU的几何约束。
- 动力学模型
与传统EKF基本一致。 - 测量模型
包括点特征和线特征,本文只描述结构线的测量模型。为了推导结构线的测量模型,首先需要计算出给定结构线在图像平面的投影。
假设结构线在第 k k k 帧图像平面的投影为 i m l = ( l 1 , l 2 , l 3 ) T ^{im}l=(l_1,l_2,l_3)^T iml=(l1,l2,l3)T 和对应观测线段的两个端点为 s a , s b s_a,s_b sa,sb (齐次坐标)。采用有符号距离来描述观测预测线投影和观测线特征的相似度如下:
(本质上是求两个观测直线端点到预测投影直线的距离所组成的向量)
通过线性化线特征参数的最新估计,残差 r k r_k rk 近似为:
通过堆叠所有可见帧的两侧,可以得到单条结构线的测量方程如下:
然后将残差 z z z 投影到 H l H_l Hl 的左零空间来得到一个新的残差 z ( 0 ) z^{(0)} z(0)如下:
写成简单形式为:
由此,将结构线本身从状态估计中解耦,从而显著减少了状态变量的维度。
注意,与点特征不同的是,结构线的测量模型有一个与给定manhattan世界方向相关联的新成分。意味着这个manhattan世界的水平方向可以通过滤波器估计,从而允许我们使用多个箱型世界来建模复杂环境。
由于解析形式过于复杂、难以计算,本文采用数值微分来计算所有的Jacobian矩阵。
- EKF更新
有两种EKF更新的触发方式:a) 一条结构线不再被跟踪到:所有的该结构线的量测将被用来进行EKF更新。考虑到偶然的跟踪失败,我们不会立即触发EKF更新,而是进行等待直到跟踪器无法恢复一定数量的帧。这种延迟更新策略显著改善了线跟踪。 b) 位姿数量超出窗口阈值。这种情况下,我们从第二旧帧开始,选择1/3的位姿状态,并且使用所有的对应点特征和线特征来构建量测方程。 - 状态管理
同样被两种事件触发:a) 新帧到来:状态增广和协方差矩阵增广。b) 检测到新的manhattan世界:状态增广和协方差矩阵增广。
注意,由于manhattan世界的方向受多种因素影响,尽管考虑到所有误差因子可以获取准确的对应协方差,但是在实践中采取简单形式: P k x ϕ = 0 P^{x\phi}_k=0 Pkxϕ=0, P k ϕ = σ ϕ 2 P^\phi_k=\sigma^2_\phi Pkϕ=σϕ2,其中 σ ϕ = 5 o \sigma_{\phi}= 5^o σϕ=5o。
关键实现
- 直线检测和跟踪
-
- LSD直线检测
-
- 3D-2D匹配来跟踪直线量测。优点是利用了相机运动去预测结构线的可能位置,减少了搜索范围,同时也方便处理偶然的跟踪丢失。对每条结构线,除了之前提到的几何参数,同样引入一个变量 r = ( r s , r e ) T r=(r_s,r_e)^T r=(rs,re)T 来表示3d空间中结构线的范围。结构线的两个端点在参数空间的坐标为: L l s = ( a , b , r s ) T ^Ll_s=(a,b,r_s)^T Lls=(a,b,rs)T 和 L l s = ( a , b , r e ) T ^Ll_s=(a,b,r_e)^T Lls=(a,b,re)T 。当新帧图片到来时,结构线将基于IMU预积分预测的相机位姿投影到图片中。下一步就是在新图像中搜索与投影相近的线段,可以通过检查位置和方向差别来实现。我们将会得到一系列的备选线段。然后,将尝试在候选线段中辨别出真实的线段关联(line correspondence)。不仅考虑线段的中点,我们将同时考虑线段上的多个点来提升跟踪的鲁棒性。
-
- 我们在结构线的两端点之间上均匀采样为: r s , r ( 1 ) , r ( 2 ) , . . . , r ( K − 2 ) , r e r_s,r^{(1)}, r^{(2)},...,r^{(K-2)},r_e rs,r(1),r(2),...,r(K−2),re ,对每个采样点,我们保存它在最新一帧的投影附近的图片块,然后通过ZNCC图像匹配方法来在候选线段中搜索它的相关点,如下图所示。最终我们选取相关点最多的线段。
- 结构线段识别
-
- 对于每个结构线段,我们试图找到它们位于哪个局部Manhattan World下,并且分类到不同的方向。
-
- 首先,计算某个Manhattan World三个方向的所有消影点。则
Z
Z
Z 轴上的消影点为:
v z = K W C R [ 0 0 1 ] T v_z=K {^C_WR}[0 \ \ 0\ \ 1]^T vz=KWCR[0 0 1]T同理,对于 X X X 和 Y Y Y 方向,计算如下:
v x ϕ i = K W C R [ cos ϕ i , -sin ϕ i , 0 ] T v y ϕ i = K W C R [ sin ϕ i , -cos ϕ i , 0 ] T \begin{aligned} v^{\phi_i}_x&=K{^C_W}R[\text{cos}\phi_i, \text{-sin}\phi_i, 0]^T\\ v^{\phi_i}_y&=K{^C_W}R[\text{sin}\phi_i, \text{-cos}\phi_i, 0]^T \end{aligned} vxϕivyϕi=KWCR[cosϕi,-sinϕi,0]T=KWCR[sinϕi,-cosϕi,0]T 注意到只有水平方向的消影点才会依赖局部Manhattan World的方向 ϕ i \phi_i ϕi 。因此,即使没有检测到Manhattan World,我们也可以识别出竖直线。
- 首先,计算某个Manhattan World三个方向的所有消影点。则
Z
Z
Z 轴上的消影点为:
-
- 为了从所有的检测线段中识别出结构线,从每个消影点引一条向线段 S S S 中点的射线,然后检测该射线和线段 S S S 的一致性。对每条线段都进行所有消影点的测试。一条线段若是与消影点其中一条一致,则判断为结构线。与该线段相关的Manhattan World就可以以对应的消影点确定下来。
-
- 当存在多个一致的消影点时选取一致性最高的那个。不属于任何消影点的直线段将从状态估计器中排除。注意即使没有检测出manhattan世界,竖直方向的线段仍然可以被识别出来。
- 结构线初始化
-
- 选择最具信息量的线段用来初始化:1)长度最长的。2)与已初始化线段距离较远的。
-
- 初始化一个新结构线 l = ( θ , ρ ) l=(\theta,\rho) l=(θ,ρ) 的关键是找到角度参数 θ \theta θ , 因为逆深度值可以设定为一个预设值。初始化第一步是为结构线建立一个起始帧。
-
- 对所有竖直方向的结构线,令起始帧与世界坐标系轴对齐,或者建立一个虚拟的Manhattan World坐标系让 ϕ 0 = 0 \phi_0=0 ϕ0=0。
-
- 对水平方向的结构线,起始帧选择为与局部Manhattan World坐标一致, ϕ i \phi_i ϕi。
-
- 角度参数
θ
\theta
θ 由相机光心到局部参数空间的XY平面上的线的方向决定。该方向可以近似为相机光心到线段
s
s
s 中点的方向。
m
m
m 为
s
s
s 中点的齐次坐标,其相机平面的逆投影射线为
K
−
1
m
K^{-1}m
K−1m,通过下式可以进一步转化到局部参数空间
L m = S L R W S R ( ϕ i ) C W R K − 1 m ^Lm={^L_SR}{^S_WR}(\phi_i){^W_CRK^{-1}m} Lm=SLRWSR(ϕi)CWRK−1m其中对于竖直方向的结构线,前两个旋转皆为单位阵,简化为 L m = C W R K − 1 m ^Lm={^W_C}RK^{-1}m Lm=CWRK−1m角度参数 θ \theta θ 就可以通过局部参数空间的水平方向来确定。让 L m = ( m x , m y , m z ) T ^Lm=(m_x,m_y,m_z)^T Lm=(mx,my,mz)T, 则角度参数为 θ 0 = atan2 ( m y , m x ) \theta_0=\text{atan2}(m_y,m_x) θ0=atan2(my,mx),逆深度都初始化为 ρ 0 \rho_0 ρ0。
- 角度参数
θ
\theta
θ 由相机光心到局部参数空间的XY平面上的线的方向决定。该方向可以近似为相机光心到线段
s
s
s 中点的方向。
m
m
m 为
s
s
s 中点的齐次坐标,其相机平面的逆投影射线为
K
−
1
m
K^{-1}m
K−1m,通过下式可以进一步转化到局部参数空间
-
- 将该初始化过程表示如下:
其不确定性可以通过如下协方差矩阵设定:
通常 σ θ 0 \sigma_{\theta_0} σθ0 是一个逆投影回来的微小值(相对于线段检测的平均误差2~4个像素)。Manhattan世界方向和当前相机位姿的不确定性可以从滤波器获得。逆深度的不确定性被人工设置为一个较大值 5 ,来覆盖从 0.2米到正无穷的范围。
- 将该初始化过程表示如下:
- 用先验知识进行结构线三角化
- 每次状态更新后都会进行三角化来更新所有线的参数,将通过最小化所有该线段的可见帧的重投影均方误差和来完成。若仅使用状态内的可见帧来进行三角化,则可能因为较小的运动视差而导致不准确的三角化结果;若使用全部可见帧来进行三角化,则会显著增加计算消耗,且滑窗外的较老的位姿估计也可能导致更大的误差。
- 为了解决上述问题,我们将为每条结构线维护一个先验分布 N ( l 0 , Σ 0 ) N(l_0, \Sigma_0) N(l0,Σ0) (线参数和协方差),来保存初始的先验分布,或者是从历史观测中推导出的先验分布。
- 总的目标函数为:
即使初始值偏离真实值较远(由于逆深度参数化引入的高线性化),该非线性最小二乘问题也可以用高斯牛顿法在3-5个迭代中求解。
- 处理长时间跟踪的线的观测丢失
- 与滑动窗口估计器类似,本文提出估计器存在的问题是:特征被跟踪的时间可能长于状态估计器中存储的图像帧的总时间。在现存的滑动窗口估计器中,这些滑窗外的观测会被简单的抛弃,因此可能导致线参数的不准确估计。在[27]中,作者将这些特征加进滤波器状态向量中(类似经典EKF-based框架)。缺点就是需要严格控制特征点的数量,以控制状态向量维度不会过高。
- 本文提出一种新方法,将这些丢失的观测转化成几何信息的先验分布,来辅助未来的更新。称这个过程为信息积累(information accumulation)。
- D D D 为将要从滑动窗口中移除的一系列位姿, l 0 o l d l^{old}_0 l0old 和 Σ 0 o l d \Sigma^{old}_0 Σ0old 为旧先验分布。当 D D D 帧集合从滑动窗口中移除的时候,先验分布将进行更新以保留被移除观测的几何约束信息。
- 新的均值
l
0
n
e
w
l^{new}_0
l0new 将通过最小化目标函数计算如下:
与前一个目标函数不同,该公式的最小化求解仅作用于丢失的观测。 使 Λ δ l = Y Λ \delta l=Y Λδl=Y 为最后一次高斯牛顿法迭代求解出的正规方程。在下一次三角化之前,将新的协方差更新为: Σ 0 n e w = Λ − 1 \Sigma^{new}_0=Λ^{-1} Σ0new=Λ−1 - 每条结构线都会锚定在状态向量中的其中一个相机位姿(起始帧)。若该起始帧从状态向量中移除,我们需要首先改变其起始帧为状态向量中存在的某帧,并更新相应的线参数和协方差。
- 检测与合并Manhattan世界
- Manhattan世界的检测要通过聚类平行线来提取消影点,然后借助消影点来确定三轴方向。若IMU可用,则这个过程可以变得更简单,因为加速度计提供了可观的垂直重力方向。
- 只要垂直方向被确定下来,就会开启Manhattan世界检测。水平地面的消影线(世界坐标系下的XY平面)计算如下:
(法向量) - 然后,使用 1-line的 ransac算法来检测可能的manhattan世界,步骤如下:
-
- 随机选取一条未被识别成结构线的线段,将其延伸与水平消影线 l ∞ l_{\infty} l∞ 相较于一个消影点 v x v_x vx,我们假设它就是该可能的manhattan世界的x轴方向的投影。因为垂直方向已知,我们可以进一步得到manhattan世界方向 ϕ \phi ϕ 和 v y v_y vy。
-
- 获取与 v x v_x vx 和 v y v_y vy 方向一致的线段的数量。
-
- 重复上述步骤直到迭代最大次数。
- 最终,若最大一致线段超出阈值,则manhattan世界被确定。
- 让 ϕ ∗ \phi^* ϕ∗ 为检测出的manhattan世界的方向,若与现存的manhattan世界方向偏差大于 5 o 5^o 5o,该manhattan世界方向 ϕ ∗ \phi^* ϕ∗ 将被加进状态向量中,协方差也如前文所述进行更新。
- 在EKF更新后,两个manhattan世界的方向偏差偶尔会小于阈值。这种情况下,我们将合并两个manhattan世界(移除最新的那个并且更新协方差,相应的结构线也会被移动到剩余的那个manhattan世界)。
- 外点排除
在EKF更新前将通过卡方检验来检测外点:
没有通过卡方检验的结构线将被排除。在EKF更新后,我们将重新三角化所有的结构线,然后进一步检查重投影误差。重投影误差大于一定阈值(本文为4个像素)的将被抛弃。
二、Leveraging Structural Information to Improve Point Line Visual-Inertial Odometry(2021)
摘要
- 和线特征一样,结构线特征同样能帮助改善基于点特征的视觉惯性里程计系统。本文提出了一个基于点和线的VIO系统,分为结构线和非结构线。
- 与 4 参数的正交表示法不同,我们只使用了两个参数来最小化结构线和非结构线的表示方法。
- 我们设计了一种基于采样点的直线匹配策略,来改善直线匹配的效率和成功率。
I. 引言
- 在一些纹理缺失和照明挑战的环境里,基于点的VIO会产生位姿估计失败。
- 线特征可以分为两类:结构线和非结构线。非结构线的处理不受环境限制,更具通用性和鲁棒性,但因为不像结构线一样具备有效的方向约束信息,导致收敛速度较慢。
- 结构性的方向取决于对应的消影点方向,可以提供高效的几何约束来改善SLAM系统性能。
- 点线VIO:
-
- PL-VIO,将直线特征加进基于点的VIO系统中。用普吕克坐标来表示直线,用最小4参数正交表示法来进行优化(因为普吕克表示法过参数化了)。
-
- Trifo- VIO 使用卡尔曼滤波来讲直线特征加入VIO。用直线的正交向量所表示的几何约束来构造误差约束。
- 结构线VIO:
-
- Camposeco 讲消影点加进VIO系统中,使用消影点的全局约束信息来修正航向误差。然而该系统仅仅把结构线作为中间结果,没有完全利用结构线的几何约束信息来修正平移。
-
- Zou 提出了一种结构线的新参数化方法并且将其加入VIO。但是并没有利用上非结构性特征。
- 本文提出了一种紧耦合的单目VIO系统 PLS-VIO (点、结构线、非结构线)。为了减少线优化过程的计算消耗,对非结构线采用了不同的参数化表示法。我们详细描述了对所有直线特征的分类、匹配和初始化。主要贡献如下:
-
- 为了更好的利用结构线信息,设计了一种线段分类策略将直线分为结构线和非结构线。
-
- 为了高效的匹配和表示线特征,设计了一种基于点特征的2d-2d和2d-3d的线特征匹配算法来改善匹配的速度和精度。而且,与4参数正交表示法相比,我们引入一种2参数表示方法来改善线估计效率。
-
- 提供了一系列实验来验证提出算法的有效性。
II. 系统概述
本文提出的方法基于VINS-Mono,我们加入了非结构线和结构线并构造了相应的约束。
- 前端为IMU和图像的预处理,包括IMU预积分、点特征检测、线特征检测和垂直线段分类等。
-后端将主要介绍非结构线和结构线的处理。 -
- 我们将垂直线(对其到重力方向的)和非垂直线加入后端,而X Y方向的线特征和非结构线特征将进一步分类。
-
- 然后,使用两种不同的线段匹配策略,为了代码处理效率和线段数据存储的简化,我们将线段特征的匹配放进了后端。
-
- 然后,初始化直线来获取3D的线段路标。
-
- 最后,IMU的位姿和地图中的3D路标将通过最小化残差(先验残差、IMU残差、点重投影误差、线重投影误差)来进行优化。
III. 结构线和非结构线
将分别介绍点和线的参数化表示法、线分类方法、匹配方法和初始化方法。
A. 路标表示
- 点:逆深度表示法
- 非结构线:如图所示
平面 π \pi π 由两个3d空间的端点 s ′ s^{\prime} s′ 和 e ′ e^{\prime} e′ 和光心 O O O 组成。3D线可以表示为 L = [ c n T , c v T ] L=[^c\textbf{n}^T, {^c\textbf{v}^T}] L=[cnT,cvT],其中 c n ^c\textbf{n} cn 为 π \pi π 的法向量, c v ^c\textbf{v} cv 为 L L L 的方向向量。由于3D线 L L L 在平面 π \pi π 上,我们可以使用 π \pi π 平面上的局部坐标系统 { P } \{P\} {P} 来表示 L L L。为了简化 L L L 的表示方法,让 { P } \{P\} {P} 的原点为 s ′ s^{\prime} s′ ,让 y 轴方向对准到从 O O O 到 s ′ s^{\prime} s′ 引出的一条射线,让z轴方向与 c n ^c\textbf{n} cn 平行。由于坐标轴的正交性,x轴与y轴和z轴垂直。 从光心 O O O 到 s ′ s^{\prime} s′ 的距离为 d d d 。
为了减少线参数的数量,我们提出了一种压缩参数化方法,仅两个参数: θ \theta θ 和 ρ = 1 / d \rho=1/d ρ=1/d。其中 θ \theta θ 为局部坐标系 { P } \{P\} {P} 中的线段方向 c v ′ ^c\textbf{v}^{\prime} cv′ 和 { P } \{P\} {P} 的x轴的方向角。其中 c v ′ = R C P c v ^c\textbf{v}^{\prime}=\textbf{R}^P_C {^c\textbf{v}} cv′=RCPcv,其中 R C P \textbf{R}^P_C RCP 为相机坐标系相对于局部坐标系 { P } \{P\} {P} 的旋转矩阵。 - 结构线:原封不动的用了和StructVIO一样的2参数表示法,同样采用了Atlanda假设。
B. 结构线和非结构线分类 (没有创新)
我们使用LSD来检测图像中线段并从中分类出结构线段,剩下的即为非结构线段。我们使用消影点来识别结构线。例如,为了分类z方向的结构线,我们从消影点
v
z
\textbf{v}_z
vz 引一条到线段
S
S
S 的中点的射线,当该射线和线段
S
S
S 的角度偏差
A
e
r
r
A_{err}
Aerr 和 距离
D
e
r
r
D_{err}
Derr 小于阈值的时候,认为属于z轴方向的结构线。
C. 2D-2D 和 2D-3D的线段匹配
为了提高线段匹配的稳定性和精度,我们联合两种跟踪策略,即 帧到帧(frame-to-frame)跟踪和 帧到地图(frame-to-map)跟踪。我们使用 帧到帧跟踪方法来跟踪新检测到的直线。若匹配到的线段数过少,将继续执行帧到地图的方法来增加匹配线。
- 帧到帧线跟踪方法:采样前一帧的所有直线,得到一系列采样点 p i ∈ { P _ s a m p l e } p_i \in \{P \_{sample}\} pi∈{P_sample},然后使用对极搜索方法(epipolar searching method)来找到当前帧对应的候选匹配点。ZMSSD(zero-mean sum of squared differences) 模版用来计算两个点的匹配分数。我们选择得分最高的候选点作为被跟踪点,得到 p i ′ ∈ { P _ t r a c k e d } p^{\prime}_i \in \{P \_{tracked}\} pi′∈{P_tracked}。若被跟踪点到当前帧直线的距离若小于阈值 m _ t h m\_th m_th (5个像素),则判定有效。当有效点的数量大于前一帧直线的采样点的0.8倍,则认为该线为最佳匹配。
- 帧到地图线跟踪方法:使用ZNCC(Zero-normalized-cross-correlation)匹配方法。对一条3D线,我们得到该条线对应的最新的历史观测帧 F i F_i Fi ,然后使用ZNCC方法匹配当前帧和 F i F_i Fi 的线段。由于相机的快速运动和场景遮挡,匹配到的线段长度可能非常不同以至于影响匹配精度。我们使用对极几何约束来确定线段的采样范围去辅助ZNCC匹配:即计算出 F i F_i Fi 和当前帧的本征矩阵 E \textbf{E} E 后, F i F_i Fi 的线段端点将被投影到当前帧,通过将当前帧的匹配线与两条极线相交来确定相应的采样范围,以此来提升ZNCC匹配的成功率。
D. 结构线和非结构线的初始化
线段初始化的稳定性和精度显著影响位姿估计。对2参数表示的结构线和非结构线,我们使用不同的初始化方法来确定合理的初始值。
- 非结构线初始化。归一化平面的线段可以用两个端点来表示,和光心一起即可确定一个平面
π
\pi
π ,给定相机帧
c
1
c_1
c1 和
c
2
c_2
c2 的两个平面
π
1
\pi_1
π1 和
π
2
\pi_2
π2 ,对偶普吕克矩阵可以计算如下:
我们从对偶普吕克矩阵中获取对应的普吕克坐标 。我们把相机坐标系下的 c v ^c\textbf{v} cv 转换到局部坐标系 { P } \{P\} {P} 下, c v ′ = R C P c v ^c\textbf{v}^{\prime}=R^P_C{^c\textbf{v}} cv′=RCPcv。因此,我们可以初始化 θ \theta θ 为局部坐标系x轴方向和 c v ′ ^c\textbf{v}^{\prime} cv′ 的方向偏差,而逆深度 ρ \rho ρ 通常初始化为 ρ 0 = 0.2 \rho_0=0.2 ρ0=0.2。
- 结构线初始化(与StructVIO略有不同)。同样要首先计算出普吕克坐标
c
L
=
[
c
n
T
,
c
v
T
]
^cL=[^c\textbf{n}^T, {^c\textbf{v}^T}]
cL=[cnT,cvT]。然后我们使用线段修建(line trimming) 来得到线段端点在世界坐标系下的3D表示
w
L
=
[
s
w
T
,
e
w
T
]
T
^wL=[\textbf{s}^{wT}, \textbf{e}^{wT}]^T
wL=[swT,ewT]T。为了获取
w
L
^wL
wL 和XY 平面在世界坐标系系统下的交点,我们将
l
p
^l\textbf{p}
lp 转换到世界坐标下:
然后我们将 w L ^wL wL 与平面 w p ^w\textbf{p} wp 相交来获取交叉点 w l p ^w\textbf{l}_p wlp, 并且将 w l p ^w\textbf{l}_p wlp 重新转化到参数空间来获取 l l p ^l\textbf{l}_p llp, 如下:
从而获得structVIO里的参数化表示方法。
- 本文里所提出的方法,是将参数空间的XY平面转换到世界坐标系下,再与世界坐标系下用普吕克表示的直线相交,求出交点坐标,然后再将该交点转换回参数空间获得参数化表示。
- 而StructVIO里所提出的方法,是直接将相机坐标系下观测到的线段中点视作参数空间直线与XY平面的交点在相机坐标系下的表示,因此将该交点直接从世界坐标系转化到参数空间。
因此可以说,本质区别就在于这个参数空间的交点在世界坐标系下坐标的计算方法不同,前者理论上更精确,但后者处理上更方便。
IV. 点线结合的VIO
A. VIO 系统组成
其中
r
f
k
c
i
r^{c_i}_{f_k}
rfkci 为点的重投影误差,
r
L
l
c
i
r^{c_i}_{L_l}
rLlci 和
r
L
s
c
i
r^{c_i}_{L_s}
rLsci 分别为非结构线和结构线的重投影误差。
ρ
(
⋅
)
\rho (\cdot)
ρ(⋅) 为cauchy 鲁棒核函数,用来抑制外点。
B. 点特征测量模型(与VINS-Mono一致)
C. 非结构线测量模型
局部坐标系下的线段方向向量可以计算如:
d
v
=
[
c
o
s
θ
,
s
i
n
θ
,
0
]
T
^d\textbf{v}=[cos\theta, sin\theta, 0]^T
dv=[cosθ,sinθ,0]T ,然后将
d
v
^d\textbf{v}
dv从局部坐标系转换到相机坐标系获得
c
v
^c\textbf{v}
cv,即
其中
R
P
C
\textbf{R}^C_P
RPC 的列由局部坐标系
{
P
}
\{P\}
{P} 的x,y,z轴组成。
相机坐标系下,线段的其中一个端点
s
\textbf{s}
s 可以计算如下:
其中
y
\textbf{y}
y 为局部坐标系的y轴。
为了获得一条直线在归一化图像坐标系的投影,需要将第一个观测关键帧下的端点
s
\textbf{s}
s 和 方向
c
v
^c\textbf{v}
cv 转换到目标帧上去:
我们就可以得到在目标归一化图像平面下的线段方程:
归一化平面的线段测量
z
L
l
m
i
\textbf{z}^{m_i}_{L_l}
zLlmi 包含两个端点的测量
s
l
m
i
\textbf{s}^{m_i}_l
slmi 和
e
l
m
i
\textbf{e}^{m_i}_l
elmi,则线的重投影误差为:
其中
d
(
s
,
l
)
d(\textbf{s},\textbf{l})
d(s,l) 是端点
s
\textbf{s}
s 到 投影直线
l
\textbf{l}
l 的距离:
D. 结构线的测量模型
结构线的残差是将参数空间中的线参数转移到起始帧再转移到目标帧,然后构建重投影误差。
其中
R
W
C
\textbf{R}^C_W
RWC 为从世界坐标系到目标坐标系的旋转矩阵,
l
l
p
=
[
a
,
b
,
0
]
T
^l\textbf{l}_p =[a,b,0]^T
llp=[a,b,0]T,其中 a,b可以用 参数
θ
\theta
θ 和
ρ
\rho
ρ 来表示。
在参数空间中的结构线的方向向量表示为
l
v
=
[
0
,
0
,
1
]
T
^l\textbf{v}=[0,0,1]^T
lv=[0,0,1]T ,将其变换到相机坐标系下:
则可以进一步得到目标帧下的直线方程:
与非结构线类似,结构线的残差同样定义为观测线段端点到投影直线的距离。