CN101290326A - 石英挠性加速度计测量组件的参数辨识标定方法 - Google Patents

Abstract

本发明提供的是一种石英挠性加速度计测量组件的参数辨识标定方法。本发明的技术方案通过控制三轴位置速率转台到不同的位置来激励出石英挠性加速度计测量组件的静态误差，通过控制三轴位置速率转台以不同的速率转动激励出石英挠性加速度计测量组件的动态误差，并应用Kalman滤波器辨识出石英挠性加速度计测量组件的误差模型参数，最后利用所辨识的参数确定石英挠性加速度计测量组件的误差模型。本发明提供的技术方案对于提高石英挠性加速度计测量组件的标定精度，从而提高装配有石英挠性加速度计测量组件的捷联惯性导航系统的精度有着积极意义。

Description

石英挠性加速度计测量组件的参数辨识标定方法

(一)技术领域

本发明涉及的是一种元器件的标定方法，具体地说涉及一种石英挠性加速度计测量组件的参数辨识标定方法。

(二)背景技术

在捷联惯性导航系统中，加速度计测量组件作为其重要的传感器部件敏感载体的比力加速度。捷联惯性导航系统的测量精度与加速度计测量组件的精度密切相关。高精度的捷联惯性导航系统要求加速度计测量组件有很高的分辨率，要能清晰而精确地反映出载体的加速度，并能给出与之相应的信号，同时还要求加速度计的零位偏差尽可能小的，且数值稳定。

石英挠性加速度计具有结构工艺简单、成本低、可靠性高等优点。因此，由石英挠性加速度计构成的石英挠性加速度计测量组件被广泛地应用于捷联式惯性导航系统中。目前，提高石英挠性加速度计测量组件的精度主要有硬件、软件两条途径：

1、硬件方面对石英挠性加速度计测量组件从物理结构及工艺上进行改进；研制开发新型的、性能更为优越的石英挠性加速度计测量组件。

2、软件方面是对石英挠性加速度计测量组件进行测试，建立误差模型方程，通过误差补偿来提高石英挠性加速度计测量组件的实际使用精度。

从硬件方面提高石英挠性加速度计测量组件的精度受到加工与制造工艺水平的限制，同时也会提高成本。因此，从软件方面设计补偿石英挠性加速度计测量组件的误差补偿方法来提高其精度成为一条可行的途径。

在捷联惯性导航系统中，石英挠性加速度计直接安装在载体上，载体的恶劣动力学环境如过载、冲击、振动以及机动等都会给石英挠性加速度计测量组件和捷联惯性导航系统带来动态误差。因此误差补偿技术是提高捷联惯性导航系统精度的关键技术之一。误差补偿的关键是必须准确的取得误差的参数值，为得到有关的石英挠性加速度计测量组件的参数，就需要做精密的标定试验。

标定试验就是建立惯性传感器和惯性传感器测量组件的模型方程，利用专门的测试设备，标定出惯性传感器测量组件的误差项，并给予补偿，提高惯性传感器测量组件的使用精度。针对石英挠性加速度计测量组件，标定工作有以下几项：

1、建立与实际工程应用环境相适应的石英挠性加速度计测量组件的数学模型或误差数学模型。

2、给石英挠性加速度计测量组件精确的已知输入量。

3、观测并记录石英挠性加速度计测量组件的输出。

4、确定石英挠性加速度计测量组件的输入与输出之间的关系，从而确定其传递函数。

对于石英挠性加速度计测量组件的标定，根据标定的场所可分为内场标定和外场标定，内场标定是外场标定的基础。内场标定是指在实验室内利用惯性测试设备标定系统的参数。外场标定则是石英挠性加速度计测量组件安装在载体后进行的标定，由于安装的关系，载体实际的姿态角与系统解算出的姿态角存在系统误差。外场标定的任务即确定此项误差。

对于石英挠性加速度计测量组件的标定，根据观测量的不同可以分为分立标定法和系统级标定法。分立标定法直接利用石英挠性加速度计的输出作为观测量，一般采用最小二乘法。系统级标定则利用石英挠性加速度计的输出进行导航解算，以导航误差(位置误差、速度误差及姿态误差)作为观测量来标定系统的误差参数。

石英挠性加速度计测量组件的标定通常采用分立标定法。常规标定方法直接利用石英挠性加速度计的输出进行标定。使用该方法时，由惯性测试转台提供多个转动位置，这样可以给每个石英挠性加速度计提供不同的重力加速度分量输入，在每一个位置上采集石英挠性加速度计的输出，以此可以标定出石英挠性加速度计的参数，通常有四位置法、八位置法、十二位置法。需要标定的参数越多，则需要的测量位置越多。常规标定方法有如下一些缺点：

1、事后处理，实时性不强。

2、数据量大，需要记录的数据多，而且随着标定参数的增加，数据量剧增，耗费时间长。

3、标定精度直接取决于测试转台的精度。

因此，探索适合于石英挠性加速度计测量组件的新型标定方法，克服常规加速度计标定方法的缺点，实时而准确地给出石英挠性加速度计测量组件的误差模型具有重要意义。

(三)发明内容

本发明的目的在于提供一种能够提高石英挠性加速度计测量组件的标定精度，缩短标定时间，从而提高捷联惯性导航系统的精度的石英挠性加速度计测量组件的参数辨识标定方法。

本发明的目的是这样实现的：

步骤1：将配备石英挠性加速度计测量组件的捷联惯性导航系统放置于三轴位置速率转台上，石英挠性加速度计测量组件的X、Y、Z轴陀螺的主轴分别与转台的内、中、外框的自转轴平行，捷联惯性导航系统进行预热，然后采集陀螺仪和加速度计输出的数据；

步骤2：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理东向，y轴指向地理北向，z轴指向地理天向。即记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角：纵摇角θ、横摇角γ和航向角ψ；

以地理坐标系n系下的重力加速度[0，0，g]^T作为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_z、标度因数误差S_z、二次误差项d_z进行参数辨识；

步骤3：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理东向，y轴指向地理北向，z轴指向地理天向，以此作为初始位置。石英挠性加速度计测量组件的y轴始终朝北，绕y轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角：纵摇角θ、横摇角γ和航向角ψ；

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_yz和τ_yx进行参数辨识；

步骤4：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理东向，y轴指向地理北向，z轴指向地理天向，以此作为初始位置。石英挠性加速度计测量组件的y轴始终朝北，绕y轴将石英挠性加速度计测量组件按正方向匀速旋转，

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ和围绕三轴的旋转速率ω_x、ω_y、ω_z，

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_z进行参数辨识；

步骤5：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理正南方向，y轴水平向下，z轴指向地理天向。记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_y、标度因数误差S_y、二次误差项d_z进行参数辨识；

步骤6：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理正南方向，y轴水平向下，z轴指向地理天向，以此作为初始位置。石英挠性加速度计测量组件的x轴始终朝南，绕x轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_xz和τ_xy进行参数辨识；

步骤7：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理正南方向，y轴水平向下，z轴指向地理天向，以此作为初始位置。石英挠性加速度计测量组件的x轴始终朝南，绕x轴将石英挠性加速度计测量组件按正方向匀速旋转，旋转角速度设置为1.6度每秒，

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ、航向角ψ和围绕三轴的旋转速率ω_x、ω_y、ω_z，

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_y进行参数辨识；

步骤8：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴水平向下，y轴指向地理正西方向，z轴指向地理北向。记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_x、标度因数误差S_x、 二次误差项d_x进行参数辨识；

步骤9：操作三轴位置速率转台使石英挠性加速度计测量组件的的x轴水平向下，y轴指向地理正西方向，z轴指向地理北向，以此作为初始位置。石英挠性加速度计测量组件的z轴始终朝北，绕z轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_zx和τ_zy进行参数辨识；

步骤10：操作三轴位置速率转台使石英挠性加速度计测量组件的的x轴水平向下，y轴指向地理正西方向，z轴指向地理北向，以此作为初始位置。石英挠性加速度计测量组件的z轴始终朝北，绕z轴将石英挠性加速度计测量组件按正方向匀速旋转，旋转角速度设置为1.6度每秒，

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ、航向角ψ和围绕三轴的旋转速率ω_x、ω_y、ω_z，

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_x进行参数辨识；

步骤11：将步骤2至步骤10中经过Kalman滤波器辨识得到的参数带入到石英挠性加速度计测量组件误差模型

$\left[\begin{array}{c}{\mathrm{\δf}}_x\\ {\mathrm{\δf}}_y\\ {\mathrm{\δf}}_z\end{array}\right]=\left(\begin{array}{ccc}{S}_x& {\mathrm{\τ}}_{\mathrm{yx}}& {\mathrm{\τ}}_{\mathrm{zx}}\\ {\mathrm{\τ}}_{\mathrm{xy}}& S_y& {\mathrm{\τ}}_{\mathrm{zy}}\\ {\mathrm{\τ}}_{\mathrm{xz}}& {\mathrm{\τ}}_{\mathrm{yz}}& S_z\end{array}\right)\left[\begin{array}{c}{f}_x\\ f_y\\ f_z\end{array}\right]+\left[\begin{array}{c}{b}_x\\ b_y\\ b_z\end{array}\right]+\left(\begin{array}{ccc}{d}_x& 0& 0\\ 0& d_y& 0\\ 0& 0& d_z\end{array}\right)\left[\begin{array}{c}{f}_x^2\\ f_y^2\\ f_z^2\end{array}\right]+\left(\begin{array}{ccc}{r}_x& 0& 0\\ 0& r_y& 0\\ 0& 0& r_z\end{array}\right)\left[\begin{array}{c}{\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2\\ {\mathrm{\ω}}_z^2+{\mathrm{\ω}}_x^2\\ {\mathrm{\ω}}_x^2+{\mathrm{\ω}}_y^2\end{array}\right]---\left(1\right)$

式中下标x，y，z表示载体坐标系的三个轴，即加速度计测量组件的三个轴；δf表示石英挠性加速度计测量组件的实际输出；f表示石英挠性加速度计测量组件的理论输出；S表示加速度计测量组件的刻度因数；τ表示加速度计测量组件的安装误差；b表示加速度计测量组件的常值偏差；d表示加速度计测量组件的二次误差项；r表示加速度计测量组件的尺寸效应误差，

确定上述模型系数，按照式(1)建立起石英挠性加速度计测量组件误差模型，完成整个标定过程。

为了提高石英挠性加速度计测量组件的标定精度，本发明的技术方案通过控制三轴位置速率转台到不同的位置来激励出石英挠性加速度计测量组件的静态误差，通过控制三轴位置速率转台以不同的速率转动激励出石英挠性加速度计测量组件的动态误差，并应用Kalman滤波器辨识出石英挠性加速度计测量组件的误差模型。

平台惯性导航系统的平台自对准自标定技术是对平台惯性导航系统来说的，由于常规的捷联惯性导航系统中没有实体平台，也没有高精度的多自由度机械转动功能，因此目前的捷联惯性导航系统不能实现自标定功能。但在实验室中，可以利用三轴位置速率转台，借鉴平台惯性导航系统的自标定技术，设计适合于在捷联惯性导航系统中石英挠性加速度计测量组件的连续参数辨识标定技术，该技术对于提高石英挠性加速度计测量组件的标定精度，从而提高捷联惯性导航系统的精度有着积极意义。

本发明适用于配备石英挠性加速度计测量组件的高精度捷联惯性导航系统。

在Matlab中利用传统的24位置常规标定方法，本发明中适合于石英挠性加速度计测量组件的参数辨识标定方法对石英挠性加速度计测量组件的误差模型参数进行标定。石英挠性加速度计测量组件的误差模型参数的设定值以及两种方法的标定结果如表1、表2、表3所示。

表1误差模型参数I

表2误差模型参数II

表3误差模型参数III

在同等的测试标定条件下对石英挠性加速度计测量组件的误差模型参数进行标定时，适合于石英挠性加速度计测量组件的参数辨识标定方法的标定精度高于传统的24位置常规标定方法。

(四)附图说明

图1为石英挠性加速度计测量组件的参数辨识标定方法的步骤2中，石英挠性加速度计测量组件所处位置，同时也是步骤3和步骤4中，石英挠性加速度计测量组件所处的初始位置。

图2为石英挠性加速度计测量组件的参数辨识标定方法的步骤5中，石英挠性加速度计测量组件所处位置，同时也是步骤6和步骤7中，石英挠性加速度计测量组件所处的初始位置。

图3为石英挠性加速度计测量组件的参数辨识标定方法的步骤8中，石英挠性加速度计测量组件所处位置，同时也是步骤9和步骤10中，石英挠性加速度计测量组件所处的初始位置。

图4为石英挠性加速度计测量组件的参数辨识标定方法的流程图。

(五)具体实施方式

下面结合附图举例对本发明做更详细地描述：

步骤1、将配备石英挠性加速度计测量组件的捷联惯性导航系统放置于三轴位置速率转台上，石英挠性加速度计测量组件的X、Y、Z轴陀螺的主轴分别与转台的内、中、外框的自转轴平行。捷联惯性导航系统进行预热，然后采集陀螺仪和加速度计输出的数据。预热时间根据具体系统设定。

步骤2、操作三轴位置速率转台使捷联惯性导航系统(即石英挠性加速度计测量组件)调整到图1所示位置。记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ和航向角ψ)。

以地理坐标系n系下的重力加速度[0，0，g]^T作为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_z、标度因数误差S_z、二次误差项d_z进行参数辨识。

步骤3、以图1所示位置作为石英挠性加速度计测量组件的初始位置，石英挠性加速度计测量组件的y轴始终朝北，绕y轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ和航向角ψ)。

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_yz和τ_yx进行参数辨识。

步骤4、以图1所示位置作为石英挠性加速度计测量组件的初始位置，石英挠性加速度计测量组件的y轴始终朝北，绕y轴将石英挠性加速度计测量组件按正方向匀速旋转。为了减小加速度计动态误差对标定精度的影响，旋转角速度设置为1.6度每秒。

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ、航向角ψ)和围绕三轴的旋转速率(ω_x、ω_y、ω_z)。

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_z进行参数辨识。

步骤5、操作三轴位置速率转台使捷联惯性导航系统(即石英挠性加速度计测量组件)调整到图2所示位置。记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ和航向角ψ)。

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_y、标度因数误差S_y、二次误差项d_z进行参数辨识。

步骤6、以图2所示位置作为石英挠性加速度计测量组件的初始位置，石英挠性加速度计测量组件的x轴始终朝南，绕x轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ和航向角ψ)。

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_xz和τ_xy进行参数辨识。

步骤7、以图2所示位置作为石英挠性加速度计测量组件的初始位置，石英挠性加速度计测量组件的x轴始终朝南，绕x轴将石英挠性加速度计测量组件按正方向匀速旋转。旋转角速度设置为1.6度每秒。

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ、航向角ψ)和围绕三轴的旋转速率(ω_x、ω_y、ω_z)。

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_y进行参数辨识。

步骤8、操作三轴位置速率转台使捷联惯性导航系统(即石英挠性加速度计测量组件)调整到图3所示位置。记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ和航向角ψ)。

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_x、标度因数误差S_x、二次误差项d_x进行参数辨识。

步骤9、以图3所示位置作为石英挠性加速度计测量组件的初始位置，石英挠性加速度计测量组件的z轴始终朝北，绕z轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ和航向角ψ)。

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_zx和τ_zy进行参数辨识。

步骤10、以图3所示位置作为石英挠性加速度计测量组件的初始位置，石英挠性加速度计测量组件的z轴始终朝北，绕z轴将石英挠性加速度计测量组件按正方向匀速旋转。旋转角速度设置为1.6度每秒。

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角(即纵摇角θ、横摇角γ、航向角ψ)和围绕三轴的旋转速率(ω_x、ω_y、ω_z)。

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_x进行参数辨识。

步骤11、将步骤2至步骤10中经过Kalman滤波器辨识得到的参数带入到石英挠性加速度计测量组件误差模型

$\left[\begin{array}{c}{\mathrm{\δf}}_x\\ {\mathrm{\δf}}_y\\ {\mathrm{\δf}}_z\end{array}\right]=\left(\begin{array}{ccc}{S}_x& {\mathrm{\τ}}_{\mathrm{yx}}& {\mathrm{\τ}}_{\mathrm{zx}}\\ {\mathrm{\τ}}_{\mathrm{xy}}& S_y& {\mathrm{\τ}}_{\mathrm{zy}}\\ {\mathrm{\τ}}_{\mathrm{xz}}& {\mathrm{\τ}}_{\mathrm{yz}}& S_z\end{array}\right)\left[\begin{array}{c}{f}_x\\ f_y\\ f_z\end{array}\right]+\left[\begin{array}{c}{b}_x\\ b_y\\ b_z\end{array}\right]+\left(\begin{array}{ccc}{d}_x& 0& 0\\ 0& d_y& 0\\ 0& 0& d_z\end{array}\right)\left[\begin{array}{c}{f}_x^2\\ f_y^2\\ f_z^2\end{array}\right]+\left(\begin{array}{ccc}{r}_x& 0& 0\\ 0& r_y& 0\\ 0& 0& r_z\end{array}\right)\left[\begin{array}{c}{\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2\\ {\mathrm{\ω}}_z^2+{\mathrm{\ω}}_x^2\\ {\mathrm{\ω}}_x^2+{\mathrm{\ω}}_y^2\end{array}\right]---\left(1\right)$

式中下标x，y，z表示载体坐标系的三个轴，即加速度计测量组件的三个轴；δf表示石英挠性加速度计测量组件的实际输出；f表示石英挠性加速度计测量组件的理论输出；S表示加速度计测量组件的刻度因数；τ表示加速度计测量组件的安装误差；b表示加速度计测量组件的常值偏差；d表示加速度计测量组件的二次误差项；r表示加速度计测量组件的尺寸效应误差。

确定上述模型系数，按照式(1)建立起石英挠性加速度计测量组件误差模型，完成整个标定过程。

在步骤2至步骤10中，所述的利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的误差系数进行参数辨识，其具体步骤如下：

步骤a：针对石英挠性加速度计测量组件，建立系统的状态方程和量测方程如下：

$\stackrel{\·}{x}\left(t\right)=W\left(t\right)$ (2)

z(t)＝H(t)x(t)+V(t)

系统状态方程的状态变量为

x＝[b_x，b_y，b_z，S_x，S_y，S_z，τ_xz，τ_xy，τ_yz，τ_yx，τ_zy，τ_zx，d_x，d_y，d_z，r_x，r_y，r_z]

系统噪声为W，系统模型的方差为Q

Q＝WW^T

观测向量为

$z=\left[\begin{array}{ccc}{C}_{11}& C_{12}& C_{13}\\ C_{21}& C_{22}& C_{23}\\ C_{31}& C_{32}& C_{33}\end{array}\right]\left[\begin{array}{c}{\mathrm{\δf}}_x\\ {\mathrm{\δf}}_y\\ {\mathrm{\δf}}_z\end{array}\right]-\left[\begin{array}{c}0\\ 0\\ g\end{array}\right]$

式中的C_ij为捷联姿态矩阵C_b ⁿ中的元素；δf表示石英挠性加速度计测量组件的实际输出；g表示重力加速度。捷联姿态矩阵通过三轴位置速率转台输出的姿态角，即纵摇角θ、横摇角γ和航向角ψ获得。

$C_b^n=\left[\begin{array}{ccc}{C}_{11}& C_{12}& C_{13}\\ C_{21}& C_{22}& C_{23}\\ C_{31}& C_{32}& C_{33}\end{array}\right]$

$=\left[\begin{array}{ccc}\cos \mathrm{\ψ}\cos \mathrm{\θ}-\sin \mathrm{\θ}\sin \mathrm{\ψ}\sin \mathrm{\γ}& -\sin \mathrm{\ψ}\cos \mathrm{\θ}& \cos \mathrm{\ψ}\sin \mathrm{\γ}+\sin \mathrm{\ψ}\sin \mathrm{\γ}\cos \mathrm{\γ}\\ \cos \mathrm{\ψ}\sin \mathrm{\θ}\sin \mathrm{\γ}+\sin \mathrm{\ψ}\cos \mathrm{\γ}& \cos \mathrm{\θ}\cos \mathrm{\ψ}& \sin \mathrm{\ψ}\sin \mathrm{\γ}-\cos \mathrm{\ψ}\sin \mathrm{\γ}\cos \mathrm{\γ}\\ -\cos \mathrm{\θ}\sin \mathrm{\γ}& \cos \mathrm{\θ}& \cos \mathrm{\γ}\cos \mathrm{\θ}\end{array}\right]$

观测矩阵为

$H=\left[\begin{array}{ccccc}{C}_b^n& H_S& H_{\mathrm{\τ}}& H_d& H_r\end{array}\right]$

式中

$H_S=\left[\begin{array}{ccc}-C_{11}{C}_{31}g& -C_{12}{C}_{32}g& -C_{13}{C}_{32}g\\ -C_{21}{C}_{31}g& {-C}_{22}{C}_{32}g& -C_{23}{C}_{33}g\\ -C_{31}{C}_{31}g& -C_{32}{C}_{32}g& -C_{33}{C}_{33}g\end{array}\right]$

$H_{\mathrm{\τ}}=\left[\begin{array}{cccccc}-C_{11}{C}_{32}g& C_{11}{C}_{33}g& C_{12}{C}_{31}g& -C_{12}{C}_{33}g& -C_{13}{C}_{31}g& C_{13}{C}_{32}g\\ -C_{21}{C}_{32}g& C_{21}{C}_{33}g& C_{22}{C}_{31}g& -C_{22}{C}_{33}g& -C_{23}{C}_{31}g& C_{23}{C}_{32}g\\ -C_{31}{C}_{32}g& C_{31}{C}_{33}g& C_{32}{C}_{31}g& -C_{32}{C}_{33}g& -C_{33}{C}_{31}g& C_{33}{C}_{32}g\end{array}\right]$

$H_d=\left[\begin{array}{ccc}{C}_{11}{C}_{31}{C}_{31}{g}^2& C_{12}{C}_{32}{C}_{32}{g}^2& C_{13}{C}_{33}{C}_{33}{g}^2\\ C_{21}{C}_{31}{C}_{31}{g}^2& C_{22}{C}_{32}{C}_{32}{g}^2& C_{23}{C}_{32}{C}_{32}{g}^2\\ C_{31}{C}_{31}{C}_{31}{g}^2& C_{32}{C}_{32}{C}_{32}{g}^2& C_{33}{C}_{33}{C}_{33}{g}^2\end{array}\right]$

$H_r=-\left[\begin{array}{ccc}{C}_{11}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{12}({\mathrm{\ω}}_x^2+{\mathrm{\ω}}_z^2)& C_{13}({\mathrm{\ω}}_x^2+{\mathrm{\ω}}_y^2)\\ C_{21}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{22}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{23}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)\\ C_{31}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{32}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{33}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)\end{array}\right]$

系统噪声为V，系统模型的方差为R

R＝VV^T

步骤b：对系统的状态变量进行估计，需要对步骤a中系统方程进行离散化。离散化采用泰勒级数展开的方法

Φ(k+1，k)＝I+…    (7)

其中T为滤波周期。

对步骤a中观测方程进行离散化

H(k)＝H(t_k)    (8)

步骤c：Kalman滤波器迭代工作，估计出相应的状态变量，即石英挠性加速度计测量组件的误差模型中的误差系数。

第k+1步的量测值为z(k+1)，则x(k+1)的卡尔曼滤波估计值

按下述方程求解：

计算系统状态的一步预测

$\hat{x}(k+1|k)=\mathrm{\Φ}(k+1,k)\hat{x}\left(k\right)---\left(9\right)$

计算系统状态的估计

$\hat{x}(k+1)=\hat{x}(k+1|k)+K(k+1)[z(k+1)-H(k+1)\hat{x}(k+1|k)]---\left(10\right)$

计算卡尔曼滤波器的增益

K(k+1)＝P(k+1|k)H^T(k+1)[H(k+1)P(k+1|k)H^T(k+1)+R(k+1)]^-1  (11)

计算预测误差的方差

P(k+1|k)＝Φ(k+1，k)P(k)Φ^T(k+1，k)+Γ(k+1，k)Q(k)Γ^T(k+1，k)  (12)

计算滤波估计误差的方差

P(k+1)＝(I-K(k+1)H(k+1))P(k+1|k)    (13)。

Claims (2)

1、一种石英挠性加速度计测量组件的参数辨识标定方法，其特征是：

步骤1：将配备石英挠性加速度计测量组件的捷联惯性导航系统放置于三轴位置速率转台上，石英挠性加速度计测量组件的X、Y、Z轴陀螺的主轴分别与转台的内、中、外框的自转轴平行，捷联惯性导航系统进行预热，然后采集陀螺仪和加速度计输出的数据；

步骤2：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理东向，y轴指向地理北向，z轴指向地理天向，记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角：纵摇角θ、横摇角γ和航向角ψ；

以地理坐标系n系下的重力加速度[0，0，g]^T作为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_z、标度因数误差S_z、二次误差项d_z进行参数辨识；

步骤3：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理东向，y轴指向地理北向，z轴指向地理天向，以此作为初始位置，石英挠性加速度计测量组件的y轴始终朝北，绕y轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角：纵摇角θ、横摇角γ和航向角ψ；

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_yz和τ_yx进行参数辨识；

步骤4：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理东向，y轴指向地理北向，z轴指向地理天向，以此作为初始位置，石英挠性加速度计测量组件的y轴始终朝北，绕y轴将石英挠性加速度计测量组件按正方向匀速旋转，

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ和围绕三轴的旋转速率ω_x、ω_y、ω_z，

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_z进行参数辨识；

步骤5：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理正南方向，y轴水平向下，z轴指向地理天向，记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_y、标度因数误差S_y、二次误差项d_z进行参数辨识；

步骤6：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理正南方向，y轴水平向下，z轴指向地理天向，以此作为初始位置，石英挠性加速度计测量组件的x轴始终朝南，绕x轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_xz和τ_xy进行参数辨识；

步骤7：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴指向地理正南方向，y轴水平向下，z轴指向地理天向，以此作为初始位置，石英挠性加速度计测量组件的x轴始终朝南，绕x轴将石英挠性加速度计测量组件按正方向匀速旋转，旋转角速度设置为1.6度每秒，

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ、航向角ψ和围绕三轴的旋转速率ω_x、ω_y、ω_z，

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_y进行参数辨识；

步骤8：操作三轴位置速率转台使石英挠性加速度计测量组件的x轴水平向下，y轴指向地理正西方向，z轴指向地理北向，记录石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：常值偏差b_x、标度因数误差S_x、二次误差项d_x进行参数辨识；

步骤9：操作三轴位置速率转台使石英挠性加速度计测量组件的的x轴水平向下，y轴指向地理正西方向，z轴指向地理北向，以此作为初始位置，石英挠性加速度计测量组件的z轴始终朝北，绕z轴将石英挠性加速度计测量组件按正方向依次转动45度，连续转动7次，记录下每个位置上石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ和航向角ψ，

以n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的静态误差系数：安装误差τ_zx和τ_zy进行参数辨识；

步骤10：操作三轴位置速率转台使石英挠性加速度计测量组件的的x轴水平向下，y轴指向地理正西方向，z轴指向地理北向，以此作为初始位置，石英挠性加速度计测量组件的z轴始终朝北，绕z轴将石英挠性加速度计测量组件按正方向匀速旋转，旋转角速度设置为1.6度每秒，

以加速度计采样频率记录下石英挠性加速度计测量组件的原始输出N_x ^b、N_y ^b和N_z ^b，以及三轴位置速率转台输出的姿态角纵摇角θ、横摇角γ、航向角ψ和围绕三轴的旋转速率ω_x、ω_y、ω_z，

以地理坐标系n系下的重力加速度为外观测量，利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的动态误差系数：尺寸效应误差r_x进行参数辨识；

步骤11：将步骤2至步骤10中经过Kalman滤波器辨识得到的参数带入到石英挠性加速度计测量组件误差模型

$\left[\begin{array}{c}\mathrm{\δ}{f}_x\\ \mathrm{\δ}{f}_y\\ \mathrm{\δ}{f}_z\end{array}\right]=\left(\begin{array}{ccc}{S}_x& {\text{\τ}}_{\mathrm{yx}}& {\mathrm{\τ}}_{\mathrm{zx}}\\ {\mathrm{\τ}}_{\mathrm{xy}}& S_y& {\mathrm{\τ}}_{\mathrm{zy}}\\ {\mathrm{\τ}}_{\mathrm{xz}}& {\mathrm{\τ}}_{\mathrm{yz}}& S_z\end{array}\right)\left[\begin{array}{c}{f}_x\\ f_y\\ f_z\end{array}\right]+\left[\begin{array}{c}{b}_x\\ b_y\\ b_z\end{array}\right]+\left(\begin{array}{ccc}{d}_x& 0& 0\\ 0& d_y& 0\\ 0& 0& d_z\end{array}\right)\left[\begin{array}{c}{f}_x^2\\ f_y^2\\ f_z^2\end{array}\right]+\left(\begin{array}{ccc}{r}_x& 0& 0\\ 0& r_y& 0\\ 0& 0& r_z\end{array}\right)\left[\begin{array}{c}{\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2\\ {\mathrm{\ω}}_z^2+{\mathrm{\ω}}_x^2\\ {\mathrm{\ω}}_x^2+{\mathrm{\ω}}_y^2\end{array}\right]---\left(1\right)$

式中下标x，y，z表示载体坐标系的三个轴，即加速度计测量组件的三个轴；δf表示石英挠性加速度计测量组件的实际输出；f表示石英挠性加速度计测量组件的理论输出；S表示加速度计测量组件的刻度因数；τ表示加速度计测量组件的安装误差；b表示加速度计测量组件的常值偏差；d表示加速度计测量组件的二次误差项；r表示加速度计测量组件的尺寸效应误差，

确定上述模型系数，按照式(1)建立起石英挠性加速度计测量组件误差模型，完成整个标定过程。

2、根据权利要求1所述的石英挠性加速度计测量组件的参数辨识标定方法，其特征是：

在步骤2至步骤10中，所述的利用Kalman滤波器对石英挠性加速度计测量组件的误差模型中的误差系数进行参数辨识，其具体步骤如下：

步骤a：针对石英挠性加速度计测量组件，建立系统的状态方程和量测方程如下：

$\stackrel{\·}{x}\left(t\right)=W\left(t\right)$                        (2)

z(t)＝H(t)x(t)+V(t)

系统状态方程的状态变量为

x＝[b_x，b_y，b_z，S_x，S_y，S_z，τ_xz，τ_xy，τ_yz，τ_yx，τ_zy，τ_zx，d_x，d_y，d_z，r_x，r_y，y_z]系统噪声为W，系统模型的方差为Q

Q＝WW^T

观测向量为

$z=\left[\begin{array}{ccc}{C}_{11}& C_{12}& C_{13}\\ C_{21}& C_{22}& C_{23}\\ C_{31}& C_{32}& C_{33}\end{array}\right]\left[\begin{array}{c}\mathrm{\δ}{f}_x\\ {\mathrm{\δf}}_y\\ \mathrm{\δ}{f}_z\end{array}\right]-\left[\begin{array}{c}0\\ 0\\ g\end{array}\right]$

式中的C_ij为捷联姿态矩阵C_b ⁿ中的元素、δf表示石英挠性加速度计测量组件的实际输出、g表示重力加速度，捷联姿态矩阵通过三轴位置速率转台输出的姿态角，即纵摇角θ、横摇角γ和航向角ψ获得，

$C_b^n=\left[\begin{array}{ccc}{C}_{11}& C_{12}& C_{13}\\ C_{21}& C_{22}& C_{23}\\ C_{31}& C_{32}& C_{33}\end{array}\right]$

$=\left[\begin{array}{ccc}\cos \mathrm{\ψ}\cos \mathrm{\θ}-\sin \mathrm{\θ}\sin \mathrm{\ψ}\sin \mathrm{\γ}& -\sin \mathrm{\ψ}\cos \mathrm{\θ}& \cos \mathrm{\ψ}\sin \mathrm{\γ}+\sin \mathrm{\ψ}\sin \mathrm{\γ}\cos \mathrm{\γ}\\ \cos \mathrm{\ψ}\sin \mathrm{\θ}\sin \mathrm{\γ}+\sin \mathrm{\ψ}\cos \mathrm{\γ}& \cos \mathrm{\θ}\cos \mathrm{\ψ}& \sin \mathrm{\ψ}\sin \mathrm{\γ}-\cos \mathrm{\ψ}\sin \mathrm{\γ}\cos \mathrm{\γ}\\ -\cos \mathrm{\θ}\sin \mathrm{\γ}& \cos \mathrm{\θ}& \cos \mathrm{\γ}\cos \mathrm{\θ}\end{array}\right]$

观测矩阵为

$H=\left[\begin{array}{ccccc}{C}_b^n& G_S& H_{\mathrm{\τ}}& H_d& H_r\end{array}\right]$

式中

$H_S=\left[\begin{array}{ccc}-C_{11}{C}_{31}g& -C_{12}{C}_{32}g& -C_{13}{C}_{32}g\\ -C_{21}{C}_{31}g& -C_{22}{C}_{32}g& -C_{23}{C}_{33}g\\ -C_{31}{C}_{31}g& -C_{32}{C}_{32}g& -C_{33}{C}_{33}g\end{array}\right]$

$H_{\mathrm{\τ}}=\left[\begin{array}{cccccc}-C_{11}{C}_{32}g& C_{11}{C}_{33}g& C_{12}{C}_{31}g& -C_{12}{C}_{33}g& -C_{13}{C}_{31}g& C_{13}{C}_{32}g\\ -C_{21}{C}_{32}g& C_{21}{C}_{33}g& C_{22}{C}_{31}g& -C_{22}{C}_{33}g& -C_{23}{C}_{31}g& C_{23}{C}_{32}g\\ -C_{31}{C}_{32}g& C_{31}{C}_{33}g& C_{32}{C}_{31}g& -C_{32}{C}_{33}g& -C_{33}{C}_{31}g& C_{33}{C}_{32}g\end{array}\right]$

$H_d=\left[\begin{array}{ccc}{C}_{11}{C}_{31}{C}_{31}{g}^2& C_{12}{C}_{32}{C}_{32}{g}^2& C_{13}{C}_{33}{C}_{33}{g}^2\\ C_{21}{C}_{31}{C}_{31}{g}^2& C_{22}{C}_{32}{C}_{32}{g}^2& C_{23}{C}_{32}{C}_{32}{g}^2\\ C_{31}{C}_{31}{C}_{31}{g}^2& C_{32}{C}_{32}{C}_{32}{g}^2& C_{33}{C}_{33}{C}_{33}{g}^2\end{array}\right]$

$H_r=\left[\begin{array}{ccc}{C}_{11}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{12}({\mathrm{\ω}}_x^2+{\mathrm{\ω}}_z^2)& C_{13}({\mathrm{\ω}}_x^2+{\mathrm{\ω}}_y^2)\\ C_{21}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{22}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{23}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)\\ C_{31}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{32}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)& C_{33}({\mathrm{\ω}}_y^2+{\mathrm{\ω}}_z^2)\end{array}\right]$

系统噪声为V，系统模型的方差为R

R＝VV^T；

步骤b：对系统的状态变量进行估计，需要对步骤a中系统方程进行离散化，离散化采用泰勒级数展开的方法

Φ(k+1，k)＝I+…    (7)

其中T为滤波周期，

对步骤a中观测方程进行离散化

H(k)＝H(t_k)；    (8)

步骤c：Kalman滤波器迭代工作，估计出相应的状态变量，即石英挠性加速度计测量组件的误差模型中的误差系数，

第k+1步的量测值为z(k+1)，则x(k+1)的卡尔曼滤波估计值按下

述方程求解：

计算系统状态的一步预测

$\hat{x}(k+1|k)=\mathrm{\Φ}(k+1,k)\hat{x}\left(k\right)---\left(9\right)$

计算系统状态的估计

$\hat{x}(k+1)=\hat{x}(k+1|k)+K(k+1)[z(k+1)-H(k+1)\hat{x}(k+1|k)]---\left(10\right)$

计算卡尔曼滤波器的增益

K(k+1)＝P(k+1|k)H^T(k+1)[H(k+1)P(k+1|k)H^T(k+1)+R(k+1)]^-1(11)

计算预测误差的方差

P(k+1|k)＝Φ(k+1，k)P(k)Φ^T(k+1，k)+Γ(k+1，k)Q(k)Γ^T(k+1，k)  (12)

计算滤波估计误差的方差

P(k+1)＝(I-K(k+1)H(k+1))P(k+1|k)    (13)。

Images