CN111141313A - 一种提高机载局部相对姿态匹配传递对准精度的方法 - Google Patents

一种提高机载局部相对姿态匹配传递对准精度的方法 Download PDF

Info

Publication number
CN111141313A
CN111141313A CN202010010721.6A CN202010010721A CN111141313A CN 111141313 A CN111141313 A CN 111141313A CN 202010010721 A CN202010010721 A CN 202010010721A CN 111141313 A CN111141313 A CN 111141313A
Authority
CN
China
Prior art keywords
measurement unit
relative attitude
equation
sub
inertial measurement
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.)
Granted
Application number
CN202010010721.6A
Other languages
English (en)
Other versions
CN111141313B (zh
Inventor
吴亚丽
杨文博
张亚崇
刘庆
江巧永
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xian University of Technology
Original Assignee
Xian University of Technology
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Xian University of Technology filed Critical Xian University of Technology
Priority to CN202010010721.6A priority Critical patent/CN111141313B/zh
Publication of CN111141313A publication Critical patent/CN111141313A/zh
Application granted granted Critical
Publication of CN111141313B publication Critical patent/CN111141313B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • 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
    • G01C25/005—Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass initial alignment, calibration or starting-up of inertial devices
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/12—Simultaneous equations, e.g. systems of linear equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13—Differential equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Algebra (AREA)
  • General Engineering & Computer Science (AREA)
  • Operations Research (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Evolutionary Biology (AREA)
  • Probability & Statistics with Applications (AREA)
  • Manufacturing & Machinery (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Computing Systems (AREA)
  • Navigation (AREA)

Abstract

本发明公开了一种提高机载局部相对姿态匹配传递对准精度的方法，包括以下步骤：建立主子惯导性测量单元系统，对于子惯性测量单元陀螺漂移进行补偿，通过对随机漂移建立自回归滑动平均模型，并对自回归模型进行了卡尔曼滤波；求解主、子惯性测量单元之间相对姿态阵
Figure DDA0002357051090000011
；建立相对姿态匹配传递对准的相对姿态误差方程、相对姿态量测方程、相对姿态状态方程，得到系统状态空间模型；利用Sage_Husa自适应滤波对系统状态空间模型中的状态方程、量测方程进行解算，完成传递对准。对子惯性测量单元陀螺漂移误差进行了补偿，子惯性测量单元输出惯导参数精确度得以提高。

Description

一种提高机载局部相对姿态匹配传递对准精度的方法
技术领域
本发明属于机载惯性导航技术领域，涉及一种提高机载局部相对姿态匹配传递对准精度的方法。
背景技术
对战斗机而言，机载武器系统能否对目标进行精确打击是评价战斗机性能最主要的参考指标之一。在战斗机对目标进行打击时，各任务传感器组件的运作都需要使用实时有效的局部姿态数据，当传感器坐标系相对机体系存在明显的时变失准角时，将直接影响战斗机对目标进行精确打击的能力。目前战斗机一般都在机身重心处安装有高精度主惯性测量单元，而在重要机载任务传感器附近安装有低成本子惯性测量单元。
现有获得战斗机实时有效的局部姿态数据方法，主要是根据系统状态空间模型，通过卡尔曼滤波进行相对姿态匹配传递对准，来实时获得惯性测量单元相对载机机体系的相对姿态。然而由于低成本子惯性测量单元陀螺漂移过大，影响子惯性测量单元导航参数输出精度；采用卡尔曼滤波必须有准确的对象数学模型和噪声统计，但在实际系统中，该条件往往无法满足，由于模型参数和噪声统计具有未知性或近似性，当模型参数选取不当时，导致状态估计值增大，甚至出现滤波发散的现象。
所以亟需一种即可对陀螺随机漂移误差进行补偿，又可对噪声统计特性进行准确估计的方法，从而提高机载局部相对姿态匹配传递对准精度。
发明内容
本发明的目的是提供一种提高机载局部相对姿态匹配传递对准精度的方法，解决了机载任务传感器相对姿态匹配对准精度不高，对相对姿态误差角估计值不准确的问题。
本发明所采用的技术方案是，一种提高机载局部相对姿态匹配传递对准精度的方法，包括以下步骤：
步骤1.建立主子惯导性测量单元系统，其中主惯性测量单元安装在机体重心位置，子惯性测量单元安装在任务传感器腔体内；
步骤2.对于子惯性测量单元陀螺漂移进行补偿，通过对随机漂移建立自回归滑动平均模型，并对自回归模型进行了卡尔曼滤波，对陀螺漂移进行误差补偿；
步骤3.求解主、子惯性测量单元之间相对姿态阵
Figure RE-GDA0002414722850000021
其初值通过“双矢量定姿法”求解；
步骤4.建立相对姿态匹配传递对准的相对姿态误差方程、相对姿态量测方程、相对姿态状态方程，得到系统状态空间模型；
步骤5.利用Sage_Husa自适应滤波对系统状态空间模型中的状态方程、量测方程进行解算，估计出机载任务传感器相对姿态误差值，完成传递对准。
步骤2部分，具体过程是：
利用时间序列分析法建立自回归滑动平均模型并进行参数估计，基于该自回归滑动平均模型进行卡尔曼滤波；
步骤2.1：建立自回归滑动平均模型(ARMA(p,q))并进行参数估计；
步骤2.2：基于自回归滑动平均模型ARMA(p,q)，进行卡尔曼滤波。
步骤2.2中卡尔曼滤波包括以下五个基本方程：
2.2.1：对状态做一步预测；
2.2.2：状态一步预测的均方误差阵；
2.2.3：滤波增益矩阵；
2.2.4：状态估计方程；
2.2.5：估计均方误差方程。
步骤3部分，具体过程是：
步骤3.1：求解相对姿态阵
Figure RE-GDA0002414722850000031
首先需求出相对姿态初值通过“双矢量定姿法”求解；
步骤3.2：当机翼产生挠曲变形时，主、子惯性测量单元之间的相对姿态矩阵
Figure RE-GDA0002414722850000032
是一个时变函数，
Figure RE-GDA0002414722850000033
的动态变化过程可用姿态矩阵的微分方程来描述。
步骤3.1具体步骤为：
相对姿态初始值的确定阶段为飞机起飞离地之前，即飞机起飞前静止阶段和滑跑加速阶段；
步骤3.1.1：在载机静止条件下，采集主、子惯性测量单元一段时间内加速度计输出的比力均值，得到主子惯性测量单元比力向量；
步骤3.1.2：接着将调整载机水平姿态角(俯仰或横滚，无需确切知道大小)后静止；
步骤3.1.3：采集子惯性测量单元中加速度计输出的比力均值，得到,该位置子惯性测量单元的比力向量；
步骤3.1.4：通过子惯惯性测量单元两种情况下不共线的比力向量，构造主子惯性测量单元共同的中间过渡坐标系V，得到子惯性测量单元与该中间系的旋转变换阵，确定子惯性测量单元相对主惯惯性测量单元之间的相对姿态。
步骤3.2具体步骤为：
3.2.1：当机翼产生挠曲变形时，主、子惯性测量单元之间的相对姿态矩阵
Figure RE-GDA0002414722850000041
是一个时变函数，
Figure RE-GDA0002414722850000042
的动态变化过程可用姿态矩阵的微分方程来描述：
Figure RE-GDA0002414722850000043
式中，
Figure RE-GDA0002414722850000044
表示由括号内向量组成的反对称矩阵；
Figure RE-GDA0002414722850000045
表示子惯性测量单元相对于主惯性测量单元的角速度在主惯性测量单元系下的投影，其值可由主、子惯性测量单元的陀螺仪输出计算得到；
3.2.2：相对姿态微分方程的最终形式可表示为：
Figure RE-GDA0002414722850000046
式中，
Figure RE-GDA0002414722850000047
分别为主、子惯性测量单元陀螺仪输出值；
记
Figure RE-GDA0002414722850000048
Figure RE-GDA0002414722850000051
故，根据(11)、(12)式可得到由9个未知数构成的线性微分方程组，在初值已知的情况下，可根据龙格-库塔法对该微分方程的进行数值运算，为了简化计算量只考虑式(11)后六个方程，再根据式(13) 式中的约束条件对其进行求解，最终确定出
Figure RE-GDA0002414722850000052
步骤4部分，具体过程是：
步骤4.1：采用二阶马尔科夫过程描述挠曲变形；
步骤4.2：描述主子惯性测量单元相对姿态误差值，即相对姿态角的计算值与实际值之间的差值；
步骤4.3：建立相对姿态量测方程；
步骤4.4：建立相对姿态状态方程；
最终得到系统状态空间模型。
步骤5部分，具体过程是：
利用Sage_Husa自适应滤波，进行相对姿态匹配传递对准； Sage_Husa自适应滤波方程由以下八个基本方程构成，分别为：
步骤5.1：计算加权系数dk；
步骤5.2：状态一步预测方程；
步骤5.3：一步预测的均方误差阵；
步骤5.4：新息序列；
步骤5.5：估计量测噪声；
步骤5.6：滤波增益矩阵；
步骤5.7：状态估计方程；
步骤5.8：估计均方误差方程；
对姿态的估计值进行校正处理，得到更优的估计值。
本发明的有益效果是：
从两个方面对相对姿态匹配过程进行优化，以此提高对载机局部相对姿态匹配传递对准精度。
第一方面，对子惯性测量单元进行陀螺漂移补偿，通过建立自回归滑动平均模型，使用卡尔曼滤波对自回归滑动平均模型进行处理，提高子惯性测量单元导航参数输出精度。
第二方面，依据系统状态空间模型，采用Sage_Husa自适应滤波，实时估计系统的噪声参数，以此降低系统状态误差估计值。
通过上述方法，将有效提高机载局部相对姿态匹配传递对准精度。
1)由于对子惯性测量单元陀螺漂移误差进行了补偿，子惯性测量单元输出惯导参数精确度得以提高。
2)在动态随机系统噪声统计特性参数未知的情况下，采用 Sage_Husa自适应滤波，确保其在进行状态估计的同时还可以通过量测输出系统的噪声参数。
附图说明
图1是本发明提高机载局部相对姿态匹配传递对准的流程框图；
图2是本发明子惯性测量单元进行陀螺随机漂移误差补偿流程框图；
图3是本发明中相关坐标系关系示意图；
具体实施方式
下面结合附图和具体实施方式对本发明进行详细说明。
参照图1，本发明一种提高机载局部相对姿态匹配传递对准精度的方法，将主惯性测量单元的导航参数作为后续步骤的参考基准，对子惯性测量单元陀螺漂移进行补偿后，结合主、子惯性测量单元输出导航数据，利用双矢量定姿法求得相对姿态初始值，并进行相对姿态解算得到相对姿态误差方程，然后建立系统状态空间模型，这一过程中需考虑挠曲变形、子惯性测量单元陀螺漂移输出噪声等影响因素。为了得到系统状态空间模型中相对姿态误差估计值，采用Sage_Husa 自适应滤波器进行相对姿态匹配传递对准，使用Sage_Husa自适应滤波器，可在进行相对姿态估计的同时，通过量测输出实时地估计噪声统计参数，对加快滤波收敛速度、防止滤波发散具有很好的效果。
基于上述原理，参照图1所示流程，具体包括以下步骤：
步骤1.将主子惯性测量单元分别安装在机身(重心)及机翼位置处，组成主子惯性测量单元传递对准结构。战斗机实际飞行过程中，需要具有较高精度的相对姿态估计结果，待估计值达到要求的相对姿态估计精度后，对实时计算的姿态阵进行修正。
步骤2.对于子惯性测量单元陀螺漂移进行补偿，通过对随机漂移建立自回归滑动平均模型，并对自回归模型进行了卡尔曼滤波，对陀螺漂移进行误差补偿。
本例中采用时间序列分析法建立自回归滑动平均模型并进行参数估计，基于该自回归滑动平均模型进行卡尔曼滤波，具体包括：
步骤2.1：建立自回归滑动平均模型，并进行模型参数估计：
自回归滑动平均模型(ARMA(p,q))形式为：
xk-φ1xk-1-…-φpxk-p＝ak-θ1ak-1-…-θqak-q (1)
其中，xk为随机时间序列；p、q分别为自回归阶数和滑动平均阶数；φii(i＝1,2,…,p)、θi(i＝1,2,…,q)分别为自回归系数和滑动平均系数；ak是均值为零方差σ2的白噪声。
为了便于使用Yule-walker法，在此将ARMA(p,q)模型参数写成矩阵的形式：
Figure RE-GDA0002414722850000081
其中，式中所涉及的参数：Y＝[xk xk+1…xn]T，
Figure RE-GDA0002414722850000082
Figure RE-GDA0002414722850000083
用Yule-walker法估计(2)式中
Figure RE-GDA0002414722850000084
和θ，可得到相对应的估计值
Figure RE-GDA0002414722850000085
根据
Figure RE-GDA0002414722850000086
的值可确定状态空间模型中状态转移矩阵和噪声驱动矩阵。
参照图2，基于选取的ARMA(p,q)模型，建立系统状态空间模型，利用卡尔曼滤波算法来抑制误差，实现对陀螺漂移误差的补偿。这里选择的陀螺漂移模型为三阶，即自回归系数和滑动平均系数分别为3 和2时的情况，其模型可表示为：
Figure RE-GDA0002414722850000091
其中选取状态、观测量分别为：Xk＝[yk yk-1…yk-2]T， Zk＝yk，据此建立出该模型的状态方程和量测方程，为卡尔曼滤波算法提供必要的动态系统模型，该动态系统模型可表示为：
Figure RE-GDA0002414722850000092
式中，A、B、C为是已知的系统结构参数，Xk是n维状态向量、Zk是m维量测向量，Wk-1和Vk-1分别是均值为零方差Qk和Rk的白噪声且互不相关。
步骤2.2：基于自回归滑动平均模型ARMA(p,q)，进行卡尔曼滤波；根据得到的系统状态方程以及量测方程，构建卡尔曼滤波算法来抑制陀螺漂移误差，卡尔曼滤波器由以下五个基本方程构成：
2.2.1)对状态做一步预测：
Figure RE-GDA0002414722850000093
其中，X为步骤2.1中的状态向量，右下标k表示在k时刻的状态向量，上标
Figure RE-GDA0002414722850000094
表示对X的状态一步估值(下同)，A为一步转移矩阵，
Figure RE-GDA0002414722850000095
是根据
Figure RE-GDA0002414722850000096
计算得到的对Xk的一步状态预测；
2.2.2)状态一步预测的均方误差阵：
Pk/k-1＝APk-1AT+BQk-1BT (6)
式(6)对预测的质量优劣做了定量描述。其中，B为系统噪声驱动阵，Qk为系统噪声序列的方差阵，为非负定阵，Pk表示k时刻的均方误差，Pk-1表示k-1时刻的均方误差，Pk/k-1为由tk-1时刻至tk时刻的均方误差，Q对应表示与P相同；
2.2.3)滤波增益矩阵：
Kk＝Pk/k-1CT(CPk/k-1CT+Rk)-1 (7)
式(7)对滤波增益矩阵进行求解。其中，C为量测阵，Rk为量测噪声序列的方差阵；
2.2.4)状态估计方程：
Figure RE-GDA0002414722850000101
其中，Kk为滤波增益矩阵，Zk＝yk为量测方程；
2.2.5)估计均方误差方程：
Pk＝(I-KkC)Pk/k-1 (9)
其中，I为单位阵；
基于上述方程，在给定初始状态向量和误差方差阵的情况下，可对陀螺各个时刻输出状态值进行估计，通过卡尔曼滤波器可有效减小随机误差，降低随机漂移，提高导航参数输出精度。
步骤3.求解相对姿态阵
Figure RE-GDA0002414722850000102
相对姿态初值通过“双矢量定姿法”求解。
步骤3.1：相对姿态初始值的确定阶段为飞机起飞离地之前，即飞机起飞前静止阶段和滑跑加速阶段。
步骤3.1.1：在载机静止条件下，采集主、子惯性测量单元一段时间内加速度计输出的比力均值，得到主子惯性测量单元比力向量；
步骤3.1.2：接着将调整载机水平姿态角(俯仰或横滚，无需确切知道大小)后静止；
步骤3.1.3：采集子惯性测量单元中加速度计输出的比力均值，得到,该位置子惯性测量单元的比力向量；
步骤3.1.4：通过子惯惯性测量单元两种情况下不共线的比力向量，构造主子惯性测量单元共同的中间过渡坐标系V，得到子惯性测量单元与该中间系的旋转变换阵，确定子惯性测量单元相对主惯惯性测量单元之间的相对姿态；
步骤3.2：当机翼产生挠曲变形时，主、子惯性测量单元之间的相对姿态矩阵
Figure RE-GDA0002414722850000111
是一个时变函数，
Figure RE-GDA0002414722850000112
的动态变化过程可用姿态矩阵的微分方程来描述：
Figure RE-GDA0002414722850000113
式中，
Figure RE-GDA0002414722850000114
表示由括号内向量组成的反对称矩阵。
Figure RE-GDA0002414722850000115
表示子惯性测量单元相对于主惯性测量单元的角速度在主惯性测量单元系下的投影，其值可由主、子惯性测量单元的陀螺仪输出计算得到。
相对姿态微分方程的最终形式可表示为：
Figure RE-GDA0002414722850000116
式中，
Figure RE-GDA0002414722850000117
分别为主、子惯性测量单元陀螺仪输出值。
记
Figure RE-GDA0002414722850000118
Figure RE-GDA0002414722850000119
故，根据(11)、(12)式可得到由9个未知数构成的线性微分方程组，在初值已知的情况下，可根据龙格-库塔法对该微分方程的进行数值运算，为了简化计算量只考虑式(12)后六个方程，再根据式(13) 式中的约束条件对其进行求解，最终确定出
Figure RE-GDA0002414722850000121
步骤4：建立“相对姿态匹配”传递对准的相对姿态误差方程、相对姿态量测方程、相对姿态状态方程，得到系统状态空间模型，为后续利用Sage_Husa自适应滤波进行相对姿态匹配对准做准备。
步骤4.1：采用二阶马尔科夫过程描述挠曲变形，表达式为：
Figure RE-GDA0002414722850000122
其中，
Figure RE-GDA0002414722850000123
Var(w)＝Q＝4β3σ2，θ,ω为变形角和角速率，σ变形角方差，τ二阶马尔科夫过程的相关时间，w,Q激励噪声及其方差。
步骤4.2：主子惯性测量单元相对姿态误差指的是其相对姿态角的计算值与实际值之间的差值。其表达式为：
Figure RE-GDA0002414722850000124
式中，η为相对姿态误差，
Figure RE-GDA0002414722850000125
为子惯性测量单元陀螺的随机常值漂移(经上述步骤2补偿后，随机漂移误差已降低)，
Figure RE-GDA0002414722850000126
为白噪声。
步骤4.3：建立相对姿态量测方程，其表达式为：
构造量测阵表达式为：
Figure RE-GDA0002414722850000127
考虑h、ξ和q为小角度，其坐标系关系示意图如图3所示，其中a系与m系之间为安装误差角关系，s系和a系之间为弹性变形角关系q，
Figure RE-GDA0002414722850000137
系表示计算相对参考系，相对姿态初值对应的相对姿态阵初值为
Figure RE-GDA0002414722850000131
相应的坐标系为
Figure RE-GDA0002414722850000138
系，它与a系之间偏差角是固定的安装角估计误差ξ。则量测方程为：
Figure RE-GDA0002414722850000132
取相对姿态量测向量为：
Figure RE-GDA0002414722850000133
可以看出，量测向量中包含相对姿态误差角h、初始相对姿态误差角ξ、挠曲变形角θ。
步骤4.4：相对姿态状态方程
选取相对姿态匹配传递对准的系统状态为：
Figure RE-GDA0002414722850000134
状态方程为：
Figure RE-GDA0002414722850000135
由上式可得系统的状态方程为：
Figure RE-GDA0002414722850000136
式中，
Figure RE-GDA0002414722850000141
Figure RE-GDA0002414722850000142
Figure RE-GDA0002414722850000143
其中，
Figure RE-GDA0002414722850000144
为m系(主惯性测量单元坐标系)相对于i系(惯性坐标系)的运动角速率，
Figure RE-GDA0002414722850000145
为
Figure RE-GDA0002414722850000146
的反对称矩阵(下同)，
Figure RE-GDA0002414722850000147
Figure RE-GDA0002414722850000148
为子惯性测量单元x,y,z三轴随机游走误差，
Figure RE-GDA0002414722850000149
为加在其上激励噪声；
系统的量测方程为：
Z＝HX+V (25)
其中，H量测矩阵、量测噪声V分别为：
Figure RE-GDA00024147228500001410
Figure RE-GDA00024147228500001411
综上，根据相对姿态误差方程、相对姿态量测方程、相对姿态状态方程，可得到系统状态空间模型。
步骤5.利用Sage_Husa自适应滤波对系统状态空间模型中的状态方程、量测方程进行解算，估计出机载任务传感器相对姿态误差值，完成传递对准。
在实际应用中，系统噪声统计特性是未知的、近似的或部分已知的，采用Sage_Husa自适应滤波通过量测输出实时地估计系统的噪声参数，可有效降低系统状态误差、抑制系统滤波发散、加快收敛速度。
根据步骤4中得到的系统相对姿态状态空间模型，利用 Sage_Husa自适应滤波，进行相对姿态匹配传递对准。Sage_Husa自适应滤波方程由以下八个基本方程构成，分别为：
步骤5.1：计算加权系数dk：
dk＝(1-b)/(1-bk) (28)
式中，b为遗忘因子，取值范围在(0,1)之间，右上标k表示k时刻；
步骤5.2：状态一步预测方程：
Figure RE-GDA0002414722850000151
其中，X为步骤3.3.4)中的状态变量，右下标k表示在k时刻的状态量，上标
Figure RE-GDA0002414722850000152
表示对X的状态一步估值(下同)，Φk/k-1为tk-1时刻至tk时刻一步转移矩阵，
Figure RE-GDA0002414722850000153
是根据
Figure RE-GDA0002414722850000154
计算得到的对
Figure RE-GDA0002414722850000155
的一步状态预测；
步骤5.3：一步预测的均方误差阵：
Figure RE-GDA0002414722850000156
其中，Γk-1为系统噪声驱动阵，Qk为系统噪声序列的方差阵，为非负定阵，Pk表示k时刻的均方误差，Pk-1表示k-1时刻的均方误差， Pk/k-1为由tk-1时刻至tk时刻的均方误差，Q对应表示与P相同。
步骤5.4：新息序列：
νk＝Zk-HkXk/k-1 (31)
其中，Zk为步骤3.3.4中求得的量测方程，Hk为量测阵；
步骤5.5：估计量测噪声：
Figure RE-GDA0002414722850000161
其中，Rk-1为由tk-1时刻量测均方误差,初值R0可根据实际情况选取。
步骤5.6：滤波增益矩阵：
Figure RE-GDA0002414722850000162
其中，Rk为估计量测噪声序列的方差阵；
步骤5.7：状态估计方程：
Xk＝Xk/k-1+Kkvk (34)
其中，Kk为滤波增益方程；
步骤5.8：估计均方误差方程：
Pk＝(I-KkHk)Pk/k-1 (35)
其中，I为单位阵；
本发明的方法，为提高机载局部相对姿态匹配传递对准精度，即提高相对姿态角估计值。在步骤2中，对低精度子惯性测量单元输出导航参数误差较大的缺陷，通过对陀螺随机漂移建立自回归滑动平均模型，并利用Yule-walker法对模型参数进行估计，然后通过卡尔曼滤波有效减少了随机误差，提高了子惯性测量单元输出导航参数精度；
在步骤3中，通过“双矢量定姿”法，求得主、子惯性测量单元之间相对姿态初值后，在步骤4中，考虑机体挠曲的影响以及导航坐标系与真实导航坐标系之间存在误差等因素的影响，得到相对姿态误差方程，并求得相应量测方程、状态方程，得到系统的状态空间模型。
在步骤5中，利用Sage_Husa自适应滤波方法完成相对姿态匹配传递对准过程。考虑到用不准确的噪声统计特性参数设计出的传统卡尔曼滤波将导致状态估计值误差增大，严重时还可能引起滤波发散的现象，利用Sage_Husa自适应滤波来提高相对姿态的估计值并对其进行校正处理，最后得到更优的估计值。

Claims (8)

1.本发明所采用的技术方案是，一种提高机载局部相对姿态匹配传递对准精度的方法，包括以下步骤：
步骤1.建立主子惯导性测量单元系统，其中主惯性测量单元安装在机体重心位置，子惯性测量单元安装在任务传感器腔体内；
步骤2.对于子惯性测量单元陀螺漂移进行补偿，通过对随机漂移建立自回归滑动平均模型，并对自回归模型进行了卡尔曼滤波，对陀螺漂移进行误差补偿；
步骤3.求解主、子惯性测量单元之间相对姿态阵
Figure RE-FDA0002414722840000011
其初值通过“双矢量定姿法”求解；
步骤4.建立相对姿态匹配传递对准的相对姿态误差方程、相对姿态量测方程、相对姿态状态方程，得到系统状态空间模型；
步骤5.利用Sage_Husa自适应滤波对系统状态空间模型中的状态方程、量测方程进行解算，估计出机载任务传感器相对姿态误差值，完成传递对准。
2.根据权利要求1所述的一种提高机载局部相对姿态匹配传递对准精度的方法，其特征在于，所述的步骤2部分，具体过程是：
利用时间序列分析法建立自回归滑动平均模型并进行参数估计，基于该自回归滑动平均模型进行卡尔曼滤波；
步骤2.1：建立自回归滑动平均模型(ARMA(p,q))并进行参数估计；
步骤2.2：基于自回归滑动平均模型ARMA(p,q)，进行卡尔曼滤波。
3.根据权利要求2所述的一种提高机载局部相对姿态匹配传递对准精度的方法，其特征在于，所述步骤2.2中卡尔曼滤波包括以下五个基本方程：
2.2.1：对状态做一步预测；
2.2.2：状态一步预测的均方误差阵；
2.2.3：滤波增益矩阵；
2.2.4：状态估计方程；
2.2.5：估计均方误差方程。
4.根据权利要求1所述的一种提高机载局部相对姿态匹配传递对准精度的方法，其特征在于，所述的步骤3部分，具体过程是：
步骤3.1：求解相对姿态阵
Figure RE-FDA0002414722840000021
首先需求出相对姿态初值通过“双矢量定姿法”求解；
步骤3.2：当机翼产生挠曲变形时，主、子惯性测量单元之间的相对姿态矩阵
Figure RE-FDA0002414722840000022
是一个时变函数，
Figure RE-FDA0002414722840000023
的动态变化过程可用姿态矩阵的微分方程来描述。
5.根据权利要求4所述的一种提高机载局部相对姿态匹配传递对准精度的方法，其特征在于，所述的步骤3.1的具体步骤为：
相对姿态初始值的确定阶段为飞机起飞离地之前，即飞机起飞前静止阶段和滑跑加速阶段；
步骤3.1.1：在载机静止条件下，采集主、子惯性测量单元一段时间内加速度计输出的比力均值，得到主子惯性测量单元比力向量；
步骤3.1.2：接着将调整载机水平姿态角(俯仰或横滚，无需确切知道大小)后静止；
步骤3.1.3：采集子惯性测量单元中加速度计输出的比力均值，得到,该位置子惯性测量单元的比力向量；
步骤3.1.4：通过子惯惯性测量单元两种情况下不共线的比力向量，构造主子惯性测量单元共同的中间过渡坐标系V，得到子惯性测量单元与该中间系的旋转变换阵，确定子惯性测量单元相对主惯惯性测量单元之间的相对姿态。
6.根据权利要求4所述的一种提高机载局部相对姿态匹配传递对准精度的方法，其特征在于，所述的步骤3.2的具体步骤为：3.2.1：当机翼产生挠曲变形时，主、子惯性测量单元之间的相对姿态矩阵
Figure RE-FDA0002414722840000031
是一个时变函数，
Figure RE-FDA0002414722840000032
的动态变化过程可用姿态矩阵的微分方程来描述：
Figure RE-FDA0002414722840000033
式中，
Figure RE-FDA0002414722840000034
表示由括号内向量组成的反对称矩阵；
Figure RE-FDA0002414722840000035
表示子惯性测量单元相对于主惯性测量单元的角速度在主惯性测量单元系下的投影，其值可由主、子惯性测量单元的陀螺仪输出计算得到；
3.2.2：相对姿态微分方程的最终形式可表示为：
Figure RE-FDA0002414722840000036
式中，
Figure RE-FDA0002414722840000037
分别为主、子惯性测量单元陀螺仪输出值；
记
Figure RE-FDA0002414722840000041
Figure RE-FDA0002414722840000042
故，根据(11)、(12)式可得到由9个未知数构成的线性微分方程组，在初值已知的情况下，可根据龙格-库塔法对该微分方程的进行数值运算，为了简化计算量只考虑式(11)后六个方程，再根据式(13)式中的约束条件对其进行求解，最终确定出
Figure RE-FDA0002414722840000043
7.根据权利要求1所述的一种提高机载局部相对姿态匹配传递对准精度的方法，其特征在于，所述的步骤4部分，具体过程是：
步骤4.1：采用二阶马尔科夫过程描述挠曲变形；
步骤4.2：描述主子惯性测量单元相对姿态误差值，即相对姿态角的计算值与实际值之间的差值；
Figure RE-FDA0002414722840000044
步骤4.3：建立相对姿态量测方程；
Figure RE-FDA0002414722840000045
步骤4.4：建立相对姿态状态方程；
Figure RE-FDA0002414722840000051
最终得到系统状态空间模型。
8.根据权利要求1所述的一种提高机载局部相对姿态匹配传递对准精度的方法，其特征在于，所述的步骤5部分，具体过程是：
利用Sage_Husa自适应滤波，进行相对姿态匹配传递对准；Sage_Husa自适应滤波方程由以下八个基本方程构成，分别为：
步骤5.1：计算加权系数dk；
步骤5.2：状态一步预测方程；
步骤5.3：一步预测的均方误差阵；
步骤5.4：新息序列；
步骤5.5：估计量测噪声；
步骤5.6：滤波增益矩阵；
步骤5.7：状态估计方程；
步骤5.8：估计均方误差方程；
对姿态的估计值进行校正处理，得到更优的估计值。
CN202010010721.6A 2020-01-06 2020-01-06 一种提高机载局部相对姿态匹配传递对准精度的方法 Active CN111141313B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010010721.6A CN111141313B (zh) 2020-01-06 2020-01-06 一种提高机载局部相对姿态匹配传递对准精度的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010010721.6A CN111141313B (zh) 2020-01-06 2020-01-06 一种提高机载局部相对姿态匹配传递对准精度的方法

Publications (2)

Publication Number Publication Date
CN111141313A true CN111141313A (zh) 2020-05-12
CN111141313B CN111141313B (zh) 2023-04-07

Family

ID=70523699

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010010721.6A Active CN111141313B (zh) 2020-01-06 2020-01-06 一种提高机载局部相对姿态匹配传递对准精度的方法

Country Status (1)

Country Link
CN (1) CN111141313B (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111707292A (zh) * 2020-07-18 2020-09-25 东南大学 一种自适应滤波的快速传递对准方法
CN113124903A (zh) * 2021-04-23 2021-07-16 中国电子科技集团公司第二十六研究所 传递对准下基于姿态匹配的最小二乘陀螺零偏快速估计
CN113341385A (zh) * 2021-03-30 2021-09-03 西南电子技术研究所（中国电子科技集团公司第十研究所） 机载平台协同综合传感器系统马尔科夫链误差传递模型
CN113834499A (zh) * 2021-08-26 2021-12-24 北京航天发射技术研究所 一种车载惯组与里程计行进间对准方法及系统
CN116026367A (zh) * 2023-03-29 2023-04-28 中国人民解放军火箭军工程大学 基于数字孪生技术的激光惯组故障诊断方法、系统及设备

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050240347A1 (en) * 2004-04-23 2005-10-27 Yun-Chun Yang Method and apparatus for adaptive filter based attitude updating
CN108731702A (zh) * 2018-07-03 2018-11-02 哈尔滨工业大学 一种基于Huber方法的大失准角传递对准方法
WO2018214227A1 (zh) * 2017-05-22 2018-11-29 深圳市靖洲科技有限公司 一种无人车实时姿态测量方法
CN109708670A (zh) * 2019-03-08 2019-05-03 哈尔滨工程大学 基于改进的Sage-Husa自适应卡尔曼滤波的极区传递对准挠曲变形的噪声补偿方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050240347A1 (en) * 2004-04-23 2005-10-27 Yun-Chun Yang Method and apparatus for adaptive filter based attitude updating
WO2018214227A1 (zh) * 2017-05-22 2018-11-29 深圳市靖洲科技有限公司 一种无人车实时姿态测量方法
CN108731702A (zh) * 2018-07-03 2018-11-02 哈尔滨工业大学 一种基于Huber方法的大失准角传递对准方法
CN109708670A (zh) * 2019-03-08 2019-05-03 哈尔滨工程大学 基于改进的Sage-Husa自适应卡尔曼滤波的极区传递对准挠曲变形的噪声补偿方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
徐英蛟;: "一种改进自适应增量Kalman滤波的传递对准算法" *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111707292A (zh) * 2020-07-18 2020-09-25 东南大学 一种自适应滤波的快速传递对准方法
CN113341385A (zh) * 2021-03-30 2021-09-03 西南电子技术研究所（中国电子科技集团公司第十研究所） 机载平台协同综合传感器系统马尔科夫链误差传递模型
CN113341385B (zh) * 2021-03-30 2023-09-05 西南电子技术研究所（中国电子科技集团公司第十研究所） 机载平台协同综合传感器系统马尔科夫链误差传递模型
CN113124903A (zh) * 2021-04-23 2021-07-16 中国电子科技集团公司第二十六研究所 传递对准下基于姿态匹配的最小二乘陀螺零偏快速估计
CN113834499A (zh) * 2021-08-26 2021-12-24 北京航天发射技术研究所 一种车载惯组与里程计行进间对准方法及系统
CN116026367A (zh) * 2023-03-29 2023-04-28 中国人民解放军火箭军工程大学 基于数字孪生技术的激光惯组故障诊断方法、系统及设备

Also Published As

Publication number Publication date
CN111141313B (zh) 2023-04-07

Similar Documents

Publication Publication Date Title
CN111141313B (zh) 一种提高机载局部相对姿态匹配传递对准精度的方法
CN109211276B (zh) 基于gpr与改进的srckf的sins初始对准方法
CN110398257B (zh) Gps辅助的sins系统快速动基座初始对准方法
CN109974714B (zh) 一种Sage-Husa自适应无迹卡尔曼滤波姿态数据融合方法
CN111156987B (zh) 基于残差补偿多速率ckf的惯性/天文组合导航方法
CN108827310B (zh) 一种船用星敏感器辅助陀螺仪在线标定方法
CN110146075B (zh) 一种增益补偿自适应滤波的sins/dvl组合定位方法
CN110887481B (zh) 基于mems惯性传感器的载体动态姿态估计方法
CN110567461B (zh) 一种考虑无陀螺仪的非合作航天器姿态和参数估计方法
CN110553641B (zh) 一种基于相关系数的提高惯性制导精度的方法
CN110440830B (zh) 动基座下车载捷联惯导系统自对准方法
CN111121773B (zh) 一种mems惯性测量组合
CN110887480A (zh) 基于mems传感器的飞行姿态估计方法及系统
CN110553642A (zh) 一种提高惯性制导精度的方法
CN112257186B (zh) 针对小型四旋翼飞行器气动参数的时域辨识方法
WO2022222938A1 (zh) 一种基于运动状态监测的自适应水平姿态测量方法
CN111238469A (zh) 一种基于惯性/数据链的无人机编队相对导航方法
CN113008272B (zh) 一种用于微小卫星的mems陀螺在轨常值漂移标定方法和系统
CN111707292B (zh) 一种自适应滤波的快速传递对准方法
CN110736459B (zh) 惯性量匹配对准的角形变测量误差评估方法
CN110375773B (zh) Mems惯导系统姿态初始化方法
CN112986977A (zh) 一种克服雷达扩展卡尔曼航迹滤波发散的方法
CN109674480B (zh) 一种基于改进互补滤波的人体运动姿态解算方法
CN112287289A (zh) 一种面向云控智能底盘的车辆非线性状态融合估计方法
CN113252029B (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