Tightly Coupled LiDAR Inertial Odometry and Mapping源码解析（三）

Tightly Coupled LiDAR Inertial Odometry and Mapping源码解析（三）

• 3. Joint optimization

• • 3.2 IMU preintegration中的数学基础

  • • 3.2.1 旋转矩阵、群与流形
    • 3.2.2 切空间与李代数
    • 3.2.3 指数映射
    • 3.2.4 叉积与对数映射

  • 3.3 IMU preintegration

3. Joint optimization

3.2 IMU preintegration中的数学基础

流形上的IMU预积分涉及到一些群、流行、切空间、指数与对数映射的概念，可能笔者理解能力没那些大佬那么强，又喜欢刨根问底(wtf，这公式对吗，怎么推出来的啊)，因此刚开始接触IMU预积分时，对预积分的推导理解的不是很透彻。相信也有一些同学可能会遇到类似的问题。因此，在介绍IMU预积分理论之前，先对涉及到的数学知识进行了梳理，力求严谨且容易理解。
对于其中可能出现的错误与不足，还希望看到的同学帮忙指正。

3.2.1 旋转矩阵、群与流形

旋转矩阵(Rotation Matrix) ：利用矩阵R表示旋转的一种方法，若为三维空间中的旋转，则旋转矩阵R为R_{3\times3}；若为二维空间中的旋转，则旋转矩阵R为R_{2\times2}。其他表示旋转的方法有：欧拉角，四元数，轴角等

群(Group) ：所谓群，其实就是满足一定规则的元素的集合。例如我们常见的SO(3)群和SO(2)群。我们以SO(3)群为例，其表达式就是：
\left[R\in R_{3\times3},R^TR=I,det(R)=1\right]
满足上式的3\times3矩阵就构成了一个SO(3)群。从上式中也可以看到旋转矩阵的一个非常重要的性质：就是旋转矩阵R的逆矩阵就是其转置，因此就不需要浪费时间计算R的旋转矩阵了，直接求转置就可以。

另外需要注意的是，群除了需要满足上式的要求以外，还应该具备以下性质：

• 首先是封闭性，即对于群内的任意两个元素a和b，和定义在该群上的操作符.，需要满足a.b也一定要在群内。群就像是一个孤岛，无论你咋操作都逃离不出这个群，你说气不气。对应于SO(3)群，即是任意两个旋转矩阵的乘积依然是一个旋转矩阵
• 其次是单位性，即存在唯一一个元素x，使群内的任意元素a，满足a.x=x.a=a。对于SO(3)群来说，其实就是单位矩阵I
• 最后是对于任意一个元素x，存在唯一一个y，使得x.y=y.x=I。对于SO(3)群来说，其实就是逆矩阵R^{-1}R=I

流形(Manifold) ：所谓流形，就是高维空间中的曲面。一组SO(3)群也构成了一个光滑流形。例如，假设现在有一个表示绕z轴旋转10度的旋转矩阵R:
R=\left[ \begin{matrix} 1 & 0 & 0\\ 0 & 0.98 & -0.17\\ 0 & 0.17 & 0.98 \end{matrix} \right]
我们可以把旋转矩阵R用向量v来表示，即：
v=\left[1,0,0,0,0.98,-0.17,0,0.17,0.98\right]^T
这样表示一下是不是可以理解为是9维空间内的一点，而SO(3)群里有无穷多个旋转矩阵，这些旋转矩阵共同构成了9维空间内的一个流形曲面，就像二维空间内一条直线是由其上的无穷多个点构成一样。
但是仔细想想好像不大对，明明三维空间内的旋转，怎么一下子搞到9维空间里了。那是因为9维空间里并不全是SO(3)，还有很多很多其他的9维向量，SO(3)只是其中的一个子空间。这么说可能有点抽象，可以参考一下下图。

3.2.2 切空间与李代数

切空间(tagent space) ：为解释什么是切空间，我们先来看SO(3)群中的一个单位阵I:
I=\left[ \begin{matrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 0 & 0 & 1 \end{matrix} \right]
现在我们想找该单位阵附近小区间的一个旋转矩阵，为实现该目标，我们在矩阵R的9个元素上分别添加一个无穷小量得到R
R=\left[ \begin{matrix} 1+a & b & c\\ d & 1+e & f\\ g & h & 1+i \end{matrix} \right]
为了确保该矩阵R还是旋转矩阵，将其带入R^TR=I得到：
\left[ \begin{matrix} 1+a & d & g\\ b & 1+e & h\\ c & f & 1+i \end{matrix} \right]\left[ \begin{matrix} 1+a & b & c\\ d & 1+e & f\\ g & h & 1+i \end{matrix} \right]=I
在忽略掉二阶无穷小项如a^2,e^2,i^2,ab等后，可以得到：
\left[ \begin{matrix} 1+2a & b+d & c+g\\ b+d & 1+2e & f+h\\ c+g & f+h & 1+2i \end{matrix} \right]=I
进一步化简得到：
a=e=i=0,b=-d,c=-g,f=-h\\ R=\left[ \begin{matrix} 1 & b & c\\ -b & 1 & f\\ -c & -f & 1 \end{matrix} \right]
上式可以进一步分解为：
R=I+\alpha_1G_1+\alpha_2G_2+\alpha_3G_3
其中：
G_1=\left[ \begin{matrix} 0 & 0 & 0\\ 0 & 0 & -1\\ 0 & 1 & 0 \end{matrix} \right], G_2=\left[ \begin{matrix} 0 & 0 & 1\\ 0 & 0 & 0\\ -1 & 0 & 0 \end{matrix} \right],G_3=\left[ \begin{matrix} 0 & -1 & 0\\ 1 & 0 & 0\\ 0 & 0 & 0 \end{matrix} \right]
看上式是不是有点熟悉，这不就是三个矩阵G_i的线性叠加吗。还记得我们线性代数中学过的基底吗，例如二维平面空间中的基底就是\vec{i}=[1\ 0]^T,\vec{j}=[0\ 1]^T，然后就可以利用这两个基底表示二维平面空间中的任意一点\vec{v}=x\vec{i}+y\vec{j}。只不过在现在的切空间中，基底不是向量，而是变成了矩阵。当三个系数\alpha_i不再是无穷小时，R就构成了一个三维空间内的平面。这个平面就叫做SO(3)群的切空间(tagent space) ，记为so(3)。对应的数学运算为李代数(Lie Algebra) 。

3.2.3 指数映射

指数映射(exponential map) ：对于切空间so( 3)内的任意一点A，即
A=\sum_{i} \alpha_iG_i=\left[ \begin{matrix} 0 & -\alpha_3 & \alpha_2\\ \alpha_3 & 0 & -\alpha_1\\ -\alpha_2 & \alpha_1 & 0 \end{matrix} \right]
我们来对这个点A做一个有趣的操作，对它求一个指数，也就是指数映射e^A了，我们来看一下会得到什么结果，首先我们可以利用泰勒展开求解：
e^x=1+x+\frac{x^2}{2}+\frac{x^3}{6}+...+\frac{x^n}{n!}\\ e^A=I+A+\frac{A^2}{2}+...+\frac{A^n}{n!}
多么熟悉的泰勒展开一阶近似，但是好像也没发现什么有趣的东西，但是换个角度可能就不一样了哦。关于指数函数，我们之前在微积分/高等数学中肯定还学习过这样一个极限公式：
e^x=\lim_{n\rightarrow+\infty}(1+\frac{x}{n})^n
把切空间中的A带入到上式就可以得到：
e^A=\lim_{n\rightarrow+\infty}(I+\frac{1}{n}A)^n
当n趋于无穷大时，(I+\frac{1}{n}A)趋近于3.2.2节中所述单位阵附近的旋转矩阵R了，我们已经在上一小节中证明了R是属于SO(3)群的，也就是说是一个旋转矩阵。又由于群是满足计算封闭性的，通俗来讲就是旋转矩阵相乘，无论怎么乘，乘多少次，最后的结果都还是一个旋转矩阵。因此，对(I+\frac{1}{n}A)连乘n次得到(I+\frac{1}{n}A)^n，结果肯定还是一个旋转矩阵。

Amazing！ 我们好像发现了一个不得了的事情：
对切空间 so(3)内的任意一点A求指数函数e^A，得到的结果e^A居然映射到了SO(3)内，也即：
Exponential\ map: so(3)\mapsto SO(3)

3.2.4 叉积与对数映射

叉积(Cross product) ：对于三维空间中的两个向量v_1,v_2，我们应该都知道可以对这两个向量求点积和叉积，点积是比较好求的，但是笔者之前每次求叉积的时候都要像在线性代数中求行列式的方法求，直到有一天发现原来还可以通过矩阵求。我们首先一起来看一下：
求两个向量的叉积其实可以化简为：
v_1\wedge v_2=\left[\begin{matrix} a\\ b\\ c \end{matrix}\right]\wedge \left[\begin{matrix} d\\ e\\ f \end{matrix}\right]=\left[\begin{matrix} bf-ce\\ cd - af\\ ae-bd \end{matrix}\right]=\left[\begin{matrix} 0 & -c & b\\ c & 0 & -a\\ -b & a & 0 \end{matrix}\right]\left[\begin{matrix} d\\ e\\ f \end{matrix}\right]
因此，我们可以把叉积看做是一种将向量v_1转为矩阵的操作符，即：
v_1\wedge=\left[ \begin{matrix} 0 & -c & b\\ c & 0 & -a\\ -b & a & 0 \end{matrix} \right]=a_1G_1+a_2G_2+a_3G_3=A
这就很有意思了，上一节我们把切空间中的一点表示为矩阵A，现在我们利用叉积运算符，可以直接再把切空间内的一点表示为向量的叉积符[v\wedge]。因此，SO(3)群内的一个旋转矩阵就可以表示为：
R=e^A=e^{[v\wedge]}
更有趣的是v和R的关系：向量 v的方向表示的就是旋转矩阵R对应的旋转轴，而向量v的模表示的则是绕旋转轴旋转的角度。Interesting！

根据上式我们就可以得到另外一个有趣的东西，对数映射：
log(R)=log(e^A)=log(e^{[v\wedge]})=[v\wedge]
如果我们再定义一个叉积符v \wedge的逆操作A\vee=v，即将切空间内的一点(反对称阵)转为一个向量，则可以进一步得到：
log(R)\vee=log(e^{v\wedge})\vee=v
即将旋转矩阵转为轴角向量v，最后总结一下，对数映射实现SO (3)群内一点到切空间so(3)的映射，即：
Logarithm\ map: SO(3)\mapsto so(3)

3.3 IMU preintegration

为避免篇幅过长，这部分的详细推导见下一篇博客。