CN115790589B - 一种发射系无误差捷联惯性导航方法 - Google Patents

一种发射系无误差捷联惯性导航方法 Download PDF

Info

Publication number
CN115790589B
CN115790589B CN202310028690.0A CN202310028690A CN115790589B CN 115790589 B CN115790589 B CN 115790589B CN 202310028690 A CN202310028690 A CN 202310028690A CN 115790589 B CN115790589 B CN 115790589B
Authority
CN
China
Prior art keywords
matrix
speed
initial
angular velocity
output
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
Application number
CN202310028690.0A
Other languages
English (en)
Other versions
CN115790589A (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.)
Northwestern Polytechnical University
Original Assignee
Northwestern Polytechnical University
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 Northwestern Polytechnical University filed Critical Northwestern Polytechnical University
Priority to CN202310028690.0A priority Critical patent/CN115790589B/zh
Publication of CN115790589A publication Critical patent/CN115790589A/zh
Application granted granted Critical
Publication of CN115790589B publication Critical patent/CN115790589B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Navigation (AREA)

Abstract

本发明公开了一种发射系无误差捷联惯性导航方法,属于捷联惯导领域,在陀螺仪输出角速度和加速度计输出比力满足多项式形式的条件下,根据发射系下捷联惯导的导航方程,采用泰勒级数展开的方法推导出的发射系捷联惯导数值更新算法。该算法不存在任何原理性误差,实现了对姿态、速度和位置不可交换误差的完美补偿,理论上算法精度远高于传统的二子样算法。

Description

一种发射系无误差捷联惯性导航方法
技术领域
本发明涉及捷联惯导领域,具体涉及一种发射系无误差捷联惯性导航方法。
背景技术
发射坐标系(发射系)能够满足航天和航空的双重飞行控制和导航需求,但对现有的发射坐标系导航而言,其捷联惯导数值更新算法的推导中,姿态更新求解算法采用的是主流的姿态更新方法,即先使用陀螺角增量的二子样采样计算等效旋转矢量,补偿转动不可交换误差,再使用等效旋转矢量计算姿态更新四元数;速度更新算法和位置更新算法的推导中也都包含了对等效旋转矢量的近似。这种传统的二子样算法是在等效旋转矢量方程(Bortz方程)进行二阶近似的基础上推导的,本质上存在原理性误差,特别是对于具有大机动特性的载体,由于算法原理性误差造成的计算误差更是不容忽略。
发明内容
针对现有技术中的上述不足,本发明提供了一种发射系无误差捷联惯性导航方法。
为了达到上述发明目的,本发明采用的技术方案为:
一种发射系无误差捷联惯性导航方法,包括如下步骤:
S1、设定飞行器初始状态,包括初始姿态矩阵、初始速度和初始位置;
S2、将飞行器惯性测量单元输出的角增量拟合成角速度多项式,同时将飞行器惯性测量单元输出的速度增量拟合成比力多项式;
S3、根据发射点地球自转角速度、陀螺仪输出角速度和初始姿态矩阵计算飞行器的姿态矩阵高阶导数,进一步根据姿态矩阵高阶导数计算 T时刻的姿态矩阵,完成姿态更新;
S4、根据发射点地球自转角速度的反对称矩阵、陀螺仪输出的比力以及初始速度计算速度的高阶导数,进一步根据速度的高阶导数计算 T时刻的速度,完成速度更新;
S5、根据初始位置、初始速度和速度的高阶导数计算 T时刻的位置,完成姿态更新。
进一步的,所述S2中以多项式形式拟合飞行器不含惯性器件误差的角速度多项式的具体方式为:
式中,
其中,为陀螺仪输出的角速度; t为时间,其上标为幂次方,下标为采样点编号; N为采样次数;为角增量矩阵,元素为角增量,其下标为采样点的编号;为拟合系数矩阵; p, n分别为采样次数编号,且
所述S2中将飞行器惯性测量单元输出的速度增量拟合成比力形式的具体方式为:
其中,为加速度计输出的比力,为速度增量矩阵;为时间段内的 N次加速度计速度增量采样,其下标为采样点编号;为采样间隔。
进一步的,所述S3中姿态矩阵的高阶导数的计算方式为:
其中,为初始姿态矩阵的 i阶导数;为飞行器陀螺仪输出的角速度;为发射点的地球自转角速度;为组合数且ij为非负整数,为初始姿态矩阵的 i-1阶导数。
进一步的,所述S3中进一步根据姿态矩阵高阶导数计算 T时刻的姿态矩阵的具体计算方式为:
其中,T时刻的姿态矩阵;为初始姿态矩阵,为飞行器陀螺仪输出的角速度,为飞行器陀螺仪输出的角速度的反对称矩阵;为发射点的地球自转角速度; i为非负整数;为初始姿态矩阵的 i阶导数, i为非负整数。
进一步的,所述S4中速度的高阶导数的具体计算方式为:
其中,为初始速度的 i阶导数;为加速度计输出的比力阶导数,为初始速度阶导数,为飞行器在发射坐标系下的初始重力加速度阶导数,为初始速度阶导数,为组合数且i、j为非负整数;为发射点地球自转角速度的反对称矩阵,为初始姿态矩阵的 j阶导数;
进一步的,所述S4中进一步根据速度的高阶导数计算 T时刻的速度的计算方式为:
其中,T时刻的速度,为初始速度,为发射点地球自转角速度的反对称矩阵;为飞行器在发射坐标系下的初始重力加速度;为初始速度的高阶导数;为弹体坐标系下飞行器加速度计输出的比力,为初始姿态矩阵, i为非负整数。
进一步的,所述S5中计算 T时刻的位置的具体计算方式为:
其中,T时刻的位置;为初始位置;为初始速度;为初始速度阶导数, i为非负整数。
本发明具有以下有益效果:
在陀螺仪输出角速度和加速度计输出比力满足多项式假设的条件下,根据数学知识,平滑的角运动或线运动总是可以用多项式无限逼近的,根据发射系下捷联惯导的导航方程,采用泰勒级数展开的方法导出的发射系下姿态、速度和位置捷联惯导数值更新算法。该算法不存在任何原理性误差,包含了对姿态、速度和位置更新不可交换误差的完美补偿,理论上算法精度远高于传统的二子样算法。
附图说明
图1为本发明一种发射系无误差捷联惯性导航方法流程示意图。
图2为本发明实施例发射坐标系示意图。
图3为本发明实施例发射坐标系和弹体坐标系关系示意图。
图4为本发明实施例发射系无误差捷联惯性导航数值更新计算方式示意图。
具体实施方式
下面对本发明的具体实施方式进行描述,以便于本技术领域的技术人员理解本发明,但应该清楚,本发明不限于具体实施方式的范围,对本技术领域的普通技术人员来讲,只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内,这些变化是显而易见的,一切利用本发明构思的发明创造均在保护之列。
为了实现上述目的,本发明具有如下基础前提:
1 坐标系及坐标系转换
捷联惯导系统常用的导航坐标系包括:地心惯性坐标系( i系)、地心地固坐标系( e系)、发射惯性坐标系( a系)、发射坐标系( g系)、当地水平坐标系( l系)和弹体坐标系( b系)。其中发射坐标系和弹体坐标系定义如下:
1) 发射坐标系( g系)
发射坐标系(Launch-Centered Earth-Fixed Frame,LCEF),原点为发射点, x轴在发射点水平面内并且指向发射瞄准方向, y轴垂直于发射点水平面并且指向上方, z轴与 x轴、 y轴构成右手直角坐标系,并且发射坐标系与地球固连。发射点的地心纬度、经度、高度和发射方位角(为过发射点的正北方向与发射系 x轴的夹角)确定了发射坐标系与地球之间的关系,如图2所示。
2)弹体坐标系( b系)
弹体坐标系(Body Frame),原点为飞行器的质心; x轴沿飞行器的纵轴,指向飞行器头部; y轴在飞行器的纵对称面内,垂直于 x轴指向上; z轴与 x轴、 y轴构成右手直角坐标系。
3)发射坐标系与弹体坐标系之间的转换关系
发射坐标系到弹体坐标系的方向余弦矩阵为,飞行器在发射坐标系的姿态角由俯仰角、偏航角和滚转角三个欧拉角描述,按照先绕 z轴俯仰、再绕 y轴偏航、最后绕 x轴滚转的3-2-1旋转顺序,发射系和弹体系的关系如图3所示,则姿态矩阵为                                 (1)
式中:分别为发射系俯仰角、偏航角和滚转角。
2 发射系捷联惯导算法编排
发射系下捷联惯导微分方程如下式(2)所示,
                           (2)
式中: 分别为发射系中载体的位置、速度和姿态矩阵,且如式(1)所示;为加速度计测量到的比力;系相对于系的旋转角速度,即发射点的地球自转角速度的反对称矩阵,即为载体在系中的重力矢量;系相对于系的旋转角速度的反对称矩阵,且有:
                              (3)
式中:为陀螺仪测量的角速度的反对称矩阵,
具体而言,如图1所示,提供一种发射系无误差捷联惯性导航方法,包括如下步骤:
S1、设定飞行器初始状态,包括初始姿态矩阵、初始速度和初始位置
S2、将飞行器惯性测量单元输出的角增量拟合成角速度多项式,同时将飞行器惯性测量单元输出的速度增量拟合成比力多项式。
捷联惯导系统中的陀螺仪通常是角增量输出形式,加速度计通常是速度增量输出形式,而本算法要以角速度和比力作为输入,因此,要先给出由角增量转换为角速度、速度增量转换为比力的转换方法。
若陀螺仪角增量的采样间隔为 h,在时间段内进行了 N次采样( p, n均为整数,),角增量分别记为,(),则可得角速度的关于时间 t的( N-1)次多项式拟合为
                               (4)
式中:
时表示当前姿态更新周期内的角增量采样,而当时表示利用了前面姿态更新周期的角增量采样。
同理,可得比力的多项式拟合为
                           (5)
式中:,为时间段内的 N次加速度计速度增量采样。
发射系无误差捷联惯导数值更新算法的推导中,首先给出两矩阵之乘积的求导公式,即:若,则有
       (6)
式中: ij为非负整数,为组合数;为适当维数的矩阵;右上标小括号内的数值表示求导阶次。
接下来根据式(2)所示的发射系捷联惯导姿态阵、速度和位置微分方程分别给出相应的数值更新算法。
S3、根据发射点地球自转角速度、陀螺仪输出角速度和初始姿态矩阵计算飞行器的姿态矩阵高阶导数,并进一步根据计算得到的姿态矩阵高阶导数计算 T时刻的姿态矩阵,完成姿态更新;
在发射系捷联惯导算法中,姿态阵微分方程为
                                     (7)
将式(7)代入式(1),有
  (8)
式中,为陀螺仪输出的 b系相对于惯性系( a系)的角速度,即陀螺仪输出的角速度;g系相对于 a系的旋转角速度,即发射点的地球自转角速度,为常值矩阵。
假设姿态更新时间间隔为,将 T时刻的姿态阵在0时刻展开成泰勒级数形式,并将式(8)代入,可得
   (9)
根据式(6)和式(8),有
      (10)
注意到发射点的地球自转角速度是一个常值,故式(10)中,当 j>0时,,所以式(10)可以写成
              (11)
式(11)说明,姿态阵的高阶导数可以用它的低阶导数表示,最终总可以表示成角速度的各阶导数以及初始姿态阵和地球自转 角速度的函数。利用式(9)和式(11)便实现了从的姿态阵更新,其中式(9)右端第三项(求和项)隐含了对姿态不可交换误差(圆锥误差)的精确补偿。
S4、根据发射点地球自转角速度的反对称矩阵、陀螺仪输出的比力以及和初始速度计算速度的高阶导数,并进一步根据计算得到的速度的高阶导数计算 T时刻的速度,完成速度更新;
在发射系捷联惯导算法中,速度微分方程为
                (12)
T时刻的速度展开成泰勒级数形式,并将式(12)代入为
      (13)
根据式(6)和式(12)有
(14)
同式(10)到式(11)过程,当 j>0时,;即使对高超声速飞行器这种快速运动的运载体,在姿态更新周期的短时间内,重力矢量变化都是很小的,所以在一个姿态更新周期内,可以认为重力矢量为时间的一次函数缓变量,因此有 i=1时,为一常数矩阵,且有
                        (15)
式中:为当前数值更新周期的重力矢量,为前一周期的重力矢量, T为数值更新周期。
i>1时,,则有
      (16)
式(13)右端第三项(求和项)也隐含了对速度不可交换误差(划桨误差)的精确补偿.
S5、根据初始位置、初始速度和速度的高阶导数计算 T时刻的位置,完成姿态更新。
在发射系捷联惯导算法中,位置微分方程为
                                   (17)
T时刻的位置展开成泰勒级数形式为
                       (18)
式(18)右端第三项(求和项)隐含了对位置不可交换误差(涡卷误差)的精确补偿。
经过上文推导,至此获得了完整的姿态阵更新、速度更新和位置更新的全套捷联惯导数值更新算法,其计算流程如图4所示,可以看出,新算法可以直接输出发射系下的姿态、速度和位置等导航参数,发射系与地球固连,其导航参数是相对于地球的,与很多地面发射飞行器知道控制系统需求的导航信息一致,因此此算法的另一大优势是可以给出通常需要的导航信息,可供工程上直接使用。
本发明中应用了具体实施例对本发明的原理及实施方式进行了阐述,以上实施例的说明只是用于帮助理解本发明的方法及其核心思想;同时,对于本领域的一般技术人员,依据本发明的思想,在具体实施方式及应用范围上均会有改变之处,综上所述,本说明书内容不应理解为对本发明的限制。
本领域的普通技术人员将会意识到,这里所述的实施例是为了帮助读者理解本发明的原理,应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。本领域的普通技术人员可以根据本发明公开的这些技术启示做出各种不脱离本发明实质的其它各种具体变形和组合,这些变形和组合仍然在本发明的保护范围内。

Claims (3)

1.一种发射系无误差捷联惯性导航方法,其特征在于,包括如下步骤:
S1、设定飞行器初始状态,包括初始姿态矩阵、初始速度和初始位置;
S2、将飞行器惯性测量单元输出的角增量和速度增量拟合成角速度和比力的多项式形式;
S3、根据发射点地球自转角速度、陀螺仪输出角速度和初始姿态矩阵计算飞行器的姿态矩阵高阶导数,进一步根据姿态矩阵高阶导数计算T时刻的姿态矩阵,完成姿态更新,其中,姿态矩阵的高阶导数的计算方式为:
为姿态矩阵的i阶导数;为飞行器陀螺仪输出的角速度;为发射点的地球自转角速度;为组合数且ij为非负整数,为姿态矩阵的i-1阶导数;为初始姿态矩阵的j阶导数,为飞行器陀螺仪输出的角速度的反对称矩阵;为发射点的地球自转角速度的反对称矩阵;
根据姿态矩阵高阶导数计算T时刻的姿态矩阵的具体计算方式为:
其中,T时刻的姿态矩阵;为初始姿态矩阵;
S4、根据发射点地球自转角速度的反对称矩阵、陀螺仪输出的比力以及初始速度计算速度的高阶导数,具体计算方式为:
其中,为初始速度i阶导数;为加速度计输出的比力阶导数,为初始速度阶导数,为飞行器在发射坐标系下的初始重力加速度阶导数,为初始速度阶导数,为组合数且i、j为非负整数;为发射点地球自转角速度的反对称矩阵;
进一步根据速度的高阶导数计算T时刻的速度,完成速度更新,具体计算方式为:
根据速度的高阶导数计算T时刻的速度的计算方式为:
其中,T时刻的速度;
S5、根据初始位置、初始速度和速度的高阶导数计算T时刻的位置,完成位置更新。
2.根据权利要求1所述的一种发射系无误差捷联惯性导航方法,其特征在于,所述S2中将飞行器惯性测量单元输出的角增量拟合成角速度形式的具体方式为:
式中,
其中,为陀螺仪输出的角速度;t为时间,其上标为幂次方,下标为采样点编号;N为采样次数;为角增量矩阵,元素为角增量,其下标为采样点的编号;为拟合系数矩阵;p,n分别为采样次数编号,且
所述S2中将飞行器惯性测量单元输出的速度增量拟合成比力形式的具体方式为:
其中,为加速度计输出的比力,为速度增量矩阵;为时间段内的N次加速度计速度增量采样,其下标为采样点编号;为采样间隔。
3.根据权利要求1所述的一种发射系无误差捷联惯性导航方法,其特征在于,所述S5中计算T时刻的位置的具体计算方式为:
其中,T时刻的位置;为初始位置。
CN202310028690.0A 2023-01-09 2023-01-09 一种发射系无误差捷联惯性导航方法 Active CN115790589B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310028690.0A CN115790589B (zh) 2023-01-09 2023-01-09 一种发射系无误差捷联惯性导航方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310028690.0A CN115790589B (zh) 2023-01-09 2023-01-09 一种发射系无误差捷联惯性导航方法

Publications (2)

Publication Number Publication Date
CN115790589A CN115790589A (zh) 2023-03-14
CN115790589B true CN115790589B (zh) 2023-05-02

Family

ID=85428869

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310028690.0A Active CN115790589B (zh) 2023-01-09 2023-01-09 一种发射系无误差捷联惯性导航方法

Country Status (1)

Country Link
CN (1) CN115790589B (zh)

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2000213906A (ja) * 1998-12-30 2000-08-04 Netmor Ltd 移動物体位置のトラッキング装置およびトラッキング方法
CN114911265A (zh) * 2022-06-13 2022-08-16 北京航空航天大学杭州创新研究院 一种四旋翼无人机编队协同机动控制方法

Family Cites Families (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7522099B2 (en) * 2005-09-08 2009-04-21 Topcon Gps, Llc Position determination using carrier phase measurements of satellite signals
US8344303B2 (en) * 2010-11-01 2013-01-01 Honeywell International Inc. Projectile 3D attitude from 3-axis magnetometer and single-axis accelerometer
FR2977313B1 (fr) * 2011-06-28 2013-08-09 Centre Nat Etd Spatiales Engin spatial muni d'un dispositif d'estimation d'un vecteur vitesse et procede d'estimation correspondant
CN108489485B (zh) * 2018-03-20 2021-07-06 西北工业大学 一种无误差的捷联惯导数值更新方法
CN109489690B (zh) * 2018-11-23 2020-10-23 北京宇航系统工程研究所 一种适用于高动态翻滚再入的助推器导航定位解算方法
CN110057382B (zh) * 2019-04-23 2021-07-09 西北工业大学 一种基于发射坐标系的捷联惯导数值更新方法
CN111721291B (zh) * 2020-07-17 2021-12-07 河北斐然科技有限公司 一种发射系下捷联惯组导航的工程算法
CN115248038B (zh) * 2022-09-21 2022-12-30 河北斐然科技有限公司 一种发射系下的sins/bds组合导航工程算法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2000213906A (ja) * 1998-12-30 2000-08-04 Netmor Ltd 移動物体位置のトラッキング装置およびトラッキング方法
CN114911265A (zh) * 2022-06-13 2022-08-16 北京航空航天大学杭州创新研究院 一种四旋翼无人机编队协同机动控制方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
捷联式惯性导航系统误差处理技术的新进展;刘晓光,谢玲,戴亚平,陈家斌,张鸿业;火力与指挥控制;26(04);全文 *

Also Published As

Publication number Publication date
CN115790589A (zh) 2023-03-14

Similar Documents

Publication Publication Date Title
US8519313B2 (en) Projectile navigation enhancement method
CN109823571A (zh) 一种遥感微纳卫星的多阶段姿态控制方法
CN109489690B (zh) 一种适用于高动态翻滚再入的助推器导航定位解算方法
CN109724624B (zh) 一种适用于机翼挠曲变形的机载自适应传递对准方法
CN109269504B (zh) 一种具有末端约束的姿态机动路径规划方法
CN108489485B (zh) 一种无误差的捷联惯导数值更新方法
CN104281150A (zh) 一种姿态机动的轨迹规划方法
CN109211230A (zh) 一种基于牛顿迭代法的炮弹姿态和加速度计常值误差估计方法
CN108583938B (zh) 一种可应用于运行于太阳同步晨昏轨道的全向天线通信卫星姿态控制系统及其方法
CN111207745A (zh) 一种适用于大机动无人机垂直陀螺仪的惯性测量方法
CN103231810A (zh) 一种利用卫星俯仰轴姿态机动卸载俯仰轴角动量的方法
CN113587925A (zh) 一种惯性导航系统及其全姿态导航解算方法与装置
CN115790589B (zh) 一种发射系无误差捷联惯性导航方法
CN105799949B (zh) 一种亚轨道卫星的压心设计方法、姿态控制方法和系统
CN103869097B (zh) 旋转弹航向角、俯仰角角速率测量方法
Brauer et al. Program to optimize simulated trajectories (POST). Volume 1: Formulation manual
CN114353784B (zh) 一种基于运动矢量的制导炮弹空中姿态辨识方法
CN114295145A (zh) 一种基于车载发射平台的捷联惯导系统轨迹发生器设计方法
CN115060256B (zh) 一种基于发射坐标系的制导炮弹空中姿态辨识方法
CN112464432A (zh) 一种基于双速数值积分结构的惯性预积分的优化方法
CN108088464B (zh) 一种闭环修正安装误差的挠曲估计方法
CN114383614B (zh) 一种弹道环境下的多矢量空中对准方法
CN116700020B (zh) 变后掠翼无人机的控制方法、系统、无人机及存储介质
CN110658862B (zh) 基于分布式角动量的柔性结构振动能量一体化控制方法
CN116642485A (zh) 惯性导航6参数坐标变换矩阵解算方法

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