一种任意基准系下的整星偏置角动量的偏航估计方法

摘要

一种任意基准系下的整星偏置角动量的偏航估计方法，在滚动角或俯仰角测量由敏感器给出，偏航角没有敏感器测量的情况下，仅通过偏航估计的方法来估算偏航角。在此基础上，整星通过角动量喷气管理的方法来实现基于任意基准系下的卫星三轴姿态稳定控制。本发明方法可根据用户需要或用滚动角测量值或用俯仰角测量作为输入值，具有很强的适应性，即仅通过调整初始的参数设计值，就可以获得需要的偏航估计结果，便于偏航估计的通用规范化方法编制。

说明书

技术领域

本发明属于卫星控制技术领域，涉及一种卫星姿态角的预估方法。

背景技术

卫星姿态预估是卫星控制中的常见问题。

现有的卫星偏航估计方法多是基于轨道坐标系的，主要是因为在轨卫星多以轨道坐标系为基准。但是在研型号和背景预研阶段的卫星，则根据有效载荷的任务需求不同，作为本体坐标系参考的基准坐标系也会有所不同。例如，卫星为了适应天线对国土形状的覆盖，要求卫星以东南坐标系为参考基准坐标系，由于东南坐标系与轨道坐标系仅相差一个周期性变化的偏航角，现有的操作方法是在基于轨道坐标系的偏航估计方法中再增加一部分标称偏航角的修正量来实现。然而，未来卫星型号的参考基准坐标系相对轨道坐标系将变得更为复杂，或是要求在不同任务阶段采用不同的参考坐标系，那么现有的基于轨道坐标系的偏航估计的修补方法将不再适用。

发明内容

本发明解决的技术问题是：克服现有技术的不足，提供了一种通用性的、可以在任意基准坐标系下预估卫星偏航角的方法。

本发明的技术解决方案是：一种任意基准系下的整星偏置角动量的偏航估计方法，包括如下步骤：

(1)在卫星的滚动轴、俯仰轴和/或偏航轴上分别安装动量轮，同时在卫星的滚动轴和俯仰轴上安装姿态敏感器；

(2)利用所述姿态敏感器，获取卫星的滚动角测量值和俯仰角测量值；

(3)利用滚动角测量值构建如下的偏航估计方程，

观测系数l₁,l₂满足特征方程的根为负数；

或者利用俯仰角测量值构建如下的偏航估计方程，

观测系数l₁,l₂满足特征方程的根为负数；

其中e₁和e₂为偏航角和滚动角的误差状态变量，ω_rx,ω_ry,ω_rz为基准系相对惯性系的三个方向的角速度在基准系下的描述；h_x,h_y,h_z分别为由动量轮组在卫星三轴方向产生的角动量，Δh_xΔh_y和Δh_z分别为三轴角动量偏差值，θ_b分别为滚动角和俯仰角的偏置值，T_dx,T_dy,T_dz为作用在星体三个方向上的干扰力矩分量，T_{dx_com},T_{dz_com}为滚动方向和偏航方向干扰力矩补偿值，分别为俯仰角估计值、滚动角估计值和偏航角估计值；

(4)采用步骤(3)中建立的偏航估计方程进行偏航角的估计。

所述步骤(4)中，为了进一步提高偏航角估计的精度，在利用建立的偏航估计方程进行偏航角的估计时，还对卫星上滚动轴、俯仰轴和偏航轴上安装的动量轮分别进行角动量卸载，其中：

利用滚动角测量值构建估计方法的卸载方式如下：

俯仰方向卸载槛值δ_y满足-δ_y≤h_synthesis_b(2)-h_y≤δ_y；

滚动方向卸载槛值为δ_x＝h_synthesis_i(1)-h0_unload_i(1)；

偏航方向卸载槛值为δ_z＝h_synthesis_i(3)-h0_unload_i(3)；

利用俯仰角测量值构建估计方法的卸载方式如下：

滚动方向卸载槛值δ_x满足-δ_x≤h_synthesis_b(1)-h_x≤δ_x；

俯仰方向卸载槛值为δ_y＝h_synthesis_i(2)-h0_unload_i(2)

偏航方向卸载槛值为δ_z＝h_synthesis_i(3)-h0_unload_i(3)

其中，为星本体角动量设计值在惯性系下的描述，为本体角动量在惯性系下的描述，h_synthesis_i＝C_ib*h_synthesis_b，C_ib为本体坐标系相对于惯性坐标系的姿态矩阵。

本发明与现有技术相比的优点在于：在卫星偏航方向缺少测量的情况下，现有技术是在轨道系下基于整星角动量垂直于轨道面的前提来设计偏航估计器的，本发明方法则扩展了上述前提，一方面整星角动量可以不垂直于轨道面，另一方面，在任意坐标系下都可以进行偏航估计器的设计。另外，本发明还在降低角动量偏差累积效应对偏航估计精度的影响方面，提出了角动量卸载管理方法。现有技术是将星体角动量在轨道系下进行管理，而本发明在设计模型构建时基于的是任意参考坐标系，给出的是惯性系相对任意参考系的相对角速度描述，在设计偏航估计方程中考虑了设计模型中由相对角速度带来的偏差部分，使得估计误差在没有外力矩干扰和残留角动量的影响下经过过渡过程时间后接近零值，即偏航估计值接近真值，因此本发明方法具有更强的通用性。

附图说明

图1为本发明方法的流程框图；

图2为本发明实施例中卫星相对于东南系的偏航估计姿态角；

图3为本发明实施例中卫星相对于东南系的滚动和俯仰测量角；

图4为本发明实施例中卫星相对于轨道系的三轴姿态角；

图5为本发明实施例中惯性系下的星本体角动量设计结果；

图6为本发明实施例中惯性系下的星本体角动量实测结果。

具体实施方式

为了适应未来更广阔更复杂的卫星任务需要，以及从算法的规范化、软件的维护性需求出发，本发明提出了一种适用于任意参考基准坐标系的偏航估计方法，使用者仅需给定当前本体坐标系需要参考的基准坐标系即可。

偏航估计方法主要是基于卫星长周期的运动控制模型，当有三个主要条件成立时才能达到估计效果：其一，基准坐标系相对于惯性坐标系要具备一定的角速度；其二，在此相对角速度方向存在一个可设计的标称角动量；其三，星上喷气等设备要具备角动量管理能力。

本发明方法还能够根据输入的基准坐标系、标称角动量、动量轮配置情况、太阳光压等干扰力矩补偿残余量来判断是否具备偏航估计的条件以及在条件具备下的估计精度。

如图1所示，为本发明方法的流程框图，采用本发明方法，要求卫星上安装有测量滚动方向或者俯仰方向的地球敏感器或太阳敏感器或其它敏感器，动量轮组和喷气装置。本发明方法的主要步骤如下：

(1)整星的动力学方程可以描述为：

$>\stackrel{·}{H}+{\omega }^{\times }H=T_d---\left(1\right)$>

其中，H＝[H_x>y>z]^T为整星的三轴角动量(角标x,y,z分别表示俯仰轴、滚动轴和偏航轴)，即星体内动量轮组产生的角动量和星体自身产生的角动量之和，可以写成如下形式：

H＝[H_x>y>z]^T＝[h_x+I_xω_x>y+I_yω_y>z+I_zω_z]^T(2)

其中，h_x,h_y,h_z分别为由动量轮组在卫星三轴方向产生的角动量，可计算求得，I_xω_x,I_yω_x,I_zω_x为星体自身产生的角动量，ω^×表示向量ω＝[ω_x,ω_y,ω_z]^T的反对称矩阵，ω为星体相对惯性坐标系的角速度，可表示为：

$>\omega =A\stackrel{·}{\alpha }+C_{br}{\omega }_{ri}^r---\left(3\right)$>

其中，为卫星的滚动角、俯仰角和偏航角，为基准坐标系相对惯性坐标系的三个方向的角速度在基准坐标系下的描述。C_br为卫星本体坐标系相对于基准坐标系的姿态矩阵，对于3-1-2转序下C_br的数学描述如下：

矩阵A为欧拉角的导数与ω的关系，对于3-1-2转序下A的数学描述如下：

其中，

$>\hat{1}=\left(\begin{array}{c}1\\ 0\\ 0\end{array}\right),\hat{2}=\left(\begin{array}{c}0\\ 1\\ 0\end{array}\right),\hat{3}=\left(\begin{array}{c}0\\ 0\\ 1\end{array}\right)$>

当卫星本体坐标系与基准坐标系之间的夹角为小角度时，A为单位矩阵。式(3)可写成如下形式：

$>\omega =\left(\begin{array}{c}\stackrel{·}{\phi }\\ \stackrel{·}{\theta }\\ \stackrel{·}{\psi }\end{array}\right)+{\left({\omega }_{ri}^r\right)}^{\times }\left(\begin{array}{c}\phi \\ \theta \\ \psi \end{array}\right)+{\omega }_{ri}^r---\left(4\right)$>

为的反对称矩阵；

将式(2)、(3)、(4)代入式(1)中，且忽略二阶小量后，经整理得

$>h^{\times }\stackrel{·}{\alpha }+h^{\times }{\left({\omega }_{ri}^r\right)}^{\times }\alpha +h^{\times }{\omega }_{ri}^r=-T_d---\left(5\right)$>

其中，h＝[h_x>y>z]^T，T_d＝[T_dx,T_dy,T_dz]^T为外界干扰力矩，如太阳光压、气动力矩等，T_dx,T_dy,T_dz为作用在星体三轴方向上的干扰力矩分量。

本发明估计方法的输入可以选择采用滚动角测量，也可以选择采用俯仰角测量。

在不考虑角动量累积的情况下，若采用滚动角测量，则设置

h＝[0；h_y；0]^T(6a)

若采用俯仰角测量，则设置

h＝[h_x；0；0]^T(6b)

其中，h_x和h_y分别为动量轮组产生的俯仰方向角动量偏置设计值和滚动方向角动量偏置设计值。

但在实际工作过程中，受到外界干扰，星体和动量轮之间通过角动量交换，使得动量轮吸收了外界的干扰角动量，三个方向的角动量就会偏离设计值(6a)或(6b)，则实际三个方向的角动量可以写成如下形式：

对应于式(6a)，实际状况下的三轴角动量为

h_real＝[Δh_x>y+Δh_yΔh_z]^T(7a)

对应于式(6b)，实际状况下的三轴角动量为

h_real＝[h_x+Δh_xΔh_yΔh_z]^T(7b)

Δh_xΔh_y和Δh_z分别为卫星三轴角动量偏差值。

(2)若采用滚动角测量，采用理论设计值(6a)，则方程(5)可以简化为

其中，为滚动角测量值。方程(8)中仅考虑偏航和滚动的耦合设计方程部分，俯仰方向有独立的测量值和控制量，因此不予在方程中给出。

若采用(7a)，则实际对象表述为

进入步骤(3)；

若采用俯仰角测量，采用理论设计值(6b)，则待设计的对象为

$>y={\theta }_m=\left(\begin{array}{cc}0& 1\end{array}\right)\left(\begin{array}{c}\psi \\ \theta \end{array}\right)$>

其中，θ_m为俯仰角测量值。方程(10)中仅考虑偏航和俯仰的耦合设计方程部分，滚动方向有独立的测量值和控制量，因此不予在方程中给出。

若采用(7b)，则实际对象表述为

进入步骤(4)；

(3)根据式(8)设计的偏航估计方程为

其中，分别为滚动角估计值和偏航角估计值。

利用式(9)和(12)，则有

$>\begin{array}{c}\left(\begin{array}{cc}{h}_y& 0\\ 0& -h_y\end{array}\right)\left(\begin{array}{c}{\stackrel{·}{e}}_1\\ {\stackrel{·}{e}}_2\end{array}\right)+\left(\begin{array}{cc}0& -h_y{\omega }_{ry}+l_1\\ -h_y{\omega }_{ry}& l_2\end{array}\right)\left(\begin{array}{c}{e}_1\\ e_2\end{array}\right)\\ =\left(\begin{array}{c}-h_y{\omega }_{rx}(\theta -{\theta }_b)-(T_{dx}-T_{dx\_com})-{\mathrm{\Delta h}}_y{\omega }_{rz}+{\mathrm{\Delta h}}_z{\omega }_{ry}\\ -h_y{\omega }_{rz}(\theta -{\theta }_b)-(T_{dz}-T_{dz\_com})-{\mathrm{\Delta h}}_x{\omega }_{ry}+{\mathrm{\Delta h}}_y{\omega }_{rx}\end{array}\right)\end{array}---\left(13\right)$>

其中，

待设计的观测系数l₁,l₂满足以下特征方程根为负数：

$>s^2-\frac{l_2}{h_y}s+{\omega }_{ry}^2-\frac{l_1{\omega }_{ry}}{h_y}=0---\left(14\right)$>

θ_b为俯仰角偏置值，该偏置值为控制任务需要预期设置值，T_{dx_com},T_{dz_com}为滚动方向和偏航方向干扰力矩补偿值。偏航估计的稳态精度为

$>{\{\hat{\psi }-\psi \}}_{(t\to \infty )}=\frac{-h_y{\omega }_{rx}(\theta -{\theta }_b)-(T_{dx}-T_{dx\_com})-{\mathrm{\Delta h}}_y{\omega }_{rz}+{\mathrm{\Delta h}}_z{\omega }_{ry}}{{\omega }_{ry}^2-l_1{\omega }_{ry}/h_y}---\left(15\right)$>

进入步骤(5)；

(4)根据式(10)设计的偏航估计方程为

利用式(11)和(16)，则有

其中，为俯仰角估计值，

$>\left(\begin{array}{c}{e}_1\\ e_2\end{array}\right)=\left(\begin{array}{c}\psi -\hat{\psi }\\ \theta -\hat{\theta }\end{array}\right)$>

待设计的观测系数l₁,l₂满足以下特征方程根为负数：

$>s^2+\frac{l_2}{h_x}s+{\omega }_{rx}^2-\frac{l_1{\omega }_{rx}}{h_x}=0$>

为滚动角偏置值，该偏置值为控制任务需要预期设置值，T_{dy_com},T_{dz_com}为俯仰方向和偏航方向干扰力矩补偿值。偏航估计的稳态精度为

进入步骤(5)；

(5)角动量卸载。由式子(15)和式子(18)可知，当星体姿态接近目标姿态且干扰力矩被补偿后，影响偏航估计误差的因素主要为Δh_xΔh_yΔh_z，即

由式(15)得到的偏航误差为：

$>\frac{-{\mathrm{\Delta h}}_y{\omega }_{rz}+{\mathrm{\Delta h}}_z{\omega }_{ry}}{{\omega }_{ry}^2-l_1{\omega }_{ry}/h_y}---\left(19\right)$>

由式(18)得到的偏航误差为：

$>\frac{-{\mathrm{\Delta h}}_z{\omega }_{rx}+{\mathrm{\Delta h}}_x{\omega }_{rz}}{{\omega }_{rx}^2-l_1{\omega }_{rx}/h_x}---\left(20\right)$>

若能消除Δh_xΔh_yΔh_z，则偏航估计值将趋近于真值。因此，需要进行角动量卸载，使得Δh_xΔh_yΔh_z维持在给定的范围内，达到偏航估计精度要求。

动量轮转速合成的三轴角动量在卫星本体坐标系下的描述记为h_synthesis_b，则其在惯性坐标系下的描述为

h_synthesis_i＝C_ib*h_synthesis_b

其中，C_ib为本体坐标系相对于惯性坐标系的姿态矩阵，h_synthesis_i为卫星本体角动量在惯性坐标系下的描述。h_synthesis_b由动量轮转速合成方法如下：

h_synthesis_b＝K_mC_MW*[r_MW1,r_MW2,…,r_MWn]^T+Iω(21)

其中，C_MW为动量轮安装矩阵，r_MWi，i＝1,2,…,n为动量轮i的转速，K_m为转速与角动量常值转换系数，根据不同类型的动量轮而不同。I＝diag[I_x,I_y,I_z]为星体三轴转动惯量，ω＝[ω_x,ω_y,ω_z]^T为陀螺测量的星体三轴角速度。

若采用滚动角测量，角动量设计值在卫星本体坐标系下的描述为h＝[0；h_y；0]^T，则有

h0_unload_i＝C_ib*[0,h_synthesis_b(2),0]^T(22)

其中，h0_unload_i为星本体角动量设计值在惯性坐标系下的描述，由于干扰角动量对卫星影响分析一般是以一个轨道周期为基础的，因此动量轮卸载基线也应该按照一个轨道周期来计算。这里，-δ_y≤h_synthesis_b(2)-h_y≤δ_y，δ_y为本体俯仰轴方向的角动量卸载槛值。滚动轴和偏航轴的角动量卸载槛值分别为：

δ_x＝h_synthesis_i(1)-h0_unload_i(1)

δ_z＝h_synthesis_i(3)-h0_unload_i(3)

其中，h_synthesis_i(k),k＝1,2,3表示向量的第k个元素，以下类同。

若采用俯仰角测量，角动量设计值在本体系下描述为h＝[h_x；0；0]^T，则有

h0_unload_i＝C_ib*[h_synthesis_b(1),0,0]^T(23)

其中，-δ_x≤h_synthesis_b(1)-h_x≤δ_x，δ_x为本体滚动轴方向的角动量卸载槛值。俯仰轴和偏航轴的角动量卸载槛值分别为：

δ_y＝h_synthesis_i(2)-h0_unload_i(2)

δ_z＝h_synthesis_i(3)-h0_unload_i(3)

角动量卸载本质就是将由干扰力矩引起的偏离标称值的那部分角动量通过其他执行机构给消除掉，例如，喷气推力器、磁力矩器等。

实施例

此处选用卫星的基准坐标系为东南坐标系，偏航估计基于东南坐标系设计。

东南坐标系OXdYdZd定义：O为卫星质心，Zd轴指向地心，Xd轴指向卫星运行方向，且平行于赤道平面，Yn轴与Xn、Zn轴右手正交。

卫星运行在地球同步轨道上，轨道倾角为5.5度。

初始时间(北京时)：2015年8月6日16时40分33秒

卫星位置(km)：-16504.192625 -38628.441000 -3752.763227

卫星速度(km/s)：2.821086 -1.215998 0.110196

东南系相对于惯性系的角速度描述为

$>{\overrightarrow{\omega }}_{di}^d=\left(\begin{array}{ccc}{\mathrm{cos\psi }}_l& {\mathrm{sin\psi }}_l& 0\\ -{\mathrm{sin\psi }}_l& {\mathrm{cos\psi }}_l& 0\\ 0& 0& 1\end{array}\right){\left(\begin{array}{ccc}0& {\omega }_0& 0\end{array}\right)}^T=\left(\begin{array}{c}{\omega }_0{\mathrm{sin\psi }}_l\\ {\omega }_0{\mathrm{cos\psi }}_l\\ 0\end{array}\right)---\left(24\right)$>

其中，ω₀为轨道角速度，约为-0.00416度/秒；ψ_l为轨道系相对于东南系的偏航角，为日周期正余弦变化，最大幅值为5.5度。

采用滚动角测量输入时，其角动量设计值为

$>h=\left(\begin{array}{c}0\\ h_y\\ 0\end{array}\right)=\left(\begin{array}{c}0\\ -92\\ 0\end{array}\right)---\left(25\right)$>

将(24)和(25)代入下式中

$>h^{\times }\stackrel{·}{\alpha }+h^{\times }{\omega }^{\times }\alpha +h^{\times }\omega =-T_d---\left(26\right)$>

经整理后得，

$>\begin{array}{c}-92\stackrel{·}{\psi }-0.00667{\mathrm{cos\psi }}_l\phi =-0.00667{\mathrm{sin\psi }}_l\theta -T_{dx}\\ 92\stackrel{·}{\phi }-0.00667{\mathrm{cos\psi }}_l\psi =0.00667{\mathrm{sin\psi }}_l-T_{dz}\end{array}---\left(27\right)$>

偏航估计器为

T_{dx_com}和T_{dz_com}为滚动轴和偏航轴方向的太阳光压力矩补偿分量。

仿真时间长度为10天，仿真结果见图2～图6。

图2为卫星相对于东南系的偏航估计姿态角。图中，横坐标为时间，纵坐标为偏航估计姿态角，该估计值为本方法估计所得。图3为卫星相对于东南系的滚动和俯仰测量角。图中，横坐标为时间，纵坐标为分别为滚动角测量值和俯仰角测量值。图4为卫星相对于轨道系的三轴姿态角，图中，横坐标为时间，纵坐标分别为滚动角、俯仰角和偏航角，这里给出该图的目的是为了能够清楚的表明偏航估计效果，即偏航估计在基于东南系下设计，在东南系的偏航角近零，而东南系相对于轨道系主要差别在于偏航角的日周期性变化，其幅值等于轨道倾角，从图中可以看出在轨道系下，偏航角呈现周期性变化，幅值为5.5。图5为在惯性系下的星本体角动量设计结果。图中，横坐标为时间，纵坐标分别为惯性系下的星本体滚动轴角动量、俯仰轴角动量和偏航轴角动量，星体俯仰方向角动量受到外界干扰力矩的作用呈现周期性变化，那么该角动量在惯性系下的描述为图5所示。图6在惯性系下的星本体角动量实测结果，图中，横坐标为时间，纵坐标分别为惯性系下星体实际角动量结果。从图中可以看出，当超过给定界限7Nms后，通过喷气方式对星体角动量进行卸载。

本发明说明书中未作详细描述的内容属本领域技术人员的公知技术。