[发明专利]一种适用于低精度无方位基准双轴转位设备的惯性测量单元标定方法

摘要

本发明公开了一种适用于低精度无方位基准双轴转位设备的惯性测量单元标定方法，属于惯性技术领域。本标定方法使用低精度无方位基准双轴转位设备，整个标定旋转编排需19个位置，然后以各个位置上的速度误差和天向姿态误差拟合出一阶中间参数Δ_g以及二阶中间参数最后依据中间参数与误差参数的关系，由最小二乘法计算出各个器件误差参数，为消除由转台引起的定位误差，将前一次迭代计算得到的误差参数和原有的惯性测量单元原始输出数据代入导航方程，再进行一次观测量、中间参数和误差参数残差的解算并进行残差补偿，依此类推，直至迭代计算得到的误差参数残差小于阈值。该标定方法可以大幅降低标定对转台精度的依赖性，具有很好的工程实用性。

主权项

一种适用于低精度无方位基准双轴转位设备的惯性测量单元标定方法，其特征在于：包含如下步骤：步骤一：将惯性测量单元安装在双轴转位设备上，所述惯性测量单元初始位置保证X轴朝上或朝下，惯性测量单元通电预热后开始采集输出的原始数据，所述惯性测量单元先在第0个位置上静止3‑5分钟后，再转动到第1个位置静止3‑5分钟，随后转动到第2个位置，依此类推，直至在第18个位置上静止3‑5分钟后停止采集惯性测量单元输出的原始数据；步骤二：利用步骤一采集的惯性测量单元原始数据，在第0位置上使用双矢量定姿法进行初始对准，进而得到第0位置上导航起始时刻的天向转角然后利用对准结果和第0位置上的采集数据进行导航解算，进而得到第0位置上导航过程中的实时速度以及实时天向转角θⁿ⁽⁰⁾，设第0位置上导航起始时刻的速度均为0，以速度和天向转角为观测结果拟合出第0位置上的和一阶中间参数所述包含和所述和分别为第0位置上的参数在x轴、y轴和z轴上投影的标量，所述包含和所述和分别为第0位置上的一阶中间参数在x轴、y轴和z轴上投影的标量；步骤三：利用步骤一采集的惯性测量单元原始数据，在第i位置上使用双矢量定姿法进行初始对准，所述i＝0，1，2……17，然后在第i个位置到第i+1个位置的转动过程中以及第i+1个位置上的静止过程中进行连续导航，通过导航获取转动到达第i+1个位置瞬间的速度和以及瞬间的天向转角以及转动完成后在第i+1个位置静止过程中的实时速度和以及实时天向转角θⁿ⁽ⁱ⁺¹⁾，$v_x^{n(i+1)}=v_{x0}^{n(i+1)}+{\mathrm{\Δ}}_{\mathrm{gx}}\mathrm{gT}+{\mathrm{\ω}}_{\mathrm{vx}}\frac{g{T}^2}{2}$$v_y^{n(i+1)}=v_{y0}^{n(i+1)}+{\mathrm{\Δ}}_{\mathrm{gx}}\mathrm{gT}$$v_z^{n(i+1)}=v_{z0}^{n(i+1)}+{\mathrm{\Δ}}_{\mathrm{gz}}\mathrm{gT}+{\mathrm{\ω}}_{\mathrm{vz}}\frac{g{T}^2}{2}$${\mathrm{\θ}}^{n(i+1)}={\mathrm{\θ}}_0^{n(i+1)}+{\mathrm{\ω}}_{\mathrm{vy}}T$式中：g是重力加速度，T是实时时间，ω_vx、ω_vy和ω_vz分别为系数ω_v在x轴、y轴和z轴上的分量，以速度和天向转角为观测，拟合出第i+1位置上的和一阶中间参数其中，i＝0，1，2……17，所述包含和包含和所述和分别为第i+1位置上的参数在x轴、y轴和z轴上投影的标量，所述和分别为第i+1位置上的一阶中间参数在x轴、y轴和z轴上投影的标量；步骤四：在惯性测量单元坐标系内，加速度计的误差模型为：：$\left[\begin{array}{c}\mathrm{\δ}{f}_x\\ \mathrm{\δ}{f}_y\\ \mathrm{\δ}{f}_z\end{array}\right]=\begin{array}{c}\left[\begin{array}{c}{B}_{\mathrm{ax}}\\ B_{\mathrm{ay}}\\ B_{\mathrm{az}}\end{array}\right]+\begin{array}{c}\left[\begin{array}{ccc}{K}_{\mathrm{axx}}& 0& 0\\ K_{\mathrm{ayx}}& K_{\mathrm{ayy}}& 0\\ K_{\mathrm{azx}}& K_{\mathrm{azy}}& K_{\mathrm{azz}}\end{array}\right]\left[\begin{array}{c}{f}_x\\ f_y\\ f_z\end{array}\right]+\end{array}\end{array}\left[\begin{array}{ccc}{K}_{\mathrm{ax}2}& 0& 0\\ 0& K_{\mathrm{ay}2}& 0\\ 0& 0& K_{\mathrm{az}2}\end{array}\right]\left[\begin{array}{c}{f_x}^2\\ {f_y}^2\\ {f_z}^2\end{array}\right]$上述误差模型的向量形式为：$\mathrm{\δ}{f}^b=B_a^b+K_a{f}^b+K_{a2}{\left(f^b\right)}^2$其中，f^b为载体坐标系下加速度计测得的比力，f_x^b、f_y^b和f_z^b分别为f^b在x轴、y轴和z轴上的投影，为载体坐标系下的加速度计零偏，K_a包括加速度计标度因数误差和加速度计失准角，K_a2为加速度计二次项系数，δf^b为载体坐标系下加速度计测得的比力误差；陀螺的误差模型为：$\left[\begin{array}{c}{\mathrm{\ϵ}}_x\\ {\mathrm{\ϵ}}_y\\ {\mathrm{\ϵ}}_z\end{array}\right]=\left[\begin{array}{c}{B}_{\mathrm{gx}}\\ B_{\mathrm{gy}}\\ B_{\mathrm{gz}}\end{array}\right]+\left[\begin{array}{ccc}{K}_{\mathrm{gxx}}& K_{\mathrm{gxy}}& K_{\mathrm{gxz}}\\ K_{\mathrm{gyx}}& K_{\mathrm{gyy}}& K_{\mathrm{gyz}}\\ K_{\mathrm{gzx}}& K_{\mathrm{gzy}}& K_{\mathrm{gzz}}\end{array}\right]\left[\begin{array}{c}{\mathrm{\ω}}_x\\ {\mathrm{\ω}}_y\\ {\mathrm{\ω}}_z\end{array}\right]$上述误差模型的向量形式为：${\mathrm{\ϵ}}^b={\mathrm{\ω}}_0^b+K_g{\mathrm{\ω}}^b$其中，ω^b为载体坐标系下陀螺测得的角速度，为载体坐标系下陀螺零偏，K_g包括陀螺标度因数误差和陀螺失准角，ε^b为载体坐标系下陀螺测得的角速度误差；然后将加速度计零偏B_ax、B_ay、B_az，加速度计标度因数K_axx、K_ayy、K_azz，加速度计失准角K_ayx、K_azx、K_azy，加速度计二次项系数K_ax2、K_ay2、K_az2，陀螺标度因数K_gxx、K_gyy、K_gzz，陀螺失准角K_gxy、K_gxz、K_gyx、K_gyz、K_gzx、K_gzy共计21个误差参数记为一阶误差参数K，其中，B_ax、B_ay、B_az分别为加速度计零偏B_a在x轴、y轴和z轴上投影的标量；在每个位置上，依据一阶中间参数Δ_g与一阶误差参数K的关系，利用步骤二得到的构建方程即：${\mathrm{\Δ}}_g^{\left(0\right)}=\frac{[\mathrm{\δ}{f}^{n\left(0\right)}{]}_V}{g}=\frac{{\left[C_{b\left(0\right)}^n(B_a^b+K_a{C}_n^{b\left(0\right)}{f}^n+K_{a2}{\left(C_n^{b\left(0\right)}{f}^n\right)}^2)\right]}_V}{g}$其中，δfⁿ⁽⁰⁾为导航坐标系(n系，naｖigation)下第0个位置的比力误差，g是重力加速度，是在第0个位置上从载体坐标系(b系，body)到导航坐标系的方向余弦矩阵，是在第0个位置上从导航坐标系到载体坐标系的方向余弦矩阵，为载体坐标系下的加速度计零偏，fⁿ为导航坐标系下的比力，K_a包括加速度计标度因数误差和加速度计失准角，K_a2为加速度计二次项系数，[X]_V表示垂直分量与X相同，水平分量为0的矢量；用步骤三中的构建方程其中，中i＝0，1，2…17，中i＝1，2…18，即：$\left\{\begin{array}{c}{\mathrm{\Δ}}_g^{\left(i\right)}=\frac{{-\left[\mathrm{\δ}{f}^{n(i-1)}\right]}_H+f^n\×\mathrm{\Δ}{\mathrm{\Φ}}^{n(i-1~i)}+\mathrm{\δ}{f}^{n\left(i\right)}}{g}\\ -{\left[\mathrm{\δ}{f}^{n(i-1)}\right]}_H=-{\left[C_{b(i-1)}^n(B_a^b+K_a{C}_n^{b(i-1)}{f}^n+K_{a2}{\left(C_n^{b(i-1)}{f}^n\right)}^2)\right]}_H\\ \mathrm{\δ}{f}^{n\left(i\right)}=C_{b\left(i\right)}^n(B_a^b+K_a{C}_n^{b\left(i\right)}{f}^n+K_{a2}{\left(C_n^{b\left(i\right)}{f}^n\right)}^2)\\ \mathrm{\Δ}{\mathrm{\Φ}}^{n(i-1~i)}=-\underset{t=0}{\overset{T}{\∫}}{\mathrm{\ϵ}}^n\mathrm{dt}=-\underset{t=0}{\overset{T}{\∫}}{C}_{b(i-1~i)}^n{\mathrm{\ϵ}}^b\mathrm{dt}=-\underset{t=0}{\overset{T}{\∫}}{C}_{b(i-1~i)}^n({\mathrm{\ω}}_0^b+K_g{\mathrm{\ω}}^b)\mathrm{dt}\end{array}\right.$其中，δf^n(i‑1)和δfⁿ⁽ⁱ⁾分别为导航坐标系下第i‑1和第i个位置下的比力误差，fⁿ为导航坐标系下的比力，ΔΦ^n(i‑1～i)为在导航坐标系下从第i‑1到第i个位置旋转过程中引入的姿态误差，g是重力加速度，和分别是在第i‑1个和第i个位置上从载体坐标系到导航坐标系的方向余弦矩阵，和分别是在第i‑1个和第i个位置上从导航坐标系到载体坐标系的方向余弦矩阵，为载体坐标系下的加速度计零偏，K_a包括加速度计标度因数误差和加速度计失准角，K_a2为加速度计二次项系数，[X]_H代表水平分量与X相同，垂直分量为0的矢量，为从第i‑1个位置到第i个位置的旋转过程中，载体坐标系到导航坐标系的方向余弦矩阵，εⁿ为陀螺测得的角速度误差在导航坐标系下的投影，ε^b为载体坐标系下陀螺测得的角速度误差，ω^b为载体坐标系下陀螺测得的角速度，为载体坐标系下陀螺零偏，K_g包括陀螺标度因数误差和陀螺失准角，T是实时时间；将以上方程联立得到如下方程：${\mathrm{\Δ}}_g={\left[\begin{array}{c}{\mathrm{\Δ}}_g^{\left(0\right)}\\ {\mathrm{\Δ}}_g^{\left(1\right)}\\ {\mathrm{\Δ}}_g^{\left(2\right)}\\ \·\\ \·\\ \·\\ {\mathrm{\Δ}}_g^{\left(18\right)}\end{array}\right]}_{57\×1}=\left[\begin{array}{c}{A}^{\left(0\right)}\\ A^{\left(1\right)}\\ A^{\left(2\right)}\\ \·\\ \·\\ \·\\ A^{\left(18\right)}\end{array}\right]K=\mathrm{AK}$最终构建如下方程：Δ_g＝AK步骤五：利用步骤四中的联立方程计算A的最小二乘逆矩阵进而通过计算出一阶误差参数K；步骤六：利用步骤五计算的K求解每个静止位置i对应的其中i＝0，1，2…18，计算方法如下：当i＝0时，的计算方法如下：$\left\{\begin{array}{c}{\mathrm{\ω}}_v^{0\left(0\right)}=\frac{{[f^n\×({\mathrm{\Φ}}_0^{0n\left(0\right)}\×{\mathrm{\ω}}_{\mathrm{ie}}^n)]}_H}{g}+{[{\mathrm{\Φ}}_0^{0n\left(0\right)}\×{\mathrm{\ω}}_{\mathrm{ie}}^n]}_V\\ {\mathrm{\Φ}}_0^{0n\left(0\right)}={\left[\begin{array}{ccc}\frac{\mathrm{\δ}{f}_z^{n\left(0\right)}}{g}& \frac{\tan L\·\mathrm{\δ}{f}_z^{n\left(0\right)}}{g}& -\frac{\mathrm{\δ}{f}_x^{n\left(0\right)}}{g}\end{array}\right]}^T\end{array}\right.$其中，为在第0个位置上的fⁿ为导航坐标系下的比力，为导航坐标系下第0个位置的姿态误差，包括第0个位置上初始对准引入的姿态误差，为中与一阶误差参数相关的项，为地球自转角速度在导航坐标系下的投影，g是重力加速度，[X]_H代表水平分量与X相同，垂直分量为0的矢量，[X]_V代表垂直分量与X相同，水平分量为0的矢量，[X]^T为矩阵X或向量X的转置，和分别为导航坐标系下第0个位置计算的比力误差在X轴和Z轴的分量，L为纬度；当i＝1，2……18时，的计算方法如下：$\left\{\begin{array}{c}{\mathrm{\ω}}_v^{0\left(i\right)}=\frac{{[f^n\×({\mathrm{\Φ}}_0^{0n\left(i\right)}\×{\mathrm{\ω}}_{\mathrm{ie}}^n)]}_H}{g}+{[{\mathrm{\Φ}}_0^{0n\left(i\right)}\×{\mathrm{\ω}}_{\mathrm{ie}}^n]}_V\\ {\mathrm{\Φ}}_0^{0n\left(i\right)}={\left[\begin{array}{ccc}\frac{\mathrm{\δ}{f}_z^{n(i-1)}}{g}& \frac{\tan L\·\mathrm{\δ}{f}_z^{n(i-1)}}{g}& -\frac{\mathrm{\δ}{f}_x^{n(i-1)}}{g}\end{array}\right]}^T+\mathrm{\Δ}{\mathrm{\Φ}}^{n(i-1~i)}\\ \mathrm{\δ}{f}^{n(i-1)}=C_{b(i-1)}^n(B_a^b+K_a{C}_n^{b(i-1)}{f}^n+K_{a2}{\left(C_n^{b(i-1)}{f}^n\right)}^2)\\ \mathrm{\Δ}{\mathrm{\Φ}}^{n(i-1~i)}=-\underset{t=0}{\overset{T}{\∫}}{\mathrm{\ϵ}}^n\mathrm{dt}=-\underset{t=0}{\overset{T}{\∫}}{C}_{b(i-1~i)}^n{\mathrm{\ϵ}}^b\mathrm{dt}=-\underset{t=0}{\overset{T}{\∫}}{C}_{b(i-1~i)}^n({\mathrm{\ω}}_0^b+K_g{\mathrm{\ω}}^b)\mathrm{dt}\end{array}\right.$其中，为在第i个位置上的fⁿ为导航坐标系下的比力，为导航坐标系下第i个位置的姿态误差，包括第i‑1个位置上初始对准引入的姿态误差及第i‑1个位置到第i个位置旋转过程中引入的姿态误差，为中与一阶误差参数相关的项，为地球自转角速度在导航坐标系下的投影，g是重力加速度，[X]_H代表水平分量与x轴相同，垂直分量为0的矢量，[X]_V代表垂直分量与x轴相同，水平分量为0的矢量，[X]^T为矩阵X或向量X的转置，和分别为导航坐标系下第i‑1个位置计算的比力误差在X轴和Z轴的分量，L为纬度，ΔΦⁿ(^i‑1～i)为在导航坐标系下从第i‑1到第i个位置旋转过程中引入的姿态误差，δf^n(i‑1)为导航坐标系下第i‑1个位置下的比力误差，是在第i‑1个位置上从载体坐标系到导航坐标系的方向余弦矩阵，是在第i‑1个位置上从导航坐标系到载体坐标系的方向余弦矩阵，为载体坐标系下的加速度计零偏，K_a包括加速度计标度因数误差和加速度计失准角，K_a2为加速度计二次项系数，εⁿ为陀螺测得的角速度误差在导航坐标系下的投影，ε^b为载体坐标系下陀螺测得的角速度误差，为从第i‑1个位置到第i个位置的旋转过程中，载体坐标系到导航坐标系的方向余弦矩阵，ω^b为载体坐标系下陀螺测得的角速度，为载体坐标系下陀螺零偏，K_g包括陀螺标度因数误差和陀螺失准角，T是实时时间；然后通过步骤二中的及步骤三中的计算出每个静止位置上的二阶中间参数其中，中i＝0，1，2……17，中i＝0，1，2，……，18；步骤七：将二阶误差参数即陀螺零偏B_gx、B_gy、B_gz记为列向量ω，依据二阶中间参数和二阶误差参数ω之间的关系，利用步骤六中构建方程其中i＝0，1，2……18，即：$\left\{\begin{array}{c}{\mathrm{\ω}}_v^{*\left(i\right)}=\frac{{[f^n\×({\mathrm{\Φ}}_0^{*n\left(i\right)}\×{\mathrm{\ω}}_{\mathrm{ie}}^n-{\mathrm{\ω}}_0^{n\left(i\right)})]}_H}{g}+{[{\mathrm{\Φ}}_0^{*n\left(i\right)}\×{\mathrm{\ω}}_{\mathrm{ie}}^n-{\mathrm{\ω}}_0^{n\left(i\right)}]}_V\\ {\mathrm{\Φ}}_0^{*n\left(i\right)}={\left[\begin{array}{ccc}0& -\frac{{\mathrm{\ω}}_{0z}^{n(i-1)}}{{\mathrm{\ω}}_{\mathrm{ie}}\cos L}& 0\end{array}\right]}^T\\ {\mathrm{\ω}}_0^{n(i-1)}=C_{b^{(i-1)}}^{n^8}{\mathrm{\ω}}_0^b\end{array}\right.$其中，表示第i个位置的二阶中间参数，fⁿ为导航坐标系下的比力，为导航坐标系下第i个位置的姿态误差，包括第i‑1个位置上初始对准引入的姿态误差及第i‑1个位置到第i个位置旋转过程中引入的姿态误差，为中与二阶误差参数相关的项，为地球自转角速度在导航坐标系下的投影，和分别为导航坐标系下第i‑1个和第i个位置的等效陀螺零偏，导航坐标系下第i‑1个位置的等效陀螺零偏在Z轴上的投影，ω_ie为地球自转角速率，是在第i‑1个位置上从载体坐标系到导航坐标系的方向余弦矩阵，为载体坐标系下的陀螺零偏；将以上方程联立得到如下方程：${\mathrm{\ω}}_v^{*}=\left[\begin{array}{c}{\mathrm{\ω}}_v^{*\left(0\right)}\\ {\mathrm{\ω}}_v^{*\left(1\right)}\\ {\mathrm{\ω}}_v^{*\left(2\right)}\\ .\\ .\\ .\\ {\mathrm{\ω}}_v^{*\left(18\right)}\end{array}\right]=\begin{array}{c}\left[\begin{array}{c}{B}^{\left(0\right)}\\ B^{\left(1\right)}\\ B^{\left(2\right)}\\ .\\ .\\ .\\ B^{\left(18\right)}\end{array}\right]\mathrm{\ω}=\mathrm{B\ω}\end{array}$步骤八：利用步骤七中的联立方程计算B的最小二乘逆矩阵进而通过计算出二阶误差参数ω；步骤九：当一阶误差参数K和二阶误差参数ω的残差大于阈值时，用一阶误差参数K和二阶误差参数ω残差补偿前次标定的误差参数，然后将得到的一阶误差参数K和二阶误差参数ω以及步骤一中采集的惯性测量单元原始数据代入到导航方程中，再进行一次一阶中间参数Δ_g、二阶中间参数一阶误差参数K和二阶误差参数ω残差的解算，然后对一阶误差参数K和二阶误差参数ω进行残差补偿，依此类推，经过多次迭代直至某一次迭代计算得到的一阶误差参数K和二阶误差参数ω残差小于阈值。