多传感器融合定位十三-基于图优化的建图方法其二
多传感器融合定位十二-基于图优化的建图方法其二
- 3.4 预积分方差计算
- 3.4.1 核心思路
- 3.4.2 连续时间下的微分方程
- 3.4.3 离散时间下的传递方程
- 3.5 预积分更新
- 4. 典型方案介绍
- 4.1 LIO-SAM介绍
- 5. 融合编码器的优化方案
- 5.1 整体思路介绍
- 5.2 预积分模型设计
Reference:
- 深蓝学院-多传感器融合
- 多传感器融合定位理论基础
文章跳转:
- 多传感器融合定位一-3D激光里程计其一:ICP
- 多传感器融合定位二-3D激光里程计其二:NDT
- 多传感器融合定位三-3D激光里程计其三:点云畸变补偿
- 多传感器融合定位四-3D激光里程计其四:点云线面特征提取
- 多传感器融合定位五-点云地图构建及定位
- 多传感器融合定位六-惯性导航原理及误差分析
- 多传感器融合定位七-惯性导航解算及误差分析其一
- 多传感器融合定位八-惯性导航解算及误差分析其二
- 多传感器融合定位九-基于滤波的融合方法Ⅰ其一
- 多传感器融合定位十-基于滤波的融合方法Ⅰ其二
- 多传感器融合定位十一-基于滤波的融合方法Ⅱ
- 多传感器融合定位十二-基于图优化的建图方法其一
- 多传感器融合定位十三-基于图优化的建图方法其二
- 多传感器融合定位十四-基于图优化的定位方法
- 多传感器融合定位十五-多传感器时空标定(综述)
3.4 预积分方差计算
3.4.1 核心思路
在融合时,需要给不同信息设置权重,而权重由方差得来,因此对于IMU积分,也要计算其方差。方差的计算形式与第四章相同,即:
Pi,k+1=FkPi,kFk⊤+BkQBk⊤\boldsymbol{P}_{i, k+1}=\mathbf{F}_k \boldsymbol{P}_{i, k} \mathbf{F}_k^{\top}+\mathbf{B}_k \boldsymbol{Q} \mathbf{B}_k^{\top} Pi,k+1=FkPi,kFk⊤+BkQBk⊤但需注意的是,此处 Fk\boldsymbol{F}_kFk 和 Gk\boldsymbol{G}_kGk 是离散时间下的状态传递方程中的矩阵,而我们一般是在连续时间下推导微分方程,再用它计算离散时间下的传递方程。
3.4.2 连续时间下的微分方程
连续时间下的微分方程一般写为如下形式:
x˙=Ftx+Btwx=[δαtbkδθtbkδβtbkδbatδbwt]w=[nanwnbanbw]\begin{gathered} \dot{\boldsymbol{x}}=\boldsymbol{F}_t \boldsymbol{x}+\boldsymbol{B}_t \boldsymbol{w} \\ \boldsymbol{x}=\left[\begin{array}{l} \delta \boldsymbol{\alpha}_t^{b_k} \\ \delta \boldsymbol{\theta}_t^{b_k} \\ \delta \boldsymbol{\beta}_t^{b_k} \\ \delta \boldsymbol{b}_{\boldsymbol{a}_t} \\ \delta \boldsymbol{b}_{w_t} \end{array}\right] \\ \boldsymbol{w}=\left[\begin{array}{l} \boldsymbol{n}_a \\ \boldsymbol{n}_w \\ \boldsymbol{n}_{b_a} \\ \boldsymbol{n}_{b_w} \end{array}\right] \end{gathered} x˙=Ftx+Btwx=δαtbkδθtbkδβtbkδbatδbwtw=nanwnbanbw以上变量排列顺序是为了与lio/vio等常见系统顺序保持一致。
推导方法已经在第6讲介绍,此处直接给出结果:
δθ˙tbk\delta \dot{\boldsymbol{\theta}}_t^{b_k}δθ˙tbk 的微分方程:
δθ˙tbk=−[ωt−bωt]×δθtbk+nω−δbωt\delta \dot{\boldsymbol{\theta}}_t^{b_k}=-\left[\boldsymbol{\omega}_t-\boldsymbol{b}_{\omega_t}\right]_{\times} \delta \boldsymbol{\theta}_t^{b_k}+\boldsymbol{n}_\omega-\delta \boldsymbol{b}_{\omega_t} δθ˙tbk=−[ωt−bωt]×δθtbk+nω−δbωtδβ˙tbk\delta \dot{\boldsymbol{\beta}}_t^{b_k}δβ˙tbk 的微分方程:
δβ˙tbk=−Rt[at−bat]×δθtbk+Rt(na−δbat)\delta \dot{\boldsymbol{\beta}}_t^{b_k}=-\boldsymbol{R}_t\left[\boldsymbol{a}_t-\boldsymbol{b}_{a_t}\right]_{\times} \delta \boldsymbol{\theta}_t^{b_k}+\boldsymbol{R}_t\left(\boldsymbol{n}_a-\delta \boldsymbol{b}_{a_t}\right) δβ˙tbk=−Rt[at−bat]×δθtbk+Rt(na−δbat)δα˙tbk\delta \dot{\boldsymbol{\alpha}}_t^{b_k}δα˙tbk 的微分方程:
δα˙tbk=δβtbk\delta \dot{\boldsymbol{\alpha}}_t^{b_k}=\delta \boldsymbol{\beta}_t^{b_k} δα˙tbk=δβtbk
3.4.3 离散时间下的传递方程
得到连续时间微分方程以后,就可以计算离散时间的递推方程了,表示为:
xk+1=Fkxk+Bkwk\boldsymbol{x}_{k+1}=\boldsymbol{F}_k \boldsymbol{x}_k+\boldsymbol{B}_k \boldsymbol{w}_k xk+1=Fkxk+Bkwk其中:
xk+1=[δαk+1δθk+1δβk+1δbak+1δbωk+1]xk=[δαkδθkδβkδbakδbωk]wk=[naknwknak+1nwk+1nbanbw]\boldsymbol{x}_{k+1}=\left[\begin{array}{l} \delta \boldsymbol{\alpha}_{k+1} \\ \delta \boldsymbol{\theta}_{k+1} \\ \delta \boldsymbol{\beta}_{k+1} \\ \delta \boldsymbol{b}_{a_{k+1}} \\ \delta \boldsymbol{b}_{\omega_{k+1}} \end{array}\right] \quad \boldsymbol{x}_k=\left[\begin{array}{c} \delta \boldsymbol{\alpha}_k \\ \delta \boldsymbol{\theta}_k \\ \delta \boldsymbol{\beta}_k \\ \delta \boldsymbol{b}_{a_k} \\ \delta \boldsymbol{b}_{\omega_k} \end{array}\right] \quad \boldsymbol{w}_k=\left[\begin{array}{c} \boldsymbol{n}_{a_k} \\ \boldsymbol{n}_{w_k} \\ \boldsymbol{n}_{a_{k+1}} \\ \boldsymbol{n}_{w_{k+1}} \\ \boldsymbol{n}_{b_a} \\ \boldsymbol{n}_{b_w} \end{array}\right] xk+1=δαk+1δθk+1δβk+1δbak+1δbωk+1xk=δαkδθkδβkδbakδbωkwk=naknwknak+1nwk+1nbanbw
-
δθk+1\delta \boldsymbol{\theta}_{k+1}δθk+1 的求解
由于连续时间下有
δθ˙=−[ωt−bωt]×δθ+nω−δbωt\delta \dot{\boldsymbol{\theta}}=-\left[\boldsymbol{\omega}_t-\boldsymbol{b}_{\omega_t}\right]_{\times} \delta \boldsymbol{\theta}+\boldsymbol{n}_\omega-\delta \boldsymbol{b}_{\omega_t} δθ˙=−[ωt−bωt]×δθ+nω−δbωt则离散时间下应该有
因此有
δθ˙k=−[ωk+ωk+12−bωt]×δθk+nωk+nωk+12−δbωkδθk+1=[I−[ωk+ωk+12−bωk]×δt]δθk+δtnωk+nωk+12−δtδbωk\delta \dot{\boldsymbol{\theta}}_k=-\left[\frac{\boldsymbol{\omega}_k+\boldsymbol{\omega}_{k+1}}{2}-\boldsymbol{b}_{\omega_t}\right]_{\times} \delta \boldsymbol{\theta}_k+\frac{\boldsymbol{n}_{\omega_k}+\boldsymbol{n}_{\omega_{k+1}}}{2}-\delta \boldsymbol{b}_{\omega_k}\\ \delta \boldsymbol{\theta}_{k+1}=\left[\boldsymbol{I}-\left[\frac{\boldsymbol{\omega}_k+\boldsymbol{\omega}_{k+1}}{2}-\boldsymbol{b}_{\omega_k}\right]_{\times} \delta t\right] \delta \boldsymbol{\theta}_k+\delta t \frac{\boldsymbol{n}_{\omega_k}+\boldsymbol{n}_{\omega_{k+1}}}{2}-\delta t \delta \boldsymbol{b}_{\omega_k} δθ˙k=−[2ωk+ωk+1−bωt]×δθk+2nωk+nωk+1−δbωkδθk+1=[I−[2ωk+ωk+1−bωk]×δt]δθk+δt2nωk+nωk+1−δtδbωk令 ω‾=ωk+ωk+12−bωk\overline{\boldsymbol{\omega}}=\frac{\omega_k+\omega_{k+1}}{2}-\boldsymbol{b}_{\omega_k}ω=2ωk+ωk+1−bωk ,则上式可以重新写为
δθk+1=[I−[ω‾]×δt]δθk+δt2nωk+δt2nωk+1−δtδbωk\delta \boldsymbol{\theta}_{k+1}=\left[\boldsymbol{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t\right] \delta \boldsymbol{\theta}_k+\frac{\delta t}{2} \boldsymbol{n}_{\omega_k}+\frac{\delta t}{2} \boldsymbol{n}_{\omega_{k+1}}-\delta t \delta \boldsymbol{b}_{\omega_k} δθk+1=[I−[ω]×δt]δθk+2δtnωk+2δtnωk+1−δtδbωk -
δβk+1\delta \boldsymbol{\beta}_{k+1}δβk+1 的求解
由于连续时间下有
δβ˙=−Rt[at−bat]×δθ+Rt(na−δbat)\delta \dot{\boldsymbol{\beta}}=-\boldsymbol{R}_t\left[\boldsymbol{a}_t-\boldsymbol{b}_{a_t}\right]_{\times} \delta \boldsymbol{\theta}+\boldsymbol{R}_t\left(\boldsymbol{n}_a-\delta \boldsymbol{b}_{a_t}\right) δβ˙=−Rt[at−bat]×δθ+Rt(na−δbat)
则离散时间下应该有
δβ˙k=−12Rk[ak−bak]×δθk−12Rk+1[ak+1−bak]×δθk+1+12Rknak+12Rk+1nak+1−12(Rk+Rk+1)δbak\begin{aligned} \delta \dot{\boldsymbol{\beta}}_k= & -\frac{1}{2} \boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times} \delta \boldsymbol{\theta}_k \\ & -\frac{1}{2} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \delta \boldsymbol{\theta}_{k+1} \\ & +\frac{1}{2} \boldsymbol{R}_k \boldsymbol{n}_{a_k} \\ & +\frac{1}{2} \boldsymbol{R}_{k+1} \boldsymbol{n}_{a_{k+1}} \\ & -\frac{1}{2}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta \boldsymbol{b}_{a_k} \end{aligned} δβ˙k=−21Rk[ak−bak]×δθk−21Rk+1[ak+1−bak]×δθk+1+21Rknak+21Rk+1nak+1−21(Rk+Rk+1)δbak把前面求得的 δθk+1\delta \boldsymbol{\theta}_{k+1}δθk+1 的表达式代入上式,可得
δβ˙k=−12Rk[ak−bak]×δθk−12Rk+1[ak+1−bak]×{[I−[ω‾]×δt]δθk+δt2nωk+δt2nωk+1−δtδbωk}+12Rknak+12Rk+1nak+1−12(Rk+Rk+1)δbak\begin{aligned} \delta \dot{\boldsymbol{\beta}}_k= & -\frac{1}{2} \boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times} \delta \boldsymbol{\theta}_k \\ & -\frac{1}{2} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times}\left\{\left[\boldsymbol{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t\right] \delta \boldsymbol{\theta}_k+\frac{\delta t}{2} \boldsymbol{n}_{\omega_k}+\frac{\delta t}{2} \boldsymbol{n}_{\omega_{k+1}}-\delta t \delta \boldsymbol{b}_{\omega_k}\right\} \\ & +\frac{1}{2} \boldsymbol{R}_k \boldsymbol{n}_{a_k} \\ & +\frac{1}{2} \boldsymbol{R}_{k+1} \boldsymbol{n}_{a_{k+1}} \\ & -\frac{1}{2}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta \boldsymbol{b}_{a_k} \end{aligned} δβ˙k=−21Rk[ak−bak]×δθk−21Rk+1[ak+1−bak]×{[I−[ω]×δt]δθk+2δtnωk+2δtnωk+1−δtδbωk}+21Rknak+21Rk+1nak+1−21(Rk+Rk+1)δbak经过一系列合并同类项以后,最终的合并结果为
δβ˙k=−12[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω‾]×δt)]δθk−δt4Rk+1[ak+1−bak]×nωk−δt4Rk+1[ak+1−bak]×nωk+1+δt2Rk+1[ak+1−bak]×δbωk+12Rknak+12Rk+1nak+1−12(Rk+Rk+1)δbak\begin{aligned} \delta \dot{\boldsymbol{\beta}}_k= & -\frac{1}{2}\left[\boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times}+\boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times}\left(\boldsymbol{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t\right)\right] \delta \boldsymbol{\theta}_k \\ & -\frac{\delta t}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \boldsymbol{n}_{\omega_k} \\ & -\frac{\delta t}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \boldsymbol{n}_{\omega_{k+1}} \\ & +\frac{\delta t}{2} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \delta \boldsymbol{b}_{\omega_k} \\ & +\frac{1}{2} \boldsymbol{R}_k \boldsymbol{n}_{a_k} \\ & +\frac{1}{2} \boldsymbol{R}_{k+1} \boldsymbol{n}_{a_{k+1}} \\ & -\frac{1}{2}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta \boldsymbol{b}_{a_k} \end{aligned} δβ˙k=−21[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω]×δt)]δθk−4δtRk+1[ak+1−bak]×nωk−4δtRk+1[ak+1−bak]×nωk+1+2δtRk+1[ak+1−bak]×δbωk+21Rknak+21Rk+1nak+1−21(Rk+Rk+1)δbak由导数形式可以得到递推形式如下
δβk+1=δβk−δt2[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω‾]×δt)]δθk−δt24Rk+1[ak+1−bak]×nωk−δt24Rk+1[ak+1−bak]×nωk+1+δt22Rk+1[ak+1−bak]×δbωk+δt2Rknak+δt2Rk+1nak+1−δt2(Rk+Rk+1)δbak\begin{aligned} \delta \boldsymbol{\beta}_{k+1}= & \delta \boldsymbol{\beta}_k \\ & -\frac{\delta t}{2}\left[\boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times}+\boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times}\left(\boldsymbol{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t\right)\right] \delta \boldsymbol{\theta}_k \\ & -\frac{\delta t^2}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \boldsymbol{n}_{\omega_k} \\ & -\frac{\delta t^2}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \boldsymbol{n}_{\omega_{k+1}} \\ & +\frac{\delta t^2}{2} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \delta \boldsymbol{b}_{\omega_k} \\ & +\frac{\delta t}{2} \boldsymbol{R}_k \boldsymbol{n}_{a_k} \\ & +\frac{\delta t}{2} \boldsymbol{R}_{k+1} \boldsymbol{n}_{a_{k+1}} \\ & -\frac{\delta t}{2}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta \boldsymbol{b}_{a_k} \end{aligned} δβk+1=δβk−2δt[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω]×δt)]δθk−4δt2Rk+1[ak+1−bak]×nωk−4δt2Rk+1[ak+1−bak]×nωk+1+2δt2Rk+1[ak+1−bak]×δbωk+2δtRknak+2δtRk+1nak+1−2δt(Rk+Rk+1)δbak -
δαk+1\delta \boldsymbol{\alpha}_{k+1}δαk+1 的求解
由于连续时间下有 δα˙t=δβt\delta \dot{\boldsymbol{\alpha}}_t=\delta \boldsymbol{\beta}_tδα˙t=δβt ,则离散时间下应该有
δα˙k=12δβk+12δβk+1=δβk−δt4[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω‾]×δt)]δθk−δt28Rk+1[ak+1−bak]×nωk−δt28Rk+1[ak+1−bak]×nωk+1+δt24Rk+1[ak+1−bak]×δbωk+δt4Rknak+δt4Rk+1nak+1−δt4(Rk+Rk+1)δbak\begin{aligned} \delta \dot{\boldsymbol{\alpha}}_k= & \frac{1}{2} \delta \boldsymbol{\beta}_k+\frac{1}{2} \delta \boldsymbol{\beta}_{k+1} \\ = & \delta \boldsymbol{\beta}_k \\ & -\frac{\delta t}{4}\left[\boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times}+\boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times}\left(\boldsymbol{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t\right)\right] \delta \boldsymbol{\theta}_k \\ & -\frac{\delta t^2}{8} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{\boldsymbol{a}_k}\right]_{\times} \boldsymbol{n}_{\omega_k} \\ & -\frac{\delta t^2}{8} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{\boldsymbol{a}_k}\right]_{\times} \boldsymbol{n}_{\omega_{k+1}} \\ & +\frac{\delta t^2}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{\boldsymbol{a}_k}\right] \times \delta \boldsymbol{b}_{\omega_k} \\ & +\frac{\delta t}{4} \boldsymbol{R}_k \boldsymbol{n}_{a_k} \\ & +\frac{\delta t}{4} \boldsymbol{R}_{k+1} \boldsymbol{n}_{a_{k+1}} \\ & -\frac{\delta t}{4}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta \boldsymbol{b}_{a_k} \end{aligned} δα˙k==21δβk+21δβk+1δβk−4δt[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω]×δt)]δθk−8δt2Rk+1[ak+1−bak]×nωk−8δt2Rk+1[ak+1−bak]×nωk+1+4δt2Rk+1[ak+1−bak]×δbωk+4δtRknak+4δtRk+1nak+1−4δt(Rk+Rk+1)δbak由导数形式可以写出递推形式
δαk+1=δαk+δtδβk−δt24[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ωˉ]×δt)]δθk−δt38Rk+1[ak+1−bak]×nωk−δt38Rk+1[ak+1−bak]×nωk+1+δt34Rk+1[ak+1−bak]×δbωk+δt24Rknak+δt24Rk+1nak+1−δt24(Rk+Rk+1)δbak\begin{aligned} \delta \boldsymbol{\alpha}_{k+1}= & \delta \boldsymbol{\alpha}_k \\ & +\delta t \delta \boldsymbol{\beta}_k \\ & -\frac{\delta t^2}{4}\left[\boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times}+\boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times}\left(\boldsymbol{I}-[\bar{\omega}]_{\times} \delta t\right)\right] \delta \boldsymbol{\theta}_k \\ & -\frac{\delta t^3}{8} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \boldsymbol{n}_{\omega_k} \\ & -\frac{\delta t^3}{8} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \boldsymbol{n}_{\omega_{k+1}} \\ & +\frac{\delta t^3}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \delta \boldsymbol{b}_{\omega_k} \\ & +\frac{\delta t^2}{4} \boldsymbol{R}_k \boldsymbol{n}_{a_k} \\ & +\frac{\delta t^2}{4} \boldsymbol{R}_{k+1} \boldsymbol{n}_{a_{k+1}} \\ & -\frac{\delta t^2}{4}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta \boldsymbol{b}_{a_k} \end{aligned} δαk+1=δαk+δtδβk−4δt2[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ωˉ]×δt)]δθk−8δt3Rk+1[ak+1−bak]×nωk−8δt3Rk+1[ak+1−bak]×nωk+1+4δt3Rk+1[ak+1−bak]×δbωk+4δt2Rknak+4δt2Rk+1nak+1−4δt2(Rk+Rk+1)δbak由以上的推导结果,便可以写出 xk+1=Fkxk+Bkwk\boldsymbol{x}_{k+1}=\boldsymbol{F}_k \boldsymbol{x}_k+\boldsymbol{B}_k \boldsymbol{w}_kxk+1=Fkxk+Bkwk 中的矩阵
Fk=[If12Iδt−14(Rk+Rk+1)δt2f150I−[ω‾]×δt00−Iδt0f32I−12(Rk+Rk+1)δtf35000I00000I]Bk=[14Rkδt2g1214Rk+1δt2g1400012Iδt012Iδt0012Rkδtg3212Rk+1δtg34000000Iδt000000Iδt]\begin{aligned} & \mathbf{F}_k= {\left[\begin{array}{cccccc} \mathbf{I} & \mathbf{f}_{12} & \mathbf{I} \delta t & -\frac{1}{4}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta t^2 & \mathbf{f}_{15} \\ \mathbf{0} & \mathbf{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t & \mathbf{0} & \mathbf{0} & -\mathbf{I} \delta t \\ \mathbf{0} & \mathbf{f}_{32} & \mathbf{I} & -\frac{1}{2}\left(\boldsymbol{R}_k+\boldsymbol{R}_{k+1}\right) \delta t & \mathbf{f}_{35} \\ \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{I} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{I} \end{array}\right] } \\ & \mathbf{B}_k=\left[\begin{array}{cccccc} \frac{1}{4} \boldsymbol{R}_k \delta t^2 & \mathbf{g}_{12} & \frac{1}{4} \boldsymbol{R}_{k+1} \delta t^2 & \mathbf{g}_{14} & \mathbf{0} & \mathbf{0} \\ \mathbf{0} & \frac{1}{2} \mathbf{I} \delta t & \mathbf{0} & \frac{1}{2} \mathbf{I} \delta t & \mathbf{0} & \mathbf{0} \\ \frac{1}{2} \boldsymbol{R}_k \delta t & \mathbf{g}_{32} & \frac{1}{2} \boldsymbol{R}_{k+1} \delta t & \mathbf{g}_{34} & \mathbf{0} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{I} \delta t & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{I} \delta t \end{array}\right] \end{aligned} Fk=I0000f12I−[ω]×δtf3200Iδt0I00−41(Rk+Rk+1)δt20−21(Rk+Rk+1)δtI0f15−Iδtf350IBk=41Rkδt2021Rkδt00g1221Iδtg320041Rk+1δt2021Rk+1δt00g1421Iδtg3400000Iδt00000Iδt上面的矩阵中,有
f12=−δt24[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω‾]×δt)]f15=δt34Rk+1[ak+1−bak]×δbωkf32=−δt2[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω‾]×δt)]f35=δt22Rk+1[ak+1−bak]×g12=−δt38Rk+1[ak+1−bak]×g14=−δt38Rk+1[ak+1−bak]×g32=−δt24Rk+1[ak+1−bak]×g34=−δt24Rk+1[ak+1−bak]×\begin{aligned} & \boldsymbol{f}_{12}=-\frac{\delta t^2}{4}\left[\boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times}+\boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times}\left(\boldsymbol{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t\right)\right] \\ & \boldsymbol{f}_{15}=\frac{\delta t^3}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \delta \boldsymbol{b}_{\omega_k} \\ & \boldsymbol{f}_{32}=-\frac{\delta t}{2}\left[\boldsymbol{R}_k\left[\boldsymbol{a}_k-\boldsymbol{b}_{a_k}\right]_{\times}+\boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times}\left(\boldsymbol{I}-[\overline{\boldsymbol{\omega}}]_{\times} \delta t\right)\right] \\ & \boldsymbol{f}_{35}=\frac{\delta t^2}{2} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \\ & \boldsymbol{g}_{12}=-\frac{\delta t^3}{8} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \\ & \boldsymbol{g}_{14}=-\frac{\delta t^3}{8} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \\ & \boldsymbol{g}_{32}=-\frac{\delta t^2}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \\ & \boldsymbol{g}_{34}=-\frac{\delta t^2}{4} \boldsymbol{R}_{k+1}\left[\boldsymbol{a}_{k+1}-\boldsymbol{b}_{a_k}\right]_{\times} \end{aligned} f12=−4δt2[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω]×δt)]f15=4δt3Rk+1[ak+1−bak]×δbωkf32=−2δt[Rk[ak−bak]×+Rk+1[ak+1−bak]×(I−[ω]×δt)]f35=2δt2Rk+1[ak+1−bak]×g12=−8δt3Rk+1[ak+1−bak]×g14=−8δt3Rk+1[ak+1−bak]×g32=−4δt2Rk+1[ak+1−bak]×g34=−4δt2Rk+1[ak+1−bak]×
3.5 预积分更新
回到bias变化时,预积分结果怎样重新计算的 问题,再次给出它的泰勒展开形式:
αbibj=α‾bibj+Jbiaαδbia+Jbigαδbigβbibj=β‾bibj+Jbiaβδbia+Jbigβδbigqbibj=q‾bibj⊗[112Jbiqqδbig]\begin{gathered} \boldsymbol{\alpha}_{b_i b_j}=\overline{\boldsymbol{\alpha}}_{b_i b_j}+\mathbf{J}_{b_i^a}^\alpha \delta \mathbf{b}_i^a+\mathbf{J}_{b_i^g}^\alpha \delta \mathbf{b}_i^g \\ \boldsymbol{\beta}_{b_i b_j}=\overline{\boldsymbol{\beta}}_{b_i b_j}+\mathbf{J}_{b_i^a}^\beta \delta \mathbf{b}_i^a+\mathbf{J}_{b_i^g}^\beta \delta \mathbf{b}_i^g \\ \mathbf{q}_{b_i b_j}=\overline{\mathbf{q}}_{b_i b_j} \otimes\left[\begin{array}{c} 1 \\ \frac{1}{2} \mathbf{J}_{b_i^q}^q \delta \mathbf{b}_i^g \end{array}\right] \end{gathered} αbibj=αbibj+Jbiaαδbia+Jbigαδbigβbibj=βbibj+Jbiaβδbia+Jbigβδbigqbibj=qbibj⊗[121Jbiqqδbig]这里,雅可比没有明确的闭式解,但是在推导 方差的更新时,我们得到了
xk+1=Fkxk+Bkwk\boldsymbol{x}_{k+1}=\boldsymbol{F}_k \boldsymbol{x}_k+\boldsymbol{B}_k \boldsymbol{w}_k xk+1=Fkxk+Bkwk通过该递推形式,可以知道
Jk+1=FkJk\mathbf{J}_{k+1}=\mathbf{F}_k \mathbf{J}_k Jk+1=FkJk即预积分结果的雅可比,可以通过这种迭代方式计算。
Jj\boldsymbol{J}_jJj 中关于 bias 的项,就是预积分泰勒展开时,各 bias 对应的雅可比。
4. 典型方案介绍
4.1 LIO-SAM介绍
论文名称:LIO-SAM: Tightly-coupled Lidar Inertial Odometry via Smoothing and Mapping
代码地址:https://github.com/TixiaoShan/LIO-SAM
最大特点:分两步完成,先通过点云特征计算出相对位姿,再利用相对位姿、IMU预积分和GPS做融合。
相比于直接一步做紧耦合,大大提高了效率,而且实测性能也很优异。
5. 融合编码器的优化方案
5.1 整体思路介绍
理论上,只要是在载体系下测量,且频率比关键帧的提取频率高,就都可以做预积分。
编码器与IMU基于预积分的融合有多种方法,此处选相对简单的一种:
把它们当做一个整体的传感器来用,IMU提供角速度,编码器提供位移增量,且不考虑编码器的误差。
此时对预积分模型,输入为
角速度: ωk=[ωxkωykωzk]\omega_k=\left[\begin{array}{l}\omega_{x k} \\ \omega_{y k} \\ \omega_{z k}\end{array}\right] \quadωk=ωxkωykωzk 位移增量: ϕk=[ϕxk00]\phi_k=\left[\begin{array}{c}\phi_{x k} \\ 0 \\ 0\end{array}\right]ϕk=ϕxk00
输出中包含位置、姿态,但不包含速度。
5.2 预积分模型设计
连续时间下,从 iii 时刻到 jjj 时刻 IMU的积分结果为
pwbj=pwbi+∫t∈[i,j]qwbtϕbtδtqwbj=∫t∈[i,j]qwbt⊗[012ωbt]δt\begin{aligned} \mathbf{p}_{w b_j} & =\mathbf{p}_{w b_i}+\int_{t \in[i, j]} \mathbf{q}_{w b_t} \boldsymbol{\phi}^{b_t} \delta t \\ \mathbf{q}_{w b_j} & =\int_{t \in[i, j]} \mathbf{q}_{w b_t} \otimes\left[\begin{array}{c} 0 \\ \frac{1}{2} \boldsymbol{\omega}^{b_t} \end{array}\right] \delta t \end{aligned} pwbjqwbj=pwbi+∫t∈[i,j]qwbtϕbtδt=∫t∈[i,j]qwbt⊗[021ωbt]δt把 qwbt=qwbi⊗qbibt\mathbf{q}_{w b_t}=\mathbf{q}_{w b_i} \otimes \mathbf{q}_{b_i b_t}qwbt=qwbi⊗qbibt 代入上式可得
pwbj=pwbi+qwbi∫t∈[i,j](qbibtϕbt)δtqwbj=qwbi[∫t∈[i,j]qbibt⊗[012ωbt]δt\mathbf{p}_{w b_j}=\mathbf{p}_{w b_i}+\mathbf{q}_{w b_i} \int_{t \in[i, j]}\left(\mathbf{q}_{b_i b_t} \boldsymbol{\phi}^{b_t}\right) \delta t \\ \mathbf{q}_{w b_j}=\mathbf{q}_{w b_i}\left[\int_{t \in[i, j]} \mathbf{q}_{b_i b_t} \otimes\left[\begin{array}{c} 0 \\ \frac{1}{2} \boldsymbol{\omega}^{b_t} \end{array}\right] \delta t\right. pwbj=pwbi+qwbi∫t∈[i,j](qbibtϕbt)δtqwbj=qwbi[∫t∈[i,j]qbibt⊗[021ωbt]δt为了整理公式,把积分相关的项用下面的式子代替
αbibj=∫t∈[i,j](qbibtϕbt)δtqbibj=∫t∈[i,j]qbibt⊗[012ωbt]δt\begin{aligned} & \boldsymbol{\alpha}_{b_i b_j}=\int_{t \in[i, j]}\left(\mathbf{q}_{b_i b_t} \boldsymbol{\phi}^{b_t}\right) \delta t \\ & \mathbf{q}_{b_i b_j}=\int_{t \in[i, j]} \mathbf{q}_{b_i b_t} \otimes\left[\begin{array}{c} 0 \\ \frac{1}{2} \boldsymbol{\omega}^{b_t} \end{array}\right] \delta t \end{aligned} αbibj=∫t∈[i,j](qbibtϕbt)δtqbibj=∫t∈[i,j]qbibt⊗[021ωbt]δt使用中值积分方法,把连续方法改成离散形式如下
ωb=12[(ωbk−bkg)+(ωbk+1−bkg)]ϕw=12(qbibkϕbk+qbibk+1ϕbk+1)\begin{aligned} \boldsymbol{\omega}^b & =\frac{1}{2}\left[\left(\boldsymbol{\omega}^{b_k}-\mathbf{b}_k^g\right)+\left(\boldsymbol{\omega}^{b_{k+1}}-\mathbf{b}_k^g\right)\right] \\ \boldsymbol{\phi}^w & =\frac{1}{2}\left(\mathbf{q}_{b_i b_k} \boldsymbol{\phi}^{b_k}+\mathbf{q}_{b_i b_{k+1}} \boldsymbol{\phi}^{b_{k+1}}\right) \end{aligned} ωbϕw=21[(ωbk−bkg)+(ωbk+1−bkg)]=21(qbibkϕbk+qbibk+1ϕbk+1)那么预积分的离散形式可以表示为
αbibk+1=αbibk+ϕwδtqbibk+1=qbibk⊗[112ωbδt]\begin{aligned} & \boldsymbol{\alpha}_{b_i b_{k+1}}=\boldsymbol{\alpha}_{b_i b_k}+\boldsymbol{\phi}^w \delta t \\ & \mathbf{q}_{b_i b_{k+1}}=\mathbf{q}_{b_i b_k} \otimes\left[\begin{array}{c} 1 \\ \frac{1}{2} \boldsymbol{\omega}^b \delta t \end{array}\right] \end{aligned} αbibk+1=αbibk+ϕwδtqbibk+1=qbibk⊗[121ωbδt]此时状态更新的公式可以整理为
[pwbjqwbjbjgj]=[pwbi+qwbiαbibjqwbiqbibjbig]\left[\begin{array}{c} \mathbf{p}_{w b_j} \\ \mathbf{q}_{w b_j} \\ \mathbf{b}_j^{g_j} \end{array}\right]=\left[\begin{array}{c} \mathbf{p}_{w b_i}+\mathbf{q}_{w b_i} \boldsymbol{\alpha}_{b_i b_j} \\ \mathbf{q}_{w b_i} \mathbf{q}_{b_i b_j} \\ \mathbf{b}_i^g \end{array}\right] pwbjqwbjbjgj=pwbi+qwbiαbibjqwbiqbibjbig残差形式可以写为
[rprqrbg]=[qwbi∗(pwbj−pwbi)−αbbbbj2[qbibj∗⊗(qwbi∗⊗qwbj)]xyzbjg−big]\left[\begin{array}{c} \mathbf{r}_p \\ \mathbf{r}_q \\ \mathbf{r}_{b g} \end{array}\right]=\left[\begin{array}{c} \mathbf{q}_{w b_i}^*\left(\mathbf{p}_{w b_j}-\mathbf{p}_{w b_i}\right)-\boldsymbol{\alpha}_{b_{b^b} b_j} \\ 2\left[\mathbf{q}_{b_i b_j}^* \otimes\left(\mathbf{q}_{w b_i}^* \otimes \mathbf{q}_{w b_j}\right)\right]_{x y z} \\ \mathbf{b}_j^g-\mathbf{b}_i^g \end{array}\right] rprqrbg=qwbi∗(pwbj−pwbi)−αbbbbj2[qbibj∗⊗(qwbi∗⊗qwbj)]xyzbjg−big其余部分(方差递推、残差对状态量雅可比、bias更新等)的推导,留作作业(仅优秀标准需要做)。