[发明专利]基于张量积多胞鲁棒H2滤波的无陀螺卫星姿态确定方法

摘要

本发明涉及一种基于张量积多胞鲁棒H2滤波的无陀螺卫星姿态确定方法，属于飞行器技术领域。本发明针对卫星姿态动力学方程和运动学方程的非线性特性，提出了基于张量积转换的多胞鲁棒H2滤波，将非线性滤波问题转化为线性滤波问题；首先建立姿态确定系统的状态方程和星敏感器测量方程，并利用雅可比线性化将非线性系统化为线性变参数误差系统，然后根据张量积模型转换建立LPV系统多胞模型描述，并结合鲁棒H2滤波获取姿态确定系统的状态估计校正量，最后利用估计校正量对EKF方法获取的姿态一步预测量进行校正，获取姿态估计值；避免了EKF方法中实时计算更新滤波增益，大大减小了滤波计算量。

主权项

1.基于张量积多胞鲁棒H2滤波的无陀螺卫星姿态确定方法，其特征在于：包括以下步骤：步骤1，建立卫星姿态的状态模型和星敏感器测量模型；选取惯性系作为参考坐标系，确定卫星本体系相对于惯性系的方位；建立卫星姿态动力学方程：Jω·+ω×Jω=T]]>其中，ω为在卫星本体系中相对于惯性系的角速度；J为卫星惯性阵，T为施加到卫星上的转动力矩，包含控制力矩和作用在卫星上的外部各种干扰力矩；采用修正罗德里格参数σ作为姿态描述参数,卫星姿态运动学方程为：σ·=M(σ)ω---(2)]]>式中，M(σ)=14[(1-σTσ)I3+2[σ×]+2σσT];]]>采用星敏感器测量得到的测量模型为：S_b=R(σ)S_i+ΔS其中，S_i为恒星在惯性系中的单位方向矢量，ΔS为星敏感器测量噪声，R(σ)为惯性系到卫星本体系的姿态转换矩阵；步骤2，采用雅可比线性化将步骤1建立的非线性系统变换为LPV系统，得到卫星姿态滤波状态误差模型和测量误差模型；设系统的状态估计x^=σ^Tω^TT,]]>其中和满足：ω^·=-J-1ω^×Jω^+J-1T^---(3)]]>σ^·=M(σ^)ω^---(4)]]>结合雅可比线性化，得到：ω·-ω^·=Δω·=FωΔω+nω]]>σ·-σ^·=Δσ·=-[ω^×]Δσ+14Δω]]>式中：Fω=J-1[(Jω^)×]-J-1[ω^×]J,]]>n_ω=J^-1ΔT，ΔT为力矩误差；选取误差状态变量Δx=Δσ^TΔω^TT,]]>则系统的滤波状态误差模型为Δx·=FΔx+Gw---(5)]]>式中：F=-[ω^×]14I3×303×3Fω,]]>G=03×3I3×3,]]>w＝n_ω星敏感器测量残差为ΔSb=Sb-S^b]]>=R(σ)Si-R(σ^)Si+ΔS]]>=[I-R(Δσ)]R(σ^)Si+ΔS]]>Δσ为小量时，忽略二阶量，R(Δσ)≈I-4[Δσ×]，则星敏感器的测量残差方程为：ΔSb=[4(R(σ^)Si)×]Δσ+ΔS]]>采用m个星敏感器对卫星姿态进行测量，得到系统测量模型为：z＝[S_b1，S_b2，...,S_bm]^T，系统的测量误差模型为Δz=H(σ^)Δx+v---(6)]]>式中，H(σ^)=4(R(σ^)Si1)×03×34(R(σ^)Si2)×03×3······4(R(σ^)Sim)×03×3]]>v=ΔS1ΔS2···ΔSm,]]>S_bm为第m个星敏感器的测量模型，S_im为第m个星敏感器在惯性系中的单位方向矢量，ΔS_m为第m个星敏感器的测量噪声；步骤3，根据步骤1中的卫星姿态动力学、运动学方程，得到变参数的上下界，代入仿射参数依赖的矩阵F，得到其多胞顶点F₁,F₂；对于非仿射参数依赖的矩阵通过张量积模型转换方法获得LPV系统的测量误差模型的多胞描述；系统滤波状态和测量误差方程的多胞描述形式为：Δx·=FΔx+Gw]]>(8)Δz=HΔx+v其中，(F,H)={Σi=12λ1iFi,Σi=1sλ2iHi|Σi=12λ1i=1,Σi=1sλ2i=1},]]>(F_i,H_i)为多胞形系统的顶点；s为多胞顶点数目，为1模式、2模式、3模式下分别保留的奇异值个数相乘；步骤4，结合鲁棒H2滤波获取姿态确定系统的状态估计校正量；将多胞描述形式离散化，得到系统滤波状态方程和测量误差方程的离散多胞描述模型：Δx_k+1=AΔx_k+Bw_kΔz_k=CΔx_k+Dw_k其中，A＝I+FT_s，B、C、D为多胞描述形式中相应参数的离散值，且(A,C)={Σi=12λ1iAi,Σi=1sλ2iCi|Σi=12λ1i=1,Σi=1sλ2i=1},]]>T_s为采样周期；Δx_k表示k时刻的状态误差，w_k为离散后的系统噪声，Δz_k表示k时刻的测量误差；根据姿态确定滤波误差系统的离散多胞描述模型，并利用基于LMI技术的鲁棒H2滤波原理，得到姿态估计校正量的计算公式为：x^mk+1=AFx^mk+BF(zk-z^k|k-1)]]>Δx^k=CFx^mk+DF(zk-z^k|k-1)]]>其中，是中间变量，是k-1时刻测量模型到k时刻的一步测量预测，是k时刻的状态估计值；滤波系数A_F、B_F、C_F、D_F计算公式为：AF=G2-1SA,BF=G2-1SB,]]>C_F=S_C、D_F=S_D其中，S_A、S_B、S_C、S_D、G₂通过求解下述优化问题得到：min trace(Z)s.t.G11+G11T-P11jG2+G21T-P12jψ1j*G2+G21T-P22jψ2j**ψ3j******]]>SA-F21TG11Bj+SBDjSA-λ2G2TG21Bj+SBDjψ4j-F11Bj-λ1SBDjP22j-λ2SA-λ2SAT-F21Bj-λ2SBDj*I>0,]]>ZI-SDCj-SC-SDDj*P11jP12j0**P22j0***I>0,j=1,2,···,2s]]>其中，λ₁,λ₂，G₁₁,G₂₁,G₂,F₁₁,F₂₁,S_A,S_B,S_C,S_D,P_11j,P_12j,P_22j为变量；ψ1j=G11Aj+SBCj-F11T,]]>ψ_2j=G₂₁A_j+S_BC_j-λ₁G₂^Tψ3j=P11j-F11Aj-λ1SBCj-AjTF11T-λ1CjTSBT]]>ψ4j=P12j-λ1SA-AjTF21T-λ2CjTSBT]]>步骤5，对EKF算法获取的系统姿态进行一步预测，得到：σ^k|k-1=σ^k-1+M(σ^k-1)ω^k-1Ts]]>ω^k|k-1=ω^k-1+(-J-1([ω^k-1×]Jω^k-1)+J-1T^)Ts]]>式中，为k-1时刻角速度估计值和MRP估计值；分别为k时刻角速度预测值和MRP预测值；步骤6，利用步骤4得到的估计校正量对步骤5得到的预测值进行校正，获得k时刻角速度和MRP状态估计为：σ^k=Δσ^k⊗σ^k|k-1]]>ω^k=ω^k|k-1+Δω^k]]>步骤7，令k=k+1，将步骤6获得的k时刻状态估计值代入步骤5，得到k+1时刻一步预测；将该预测值代入步骤4，得到k+1时刻的姿态估计校正量；再进行步骤6，得到k+1时刻的角速度和MRP状态估计；在无陀螺卫星运行过程中，按照上述方法，实现卫星姿态的实时获取。