基于残差正交性判别的定金磁测自主导航

基于残差正交性判别的定金磁测自主导航

0uf滤波算法独立导航是指卫星能够独立确定位置和速度的能力，不依赖地面。卫星独立导航技术是当前卫星控制技术的发展趋势。近地轨道由于有丰富的地磁场资源,且具有相当好的地磁场模型,可以采用高斯球谐函数来描述地磁场模型,这样磁场的强度和方向是位置的函数,因此能采用磁测的方法对卫星进行定轨导航。使用一台三轴磁强计作为敏感器,利用扩展卡尔曼滤波器即可得到卫星的位置和速度估计值,进而得到卫星的轨道根数。然而在状态方程的非线性严重时,扩展Kalman滤波的状态估计值中会产生较大偏差,滤波结果不能满足精度要求,这就需要对滤波算法做出改进和完善,本文将UKF滤波算法用于磁测自主导航系统中,用以提高滤波精度,改善滤波性能。但在其导航测量值中可能会存在野值时,这会使滤波收敛变慢甚至发散,滤波的稳定性也受到一定的影响,构造一种能够增强滤波收敛性、修正观测信息中野值的不良影响并具有故障检测与隔离功能的高鲁棒性容错滤波器,具有较高的理论和工程应用价值。本文对利用三轴磁强计,以地磁场强度向量为测量量的卫星自主导航的方法进行了理论和仿真研究,在此基础上将基于残差正交性判别的UKF容错滤波算法应用于磁测自主导航方法中,该算法在导航滤波收敛性、观测信息的野值剔除、滤波器故障检测等方面具有较为良好的性能,数学仿真结果证明了该方法的有效性。1动态模型和观测方程1.1j2004轨道动力学方程的建立对于近地卫星的轨道动力学模型采用二体运动模型,在系统状态方程中考虑了地球非球形摄动和大气阻力的摄动,将其它摄动因素等效为高斯白噪声,在J2000地心赤道惯性坐标系中建立卫星轨道动力学方程,即˙X=f(X,t)+w(t)(1)X˙=f(X,t)+w(t)(1)其中X=[xyzvxvyvz]T为卫星在惯性坐标系下的状态向量,w为过程噪声,假定为零均值的高斯白噪声。离散化得到如下方程:ΔXk=∫tk+1tkf(X,t)dt+wk(2)ΔXk=∫tk+1tkf(X,t)dt+wk(2)wk为该采样周期内的等效过程噪声。1.2地磁场向量与卫星位置的关系地球磁场是个向量场,对于近地卫星可略去地球外部磁场影响,这样高斯球谐函数描述的地磁场模型如下:V(r,θ,ϕ)=Re∞∑n=1n∑m=0(Rer)n+1[gmncos(mϕ)+hmnsin(mϕ)]Ρmn(cosθ)(3)式中,Pmn(cosθ)是Legendre函数;Re为地球参考球赤道半径;ϕ为地理经度;θ为地心余纬;gmn,hmn是高斯系数,虽地球内部运动变化而缓慢变化的需要根据测量数据校正。本文采用IGRF2000系数表。地磁场强度向量可以表示为地磁场势函数的负梯度B=-∇V,由此可得到卫星位置与地磁场向量在惯性坐标系下的关系。文中采用低阶的地磁场模型,以三轴磁强计测量得到的地磁向量作为测量量,这需要考虑卫星的姿态,这里假设卫星姿态由姿态测量系统提供,经过坐标变换,得到地磁场向量在体坐标系下的表达式。观测方程可写为z(tk)=Bb(X(tk),tk)+υ(tk)(4)υ(tk)为磁强计测量噪声,假定为的高斯白噪声,其方差为R。2干预ukf容错滤波方法2.1基于利益衡平算法的工性模式UKF是一种新的非线性滤波方法,它不需要对状态方程和量测方程线性化,因此也就没有对高阶项的截断误差。设存在野值或故障的系统模型如下:xk=f(xk-1,uk-1)+wk-1zk=h(xk-1)+vk-1+γρ(φ)(5)其中γ是随机向量,表示故障的大小;存在野值时ρ(φ)=1,不存在野值时ρ(φ)=0,φ表示故障发生的时间。则,UKF算法如下:Xk-1=[ˉxk-1ˉxk-1+(α√Ρx,k)iˉxk-1-(α√Ρx,k)i]Xk/k-1=f(Xk-1,uk-1)ˆx-k=2l∑i=0wmiXi,k/k-1Ρ-x,k=2l∑i=0wci(Xi,k/k-1-ˆx-k)(Xi,k/k-1-ˆx-k)Τ+QkΖk/k-1=h(Xk-1)ˆz-k=2l∑i=0wmiΖi,k/k-1Ρz,k=2l∑i=0wci(Ζi,k/k-1-ˆz-k)(Ζi,k/k-1-ˆz-k)Τ+RkΡxy,k=2l∑i=0wci(Xi,k/k-1-ˆx-k)(Ζi,k/k-1-ˆz-k)Τˆxk=ˆx-k+Ρxy,kΡ-1z,k(zk-ˆz-k)Ρx,k=Ρ-x,k-(Ρxy,kΡ-1z,k)Ρy,k(Ρxy,kΡ-1z,k)Τ式中,α=√l+λ,λ为sigma点的尺度分布参数,λ=β2(l+k)-l(其中β和k是调整参数,l为状态向量的维数)。(√Ρx,k)i为协方差矩阵均方根矩阵的第i列,可通过Cholesky分解(或者奇异值分解的方法)求得。Qk和Rk分别为系统噪声和量测噪声方差阵。2.2滤波器的残差正交性判别方法系统滤波器的残差为ˆek=zk-ˆz-k,可以看作zk和ˆz-k线性组合,其残差序列保持正交性。当故障未发生时,残差ˆek是零均值高斯白噪声,其方差为σk=2l∑i=0wci(Ζi,k/k-1-ˆz-k)(Ζi,k/k-1-ˆz-k)Τ+Rk当故障发生时,残差ˆek的均值不再为零。根据极大似然比检验原理,构造故障检测函数,便可得到故障检测函数:γk=ˆeΤkσ-1kˆek(6)由于残差ˆek是高斯随机向量,因此γk服从自由度为m(观测向量维数)的χ2分布,即γk～χ2(m),根据预先设定的门限值ε,以γk是否超过门限值ε来判定是否发生故障。根据Neyman-pearson准则,当误警率Pfa=α时,由Ρfa=Ρ[γk≥ε|Η0]=α(7)解出的门限值ε可使漏警率Pm=P[γk+1<ε|H1]达到最小。因此门限值可由误警率根据式(7)确定,α又称为显著水平。滤波器收敛性判据:当滤波器发散时,误差协方差阵无界,这时实际的估计误差往往比理论预计的误差大很多倍。可以根据这一特点,利用残差序列的性质来构造滤波器的收敛性判据。由于ˆeTˆe为残差序列的平方和,包含了实际估计误差的信息,而理论预计误差的信息可通过残差序列的方差阵E[ˆeˆeT]来描述。因此可用下式作为判断滤波器是否收敛的依据:ˆeΤkˆek≤γΤr[2l∑i=0wci(Ζi,k/k-1-ˆz-k)(Ζi,k/k-1-ˆz-k)Τ+Rk](8)式中,矩阵的迹被作为一种二次性能指针来代表估计误差的正常水平,γ≥1为可调系数。当上述不等式成立时,滤波器处于正常工作状态;如不满足,说明滤波器的实际误差将超过理论预计值γ倍,此时滤波器处于发散状态。2.3i=0wcii,k/k-1-z-k的滤波状态最优估计算子通过在UKF新息正交性抗野值方法的分析结果,可见在残差正交性抗野值方法中,在对野值进行判断的同时,直接运用了滤波器收敛性判据,实现了对滤波收敛性的判定。修正算法中的统计标量构造如下:d=ˆeΤk+1[2l∑i=0wci(Ζi,k/k-1-ˆz-k)(Ζi,k/k-1-ˆz-k)Τ+Rk]-1ˆek+1(9)则将前面的UKF滤波状态最优估计表达式修改如下:ˆxk=ˆx-k+f(r)Κk(zk-ˆz-k)(10)式中,f(r)为加权因子,当d≤ε时,f(r)=1,当d>ε时f(r)为d的倒数。f(r)与当前的观测值Zk有关,但其对滤波的综合影响,体现在增益矩阵Kk上,具有根据观测值自动修正增益矩阵的作用。此处d符合三自由度χ2分布,取显著水平α为0.01,通过χ2查表,野值检测门限值ε设为11.345。综上,构造了一种残差正交性判别的UKF容错滤波方法,在观测值存在野值的情况下,能够消除其影响,使滤波估计趋于收敛。3模拟仿真lm在仿真运算中,磁强计测量误差方差取R=(180nT)2,系统误差取10-6量级,卫星姿态精度为0.001°。按照上述的UKF容错滤波算法进行数学仿真。UKF滤波中的参数选取为:α=1,β=0,k=3。卫星的初始轨道参数如下:半长轴a=6799.4km,偏心率e=0.0134,轨道倾角i=65°,初始位置误差为[60,50,70]km,初始速度误差[0.4,0.2,0.3]km,采样周期T=1s,仿真时间15000s,在没有加入野值的情况下,UKF算法位置和速度估计误差如图1所示。在3000s-3030s连续30点处加入1800nT(x,y,z三方向)野值,UKF算法估计误差曲线如图2;容错滤波算法仿真曲线如图3。在2501-2511s点+1500nT(x,y,z三方向)野值,3001-3005s点+1500nT(x和y方向)野值,4001s-4005s点+1000nT(y和z方向)野值,5001s和5002s+1500nT(x和z方向)野值,5501s-5504s点+2000nT(x,y,z三方向)野值,6001s-6010s点+1500nT(y方向)野值,7001s-7005s点+1500nT(x,y,z三方向)野值,8001s点+1500nT(x,y,z三方向)野值,UKF仿真曲线如图4,容错滤波仿真曲线如图5。通过仿真结果可以看出,容错滤波器具有较强的野值修正能力及较强的野值检测灵敏度,在出现成片野值时滤波仍有很好的收敛性。另外通过仿真比