[0049] 当天文导航子系统Unscented卡尔曼滤波的估计均方误差阵Pk′大于已有目标天体星历误差均方误差阵时,即Pk′>Peph,则不进行星历校正,直接进行第9步信息融合;当Pk′
[0050] 8.对目标天体星历误差进行建模、估计并反馈校正
[0051] A.建立目标天体星历误差状态模型
[0052] 将捕获段内其误差特性考虑为常值误差,在日心惯性坐标系中建立目标天体星历误差状态模型为:
[0053]
[0054] 式中, 为日心惯性坐标系中目标天体星历的三轴位置误差的微分,可离散化后简写为:
[0055] Xerr(k)=Ferr(Xerr(k-1),k-1)+Werr(k-1) (13)
[0056] 式中,状态转移函数Ferr(Xerr(k-1),k-1)=Φerr,k,k-1Xerr,k-1,Φerr,k,k-1为第k-1时刻到第k时刻的状态转移矩阵,Xerr(k)为第k时刻目标天体星历误差状态向量,且Xerr(k)=Xerr,k,Werr(k-1)为第k-1时刻目标天体星历误差状态模型误差。
[0057] B.建立目标天体星历误差量测模型
[0058] 目标天体星历误差的量测模型可以表示为:
[0059] Zerr=H3(Xerr(k),k)+V3 (14)
[0060] 式中,H3(Xerr(k),k)为k时刻的量测函数,V3为目标天体星历误差量测噪声。
[0061] C.获取目标天体星历误差量测量
[0062] 目标天体星历误差量测量Zerr可以表示为:
[0063]
[0064] 式中, 为无线电导航子系统获得的相对于太阳的位置和速度, 为天文导航子系统获得的相对于目标天体的位置和速度, 为目标天体相对于太阳的位置和
速度,从天体星历数据库中获得。
[0065] D.对目标天体星历误差进行卡尔曼滤波估计
[0066] 根据步骤A和步骤B建立的目标天体星历误差状态模型和量测模型,以及步骤C所获取的目标天体星历误差量测量,利用卡尔曼滤波方法,对目标天体星历误差进行估计,获得目标天体星历误差估计状态向量与估计均方误差阵,具体如下:
[0067] 目标天体星历误差状态向量的一步预测
[0068]
[0069] 式中, 为k-1时刻火星星历误差一步预测状态向量。
[0070] 目标天体星历误差状态向量的估计均方误差阵一步预测
[0071] Perr,k/k-1=Φerr,k,k-1Perr,k-1Φerr,k,k-1T+Qerr,k (17)
[0072] 式中,Perr,k-1为k-1时刻目标天体星历误差状态向量的估计均方误差阵,Qerr,k为k时刻目标天体星历误差状态模型误差协方差阵。
[0073] 卡尔曼滤波增益
[0074] Kerr,k=Perr,k/k-1Herr,kT(Herr,kPerr,k/k-1Herr,kT+Rerr,k)-1 (18)[0075] 式中,Herr,k为k时刻目标天体星历误差量测矩阵,Herr,kXerr,k=H3(Xerr,k),Rerr,k为k时刻目标天体星历误差量测模型误差协方差阵。
[0076] 目标天体星历误差估计状态向量
[0077]
[0078] 式中,Zerr,k为k时刻目标天体星历误差量测量。
[0079] 目标天体星历误差估计均方误差阵
[0080] Perr,k=(I-Kerr,kHerr,k)Perr,k/k-1 (20)
[0081] 式中,I为单位阵。
[0082] E.对目标天体星历误差进行反馈校正
[0083] 将步骤D中获得的目标天体星历误差 和目标天体星历误差估计均方误差阵Perr,k反馈回深空探测器的第一状态模型和第二状态模型中,并重新确定第一状态模型和第二状态模型的模型误差协方差阵,最后将校正后的模型误差协方差阵输入至天文导航子系统Unscented卡尔曼滤波和无线电导航子系统Unscented卡尔曼滤波中,修正下一时刻的导航结果。
[0084] 9.信息融合
[0085] 利用星历修正系统获得的星历误差,将无线电导航系统的导航信息转换至以目标天体为中心的惯性坐标系中,与天文导航系统的导航信息进行融合。当探测器不在无线电测控范围内时,无线电导航子系统没有输入的无线电导航量测量,此时对天文导航子系统进行Unscented卡尔曼滤波,无线电导航子系统只利用第二状态模型进行时间更新;当探测器处于无线电测控范围内时,无线电导航子系统有输入的无线电导航量测量,对两个导航子系统同时进行Unscented卡尔曼滤波;
[0086] 最终输出第k时刻在目标天体为中心的惯性坐标系中表示探测器位置和速度的估计状态向量和估计均方误差阵,并根据修正后的目标行星星历,将该结果转换至日心惯性坐标系中,输出在日心惯性坐标系中表示探测器位置和速度的估计状态向量和估计均方误差阵,将这些导航信息分别返回天文导航子系统和无线电导航子系统中,用于k+1时刻的位置、速度导航信息的估计,k=1,2,...;
[0087] 本发明的原理是:首先选择基于太阳和八大行星引力的轨道动力学模型作为系统状态模型,根据天文导航和无线电导航测量的不同特点,分别建立在目标天体为中心的惯性坐标系和日心惯性坐标系中的两个状态模型;之后根据天文导航敏感器和无线电导航接收装置的工作原理,获取并处理的天文导航敏感器和无线电导航接收装置所测量的量测量;然后,建立天文导航子系统和无线电导航子系统的量测模型;此后,由于探测器的状态模型和量测模型都呈现非线性特性,且系统噪声为非高斯噪声,因此采用Unscented卡尔曼滤波方法,对两个子系统的探测器导航信息进行估计;其次,由于天文导航可以测量高精度的相对目标天体导航信息,而无线电导航可以测量高精度的相对于太阳导航信息;结合目标天体、探测器、太阳的几何关系,可以确定目标天体的星历,与原有的目标天体星历相比较,即可得到目标天体星历误差量测量;为了获得更为准确的目标天体星历误差,结合捕获段持续时间短、目标天体星历误差在捕获段变化小的特点,建立捕获段目标天体星历误差的状态模型和量测模型,并根据目标天体星历误差状态模型和量测模型都为线性模型的特点,采用卡尔曼滤波方法,估计目标天体星历误差,并将所估计的目标天体星历误差反馈至第一状态模型和第二状态模型中,修正目标天体星历数据,提高状态模型的模型精度;最后通过信息融合方法,有效利用天文导航子系统和无线电子系统的测量信息,提高对相对于目标天体和相对于日心的探测器导航信息的估计精度。
[0088] 本发明与
现有技术相比的优点在于:利用星历校正过程所估计的目标天体星历误差,一方面实现了对目标天体星历误差的高精度估计,同时获得了相对目标天体和相对日心的探测器高精度导航信息,另一方面修正了导航系统的状态模型,减小了目标天体星历误差对状态模型精度的影响,进一步提高了深空探测器的导航精度。
附图说明
[0089] 图1为本发明基于星历修正的天文无线电组合导航
流程图。
[0090] 图2为本发明天文导航子系统所用天文导航敏感器的成像原理图。
具体实施方式
[0091] 如图1所示,前述技术方案中所涉及的目标天体可以为火星、金星、木星、土星等太阳系内的天体,以下以火星作为
实施例,说明本发明的具体实施过程:
[0092] 1.建立深空探测器基于太阳和八大行星引力轨道动力学的状态模型
[0093] 首先初始化位置、速度,设X=[x y z vx vy vz]T为在日心惯性坐标系中的状态向量,x,y,z,vx,vy,vz分别为探测器在日心惯性坐标系中三轴的位置和速度,X′=[x′ y′ Tz′ vx′ vy′ vz′]为在火心惯性坐标系中的状态向量,x′,y′,z′,vx′,v′y,vz′分别为探测器在火心惯性坐标系中三轴的位置和速度,上述变量都是与t有关的函数,根据探测器的轨道设计,选取探测器的位置和速度初值为X(0)和X′(0);其次初始化火星星T
历误差初值为Xerr(0)=[xerr yerr zerr vxerr vyerr vzerr],xerr,yerr,zerr,vxerr,vyerr,vzerr分别为日心坐标系中火星星历三轴的位置误差和速度误差。
[0094] A.在火心惯性坐标系中建立深空探测器基于太阳和八大行星引力轨道动力学的第一状态模型,即天文导航子系统的状态模型;
[0095] 考虑太阳和火星、地球等太阳和八大行星对探测器的引力作用,选取火心惯性坐标系,可得深空探测器在火心惯性坐标系中第一状态模型,即天文导航子系统的状态模型为:
[0096]
[0097] 式中,x′,y′,z′为探测器在火心惯性坐标系中三轴位置,vx′,v′y,vz′为探测器在火心惯性坐标系中三轴速度, 为探测器在火心惯性坐标系中三轴位置的微分, 为探测器在火心惯性坐标系中三轴速度的微分,μs、μm和μic分别为太阳、火星和第ic颗行星的引力常数;r′ ps为日心到探测器的距离;r′ pm为火星到探测器的距离;r′ms为日心到火心的距离; 为第ic颗行星到探测器的距离;r′ mi为第ic颗行星质心到火心的距离;(xs′,ys′,zs′), 分别为太阳、第ic颗行星在火心惯性坐标
系中的三轴位置坐标,可根据时间由行星星历表获得,wx′,wy′,wz′分别为第一状态模型中探测器三轴的状态模型误差;ic表示太阳和八大行星中从内至外的第ic颗行星,ic=
1,2,3...,N(ic≠4),N=8,由于ic=4表示第4颗行星(火星),因此不包含在求和表达式中。
[0098] 式(1)中的各变量都是与时间t有关的变量,可简写为
[0099]
[0100] 式中,X′(t)=[x′,y′,z′,vx′,v′y,vz′]T为第一状态模型的状态向量,为X′(t)的微分,f1(X′(t),t)为第一状态模型的系统非线性连续状态转移函数,T
w′(t)=[wx′,w′y,wz′]为第一状态模型的状态模型误差。
[0101] B.在日心惯性坐标系中建立深空探测器基于太阳和八大行星引力轨道动力学的第二状态模型,即无线电导航子系统的状态模型;
[0102] 考虑太阳和火星、地球等太阳和八大行星对探测器的引力作用,选取日心惯性坐标系,可得深空探测器在日心惯性坐标系中展开为分量形式的第二状态模型,即无线电导航子系统的状态模型为:
[0103]
[0104] 式中,x,y,z为探测器在日心惯性坐标系中三轴位置,vx,vy,vz为探测器在日心惯性坐标系中三轴速度, 为探测器在日心惯性坐标系中三轴位置的微分, 为探测器在日心惯性坐标系中三轴速度的微分,μs为太阳引力常数, 为第ic个行星的引力常数;rps为日心到探测器的距离; 为第ic个行星到探测器的距离; 为第ic个行星质心到日心的距离; 分别为第ic个行星在日心惯性坐标系中的坐标,可根据时间由行星
星历表获得,wx,wy,wz分别为第二状态模型中探测器三轴的状态模型误差;
[0105] 式(3)中的各变量都是与时间t有关的变量,可简写为:
[0106]T
[0107] 式中,X(t)=[x,y,z,vx,vy,vz]为第二状态模型的状态向量, 为X(t)的微T
分,f2(X(t),t)为第二状态模型系统非线性连续状态转移函数,w(t)=[wx,wy,wz]为第二状态模型的状态模型误差。
[0108] 2.分别建立天文导航子系统和无线电导航子系统量测模型;
[0109] (1)天文导航子系统量测模型
[0110] 火星与第i颗背景恒星之间角度信息θmi的大小在不同坐标系中是相同的,因此其表达式为:
[0111]
[0112] 式中, 为在火星敏感器测量坐标系中从火星至探测器的单位矢量, 为在惯性坐标系中火星到探测器的单位矢量,可表示为:
[0113]
[0114] 式中,(x′,y′,z′)为探测器在火心惯性坐标系中三轴位置坐标, 为火星图像中第i颗背景恒星在火星敏感器测量坐标系中的单位矢量, 为在惯性坐标系中第i颗背景恒星星光的单位矢量,i=1,2,3,可由恒星星历数据库直接得到,
[0115] 同理可得火卫一和火卫二与其第i颗背景恒星之间角度信息θpi和θdi的表达式为:
[0116]
[0117]
[0118] 设天文导航子系统量测量 Z1= [θ m1,θm2,θm3,θp1,θp2,θp3,θd1,θd2,θd3]T, 天文导航子系统量测噪声
分别为测
量θm1,θm2,θm3,θp1,θp2,θp3,θd1,θd2,θd3的观测误差,由于各变量都是与时间t有关的变量,则可根据式(6)~式(8)建立天文导航子系统量测模型的表达式为
[0119] Z1(t)=h1[X′(t),t]+v1(t) (9)
[0120] 式中,h1[X′(t),t]为天文导航子系统非线性连续量测函数。
[0121] (2)无线电导航子系统量测模型
[0122] 无线电测控站可根据无线电信号的时延与多普勒频移测量探测器到测控站的距离与距离变化率,一般主测控站附近还配有多个副站,根据多个测控站测量的距离与距离变化率可以获得探测器相对于地心的位置与速度,无线电测控导航所获得的距离与距离变化率表达式为:
[0123]
[0124]
[0125] 式中,Ri为第i个测控站的测距量测信息(探测器与测控站之间的距离), 为第i个测控站的距离变化率量测信息(探测器与测控站之间的距离变化率),i=1,2,3,(xf,yf,zf)为地面测控站在地心惯性坐标系中的位置坐标。
[0126] 设无线电导航子系统量测量 无线电导航子系统量测噪声 分别为三个测控站所测量得到的探测器
与测控站之间的距离和距离变化率, 分别为测量 的
观测误差,由于各变量都是与时间t有关的变量,则可根据式(10)和式(11)建立无线电导航子系统量测模型的表达式为:
[0127] Z2(t)=h2[X(t),t]+v2(t) (12)
[0128] 式中,h2[X(t),t]为无线电导航子系统非线性连续量测函数。
[0129] 3.对步骤1和步骤2中的状态模型和量测模型进行离散化
[0130] X′(k)=F1(X′(k-1),k-1)+W′(k-1) (13)
[0131] X(k)=F2(X(k-1),k-1)+W(k-1) (14)
[0132] Z′(k)=H1(X′(k),k)+V1(k) (15)
[0133] Z(k)=H2(X(k),k)+V2(k) (16)
[0134] 式中,k=1,2,…,F1(X′(k-1),k-1)和F2(X(k-1),k-1)分别为f1(X′(t),t)和f2(X(t),t)离散后从第k-1时刻到第k时刻的非线性状态转移函数,H1(X′(k),k)和H2(X(k),k)分别为h1(X′(t),t)和h2(X(t),t)离散后第k时刻的非线性量测函数,
W′(k),W(k),V1(k),V2(k)分别为w′(t),w(t),v1(t)和v2(t)离散后第k时刻的等效噪声,且W′(k)和V1(k)、W(k)和V2(k)互不相关。
[0135] 4.天文导航和无线电导航量测量的获取及处理
[0136] (1)天文导航子系统量测量的获取及处理
[0137] 图2给出了天文导航子系统所用火星敏感器的成像原理图,火卫一敏感器、火卫二敏感器成像过程与之类似。火星敏感器主要由
光学透镜和二维成像面阵组成,在敏感器测量坐标系O′XcYcZc中沿火星到探测器的矢量方向 火星反射太阳光线射向火星敏感器,此时,火星在火星敏感器测量坐标系中的坐标为(xc,yc,zc);火星敏感器的光学透镜以焦距f将火星的光线折射后成像在二维成像面阵上,二维成像面阵将照在每个成像单元上的图像
亮度信号储存。根据敏感器的成像原理,天文导航子系统量测量的处理过程如下所示:
[0138] A.天体图像的处理
[0139] 由于火星在二维成像面阵上的图像并不是一个点,而是一个圆,通过质心识别等图像处理技术确定火星图像在二维成像平面坐标系OX2dY2d的质心(x2d,y2d),这一中心可以用像元像线坐标系OplXplYpl中的像元像线(p,l)表示。
[0140] B.质心坐标从像元像线坐标系转换至二维成像平面坐标系的转换
[0141] 根据像元像线坐标系与二维成像平面坐标系之间的关系,将火星质心坐标从像元像线坐标系转换至二维成像平面坐标系:
[0142]
[0143] 式中,(x2d,y2d)为火星在二维成像平面OX2dY2d中的坐标,p和l分别为火星在火星敏感器二维成像平面的像元和像线, 为由毫米转为
像素的敏感器转换矩阵,p0和l0分别为火星敏感器中心在像元像线坐标系OXplYpl中的像元和像线。
[0144] C.质心坐标从二维成像平面坐标系转换至敏感器测量坐标系的转换
[0145] 根据透镜成像原理,将火星质心坐标从二维成像平面坐标系转换至敏感器测量坐标系:
[0146]
[0147] 式中,f为火星敏感器的焦距, 为在火星敏感器测量坐标系中从火星至探测器的单位矢量。
[0148] 同理可获得火星图像中第i颗背景恒星在火星敏感器测量坐标系中的单位矢量i=1,2,3。
[0149] D.计算星光角距
[0150] 根据在火星敏感器测量坐标系中火星至探测器的单位矢量 与火星图像中第i颗背景恒星的单位矢量 计算两个矢量之间的星光角距θmi
[0151]
[0152] 同理可获得火卫一与其背景恒星、火卫二与其背景恒星之间的星光角距θpi,θdi。
[0153] (2)无线电导航子系统量测量的获取及处理
[0154] 选取探测器与测控站之间的距离变化率 和探测器与测控站之间的距离R作为无线电导航子系统的量测量。
[0155] 探测器与测控站之间的距离变化率 是通过测量探测器接收的载波信号频率f′rec与地面站发射的载波信号频率f0相比的频率偏移来获得的,具体可表示为:
[0156]
[0157] 式中,为探测器与地面站间的相对距离变化率,c为光在
真空中的传播速度,f0为地面站发射的无线电信号的固有频率,f′rec为探测器上的接收机接收到的无线电信号的频率,δfatm为大气层对信号的时延,δf0为由信号源
本振频率的漂移引起的误差,由于目前地面站多采用USO(Ultra Stable Oscillators),该误差的量级很小。
[0158] 探测器与测控站之间的距离R是通过测量无线电信号由地面测控站发射至探测器再返回的时间延迟而获得的,具体可表示为:
[0159] R=cτ/2+cδtR-cδtT+δρatm (21)
[0160] 式中,c为光在真空中的传播速度,τ为无线电信号由地面测控站发射至探测器再返回的时间延迟,δtR为接收装置时钟同步误差,δtT为发射装置时钟同步误差,δρatm为由
原子钟差引起的距离测量误差。
[0161] 5.对天文导航子系统进行Unscented卡尔曼滤波
[0162] 根据第一状态模型式(13)、天文导航子系统量测模型式(15)、由天文导航敏感器获得的量测量式(19),进行天文导航子系统Unscented卡尔曼滤波。具体步骤如下:
[0163] A.初始化
[0164]
[0165] 式中, 为第0时刻(初始时刻)在火心惯性坐标系中探测器的三轴位置和速度估计值,x0′为第0时刻(初始时刻)在火心惯性坐标系中探测器的三轴位置和速度真实值,P0′为状态向量的初始均方误差阵。
[0167] 在天文导航子系统第k-1时刻状态向量 附近选取一系列样本点,这些样本点的均值和均方误差阵分别为 和 设状态向量为6×1维,那么13个天文导航子系统
的样本点χ′0,k,...,χ′12,k及其权重W0′…W′12分别如下:
[0168]
[0169]
[0170]
[0171] 式中,当P′k-1=A′ TA′时, 取A′的第j行,当P′k-1=A′A′ T时,取A′的第j列,得第k-1时刻采样点χ′k-1的统一表达式为:
[0172]
[0173] C.时间更新
[0174] 天文导航子系统状态向量的一步预测χ′k|k-1为:
[0175] χ′k|k-1=F1(χ′k-1,k-1) (25)
[0176] 天文导航子系统所有采样点状态向量的一步预测加权后结果 为:
[0177]
[0178] 式中,Wj′为第j个采样点的权值;
[0179] 天文导航子系统状态向量的估计均方误差阵一步预测Pk′-为:
[0180]
[0181] 式中,Qk′为k时刻天文导航子系统的状态模型误差协方差阵;
[0182] 天文导航子系统采样点对应的量测估计向量Z′k|k-1
[0183] Z′k|k-1=H1(χ′k|k-1,k) (28)
[0184] 天文导航子系统所有采样点量测估计加权向量
[0185]
[0186] D.量测更新
[0187] 天文导航子系统量测均方误差阵 为:
[0188]
[0189] 式中,Rk′为k时刻天文导航子系统的量测噪声协方差阵;
[0190] 天文导航子系统状态向量量测量均方误差阵
[0191]
[0192] 天文导航子系统滤波增益Kk′为:
[0193]
[0194] 天文导航子系统的估计状态向量 和估计均方误差阵Pk′为
[0195]
[0196]
[0197] 6.对无线电导航子系统进行Unscented卡尔曼滤波
[0198] 根据第二状态模型式(14)、无线电导航子系统量测模型式(16)、由无线电接收装置获得的量测量式(20)和式(21),进行无线电导航子系统Unscented卡尔曼滤波。具体步骤如下:
[0199] A.初始化
[0200]
[0201] 式中, 为第0时刻(初始时刻)无线电导航子系统在日心惯性坐标系中探测器的三轴位置和速度估计值,x0为第0时刻(初始时刻)无线电导航子系统在日心惯性坐标系中探测器的三轴位置和速度真实值,P0为无线电导航子系统状态向量的初始均方误差阵。
[0202] B.计算采样点
[0203] 在无线电导航子系统第k-1时刻状态向量 附近选取一系列样本点,这些样本点的均值和均方误差阵分别为 和Pk-1。设状态向量为6×1维,那么13个无线电导航子
系统的样本点χ0,k,...,χ12,k及其权重W0…W12分别如下:
[0204]
[0205]
[0206]T T
[0207] 式中,当Pk-1=A A时, 取A的第j行,当Pk-1=AA 时, 取A的第j列,得第k-1时刻采样点χk-1的统一表达式为:
[0208]
[0209] C.时间更新
[0210] 无线电导航子系统状态向量的一步预测χk|k-1为:
[0211] χk|k-1=F2(χk-1,k-1) (38)
[0212] 无线电导航子系统所有采样点状态向量的一步预测加权后结果 为:
[0213]
[0214] 式中,Wj为第j个采样点的权值;
[0215] 无线电导航子系统状态向量的估计均方误差阵一步预测 为:
[0216]
[0217] 式中,Qk为k时刻无线电导航子系统的状态模型误差协方差阵;
[0218] 无线电导航子系统采样点对应的量测估计向量Zk|k-1
[0219] Zk|k-1=H2(χk|k-1,k) (41)
[0220] 无线电导航子系统所有采样点量测估计加权向量
[0221]
[0222] D.量测更新
[0223] 无线电导航子系统量测均方误差阵 为:
[0224]
[0225] 式中,Rk为无线电导航子系统量测噪声协方差阵;
[0226] 无线电导航子系统状态向量量测量均方误差阵 为:
[0227]
[0228] 无线电导航子系统滤波增益Kk为:
[0229]
[0230] 无线电导航子系统估计状态向量 和估计均方误差阵Pk为:
[0231]
[0232]
[0233] 7.判定是否需要进行火星星历校正
[0234] 当天文导航子系统Unscented卡尔曼滤波的估计均方误差阵Pk′大于已有火星星历误差均方误差阵时,即Pk′>Peph,则不进行星历校正,直接进行第9步信息融合;当Pk′
[0235] 8.对火星星历误差进行建模、估计并反馈校正
[0236] 无线电测控导航主要的观测量为探测器相对于地面站的距离与距离变化率等,这种导航方式可以直接测量相对于太阳的导航信息,天文导航可以直接测量相对于目标天体(如火星)的导航信息。由于目标天体的星历存在误差(200m~100km),而无线电导航方法仅能获取相对太阳的高精度导航信息,因此由该方法获得的相对目标天体导航信息精度低;虽然天文导航方法可以获得较为准确的相对目标天体导航信息,因此由该方法无法获得相对太阳高精度的导航信息。因此需要对目标天体星历误差进行估计,并反馈校正由目标天体星历误差引起的导航误差。
[0237] 天文导航子系统获得的相对于火星的位置和速度为 天文导航子系统获得的相对于太阳的位置和速度为 为火星相对于太阳的位置
和速度,可从天体星历数据库获得);无线电导航子系统获得的相对于太阳的位置和速度为 无线电导航子系统获得的相对于火星的位置和速度为 由
此可以看出,单一导航系统会受到火星星历误差的影响,无法同时满足探测器相对太阳和火星高精度导航的需求。因此可利用天文导航子系统Unscented卡尔曼滤波的结果和无线电导航子系统Unscented卡尔曼滤波的结果对火星星历误差进行估计,具体步骤如下:
[0238] A.建立火星星历误差状态模型
[0239] 考虑到捕获段持续时间短(约40天)这一特点,火星星历误差变化较小,将捕获段内火星的星历误差特性考虑为常值误差,在日心惯性坐标系中建立火星星历误差状态模型为:
[0240]
[0241] 式中, 为日心惯性坐标系中火星星历的三轴位置误差的微分,离散化后简写为:
[0242] Xerr(k)=Ferr(Xerr(k-1),k-1)+Werr(k-1) (49)
[0243] 式中,状态转移函数Ferr(Xerr(k-1),k-1)=Φerr,k,k-1Xerr,k-1,Φerr,k,k-1为第k-1时刻到第k时刻的状态转移矩阵,Xerr(k)为第k时刻火星星历误差状态向量,且Xerr(k)=Xerr,k,Werr(k-1)为第k-1时刻火星误差状态模型误差。
[0244] B.建立火星星历误差量测模型
[0245] 因此火星星历误差的量测模型可以表示为:
[0246] Zerr=H3(Xerr(k),k)+V3 (50)
[0247] 式中,H3(Xerr(k),k)为k时刻的量测函数,V3为火星星历误差量测噪声。
[0248] C.获取火星星历误差量测量
[0249] 火星星历误差量测量Zerr可以表示为:
[0250]
[0251] 式中, 为无线电导航子系统获得的相对于太阳的位置和速度, 为天文导航子系统获得的相对于火星的位置和速度, 为火星相对于太阳的位置和速度,从天体星历数据库中获得。
[0252] D.对火星星历误差进行卡尔曼滤波估计
[0253] 根据步骤A和步骤B建立的火星星历误差状态模型式(49)和量测模型式(50),以及步骤C所获取的火星星历误差量测量式(51),利用卡尔曼滤波方法,对火星星历误差进行估计,获得火星星历误差估计状态向量与估计均方误差阵,具体如下:
[0254] 火星星历误差状态向量的一步预测
[0255]
[0256] 式中, 为k-1时刻火星星历误差一步预测状态向量。
[0257] 火星星历误差状态向量的估计均方误差阵一步预测T
[0258] Perr,k/k-1=Φerr,k,k-1Perr,k-1Φerr,k,k-1+Qerr,k (53)
[0259] 式中,Perr,k-1为k-1时刻火星星历误差状态向量的估计均方误差阵,Qerr,k为k时刻火星星历误差状态模型误差均方误差阵。
[0260] 卡尔曼滤波增益T T -1
[0261] Kerr,k=Perr,k/k-1Herr,k(Herr,kPerr,k/k-1Herr,k+Rerr,k) (54)
[0262] 式中,Herr,k为k时刻火星星历误差量测矩阵,Herr,kXerr,k=H3(Xerr,k),Rerr,k为k时刻火星星历误差量测模型误差协方差阵。
[0263] 火星星历误差估计状态向量
[0264]
[0265] 式中,Zerr,k为k时刻火星星历误差量测量。
[0266] 火星星历误差估计均方误差阵
[0267] Perr,k=(I-Kerr,kHerr,k)Perr,k/k-1 (56)
[0268] 式中,I为单位阵。
[0269] E.对火星星历误差进行反馈校正
[0270] 将步骤D中获得的火星星历误差 和火星星历估计均方误差阵Perr,k反馈回深空探测器的第一状态模型和第二状态模型中,并重新确定第一状态模型和第二状态模型的模型误差协方差阵Qk′和Qk,最后将校正后的模型误差协方差阵Qk′和Qk输入至天文导航子系统Unscented卡尔曼滤波和无线电导航子系统Unscented卡尔曼滤波中,修正下一时刻的导航结果。
[0271] 9.信息融合
[0272] 当探测器不在无线电测控范围内时,无线电导航子系统没有输入的无线电导航量测量,此时对天文导航子系统进行Unscented卡尔曼滤波,包括时间更新和量测更新(第5步的步骤C和步骤D)等步骤,无线电导航子系统只进行时间更新(第6步的步骤C);当探测器处于无线电测控范围内时,无线电导航子系统有输入的无线电导航量测量,对两个子系统同时进行Unscented卡尔曼滤波,都进行时间更新和量测更新。
[0273] 在滤波过程中得到的两个局部估计状态向量两个估计均方误差阵 按下式进
行融合,得到全局的估计状态向量和全局估计均方误差阵分别为:
[0274]
[0275]
[0276] 将全局估计结果反馈给两个导航子系统,作为k时刻两个导航子系统的估计结果:
[0277]
[0278]
[0279]
[0280]
[0281]
[0282] 式中, 为信息分配因子。
[0283] 信息分配因子选择的基本原则是在满足信息守恒公式的前提下与局部导航滤波精度成正比,为了使导航系统具有更强的自适应能力和容错能力,使用基于估计均方误差阵范数的动态分配信息因子的
算法。
[0284] 令
[0285]
[0286] 式中,||·||F为Frobenius范数,即对于任意矩阵A,
[0287] 最终将式(59)和式(60)获得的第k时刻在火心惯性坐标系中和在日心惯性坐标系中的估计状态向量 和估计均方误差阵P1(k),P2(k)输出,估计状态向
量 分别表示在火心惯性坐标系中和日心惯性坐标系中探测器的位置、速度信
息,输出的估计均方误差阵P1(k),P2(k)表示了滤波估计的性能,并将这些导航信息分别返回天文导航子系统和无线电导航子系统中,用于k+1时刻的位置、速度导航信息,k=
1,2,...。
[0288] 本发明
说明书中未作详细描述的内容属于本领域专业技术人员公知的现有技术。