[发明专利]基于扩展卡尔曼滤波的SAR卫星自主定轨方法

摘要

基于扩展卡尔曼滤波的SAR卫星自主定轨方法，属于卫星自主定轨技术领域。现有卫星自主定轨技术主要采用GPS，由于GPS是一种半自主定轨方式，因此存在自主性和实时性差的问题。本发明中的定轨方法就可以解决这些问题。本发明的方法是：步骤一：建立基于轨道动力学的卫星运动方程，进而得到状态方程；步骤二：建立以SAR到地面标识点之间的距离和SAR与地面标识点之间的多普勒频移为观测量的观测方程；步骤三：建立扩展卡尔曼滤波的递推方程；步骤四：得到卫星的状态信息，即得到卫星位置和速度矢量。本发明的方法提高了卫星的定轨精度，具有自主性和实时性的特点，能够实现低轨卫星的高精度实时自主定轨。

主权项

1.基于扩展卡尔曼滤波的SAR卫星自主定轨方法，其特征在于：所述方法包括以下步骤：步骤一：建立基于轨道动力学的卫星运动方程，进而得到状态方程；建立包含J₂项的卫星轨道动力学方程，在J2000惯性坐标系下，卫星的轨道动力学方程为：X·=F[X(t),t]---(1)]]>写成状态方程，即：dxdt=vxdydt=vydzdt=vzdvxdt=-μxr3[1-J2(Rer)2(7.5z2r2-1.5)]+ΔFxdvydt=-μyr3[1-J2(Rer)2(7.5z2r2-1.5)]+ΔFydvzdt=-μzr3[1-J2(Rer)2(7.5z2r2-4.5)]+ΔFz---(2)]]>式中，上述式(2)简写为状态矢量X＝[x，y，z，v_x，v_y，v_z]^T，x，y，z，v_x，v_y，v_z分别为卫星在惯性坐标系的X，Y，Z方向上的位置和速度；μ为地心引力常数；r为卫星位置参数矢量；J₂为摄动系数；ΔF_x，ΔF_y，ΔF_z分别为X，Y，Z方向上的非球形摄动的高阶项及日、月引力摄动，太阳光压摄动以及大气阻力摄动未建模摄动的影响；上角标T为转置；R_e为地球半径；在所述简化模型中，上述这些摄动的影响用系统噪声w(t)来表示，上述式(2)为连续方程，将其离散化，得到状态方程，即：X(j+1)＝Φ_k+1，kX(k)+W(k)    k＝1，2，3...(3)式中，Φ_k+1，k为k至k+1时刻的一步状态转移矩阵；X(k)为第k时刻的状态；W(k)为系统噪声矩阵，对于白噪声有：E[W(k)]＝0，E[W(k)W(j)]＝Q_kδ_kj    j＝1，2，3...(4)式中，Q_k为系统噪声的方差矩阵；δ：表示δ函数；步骤二：建立以SAR到地面标识点之间的距离和SAR与地面标识点之间的多普勒频移为观测量的观测方程；当SAR的成像视场中出现地面标识点时，假设L_i为第i个地面标识点的经度，λ_i为第i个地面标识点的纬度；当i大于等于3时，就可以确定卫星的轨道；所述第i个地面标识点在地固坐标系下的坐标为：xdi=RecosλicosLiydi=RecosλisinLizdi=Resinλi---(5)]]>式中，R_e为地球半径；根据地固坐标系和惯性坐标系的转换关系，得到第i个地面标识点在惯性坐标系下的坐标，即：xmiymizmi=cosΩ-sinΩ0sinΩcosΩ0001xdiydizdi---(6)]]>第i个地面标识点在惯性坐标系下的速度，即：vmxivmyivmzi=0-ωe0ωe00000xmiymizmi---(7)]]>式中：ω_e为地球自转角速度；Ω＝Ω₀+ω_e(t-t_o)为t时刻的恒星时；Ω_o为t₀时刻的恒星时；根据上述式(6)和(7)，得到卫星到第i个地面标识点的距离方程，即：ρi=(x-xmi)2+(y-ymi)2+(z-zmi)2+V1,i=1,2,3---(8)]]>利用上述式(8)得到多普勒频移方程，即：fi=-2λρi[(vx-vmxi)(x-xmi)+(vy-vmyi)(y-ymi)+(vz-vmzi)(z-zmi)]+V2---(9)]]>式中，i＝1，2，3利用上述式(8)和(9)得到观测方程，即：G=ρf=h[X(t),t]+V(t)---(10)]]>式(10)中的ρ＝[ρ₁，ρ₂，ρ₃]^T，f＝[f₁，f₂，f₃]^T，V＝[V₁，V₂]^T式中：V为测量噪声；G为观测量；步骤三：建立扩展卡尔曼滤波的递推方程；设定k时刻的滤波值为将步骤一中的式(3)中的非线性向量函数F围绕滤波值展开成泰勒级数，并略去二次以上项，得到线性化方程：X(k+1)=f[X^(k|k),k]+∂F[X(k|k),k]∂X(k)|X(t)=X^(k|k)[X(k)-X^(k|k)]+Γ(X(k),k)W(k)---(11)]]>式中，Γ为系数矩阵；将步骤二中的式(10)中的非线性向量函数h离散化后，并围绕滤波值展开成泰勒级数，并略去二次以上项，得到离散化的观测方程，即：G(k+1)=h[X^(k+1|k),k+1]+∂h[X^(k+1),k+1]∂X(k+1)|X(t)=X^(k+1|k)[X(k+1)-X^(k+1|k)]+V(k+1)---(12)]]>将上述式(11)和(12)代入标准的卡尔曼滤波方程中，得到扩展卡尔曼滤波的递推方程，即：X^(k+1|k+1)=F[X^(k|k),k]+K(k+1){G(k+1)-h[X^(k+1|k),k+1]}---(13)]]>其中：K(k+1)＝P(k+1|k)H^T(k+1)[H(k+1)P(k+1|k)H^T(k+1)+R(k+1)]^-1    (14)P(k+1|k)＝Φ(k+1，k)P(k|k)Φ^T(k+1，k)+Q(k)                 (15)P(k+1|k+1)＝[I-K(k+1)H(k+1)]P(k+1|k)                       (16)滤波初始值为：X^(0|0)=E[X(0)]]]>P(0|0)＝P_x(0)上述式(15)中一步状态转移矩阵Φ(k+1,k)=∂F[X(k),k∂X(k)|X(k)=X^(k|k)]]>上述式(16)中的观测矩阵H(k+1)=∂h[X(k+1,k+1)∂X(k+1)|X(k)=X^(k+1|k)---(17)]]>式中，K为增益矩阵；P(k|k)为k时刻的协方差矩阵；P(k+1|k)为协方差矩阵的一步预测；Q为系统噪声的协方差矩阵；R为测量噪声的协方差矩阵；I为单位矩阵；步骤四：得到卫星的状态信息；利用步骤三建立的扩展卡尔曼滤波的递推方程得到卫星的状态信息，所述卫星的状态信息即为卫星的状态估计所述卫星的状态估计包括卫星位置和速度矢量