CN110186479B - 一种惯性器件误差系数确定方法 - Google Patents
一种惯性器件误差系数确定方法 Download PDFInfo
- Publication number
- CN110186479B CN110186479B CN201910464192.4A CN201910464192A CN110186479B CN 110186479 B CN110186479 B CN 110186479B CN 201910464192 A CN201910464192 A CN 201910464192A CN 110186479 B CN110186479 B CN 110186479B
- Authority
- CN
- China
- Prior art keywords
- equation
- state
- matrix
- function
- inertial device
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/10—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration
- G01C21/12—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning
- G01C21/16—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C25/00—Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass
Abstract
本发明涉及一种惯性器件误差系数确定方法,通过对非线性状态方程和观测方程的线性化处理,并利用扩展递推最小二乘法计算惯性器件误差系数,克服了现有卡尔曼滤波器估计值不唯一的不足,该方法能够适应状态方程和观测方程为非线性函数时的情况,也可满足状态方程和观测方程为线性函数时的情况。该方法能够实时计算,容易实现,具有工程实际应用的价值。
Description
技术领域
本发明涉及一种惯性器件误差系数确定方法,属于航空航天技术领域。
背景技术
当前航天飞行器的惯性导航主要采用陀螺仪和加速度计构成的捷联系统或平台系统。在实弹飞行前,需要在地面对陀螺仪和加速度计的误差系数进行标定,根据标定的结果通过误差补偿可有效提高惯性导航的使用精度。目前,经过地面标定的惯性器件,在实际飞行导航试验中,根据遥测数据计算的速度和位置的理论值仍与外测获得的真实飞行速度和位置值之间存在较大的偏差,出现所谓的“天地不一致”的情况。经分析,出现“天地不一致”的原因是地面标定方法和数据处理方法的精度不足,造成实际飞行过程中误差积累,导致飞行精度变差,因此需要对地面标定时的误差模型和数据处理方法进行修正。
现有技术中,只考虑到惯性器件的低阶线性误差模型,可采用最小二乘法或线性求解的方法分离误差系数。但该方法的缺点是对非线性方程不能适用。虽然可采用扩展卡尔曼滤波器等滤波方法估计惯性器件的参数,但受初始值取值不同的影响,会造成不同的估计收敛结果,且不同的估计收敛结果之间的差异较大。为此,需要寻找一种新的惯性器件误差系数确定方法。
发明内容
本发明要解决的技术问题是:克服现有技术的不足,提供了一种惯性器件误差系数确定方法,克服现有卡尔曼滤波器容易受到动态噪声w(t)和量测噪声v(t)的协方差矩阵Qn和Rn不准导致估值有偏差的不足。本发明方法能够非常精确的估计出各项状态变量,且相对于卡尔曼滤波器具有算法简单的特点。
本发明目的通过以下技术方案予以实现:
一种惯性器件误差系数确定方法,包括如下步骤:
S1、以惯性器件误差系数作为状态变量,建立状态方程和观测方程;
S2、设定状态变量的初始值、状态变量的协方差矩阵、采样时间,利用扩展递推最小二乘法计算惯性器件的误差系数的估计值。
上述惯性器件误差系数确定方法,所述惯性器件误差系数的状态方程和观测方程为:
式中,t为时间,x(t)为由惯性器件误差系数构成的状态向量,为x(t)相对时间t的微分,f(x(t))为由x(t)构成的状态方程函数,w(t)为动态噪声,y为观测量,h(x)为由x(t)构成的观测方程函数,v(t)为量测噪声。
上述惯性器件误差系数确定方法,在所述S1和S2之间,对S1中所述的状态方程和观测方程进行线性化,对所述状态方程函数f(x(t))和观测方程函数h(x(t))进行线性化后的函数矩阵为:
其中,为状态方程函数f(x(t))对x(t)的偏微分;为状态方程函数h(x(t))对x(t)的偏微分;A(x)为状态方程函数f(x(t))对应的线性化矩阵;C(x)为观测方程函数h(x(t))对应的线性化矩阵。
上述惯性器件误差系数确定方法,利用扩展递推最小二乘法计算惯性器件误差系数包括如下步骤:
S21、令序数n为0,对应的时刻为tn;
S23、测量得到tn+1时刻的观测量yn+1;
其中,
式中,φn+1为转移矩阵;I为单位矩阵;Δt为采样时间;
式中,cn+1为观测矩阵;
S26、计算增益矩阵Kn+1:
其中,
式中,Pn+1,n为协方差矩阵的估计值;
S28、计算状态变量的协方差矩阵Pn+1:
Pn+1=Pn+1,n-Kn+1cn+1Pn+1,n
上述惯性器件误差系数确定方法,所述惯性器件包括但不限于加速度计或陀螺仪或组合导航系统。
一种计算机可读存储介质,其上存储有计算机程序,该程序被处理器执行时,实现上述惯性器件误差系数确定方法的步骤。
本发明相比于现有技术具有如下有益效果:
(1)本发明给出了一种惯性器件误差系数确定方法,克服了传统的递推最小二乘法不能处理非线性方程的缺点,本发明既能适应线性系统,又能适应非线性系统,具有更广泛的适用性;
(2)本发明给出了一种惯性器件误差系数确定方法,相对于卡尔曼滤波器具有不依赖于初值的特点,鲁棒稳定性较好,且可实现高精度估计;
(3)本发明给出了一种惯性器件误差系数确定方法不受特定对象的限制,既适用于陀螺仪和加速度计等仪表级的测试,也适用于组合导航等系统级的误差标定;
(4)本发明给出了一种惯性器件误差系数确定方法,既可适用于地面误差系数分离,也可应用于飞行试验中遥外测误差系数分离;
(5)本发明给出了一种惯性器件误差系数确定方法,相对于最小二乘法需要对大维数矩阵的求逆运算而消耗存储资源的缺点,具有存储量小、可实时计算的优点;
(6)本发明给出了一种惯性器件误差系数确定方法,具有实现简单,易于编程的特点,具有较好的应用范围和工程价值。
附图说明
图1为本发明实施例3的流程图;
图2为实施例3中根据系数真值给出的加速度计输出误差序列值;
图3为实施例3中采用本发明给出的迭代计算过程;其中图(3a)为加速度计零偏误差的计算过程,图(3b)为加速度计线性度误差的计算过程,图(3c)为加速度计第一安装角误差的计算过程,图(3d)为加速度计第二安装角误差的计算过程;
图4为实施例3中卡尔曼滤波器选取初值P0=104、Qn=10-2、Rn=10-5时的迭代过程;其中图(4a)为加速度计零偏误差的计算过程,图(4b)为加速度计线性度误差的计算过程,图(4c)为加速度计第一安装角误差的计算过程,图(4d)为加速度计第二安装角误差的计算过程;
图5为实施例3中卡尔曼滤波器选取初值P0=104、Qn=10-5、Rn=10-5时的迭代过程;其中图(5a)为加速度计零偏误差的计算过程,图(5b)为加速度计线性度误差的计算过程,图(5c)为加速度计第一安装角误差的计算过程,图(5d)为加速度计第二安装角误差的计算过程。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚,下面将结合附图对本发明的实施方式作进一步详细描述。
实施例1:
一种惯性器件误差系数确定方法,所述惯性器件包括但不限于加速度计或陀螺仪或组合导航系统;包括如下步骤:
S1、以惯性器件误差系数作为状态变量,建立状态方程和观测方程;所述惯性器件误差系数的状态方程和观测方程为:
式中,t为时间,x(t)为由惯性器件误差系数构成的状态向量,为x(t)相对时间t的微分,f(x(t))为由x(t)构成的状态方程函数,w(t)为动态噪声,y为观测量,h(x)为由x(t)构成的观测方程函数,v(t)为量测噪声。
S2、对S1中所述的状态方程和观测方程进行线性化,对所述状态方程函数f(x(t))和观测方程函数h(x(t))进行线性化后的函数矩阵为:
其中,为状态方程函数f(x(t))对x(t)的偏微分;为状态方程函数h(x(t))对x(t)的偏微分;A(x)为状态方程函数f(x(t))对应的线性化矩阵;C(x)为观测方程函数h(x(t))对应的线性化矩阵。
S2、设定状态变量的初始值、状态变量的协方差矩阵、采样时间,利用扩展递推最小二乘法计算惯性器件的误差系数的估计值,具体包括如下步骤:
S21、令序数n为0,对应的时刻为tn;
S23、测量得到tn+1时刻的观测量yn+1;
其中,
式中,φn+1为转移矩阵;I为单位矩阵;Δt为采样时间;
式中,cn+1为观测矩阵;
S26、计算增益矩阵Kn+1:
其中,
式中,Pn+1,n为协方差矩阵的估计值;
S28、计算状态变量的协方差矩阵Pn+1:
Pn+1=Pn+1,n-Kn+1cn+1Pn+1,n
当惯性器件误差系数的状态方程和观测方程均为非线性方程时,线性化矩阵A(x)和C(x)的值取决于状态向量x(t)。当惯性器件误差系数的状态方程和观测方程均为线性方程时,线性化矩阵A(x)和C(x)的值为常值矩阵。
实施例2:
一种计算机可读存储介质,其上存储有计算机程序,该程序被处理器执行时,实现上述器件误差系数确定方法的步骤。
实施例3
以加速度计误差标定为例,具体流程如图1所示,设加速度计输出误差方程为
Ax=k0x+δkxax+kyxay+kzxaz
式中,ax、ay、az为重力加速度在加速度计壳体坐标系OXYZ三个轴OX、OY、OZ分别对应的投影分量;k0x为加速度计的零偏;δkx为加速度计的线性度;kyx和kzx分别为加速度计输入轴相对OY和OZ的安装误差角;Ax为加速度计的测量误差。
取状态变量为
由于上述4个误差系数都为常值,与时间无关,则按照本发明所述的状态方程和观测方程为:
即状态方程函数和观测方程函数分别为
h(x(t))=x1+axx2+ayx3+azx4
对状态方程函数和观测方程函数分别进行线性化,有
采用六位置分离误差系数时,6个位置如表1所示。
表1
序号 | x、y、z | a<sub>x</sub> | a<sub>y</sub> | a<sub>z</sub> |
1 | 天南东 | 1 | 0 | 0 |
2 | 南东天 | 0 | 0 | 1 |
3 | 东天南 | 0 | 1 | 0 |
4 | 北地西 | 0 | -1 | 0 |
5 | 西北地 | 0 | 0 | -1 |
6 | 地西北 | -1 | 0 | 0 |
设采样周期为1秒,在每个位置测试时间为100秒,6个位置共计600个测试数据y1、y2、…、y600。
结构矩阵为
在第一个位置,有
ci=[1 1 0 0],i=1,2,…,100
在第二个位置,有
ci=[1 0 0 1],i=101,102,…,200
在第三个位置,有
ci=[1 0 1 0],i=201,202,…,300
在第四个位置,有
ci=[1 0 -1 0],i=301,302,…,400
在第五个位置,有
ci=[1 0 0 -1],i=401,402,…,500
在第六个位置,有
ci=[1 -1 0 0],i=501,502,…,600
在进行仿真时,设加速度计误差系数的真值为k0x=1.1×10-4、δkx=-1.1×10-4、kyx=0.9×10-4、kzx=-0.9×10-4,代入加速度计误差模型后计算的加速度计输出误差值如图2所示。给定初值采用本发明的算法就可精确估计出各项误差系数,如图3所示,其中图(3a)为加速度计零偏误差的计算过程,图(3b)为加速度计线性度误差的计算过程,图(3c)为加速度计第一安装角误差的计算过程,图(3d)为加速度计第二安装角误差的计算过程。其中,4个系数的相对误差最大值为1.5×10-16,可以认为已非常精确地估计出真实值。
而在采用卡尔曼滤波器时,如果协方差矩阵Qn和Rn的值选择不合适,则会导致估计值与真实值有偏离。比如,在本例中取P0=104、Qn=10-2、Rn=10-5,则估计值为k0x=1.59×10-4、δkx=-0.606×10-4、kyx=0.9×10-4、kzx=-0.406×10-4,与真实值最大相对误差可达44.8%,估计过程如图4所示,其中图(4a)为加速度计零偏误差的计算过程,图(4b)为加速度计线性度误差的计算过程,图(4c)为加速度计第一安装角误差的计算过程,图(4d)为加速度计第二安装角误差的计算过程。而取P0=104、Qn=10-5、Rn=10-5,则可估计各项误差系数的真值,与真实值最大相对误差为4.88×10-8,估计过程如图5所示,其中图(5a)为加速度计零偏误差的计算过程,图(5b)为加速度计线性度误差的计算过程,图(5c)为加速度计第一安装角误差的计算过程,图(5d)为加速度计第二安装角误差的计算过程。可以看出,卡尔曼滤波器估计的结果与初值的选取直接相关,其收敛后的稳态值也不唯一。但采用本发明给出的算法,收敛后的稳态值唯一,且具有较好的精度。
以上所述,仅为本发明一个具体实施方式,但本发明的保护范围并不局限于此,任何熟悉本技术领域的技术人员在本发明揭露的技术范围内,可轻易想到的变化或替换,都应涵盖在本发明的保护范围之内。
本发明说明书中未作详细描述的内容属本领域技术人员的公知技术。
Claims (3)
1.一种惯性器件误差系数确定方法,其特征在于,包括如下步骤:
S1、以惯性器件误差系数作为状态变量,建立状态方程和观测方程;
S2、设定状态变量的初始值、状态变量的协方差矩阵、采样时间,利用扩展递推最小二乘法计算惯性器件的误差系数的估计值;
所述惯性器件误差系数的状态方程和观测方程为:
式中,t为时间,x(t)为由惯性器件误差系数构成的状态向量,为x(t)相对时间t的微分,f(x(t))为由x(t)构成的状态方程函数,w(t)为动态噪声,y为观测量,h(x)为由x(t)构成的观测方程函数,v(t)为量测噪声;
在所述S1和S2之间,对S1中所述的状态方程和观测方程进行线性化,对所述状态方程函数f(x(t))和观测方程函数h(x(t))进行线性化后的函数矩阵为:
其中,为状态方程函数f(x(t))对x(t)的偏微分;为状态方程函数h(x(t))对x(t)的偏微分;A(x)为状态方程函数f(x(t))对应的线性化矩阵;C(x)为观测方程函数h(x(t))对应的线性化矩阵;
利用扩展递推最小二乘法计算惯性器件误差系数包括如下步骤:
S21、令序数n为0,对应的时刻为tn;
S23、测量得到tn+1时刻的观测量yn+1;
其中,
式中,φn+1为转移矩阵;I为单位矩阵;Δt为采样时间;
式中,cn+1为观测矩阵;
S26、计算增益矩阵Kn+1:
其中,
式中,Pn+1,n为协方差矩阵的估计值;
S28、计算状态变量的协方差矩阵Pn+1:
Pn+1=Pn+1,n-Kn+1cn+1Pn+1,n
2.根据权利要求1所述的一种惯性器件误差系数确定方法,其特征在于,所述惯性器件为加速度计或陀螺仪或组合导航系统。
3.一种计算机存储介质,其上存储有计算机程序,其特征在于,该程序被处理器执行时,实现权利要求1所述方法的步骤。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910464192.4A CN110186479B (zh) | 2019-05-30 | 2019-05-30 | 一种惯性器件误差系数确定方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910464192.4A CN110186479B (zh) | 2019-05-30 | 2019-05-30 | 一种惯性器件误差系数确定方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110186479A CN110186479A (zh) | 2019-08-30 |
CN110186479B true CN110186479B (zh) | 2021-04-13 |
Family
ID=67719059
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910464192.4A Active CN110186479B (zh) | 2019-05-30 | 2019-05-30 | 一种惯性器件误差系数确定方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110186479B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111623770B (zh) * | 2020-04-28 | 2022-06-03 | 北京航天控制仪器研究所 | 一种基于速度误差开环修正的提高惯性制导精度的方法 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104019828A (zh) * | 2014-05-12 | 2014-09-03 | 南京航空航天大学 | 高动态环境下惯性导航系统杆臂效应误差在线标定方法 |
CN106643806A (zh) * | 2016-12-30 | 2017-05-10 | 深圳友铂科技有限公司 | 一种惯导系统对准精度评估方法 |
Family Cites Families (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
RU2399065C2 (ru) * | 2008-04-11 | 2010-09-10 | Корпорация Самсунг Электроникс Ко., Лтд. | Способ определения местоположения подвижного объекта посредством гибридной навигационной системы (варианты) |
CN103852085B (zh) * | 2014-03-26 | 2016-09-21 | 北京航空航天大学 | 一种基于最小二乘拟合的光纤捷联惯导系统现场标定方法 |
CN103941042B (zh) * | 2014-04-28 | 2016-03-30 | 北京航天控制仪器研究所 | 一种陀螺加速度计多位置误差系数标定方法 |
CN107273336A (zh) * | 2017-05-10 | 2017-10-20 | 北京航天控制仪器研究所 | 一种基于最小二乘法的曲线拟合方法 |
CN108663068A (zh) * | 2018-03-20 | 2018-10-16 | 东南大学 | 一种应用在初始对准中的svm自适应卡尔曼滤波方法 |
CN108593965B (zh) * | 2018-05-02 | 2019-09-13 | 福州大学 | 一种基于比力模和惯性稳定的加速度计系泊状态标定方法 |
CN109443385B (zh) * | 2018-11-13 | 2022-07-29 | 中国兵器装备集团自动化研究所有限公司 | 一种动中通天线的惯导安装误差自动标定方法 |
-
2019
- 2019-05-30 CN CN201910464192.4A patent/CN110186479B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104019828A (zh) * | 2014-05-12 | 2014-09-03 | 南京航空航天大学 | 高动态环境下惯性导航系统杆臂效应误差在线标定方法 |
CN106643806A (zh) * | 2016-12-30 | 2017-05-10 | 深圳友铂科技有限公司 | 一种惯导系统对准精度评估方法 |
Also Published As
Publication number | Publication date |
---|---|
CN110186479A (zh) | 2019-08-30 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN108507568B (zh) | 补偿温度漂移误差的方法、装置和组合导航系统 | |
CN110553641B (zh) | 一种基于相关系数的提高惯性制导精度的方法 | |
CN105387859B (zh) | Mems传感器组合温度漂移误差补偿方法 | |
CN110553642B (zh) | 一种提高惯性制导精度的方法 | |
CN111141313B (zh) | 一种提高机载局部相对姿态匹配传递对准精度的方法 | |
CN107870001A (zh) | 一种基于椭球拟合的磁力计校正方法 | |
KR101106048B1 (ko) | 센서오차의 작동 중 자동교정 방법과 이를 이용한 관성항법장치 | |
CN110579730A (zh) | 一种三轴磁传感器的转向差校正方法 | |
CN116817896B (zh) | 一种基于扩展卡尔曼滤波的姿态解算方法 | |
CN111623779A (zh) | 一种适用于噪声特性未知的时变系统自适应级联滤波方法 | |
CN112945271A (zh) | 磁力计信息辅助的mems陀螺仪标定方法及标定系统 | |
CN103913179B (zh) | 舰载单轴旋转调制捷联惯导姿态角修正方法 | |
CN110532621A (zh) | 一种飞行器气动参数在线辨识方法 | |
CN110186479B (zh) | 一种惯性器件误差系数确定方法 | |
CN110186482B (zh) | 一种提高惯性制导航天器的落点精度的方法 | |
CN110186483B (zh) | 提高惯性制导航天器落点精度的方法 | |
CN108692727B (zh) | 一种带有非线性补偿滤波器的捷联惯导系统 | |
CN112632454A (zh) | 一种基于自适应卡尔曼滤波算法的mems陀螺滤波方法 | |
CN110058324B (zh) | 利用重力场模型的捷联式重力仪水平分量误差修正方法 | |
CN111207734A (zh) | 一种基于ekf的无人机组合导航方法 | |
CN107228672B (zh) | 一种适用于姿态机动工况下的星敏和陀螺数据融合方法 | |
CN111189472A (zh) | 一种mems陀螺组合标定方法 | |
CN110186484B (zh) | 一种提高惯性制导航天器落点精度的方法 | |
CN114459478A (zh) | 一种基于姿态运动学模型的惯性测量单元数据融合方法 | |
CN107356786B (zh) | 加速度计的校准方法和装置、计算机可读存储介质 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |