CN103438892A - 一种改进的基于ekf的天文自主定轨算法 - Google Patents

一种改进的基于ekf的天文自主定轨算法 Download PDF

Info

Publication number
CN103438892A
CN103438892A CN2013104190759A CN201310419075A CN103438892A CN 103438892 A CN103438892 A CN 103438892A CN 2013104190759 A CN2013104190759 A CN 2013104190759A CN 201310419075 A CN201310419075 A CN 201310419075A CN 103438892 A CN103438892 A CN 103438892A
Authority
CN
China
Prior art keywords
coordinate system
filtering
ekf
orbit determination
filter
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
CN2013104190759A
Other languages
English (en)
Other versions
CN103438892B (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.)
Harbin Institute of Technology (Beijing) military and Civilian Integration Innovation Research Institute Co., Ltd.
Original Assignee
Harbin Institute 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 Harbin Institute of Technology filed Critical Harbin Institute of Technology
Priority to CN201310419075.9A priority Critical patent/CN103438892B/zh
Publication of CN103438892A publication Critical patent/CN103438892A/zh
Application granted granted Critical
Publication of CN103438892B publication Critical patent/CN103438892B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Abstract

一种改进的基于EKF的天文自主定轨算法，涉及一种控制领域的导航方法。所述方法步骤如下：一、EKF滤波器参数初始化；二、参考坐标系由地心惯性坐标系转为滤波坐标系；三、在滤波坐标系中线性化和离散化滤波器状态方程；四、在滤波坐标系中进行EKF滤波；五、参考坐标系由滤波坐标系转回为地心惯性坐标系。本发明在传统的基于EKF天文自主定轨算法中，提出并引入滤波坐标系，同时增加坐标变换和坐标反转换等步骤使得传统的基于EKF天文自主定轨算法在线性化和离散化的时候，航天器的位置矢量在滤波坐标系中的各分量能够相等，从而使得状态方程线性化的展开点处于一个线性度较好的位置，进而使得定轨算法的稳定性有所提高。

Description

一种改进的基于EKF的天文自主定轨算法
技术领域
本发明涉及一种控制领域的导航方法，具体涉及一种改进的基于EKF的天文自主定轨算法。
背景技术
自20世纪60年代以来，航天事业的发展迎来了天文导航技术的辉煌，天文导航技术在美国阿波罗载人登月计划中得到了成功应用。近年来，随着新一轮月球和火星探测等一系列深空探测活动的开展，天文导航以其自主性强、精度高、成本低廉等特点在航天领域也得到了越来越广泛的应用，成为一种航天器的重要自主导航方法。天文导航技术目前已成为卫星、深空探测和载人航天中必不可少的关键技术，在未来人类探索宇宙的星际航行中也必将发挥重要的作用。
目前导航系统中应用最为广泛的滤波方法为扩展卡尔曼滤波（Extended Kalman Filter，EKF）、Unscented卡尔曼滤波（UKF）、Unscented粒子滤波（UPF）3种滤波方法。其中EKF计算量较小，但是导航精度相对也比较低，同时因为EKF算法在工程上会出现滤波发散的问题。美国著名的深空1号卫星是一颗完全自主性的深空探测器，在自主天文定轨方案上也因为EKF算法的可能发散而没有采纳，而是选取了批处理最小二乘方法。因此有必要对基于EKF的天文学定轨算法进行改进以提高其稳定性。
EKF算法是一种近似方法，它将非线性模型在状态估计值附近进行泰勒级数展开，并取一阶截断，用得到的一阶近似项作为原状态方程和测量方程的近似表达形式，从而实现线性化同时假设线性化后的状态依然服从高斯分布，然后对线性化后的系统采用标准卡尔曼滤波获得状态估计。由于采用了局部线性化技术，EKF能得到局部最优解，但它能否收敛于全局最优解，取决于函数的非线性强度以及展开点（即泰勒级数展开点）的选择。
基于EFK的天文学自主定轨算法在展开点处对滤波器的非线性状态方程和非线性观测方程进行了线性化和离散化。为了选取一个好的展开点，常用的思想是通过施加一个偏置量，使得当前的状态偏置到一个线性度较好的地方，使得线性化后的近似方程与原方程之间的误差相对较小。
发明内容
本发明的目的在于对传统的基于EFK的天文学自主定轨算法，通过坐标变换的方法，提出一种改进的基于EKF的天文自主定轨算法，使得展开点处具有较好的线性度，从而达到提高传统的基于EKF的天文自主定轨算法稳定性目的。
本发明的目的是通过以下技术方案实现的：
在基于EKF的天文自主定轨算法进行到线性化和离散化步骤之前，先将滤波器状态变量的参考坐标系由地心惯性坐标系变为滤波坐标系，滤波坐标系是本发明为了使滤波器状态方程和滤波器观测方程在线性化时具有较好的展开点而设置的一个坐标系，凡是能使位置矢量在该坐标系的三个轴上的投影分量大小相等的坐标系都称为滤波坐标系。从而使得滤波器状态方程和滤波器观测方程在线性化的过程中，不会因为位置矢量各个分量的数值在数量级上不一致而导致滤波器发散。如图1所示，共分为五个步骤，具体步骤如下：
步骤一、EKF滤波器参数初始化：
1）给滤波迭代次数                                                
Figure 876407DEST_PATH_IMAGE001
赋初值为
Figure 698870DEST_PATH_IMAGE002
；
2）设滤波器的状态变量
Figure 965903DEST_PATH_IMAGE003
，其中
Figure 102486DEST_PATH_IMAGE004
为地心到航天器的位置矢量，
Figure 771365DEST_PATH_IMAGE005
为航天器的速度矢量。给初始时刻的状态变量
Figure 764729DEST_PATH_IMAGE006
估计值
Figure 456741DEST_PATH_IMAGE007
赋初值，其中上角标
Figure 459332DEST_PATH_IMAGE008
说明该变量为地心惯性坐标系；
 3）滤波时间赋初值
Figure 982717DEST_PATH_IMAGE009
；
4）给状态模型噪声的协方差矩阵赋值；
5）给观测模型噪声的协方差矩阵
Figure 821897DEST_PATH_IMAGE011
赋值；
6）给预估状态协方差
Figure 362599DEST_PATH_IMAGE012
赋初值。
 步骤二、参考坐标系由地心惯性坐标系转为滤波坐标系：
不同于传统的基于EKF的天文自主定轨算法，本发明在自主定轨算法的线性化和离散化之前，需要将参考坐标系由地心惯性坐标系转换为滤波坐标系。图2为滤波坐标系示意图。
滤波坐标系是本发明为了使得滤波器的状态方程和滤波器观测方程在进行线性化时具有较好的展开点而设置的一个坐标系，凡是能使位置矢量在该坐标系的三个轴上的投影分量大小相等的坐标系都称为滤波坐标系。本发明中用
Figure 278920DEST_PATH_IMAGE014
表示滤波坐标系。
 1）获取上一时刻航天器位置矢量的估计值
Figure 742262DEST_PATH_IMAGE015
以及这一时刻的观测信息矢量
Figure 86656DEST_PATH_IMAGE016
，其中
Figure 256737DEST_PATH_IMAGE017
表示上一时刻滤波器状态变量估计向量的前三个元素组成的列向量；
 2）通过下式计算
Figure 28384DEST_PATH_IMAGE018
与其期望的参考坐标表示
Figure 979023DEST_PATH_IMAGE019
之间的夹角
Figure 799211DEST_PATH_IMAGE020
：
Figure 151695DEST_PATH_IMAGE021
；
其中
Figure 94243DEST_PATH_IMAGE022
表示矩阵或向量的转置；
 3）通过下式计算
Figure 968396DEST_PATH_IMAGE023
和
Figure 654592DEST_PATH_IMAGE019
两个矢量所在平面的单位法方向量
Figure 861583DEST_PATH_IMAGE024
：
Figure 912715DEST_PATH_IMAGE025
；
4）通过下式计算地心惯性坐标系
Figure 837946DEST_PATH_IMAGE026
到滤波坐标系
Figure 62254DEST_PATH_IMAGE027
的转动四元数
Figure 327013DEST_PATH_IMAGE028
：
Figure 345785DEST_PATH_IMAGE029
；
5）通过下式计算地心惯性坐标系
Figure 758312DEST_PATH_IMAGE026
到滤波坐标系
Figure 723994DEST_PATH_IMAGE027
的余弦转换矩阵
Figure 905576DEST_PATH_IMAGE030
：
 
Figure 360828DEST_PATH_IMAGE031
其中
Figure 434220DEST_PATH_IMAGE032
为转动四元数的第
Figure 301999DEST_PATH_IMAGE034
个分量；
6）通过下式计算
Figure 865835DEST_PATH_IMAGE035
在滤波坐标系
Figure 987375DEST_PATH_IMAGE027
中的坐标表示
Figure 357177DEST_PATH_IMAGE036
和地心惯性坐标系
Figure 451034DEST_PATH_IMAGE026
中的速度矢量坐标表示
Figure 248089DEST_PATH_IMAGE037
在滤波坐标系中的坐标表示
Figure 702521DEST_PATH_IMAGE038
：
 组合滤波坐标系
Figure 681159DEST_PATH_IMAGE040
中的位置矢量和速度矢量，得到滤波器状态变量在滤波坐标系
Figure 777291DEST_PATH_IMAGE040
中的坐标表示为：
Figure 925113DEST_PATH_IMAGE041
，这时候
Figure 790301DEST_PATH_IMAGE042
的前三个分量相等，即满足以下等式：
7）通过下式计算观测信息矢量
Figure 450269DEST_PATH_IMAGE044
在滤波坐标系
Figure 965564DEST_PATH_IMAGE027
中的坐标表示
Figure 685259DEST_PATH_IMAGE045
：
Figure 932700DEST_PATH_IMAGE046
。
 步骤三、在滤波坐标系中线性化和离散化滤波器状态方程：
此步骤与传统的基于EKF天文自主定轨算法的不同之处在于，滤波器状态方程的雅克比矩阵和状态转移矩阵的求解都是在滤波坐标系中进行的。在以位置和速度为滤波器状态变量的基于EKF的天文自主定轨算法中，求解的雅克比矩阵中各个元素与位置有关而与速度无关；而且在步骤二，已近将位置矢量的三个分量调整为大小相等，所以位置矢量的三个分量对雅克比矩阵的贡献可以保持相近，这有助于提高EKF滤波器的稳定性。
1）计算滤波器状态方程
Figure 3424DEST_PATH_IMAGE047
的雅克比矩阵
Figure 56831DEST_PATH_IMAGE048
，其中
Figure 834294DEST_PATH_IMAGE049
为
Figure 314954DEST_PATH_IMAGE050
对时间的导数；
 2）计算
Figure 872974DEST_PATH_IMAGE051
时刻的状态转移矩阵
Figure 464493DEST_PATH_IMAGE052
，其中
Figure 592068DEST_PATH_IMAGE053
为单位矩阵，为滤波时间步长。
    步骤四、在滤波坐标系中进行EKF滤波：
此步骤与传统的基于EKF天文自主定轨算法的不同之处在于，状态预估和状态更新两个步骤都是在滤波坐标系下进行的。
1）在滤波坐标系中进行状态预估：
 
其中表示滤波器的状态方程，
Figure 436210DEST_PATH_IMAGE057
表示的是
Figure 258673DEST_PATH_IMAGE058
时刻滤波器状态变量的初步计算值，
Figure 463389DEST_PATH_IMAGE059
为
Figure 662289DEST_PATH_IMAGE060
时刻的预估状态协方差，
Figure 331168DEST_PATH_IMAGE061
为
Figure 262215DEST_PATH_IMAGE062
时刻滤波器的预估状态协方差的初步计算值；
2）在滤波坐标系中进行状态更新：
 
Figure 16544DEST_PATH_IMAGE063
其中
Figure 753556DEST_PATH_IMAGE064
表示滤波器的观测方程，
Figure 713160DEST_PATH_IMAGE065
为
Figure 143004DEST_PATH_IMAGE066
时刻的观测方程的雅克比矩阵，
Figure 322312DEST_PATH_IMAGE067
为增益矩阵，
Figure 597436DEST_PATH_IMAGE068
为单位矩阵，
Figure 240907DEST_PATH_IMAGE069
为
Figure 513756DEST_PATH_IMAGE066
时刻的预估状态协方差，为
Figure 259176DEST_PATH_IMAGE066
时刻滤波器状态变量的估计值在滤波坐标系下的表示。
步骤五、参考坐标系由滤波坐标系转回为地心惯性坐标系：
此步骤与传统的基于EKF天文自主定轨算法的不同之处在于，传统的基于EKF天文自主定轨算法没有这个步骤。
改进的基于EKF天文自主定轨算法通过步骤四滤波得出的
Figure 757153DEST_PATH_IMAGE058
时刻的滤波器状态变量的估计值
Figure 528800DEST_PATH_IMAGE071
是在滤波坐标系下的坐标表示，为了得到其在原来的地心惯性坐标系中的坐标表示，需要对结果进行坐标反变换：
  
Figure 918586DEST_PATH_IMAGE072
。
 通过
Figure 801092DEST_PATH_IMAGE073
可得到航天器的位置矢量的估计值和速度矢量的估计值；递增滤波迭代次数
Figure 153576DEST_PATH_IMAGE074
，然后跳回至步骤二，通过以上步骤的周期运行，就可以从滤波器状态变量的估计值
Figure 33807DEST_PATH_IMAGE075
中获得航天器当前时刻位置矢量的估计值
Figure 471741DEST_PATH_IMAGE076
以及速度矢量的估计值。
本发明与现有技术相比具有如下优点：
在传统的基于EKF天文自主定轨算法中，提出并引入滤波坐标系，同时增加坐标变换和坐标反转换等步骤使得传统的基于EKF天文自主定轨算法在线性化和离散化的时候，航天器的位置矢量在滤波坐标系中的各分量能够相等，从而使得状态方程线性化的展开点处于一个线性度较好的位置，进而使得定轨算法的稳定性有所提高，适用于深空探测器等需要天文自主导航功能的航天器。
附图说明
图1为本发明的流程图；
图2为滤波坐标系示意图；
图3为探月卫星的位置时间曲线；
图4为传统的基于EKF天文自主定轨算法的定位误差时间曲线；
图5为传统的基于EKF天文自主定轨算法的定速误差时间曲线；
图6为改进的基于EKF天文自主定轨算法的定位误差时间曲线；
图7为改进的基于EKF天文自主定轨算法的定速误差时间曲线。
具体实施方式
下面结合实施例对本发明的技术方案作进一步的说明，但并不局限如此，凡是对本发明技术方案进行修改或者等同替换，而不脱离本发明技术方案的精神和范围，均应涵盖在本发明的保护范围中。
本实施例将改进的基于EKF天文自主定轨算法运用到探月卫星的天文学定轨过程中，其中滤波器的观测值为地球视矢量、月球视矢量。
先给出基于日地月信息和EKF的天文学定轨算法的基本描述：
在地球惯性坐标系中，航天器的位置矢量为
Figure 302611DEST_PATH_IMAGE078
，速度矢量为
Figure 416061DEST_PATH_IMAGE079
。设状态变量
Figure 341291DEST_PATH_IMAGE080
，阳质心的坐标为
Figure 503283DEST_PATH_IMAGE081
，月球质心的坐标为
Figure 564779DEST_PATH_IMAGE082
。在只考虑日月摄动的情况下，滤波器的状态方程为：
Figure 849130DEST_PATH_IMAGE083
；
其中
Figure 697875DEST_PATH_IMAGE084
、
Figure 725874DEST_PATH_IMAGE085
、
Figure 845140DEST_PATH_IMAGE086
、
Figure 34813DEST_PATH_IMAGE087
、
Figure 934636DEST_PATH_IMAGE088
、
Figure 438429DEST_PATH_IMAGE089
分别表示的是变量
Figure 740098DEST_PATH_IMAGE090
、
Figure 100672DEST_PATH_IMAGE091
、
Figure 425474DEST_PATH_IMAGE092
、
Figure 795275DEST_PATH_IMAGE093
、、
Figure 916214DEST_PATH_IMAGE095
对时间的导数，
Figure 525050DEST_PATH_IMAGE096
为地球引力常数，
Figure 698543DEST_PATH_IMAGE097
为月球引力常数，
Figure 646907DEST_PATH_IMAGE098
为太阳引力常数。
Figure 614863DEST_PATH_IMAGE099
为航天器到地球球心的距离，
Figure 710995DEST_PATH_IMAGE100
为航天器到太阳质心的距离，
Figure 422599DEST_PATH_IMAGE101
为航天器到太阳质心的距离。
基于地球视矢量信息的观测信息矢量
Figure 225470DEST_PATH_IMAGE102
。
基于月球视矢量信息的观测信息矢量
Figure 364327DEST_PATH_IMAGE103
。
滤波器的状态方程简写为
Figure 947755DEST_PATH_IMAGE104
。
滤波器的观测方程为
Figure 135154DEST_PATH_IMAGE105
，简写为
Figure 120428DEST_PATH_IMAGE106
。
执行步骤一：改进的基于EKF天文自主定轨算法需要首先设定一些初始值。在实际工程中，由于探月卫星在近地环绕阶段可以使用GPS或者通过地面遥测来获得初始的轨道信息，因此
Figure 430186DEST_PATH_IMAGE107
的初始值可以在近地环绕阶段使用其他定轨方式来确定。状态模型噪声的协方差矩阵
Figure 937129DEST_PATH_IMAGE108
赋值取决于模型的精度，在仅考虑日月摄动的情况下可以取值为
Figure 990535DEST_PATH_IMAGE109
，预估状态协方差的初始值可以取为
Figure 310975DEST_PATH_IMAGE111
，观测模型噪声的协方差矩阵
Figure 541100DEST_PATH_IMAGE112
应该根据敏感器的噪声水平来确定，对于敏感器的信噪比在60dB左右的情况，可以设置
Figure 398197DEST_PATH_IMAGE113
，滤波时间赋初值
Figure 92484DEST_PATH_IMAGE114
，滤波迭代次数赋初值为
Figure 727044DEST_PATH_IMAGE116
。
 执行步骤二：根据步骤一计算出地心惯性坐标系
Figure 122254DEST_PATH_IMAGE117
到滤波坐标系
Figure 608730DEST_PATH_IMAGE118
的余弦转换矩阵，然后计算滤波器状态变量在滤波坐标系中的坐标表示，计算观测信息矢量在滤波坐标系中的坐标表示
Figure 963805DEST_PATH_IMAGE120
。
执行步骤三：计算滤波器状态方程的雅克比矩阵
 
Figure 162705DEST_PATH_IMAGE121
 其中：
Figure 5152DEST_PATH_IMAGE122
。
注意以上计算
Figure 264095DEST_PATH_IMAGE123
到的各个矢量的坐标表示是在滤波参考坐标系
Figure 693120DEST_PATH_IMAGE125
下的，继续按照步骤三计算时刻的状态转移矩阵
Figure 646349DEST_PATH_IMAGE127
。
执行步骤四：在滤波坐标系中进行EKF滤波得到
Figure 560079DEST_PATH_IMAGE128
时刻滤波器状态变量的估计值
Figure 100781DEST_PATH_IMAGE129
，在这个过程中观测方程的雅克比矩阵为：
 
Figure 744252DEST_PATH_IMAGE130
Figure 17102DEST_PATH_IMAGE131
。
执行步骤五：进行坐标反变换，计算出第
Figure 746023DEST_PATH_IMAGE128
步滤波器的状态变量的估计值在地心惯性坐标系中的表示形式
Figure 824838DEST_PATH_IMAGE132
；递增
Figure 322815DEST_PATH_IMAGE128
，然后返回执行步骤二。
通过以上步骤的周期运行，就可以从滤波器状态变量的估计值
Figure 265101DEST_PATH_IMAGE133
中获得探月卫星当前时刻位置矢量的估计值
Figure 481319DEST_PATH_IMAGE134
以及速度矢量的估计值
Figure 363824DEST_PATH_IMAGE135
。
图3给出了传统的基于EKF的天文自主定轨算法仿真运行2万秒过程中，探月卫星的位置时间曲线，6000s左右的时候探月卫星从环地轨道过渡到奔月轨道，并开始使用传统的基于EKF天文自主定轨算法来确定轨道，随着探月卫星不断地飞向月球的过程中，卫星位置矢量的x轴分量越来越大，这是奔月轨迹此阶段近似平行于地心惯性坐标系x轴的结果。从图4可以看出，传统的基于EKF的天文自主定轨算法随着时间的增加逐渐发散，这是由于卫星位置矢量的x轴分量过大，逐渐远远超过了位置矢量在y轴和z轴上的分量，从而导致了计算雅克比矩阵时线性化误差不断增大，并最终导致了定轨算法的发散。从图5所示的传统的基于EKF的天文自主定轨算法的定速误差结果可以看出，x轴向定速误差的均值明显非零。
在上述同等条件下仿真改进的基于EKF的天文自主定轨算法，图6和图7给出了改进的基于EKF的天文自主定轨算法运行2万秒之后的定轨结果。从图6可以看出滤波器在这2万秒之内都是稳定的，没有发散的迹象。从图7中也可以看出x轴向的定速误差的均值接近零，这使得图6中的定轨误差不会随着时间的推移越变越大。可见，本发明所提出的方法确实能够提高基于EKF的天文自主定轨算法的稳定性，具有很强的理论和工程实际意义。   

Claims (7)

1.一种改进的基于EKF的天文自主定轨算法，其特征在于所述算法步骤如下：
一、参考坐标系由地心惯性坐标系转为滤波坐标系；
二、在滤波坐标系中线性化和离散化滤波器状态方程；
三、在滤波坐标系中进行EKF滤波；
四、参考坐标系由滤波坐标系转回为地心惯性坐标系。
2.根据权利要求1所述的改进的基于EKF的天文自主定轨算法，其特征在于所述步骤一中，将参考坐标系由地心惯性坐标系转为滤波坐标系之前，需要对EKF滤波器参数初始化，具体步骤如下：
1）给滤波迭代次数                                                
Figure 290708DEST_PATH_IMAGE001
赋初值为
Figure 720553DEST_PATH_IMAGE002
；
2）设滤波器的状态变量
Figure 758916DEST_PATH_IMAGE003
，其中
Figure 299619DEST_PATH_IMAGE004
为地心到航天器的位置矢量，
Figure 943090DEST_PATH_IMAGE005
为航天器的速度矢量；给初始时刻的状态变量
Figure 340573DEST_PATH_IMAGE006
估计值
Figure 69495DEST_PATH_IMAGE007
赋初值，其中上角标
Figure 213556DEST_PATH_IMAGE008
说明该变量为地心惯性坐标系；
3）滤波时间赋初值；
4）给状态模型噪声的协方差矩阵
Figure 217601DEST_PATH_IMAGE010
赋值；
5）给观测模型噪声的协方差矩阵
Figure 496135DEST_PATH_IMAGE011
赋值；
6）给预估状态协方差
Figure 378641DEST_PATH_IMAGE012
赋初值。
3.根据权利要求1或2所述的改进的基于EKF的天文自主定轨算法，其特征在于所述步骤一的具体步骤如下：
1）获取上一时刻航天器位置矢量的估计值
Figure 793441DEST_PATH_IMAGE013
以及这一时刻的观测信息矢量
Figure 735990DEST_PATH_IMAGE014
，其中
Figure 173924DEST_PATH_IMAGE015
表示上一时刻滤波器状态变量估计向量的前三个元素组成的列向量；
2）通过下式计算
Figure 656858DEST_PATH_IMAGE016
与其期望的参考坐标表示之间的夹角
Figure 977298DEST_PATH_IMAGE018
：
Figure 964846DEST_PATH_IMAGE019
其中
Figure 189154DEST_PATH_IMAGE020
表示矩阵或向量的转置；
 3）通过下式计算
Figure 516230DEST_PATH_IMAGE021
和
Figure 862898DEST_PATH_IMAGE022
两个矢量所在平面的单位法方向量
Figure 275424DEST_PATH_IMAGE023
：
Figure 37844DEST_PATH_IMAGE024
 4）通过下式计算地心惯性坐标系
Figure 281743DEST_PATH_IMAGE025
到滤波坐标系
Figure 736996DEST_PATH_IMAGE026
的转动四元数
Figure 699135DEST_PATH_IMAGE027
：
Figure 265246DEST_PATH_IMAGE028
5）通过下式计算地心惯性坐标系
Figure 301335DEST_PATH_IMAGE025
到滤波坐标系的余弦转换矩阵
Figure 373995DEST_PATH_IMAGE029
：
Figure 478217DEST_PATH_IMAGE031
其中
Figure 696709DEST_PATH_IMAGE032
为转动四元数的第
Figure 102599DEST_PATH_IMAGE034
个分量；
6）通过下式计算在滤波坐标系
Figure 349090DEST_PATH_IMAGE026
中的坐标表示
Figure 317046DEST_PATH_IMAGE036
和地心惯性坐标系
Figure 475495DEST_PATH_IMAGE025
中的速度矢量坐标表示
Figure 921520DEST_PATH_IMAGE037
在滤波坐标系中的坐标表示
Figure 191144DEST_PATH_IMAGE039
：
组合滤波坐标系
Figure 86605DEST_PATH_IMAGE026
中的位置矢量和速度矢量，得到滤波器状态变量在滤波坐标系
Figure 71878DEST_PATH_IMAGE026
中的坐标表示为：
Figure 443954DEST_PATH_IMAGE041
，这时
Figure 514678DEST_PATH_IMAGE042
的前三个分量相等，即满足以下等式：
 7）通过下式计算观测信息矢量
Figure 470181DEST_PATH_IMAGE044
在滤波坐标系中的坐标表示
Figure 243282DEST_PATH_IMAGE045
：
Figure 162697DEST_PATH_IMAGE046
。
4.根据权利要求1所述的改进的基于EKF的天文自主定轨算法，其特征在于所述步骤二的具体步骤如下：
1）计算滤波器状态方程的雅克比矩阵
Figure 508545DEST_PATH_IMAGE048
，其中
Figure 619108DEST_PATH_IMAGE049
为
Figure 14317DEST_PATH_IMAGE050
对时间的导数；
2）计算
Figure 563110DEST_PATH_IMAGE051
时刻的状态转移矩阵
Figure 713468DEST_PATH_IMAGE052
，其中为单位矩阵，
Figure 913823DEST_PATH_IMAGE054
为滤波时间步长。
5.根据权利要求1所述的改进的基于EKF的天文自主定轨算法，其特征在于所述步骤三的具体步骤如下：
1）在滤波坐标系中进行状态预估：
其中
Figure 903961DEST_PATH_IMAGE057
表示滤波器的状态方程，
Figure 392711DEST_PATH_IMAGE058
表示的是
Figure 457619DEST_PATH_IMAGE059
时刻滤波器状态变量的初步计算值，为
Figure 145270DEST_PATH_IMAGE061
时刻的预估状态协方差，
Figure 449212DEST_PATH_IMAGE062
为
Figure 989915DEST_PATH_IMAGE059
时刻滤波器的预估状态协方差的初步计算值；
2）在滤波坐标系中进行状态更新：
 
Figure 367807DEST_PATH_IMAGE064
其中
Figure 30869DEST_PATH_IMAGE065
表示滤波器的观测方程，
Figure 494211DEST_PATH_IMAGE066
为
Figure 838605DEST_PATH_IMAGE059
时刻的观测方程的雅克比矩阵，
Figure 398899DEST_PATH_IMAGE067
为增益矩阵，为单位矩阵，
Figure 121185DEST_PATH_IMAGE069
为时刻的预估状态协方差，
Figure 418491DEST_PATH_IMAGE070
为
Figure 95460DEST_PATH_IMAGE059
时刻滤波器状态变量的估计值在滤波坐标系下的表示。
6.根据权利要求1所述的改进的基于EKF的天文自主定轨算法，其特征在于所述步骤四的具体步骤如下：
通过下式对步骤三的结果进行坐标反变换：
Figure 592782DEST_PATH_IMAGE071
通过
Figure 278978DEST_PATH_IMAGE072
可得到航天器的位置矢量的估计值和速度矢量的估计值。
7.根据权利要求6所述的改进的基于EKF的天文自主定轨算法，其特征在于递增滤波迭代次数
Figure 485969DEST_PATH_IMAGE073
，然后跳回至步骤一，通过以上步骤的周期运行，可从滤波器状态变量的估计值
Figure 661735DEST_PATH_IMAGE074
中获得航天器当前时刻位置矢量的估计值和速度矢量的估计值。
CN201310419075.9A 2013-09-16 2013-09-16 一种改进的基于ekf的天文自主定轨算法 Expired - Fee Related CN103438892B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310419075.9A CN103438892B (zh) 2013-09-16 2013-09-16 一种改进的基于ekf的天文自主定轨算法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310419075.9A CN103438892B (zh) 2013-09-16 2013-09-16 一种改进的基于ekf的天文自主定轨算法

Publications (2)

Publication Number Publication Date
CN103438892A true CN103438892A (zh) 2013-12-11
CN103438892B CN103438892B (zh) 2015-09-30

Family

ID=49692588

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310419075.9A Expired - Fee Related CN103438892B (zh) 2013-09-16 2013-09-16 一种改进的基于ekf的天文自主定轨算法

Country Status (1)

Country Link
CN (1) CN103438892B (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104048664A (zh) * 2014-07-01 2014-09-17 南京航空航天大学 一种导航卫星星座自主定轨的方法
CN108872973A (zh) * 2018-08-30 2018-11-23 北京电子工程总体研究所 一种弹道导弹目标定轨的ekf滤波方法
CN109752005A (zh) * 2018-11-23 2019-05-14 中国西安卫星测控中心 一种基于精确轨道模型的航天器初轨确定方法
CN109752006A (zh) * 2018-11-23 2019-05-14 中国西安卫星测控中心 一种非完备外测数据在实时滤波中的使用方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101782391A (zh) * 2009-06-22 2010-07-21 北京航空航天大学 机动加速度辅助的扩展卡尔曼滤波航姿系统姿态估计方法
CN102305630A (zh) * 2011-05-17 2012-01-04 哈尔滨工业大学 基于扩展卡尔曼滤波的sar卫星自主定轨方法
CN102519463A (zh) * 2011-12-13 2012-06-27 华南理工大学 一种基于扩展卡尔曼滤波的导航方法及装置
CN102830414A (zh) * 2012-07-13 2012-12-19 北京理工大学 一种基于sins/gps的组合导航方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101782391A (zh) * 2009-06-22 2010-07-21 北京航空航天大学 机动加速度辅助的扩展卡尔曼滤波航姿系统姿态估计方法
CN102305630A (zh) * 2011-05-17 2012-01-04 哈尔滨工业大学 基于扩展卡尔曼滤波的sar卫星自主定轨方法
CN102519463A (zh) * 2011-12-13 2012-06-27 华南理工大学 一种基于扩展卡尔曼滤波的导航方法及装置
CN102830414A (zh) * 2012-07-13 2012-12-19 北京理工大学 一种基于sins/gps的组合导航方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
XIAOLIN NING. ET AL: "Analysis of Filtering Methods for Satellite Autonomous Orbit Determination Using Celestial and Geomagnetic Measurement", 《MATHEMATICAL PROBLEMS IN ENGINEERING》 *
李钧: "地球同步轨道卫星容错定姿定轨研究", 《中国优秀硕士学位论文全文数据库工程科技II辑》 *
龚享铱: "一种基于极坐标的无源定位与跟踪稳定算法研究", 《信号处理》 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104048664A (zh) * 2014-07-01 2014-09-17 南京航空航天大学 一种导航卫星星座自主定轨的方法
CN108872973A (zh) * 2018-08-30 2018-11-23 北京电子工程总体研究所 一种弹道导弹目标定轨的ekf滤波方法
CN108872973B (zh) * 2018-08-30 2022-07-29 北京电子工程总体研究所 一种弹道导弹目标定轨的ekf滤波方法
CN109752005A (zh) * 2018-11-23 2019-05-14 中国西安卫星测控中心 一种基于精确轨道模型的航天器初轨确定方法
CN109752006A (zh) * 2018-11-23 2019-05-14 中国西安卫星测控中心 一种非完备外测数据在实时滤波中的使用方法

Also Published As

Publication number Publication date
CN103438892B (zh) 2015-09-30

Similar Documents

Publication Publication Date Title
Xiangyu et al. An autonomous optical navigation and guidance for soft landing on asteroids
Martin et al. Design and implementation of a low-cost observer-based attitude and heading reference system
Xiong et al. Modified unscented Kalman filtering and its application in autonomous satellite navigation
CN106873611A (zh) 一种多通道线性自抗扰控制器的设计方法
CN109931955B (zh) 基于状态相关李群滤波的捷联惯性导航系统初始对准方法
Vandersteen et al. Spacecraft attitude estimation and sensor calibration using moving horizon estimation
CN104374388A (zh) 一种基于偏振光传感器的航姿测定方法
CN102175260A (zh) 一种自主导航系统误差校正方法
Abdelrahman et al. Sigma-point Kalman filtering for spacecraft attitude and rate estimation using magnetometer measurements
CN113361013B (zh) 一种基于时间同步稳定的航天器姿态鲁棒控制方法
CN103542853B (zh) 一种估计加速度计漂移的绝对导航滤波方法
Mashtakov et al. Study of the disturbances effect on small satellite route tracking accuracy
CN105737834B (zh) 一种基于轨道平根数的相对导航鲁棒滤波方法
CN103438892A (zh) 一种改进的基于ekf的天文自主定轨算法
CN104215244A (zh) 基于发射惯性坐标系的空天飞行器组合导航鲁棒滤波方法
Lee et al. Vision-based relative state estimation using the unscented Kalman filter
Quan et al. Interlaced optimal-REQUEST and unscented Kalman filtering for attitude determination
Nguyen et al. Lateral stability control of electric vehicle based on disturbance accommodating kalman filter using the integration of single antenna GPS receiver and yaw rate sensor
Yan et al. High-precision simulator for strapdown inertial navigation systems based on real dynamics from GNSS and IMU integration
Hajiyev et al. Attitude determination and control system design of the ITU-UUBF LEO1 satellite
Qi et al. A novel maximum correntropy adaptive extended Kalman filter for vehicle state estimation under non-Gaussian noise
CN102866633B (zh) 微陀螺仪的动态滑模控制系统
CN102582850A (zh) 提高卫星磁控精度的方法
Xiang et al. In-motion initial alignment method for a laser Doppler velocimeter-aided strapdown inertial navigation system based on an adaptive unscented quaternion H-infinite filter
Nordkvist et al. Attitude feedback tracking with optimal attitude state estimation

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
TR01 Transfer of patent right

Effective date of registration: 20180628

Address after: 150001 No. 434 postal street, Nangang District, Harbin, Heilongjiang.

Patentee after: Harbin Institute of Technology Asset Investment Management Co., Ltd.

Address before: 150000 No. 92, West Da Zhi street, Nangang District, Harbin, Heilongjiang.

Patentee before: Harbin Institute of Technology

TR01 Transfer of patent right
TR01 Transfer of patent right

Effective date of registration: 20180723

Address after: 101300 Shunyi District, Beijing, the core area of the airport economy, 6 yuan, 2-35

Patentee after: Harbin Institute of Technology (Beijing) military and Civilian Integration Innovation Research Institute Co., Ltd.

Address before: 150001 No. 434 postal street, Nangang District, Harbin, Heilongjiang.

Patentee before: Harbin Institute of Technology Asset Investment Management Co., Ltd.

TR01 Transfer of patent right
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20150930

Termination date: 20200916

CF01 Termination of patent right due to non-payment of annual fee