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

权利要求书

1.一种适用于低精度无方位基准单轴转位设备的惯性测量单元标定方法，其特征在于：包含如下步骤：

步骤一：将惯性测量单元安装在单轴转位设备上，所述惯性测量单元初始位置保证X轴朝上或朝下，惯性测量单元通电预热后开始采集输出的原始数据，所述惯性测量单元先在第0个位置上静止3-5分钟后，再转动到第1个位置静止3-5分钟，随后转动到第2个位置，依此类推，直至在第4个位置上静止3-5分钟后停止采集惯性测量单元输出的原始数据；

步骤二：利用步骤一采集的惯性测量单元原始数据，在第0位置上使用双矢量定姿法进行初始对准，进而得到第0位置上导航起始时刻的天向转角然后利用对准结果和第0位置上的采集数据进行导航解算，进而得到第0位置上导航过程中的实时速度以及实时天向转角θⁿ⁽⁰⁾，设第0位置上导航起始时刻的速度均为0，以速度和天向转角为观测结果拟合出第0位置上的和一阶中间参数所述包含和所述和分别为第0位置上的参数在x轴、y轴和z轴上投影的标量，所述包含和所述和分别为第0位置上的一阶中间参数在x轴、y轴和z轴上投影的标量；

步骤三：利用步骤一采集的惯性测量单元原始数据，在第i位置上使用双矢量定姿法进行初始对准，i＝0，1，2，3，然后在第i个位置到第i+1个位置的转动过程中以及第i+1个位置上的静止过程中进行连续导航，通过导航获取转动到达第i+1个位置瞬间的速度和和瞬间的天向转角以及转动完成后在第i+1个位置静止过程中的实时速度和和实时的天向转角θⁿ⁽ⁱ⁺¹⁾，

vxn(i+1)=vx0n(i+1)+ΔgxgT+ωvxgT22]]>

vyn(i+1)=vy0n(i+1)+ΔgygT]]>

vzn(i+1)=vz0n(i+1)+ΔgzgT+ωvzgT22]]>

θn(i+1)=θ0n(i+1)+ωvyT]]>

式中：g是重力加速度，T是实时时间，

ω_vx、ω_vy和ω_vz分别为系数ω_v在x轴、y轴和z轴上的分量，

以速度和天向转角为观测，拟合出第i+1位置上的和一阶中间参数其中，i＝0，1，2，3，所述包含和包含和所述和分别为第i+1位置上的参数在x轴、y轴和z轴上投影的标量，所述和分别为第i+1位置上的一阶中间参数在x轴、y轴和z轴上投影的标量；

步骤四：在惯性测量单元坐标系内，加速度计的误差模型为：

δfxδfyδfz=BaxBayBaz+Kaxx00KayxKayy0KazxKazyKazzfxfyfz+Kax2000Kay2000Kaz2fx2fy2fz2]]>

上述误差模型的向量形式为：

δfb=Bab+Kafb+Ka2(fb)2]]>

其中，

f^b为载体坐标系下加速度计测得的比力，

f_x^b、f_y^b和f_z^b分别为f^b在x轴、y轴和z轴上的投影，

为载体坐标系下的加速度计零偏，

K_a包括加速度计标度因数误差和加速度计失准角，

K_a2为加速度计二次项系数，

δf^b为载体坐标系下加速度计测得的比力误差；

陀螺的误差模型为：

ϵxϵyϵz=BgxBgyBgz+KgxxKgxyKgxzKgyxKgyyKgyzKgzxKgzyKgzzωxωyωz]]>

上述误差模型的向量形式为：

ϵb=ω0b+Kgω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的关系，利用步骤二得到的构建方程即：

Δg(0)=[δfn(0)]Vg=[Cb(0)n(Bab+KaCnb(0)fn+Ka2(Cnb(0)fn)2)]Vg]]>

其中，

δfⁿ⁽⁰⁾为导航坐标系(n系，navigation)下第0个位置的比力误差，

g是重力加速度，

是在第0个位置上从载体坐标系(b系，body)到导航坐标系的方向余弦矩阵，

是在第0个位置上从导航坐标系到载体坐标系的方向余弦矩阵，

为载体坐标系下的加速度计零偏，

fⁿ为导航坐标系下的比力，

K_a包括加速度计标度因数误差和加速度计失准角，

K_a2为加速度计二次项系数，

[X]_v表示垂直分量与X相同，水平分量为0的矢量；

用步骤三得到的构建方程其中中i＝0，1，2，3，中i＝1，2，3，4，即：

Δg(i)=-[δfn(i-1)]H+fn×ΔΦn(i-1~i)+δfn(i)g-[δfn(i-1)]H=-[Cb(i-1)n(Bab+KaCnb(i-1)fn+Ka2(Cnb(i-1)fn)2)]Hδfn(i)=Cb(i)n+(Bab+KaCnb(i)fn+Ka2(Cnb(i)fn)2)ΔΦn(i-1~i)=-∫t=0Tϵndt=-∫t=0TCb(i-1~i)nϵbdt=-∫t=0TC0(i-1~i)(ω0b+Kgωb)dtb]]>

其中，

δ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是实时时间；

将以上方程联立得到如下方程：

ΔgΔg(0)Δg(1)Δg(2)Δg(3)Δg(4)15×1=A(0)A(1)A(2)A(3)A(4)K=AK]]>

最终构建如下方程：

Δ_g＝AK

步骤五：根据一阶中间参数与一阶误差参数的关系式，设A＝[a₁ a₂ … a_n]、K＝[k₁ k₂ … k_n]^T，其中，a_i为矩阵A的列向量，i＝1，2…n，k_i为向量K的各元素，i＝1，2…n，[k₁ k₂ … k_n]^T为行向量[k₁ k₂ … k_n]的转置，

K中包含B_ay、K_ayy、K_azy、K_ax2、K_ay2、K_az2、K_gxx、K_gzz、K_gxy、K_gxz、K_gyx、K_gyz、K_gzx的13个误差参数直接由最近一次的完整标定结果给出，这些由外部给出的误差参数在向量K中的序号为e_l，l＝1，2，...，n_e，其余参与标定计算的误差参数在向量K中的序号为c_j，j＝1，2，...，n_c，且n_c+n_e＝n；然后将矩阵A中所有序号非c_j的列置为零后形成的矩阵设为A_cal，j＝1，2，...，n_c，将矩阵A中所有序号非e_l的列置为零后形成的矩阵设为A_ext，l＝1，2，...，n_e，将向量K中所有序号非c_j的元素置为零后形成的向量设为K_cal，j＝1，2，...，n_c，将向量K中所有序号非e_l的元素置为零后形成的向量设为K_ext，l＝1，2，...，n_e，则Δ_g＝AK可以记为：

Δg=Σi=1nkiai=Σj=1nckcjacj+Σl=1nekelael=AcalKcal+AextKext]]>

然后将A_cal去掉元素全为零的行和列得到A_ca ln z，并记下这些行和列的序号，接着求A_ca ln z的最小二乘逆矩阵并对扩充元素全为零的行和列得到其中，扩充的行序号与A_cal中全为零的列序号相同，列序号与A_cal中全为零的行序号相同，

由下式求解参与标定的误差参数：

计算得到K_cal，其中，参与标定计算的误差参数序号上元素为相应的标定计算值，外部给出的误差参数序号上元素全为零，将K_cal与外部给出的K_ext相加得到一阶误差参数K；

步骤六：利用步骤五计算的K求解每个静止位置i对应的补偿分量其中i＝0，1，2，3，4，计算方法如下：

当i＝0时，的计算方法如下：

ωv0(0)=[fn×(Φ00n(0)×ωien)]Hg+[Φ00n(0)×ωien]VΦ00n(0)=δfzn(0)gtanL·δfzn(0)g-δfxn(0)gT]]>

其中，

为在第0个位置上的

fⁿ为导航坐标系下的比力，

为导航坐标系下第0个位置的姿态误差，包括第0个位置上初始对准引入的姿态误差，

为中与一阶误差参数相关的项，

为地球自转角速度在导航坐标系下的投影，

g是重力加速度，

[X]_H代表水平分量与X相同，垂直分量为0的矢量，

[X]_V代表垂直分量与X相同，水平分量为0的矢量，

[X]^T为矩阵X或向量X的转置，

和分别为导航坐标系下第0个位置计算的比力误差在X轴和Z轴的分量，

L为纬度；

当i＝1，2，3，4时，的计算方法如下：

ωv0(i)=[fn×(Φ00n(i)×ωien)]Hg+[Φ00n(i)×ωien]VΦ00n(i)=δfzn(i-1)gtanL·δfzn(i-1)g-δfxn(i-1)gT+ΔΦn(i-1~i)δfn(i-1)=Cb(i-1)n(Bab+KaCnb(i-1)fn+Ka2(Cnb(i-1)fn)2)ΔΦn(i-1~i)=-∫t=0Tϵndt=-∫t=0TCb(i-1~i)nϵbdt=-∫t=0TCb(i-1~i)n(ω0b+Kgωb)dt]]>

其中，

为在第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为纬度，

ΔФ^n(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，3，中i＝0，1，2，3，4；

步骤七：将二阶误差参数即陀螺零偏B_gx、B_gy、B_gz记为列向量ω，依据二阶中间参数和二阶误差参数ω之间的关系，利用步骤六中构建方程其中i＝0，1，2，3，4，即：

ωv*(i)=[fn×(Φ0*n(i)×ωien-ω0n(i))]Hg+[Φ0*n(i)×ωien-ω0n(i)]VΦ0*n(i)=0-ω0zn(i-1)ωiecosL0Tω0n(i-1)=Cb(i-1)nω0b]]>

其中，

表示第i个位置的二阶中间参数，

fⁿ为导航坐标系下的比力，

为导航坐标系下第i个位置的姿态误差，包括第i-1个位置上初始对准引入的姿态误差及第i-1个位置到第i个位置旋转过程中引入的姿态误差，

为中与二阶误差参数相关的项，

为地球自转角速度在导航坐标系下的投影，

和分别为导航坐标系下第i-1个和第i个位置的等效陀螺零偏，

导航坐标系下第i-1个位置的等效陀螺零偏在Z轴上的投影，

ω_ie为地球自转角速率，

是在第i-1个位置上从载体坐标系到导航坐标系的方向余弦矩阵，

为载体坐标系下的陀螺零偏；

将以上方程联立得到如下方程：

ωv*=ωv*(0)ωv*(1)ωv*(2)ωv*(3)ωv*(4)=B(0)B(1)B(2)B(3)B(4)ω=Bω]]>

步骤八：利用步骤七中的联立方程计算B的最小二乘逆矩阵进而通过计算出二阶误差参数ω；

步骤九：当一阶误差参数K和二阶误差参数ω的残差大于阈值时，用一阶误差参数K和二阶误差参数ω残差补偿前次标定的误差参数，然后将得到的一阶误差参数K和二阶误差参数ω以及步骤一中采集的惯性测量单元原始数据代入到导航方程中，再进行一次一阶中间参数Δ_g、二阶中间参数一阶误差参数K和二阶误差参数ω残差的解算，然后对一阶误差参数K和二阶误差参数ω进行残差补偿，依此类推，经过多次迭代直至某一次迭代计算得到的一阶误差参数K和二阶误差参数ω残差小于阈值。