MEMS惯性测量单元误差分析与分立标定方法

由于传感器本身的物理模型和工作原理，以及加工工艺的缺陷，MEMS IMU本身就会引入一些误差。再加上不确定的工作环境的影响，MEMS IMU就是一种精度较低的传感器。为了提高传感器的精度，人们分析了它的误差来源，建立了一些误差模型来补偿这些误差。获取模型参数的过程就是所谓的标定。

一般而言，传感器的误差可以分为两类。一类是，因为传感器自身物理特性所带来的确定性误差，这类误差是可以通过一些标定的手段测量出来并补偿的。另一类是因为工作环境的影响等一些不确定性因素导致的随机误差，这类误差一般都用均值为0方差有限的高斯白噪声来表示。

本文中，我们关注产生确定性误差的来源，并对它们建模。最后介绍传统的分立标定方法。

1. 误差源

理论上来说误差模型越精确，补偿的效果就越好，但付出的代价就是标定和补偿算法的复杂度和实现难度。一般我们认为传感器的确定性误差主要有三个方面：温漂和零偏、刻度因子、和安装误差。

温漂和零偏: 零偏是指当输入信号为零时，传感器的测量值。该误差并不是一个固定值，常伴随着一些随机信号在其中，我们只能通过统计的方式来估计。一般将静止一段时间所采集的样本均值作为零偏，方差来衡量零偏稳定度。温漂一般是指传感器的零偏值随着温度的不同而有变化的现象。

实际上零偏对于我们而言并不是一个致命的误差，因为如果零偏是一个固定的值，我们在做积分运算的时候将着一个固定偏差减去就可以。实际上零偏的稳定性对于位姿的积分运算准确度起着至关重要的作用，因为其中存在很多不可能建模或者很难建模的未知因素存在。

而温度是影响零偏稳定性的一个主要因素，一般可以通过标定来解决。但是温漂的标定过程比较麻烦，有时人们会为IMU建立一个恒温的环境，这样可以从根源上干掉温度的影响，相应的也会付出成本和系统复杂度上的代价。

刻度因子: 刻度因子描述的是传感器输出与系统输入的比例关系。如果系统输入和传感器输出的数据都是在一个量纲下，那么理想情况下刻度因子应该是1.0。但是因为各种不确定性因素的存在，该值也会有一定的误差。其误差主要体现在偏移量和非线性度两个方面。

所谓的偏移量，与零偏类似，实际的刻度因子是0.999几或者1.00几的一个值。一般我们认为该偏移量是一个固定值，可以通过标定的手段测得。我们一般假设系统的输入和输出是线性关系，即可以通过$y = ax + b$的形式描述。其中$y$是传感器的输出，$x$是系统输入，$a$就是这里的刻度因子，$b$则描述了零偏。但因为各种不确定性的影响，系统可能存在一定的非线性特征，可以通过样本$x,y$的相关性来评价。

偏移量可以认为是一种确定性的误差，通过标定来修正。而非线性度描述的是一种非确定性的因素，只能通过一些统计学的手段予以估计。此外刻度因子也会受到温度等环境因素的影响。非线性度、零偏的温度稳定性、刻度因子的温度稳定性，三者是IMU选型的重要参考。

安装误差: 为了计算方便，我们一般认为传感器的测量坐标系与系统的坐标系是完全重合的，而且各个敏感轴之间不应当存在耦合误差。但实际的生产工艺是不能做到这点的，既不可能保证坐标系轴线完全对准，也不可能完全保证各轴之间的正交关系。如果一个IMU中封装了多个传感器，这些传感器之间也存在相对的偏差。

安装的细微误差都会随着系统运行时间不断的累加，而且在运动过程中还会起到放大系统不确定性因素的作用。好在这种安装误差，可以认为是一种确定的系统误差，通过标定能够一定程度上测量得到。

关于误差来源还有很多，比如说系统的振动，工作环境的湿度，甚至在地球上的不同位置都会给传感器带来不同的影响。我们不可能面面俱到的考虑这些奇奇怪怪的因素，通常把那些不容易测量或者不能测量的误差都归结为随机误差，通过一些统计学的方法来估计。

2. 误差模型

虽然温度对于MEMS IMU而言，有着很重要的影响，但是其标定过程比较复杂，有人通过提供恒温的环境来规避这一问题。我们这里的误差模型不考虑温度的影响，假设传感器的工作过程中温度是不变的。在这种假设下，我们认为MEMS IMU的误差补偿项主要有零偏、刻度因子、安装误差三个部分。

2.1 陀螺误差模型

设陀螺的输入角速度为$\vec{\omega} = \begin{bmatrix} \omega_x & \omega_y & \omega_z \end{bmatrix}^T$，输出测量值为$\vec{\omega_g} = \begin{bmatrix} \omega_{gx} & \omega_{gy} & \omega_{gz} \end{bmatrix}^T$。那么在忽略传感器的随机噪声的情况下，陀螺的误差模型可以写作$\vec{\omega_g} = \mathbf{G} \vec{\omega} + \vec{\omega_0}$:

$$ \begin{equation} \label{f1} \begin{bmatrix} \omega_{gx} \\ \omega_{gy} \\ \omega_{gz} \end{bmatrix} = \begin{bmatrix} G_x & G_{xy} & G_{xz} \\ G_{yx} & G_y & G_{yz} \\ G_{zx} & G_{zy} & G_{z} \end{bmatrix} \begin{bmatrix} \omega_x \\ \omega_y \\ \omega_z\end{bmatrix} + \begin{bmatrix} \omega_{x0} \\ \omega_{y0} \\ \omega_{z0}\end{bmatrix} \end{equation} $$

其中，$\vec{\omega_0} = \begin{bmatrix} \omega_{x0} & \omega_{y0} & \omega_{z0}\end{bmatrix}^T$是陀螺的零偏。$G_x, G_y, G_z$分别是陀螺三轴的刻度因子， $G_{ij}, i,j \in {x,y,z}$则描述了陀螺的安装误差。通过标定，我们可以测量得到系数矩阵$\mathbf{G}$和零偏向量$\vec{\omega_0}$。在实际使用陀螺输出的时候还需要通过下式进行补偿：

$$ \begin{equation} \label{f2} \vec{\omega} = \mathbf{G}^{-1} (\vec{\omega_g} - \vec{\omega_0}) \end{equation} $$

2.2 加表误差模型

类似的在忽略随机噪声的情况下，设加表的输入加速度为$\vec{a} = \begin{bmatrix} a_x & a_y & a_z \end{bmatrix}^T$，输出测量值为$\vec{a_a} = \begin{bmatrix} a_{ax} & a_{ay} & a_{az} \end{bmatrix}^T$，那么加速度计的误差模型可以写作$\vec{a_a} = \mathbf{A}\vec{a} + \vec{a_0}$:

$$ \begin{equation} \label{f3} \begin{bmatrix} a_{ax} \\ a_{ay} \\ a_{az} \end{bmatrix} = \begin{bmatrix} A_x & A_{xy} & A_{xz} \\ A_{yx} & A_y & A_{yz} \\ A_{zx} & A_{zy} & A_{z} \end{bmatrix} \begin{bmatrix} a_x \\ a_y \\ a_z\end{bmatrix} + \begin{bmatrix} a_{x0} \\ a_{y0} \\ a_{z0}\end{bmatrix} \end{equation} $$

其中，$\vec{a_0} = \begin{bmatrix} a_{x0} & a_{y0} & a_{z0}\end{bmatrix}^T$是加表的零偏。$A_x, A_y, A_z$分别是加表三轴的刻度因子， $A_{ij}, i,j \in {x,y,z}$则描述了加表的安装误差。通过标定测得系数矩阵$\mathbf{A}$和零偏向量$\vec{a_0}$。可以通过下式进行补偿：

$$ \begin{equation} \label{f4} \vec{a} = \mathbf{A}^{-1} (\vec{a_a} - \vec{a_0}) \end{equation} $$

3. 分立标定方法

分立标定是一种传统的标定方法。它需要借助转台严格控制传感器的输入量，通过直接测量传感器的输出量来完成标定。分立标定有角速率标定和多位置静态标定两种主要形式，其中角速率标定一般用于标定陀螺的刻度因子和安装误差，而多位置静态标定一般用于标定加表的刻度因子、安装误差和零偏，以及陀螺的零偏。

3.1 加表标定方法

在地球的环境下，加表总是会受到重力加速度的影响，其输出结果中也因此总是包含有重力加速度的量。利用这一信息，我们可以将加表静置在不同的姿态下，通过加表的姿态和重力加速度之间的关系，可以计算出各个误差系数。将上式$\ref{f3}$改写成如下形式，可以通过最小二乘的形式估计各个标定参数。

$$ \begin{equation}\label{f5} \begin{bmatrix} a_{ax} \\ a_{ay} \\ a_{az} \end{bmatrix} = \begin{bmatrix} A_x & A_{xy} & A_{xz} & a_{x0} \\ A_{yx} & A_y & A_{yz} & a_{y0} \\ A_{zx} & A_{zy} & A_z & a_{z0} \end{bmatrix} \begin{bmatrix} a_x \\ a_y \\ a_z \\ 1 \end{bmatrix} \end{equation} $$

假设加表正常静置时，测得重力加速度为$\vec{g} = \begin{bmatrix} 0 \\ 0 \\ g \end{bmatrix} $。通过转台调整加表位姿并静置后，如果加表的旋转矩阵为$R$，那么此时加表的输入激励应当是$\vec{a} = R\vec{g}$。通过改变加表的位姿，在多个不同的旋转矩阵的加持下，我们可以得到一组样本$(Y, X)$。其中$Y = \begin{bmatrix} \vec{a_{a1}} & \vec{a_{a2}} & \cdots & \vec{a_{an}} \\ 1 & 1 & \cdots & 1 \end{bmatrix}$，$\vec{a_{an}}$是第$n$个样本的加表输出向量。 $X = \begin{bmatrix} R_1\vec{g} & R_2\vec{g} & \cdots & R_n\vec{g} \\ 1 & 1 & \cdots & 1 \end{bmatrix}$，$R_n$是第$n$个样本的加表旋转矩阵。记$K = \begin{bmatrix} A_x & A_{xy} & A_{xz} & a_{x0} \\ A_{yx} & A_y & A_{yz} & a_{y0} \\ A_{zx} & A_{zy} & A_z & a_{z0} \end{bmatrix}$，根据最小二乘法，有$Y = \hat{K} X + \Delta$，其中$\Delta$是测量噪声，$\hat{K}$是$K$的估计：

$$ \begin{equation}\label{f6} \hat{K} = Y \bullet X^T \bullet \left(X X^T\right)^{-1} \end{equation} $$

3.2 陀螺标定方法

我们知道地球是有自转的，所以理论上也可以像加表那样，通过静置在不同的姿态下，利用地球自转特性来标定陀螺的参数。但是，MEMS陀螺的精度很差，地球自转的微弱信号通常都湮没在各种噪声中了。所以需要利用转台来输入不同的角速度，进而完成标定。

转台的转动是相对于地面的，所以安装在转台上的陀螺输入激励有两个部分：地球在惯性空间下的自转角速率记为$\vec{\omega_{e}}$，转台相对于地球的转动角速率$\vec{\omega_{p}}$。即：

$$ \begin{equation}\label{f7} \vec{\omega} = \vec{\omega_{e}} + \vec{\omega_{p}} \end{equation} $$

这里我们忽略地球自转的影响，进行标定。但仍然要求转台尽可能的水平放置，来减少其对标定结果的影响。标定过程中，每次都绕着一个轴，先以一个角速率逆时针转动，再以相同的角速率顺时针转动。记$\omega^+, \omega_g^+$分别为逆时针转动输入角速度矢量和陀螺测量矢量，$\omega^-, \omega_g^-$分别为顺时针转动输入角速度矢量和陀螺测量矢量。有$\omega^- = -\omega^+$。根据式($\ref{f1}$)，有：

$$ \begin{equation} \label{f8} \begin{bmatrix} \omega_{gx}^+ - \omega_{gx}^- \\ \omega_{gy}^+ - \omega_{gy}^- \\ \omega_{gz}^+ - \omega_{gz}^- \end{bmatrix} = 2 \begin{bmatrix} G_x & G_{xy} & G_{xz} \\ G_{yx} & G_y & G_{yz} \\ G_{zx} & G_{zy} & G_{z} \end{bmatrix} \begin{bmatrix} \omega_x^+ \\ \omega_y^+ \\ \omega_z^+\end{bmatrix} \end{equation} $$

令$y = \omega_g^+ - \omega_g^-$，$x = 2\omega^+$，我们有样本集合$(Y, X)$，其中$Y = \begin{bmatrix} y_1 & \cdots & y_n \end{bmatrix}$， $X = \begin{bmatrix} x_1 & \cdots & x_n \end{bmatrix}$。记$K = \begin{bmatrix} G_x & G_{xy} & G_{xz} \\ G_{yx} & G_y & G_{yz} \\ G_{zx} & G_{zy} & G_{z} \end{bmatrix}$，根据最小二乘法可以得到其估计$\hat{K} = Y \bullet X^T \bullet \left(X X^T\right)^{-1}$。

理论上，我们也可以像加表的标定那样，对参数矩阵进行扩展，采集不同角速度下的陀螺输出，通过最小二乘法直接解算出零偏、刻度因子和安装误差。地球自转是客观存在的，虽然我们认为它在陀螺的输出中表现很微弱，可以忽略。但有时我们还是希望标定过程应当尽可能严谨一点，通过采集相同加速率正反转的陀螺数据，并将两者相减可以约掉零偏项以及地球自转的影响。

在不考虑地球自转的影响下，我们可以将陀螺静置，采集一段时间的数据后，计算均值作为零偏$\vec{\omega_0}$。如果考虑地球的自转，人们还研究了二十四位置法、十二位置法、六位置法，等多位置标定方法。这些标定方法的本质都是借助转台精准对北，严格控制陀螺姿态，对称的摆放传感器，达到消除自转分量的目的。这里不再细述，详细可以参考2015 - MEMS陀螺接连惯导系统标定方法研究.pdf

4. 总结

本文中我们了解到MEMS惯性器件的误差源主要有零偏、刻度因子、安装误差，可以用一个系数矩阵和零偏向量来表示。分立标定方法是一种传统的标定方法，它需要借助高精度的转台，精确的控制待标传感器的姿态。通过多位置静置来标定加表的系数矩阵和零偏向量，以及陀螺的零偏。通过逐次沿各轴按照一定角速率正反转的角速率标定来获得陀螺的刻度因子和安装误差。

但是，分立标定这种传统的方法，严重依赖外部设备的精度。而且标定过程复杂，成本高，不易大批量标定。对于需要批量生产的消费级产品而言，这往往是一个致命的缺陷。