一种基于最大角加速度和最大角速度卫星规划轨迹方法

摘要

一种基于最大角加速度和最大角速度卫星规划轨迹方法，属于卫星机动轨迹规划领域。现有的规划轨迹确定方法不能充分利用执行机构的机动能力，且不能保证机动时间最短的问题。一种基于最大角加速度和最大角速度卫星规划轨迹方法，设定与目标姿态对应的目标坐标系，计算卫星由初始姿态机动至目标姿态的欧拉轴e

说明书

技术领域

本发明涉及一种基于最大角加速度和最大角速度卫星规划轨迹方法。

背景技术

对地观测卫星特别是具有快速响应能力的成像卫星，在轨运行时需要迅速获取地面不 同目标信息或者地面同一目标连续/间断信息，因此，这类卫星都必须具备快速大角度姿 态机动的能力。

卫星实现大角度姿态机动的前提是机动轨迹规划。在以往的研究中，研究人员往往人 为设定规划轨迹的最大角加速度和最大角速度，通过试凑的方法，使设计出的规划轨迹在 执行机构的能力范围之内，以实现大角度机动。这类方法给出的规划轨迹往往不能充分利 用执行机构的机动能力，而且不能保证机动时间最短。因此，有必要给出规划轨迹的最大 角加速度和最大角速度的设计方法，使执行机构的能力得到充分利用，从而最大程度地缩 短机动时间。

发明内容

本发明的目的是为了解决现有的规划轨迹确定方法不能充分利用执行机构的机动能 力，而且不能保证机动时间最短的问题，而提出一种基于最大角加速度和最大角速度卫星 规划轨迹方法。

一种基于最大角加速度和最大角速度卫星规划轨迹方法，所述方法通过以下步骤实 现：

步骤一、设定与目标姿态对应的目标坐标系o_bx_ty_tz_t，根据：

$e_m=\frac{1}{sin\frac{{\Phi }_m}{2}}{q}_m---\left(1\right),$

Φ_m＝2arccosq_m0                                    (2)，

计算卫星由初始姿态机动至目标姿态的欧拉轴e_m和转角Φ_m；式中，q_m0表示卫星由初始 姿态机动至目标姿态的机动四元数Q_m的标部，q_m表示卫星由初始姿态机动至目标姿态的 机动四元数Q_m的矢部，Q_m为四维矢量Q_m＝[q_m0 q_m]^T，且

式中，Q_i表示初始姿态四元数，Q_t为表示目标姿态四元数；

步骤二、获得规划轨迹的最大角加速度和最大角速度的约束方程：

${\stackrel{··}{\Phi }}_{max}{N}_{\stackrel{··}{\Phi }}+{\stackrel{·}{\Phi }}_{max}{N}_{\stackrel{·}{\Phi }}\approx T_{wmax}---\left(4\right),$

其中，$N_{\stackrel{··}{\Phi }}=-{\mathrm{Ie}}_m---\left(5\right),$

$N_{\stackrel{·}{\Phi }}={{\mathrm{Ie}}_m}^{\times }{\omega }_{oI}-{e_m}^{\times }{H}_0---\left(6\right),$

式中，表示规划轨迹的最大角加速度，表示规划轨迹的最大角速度，表 示角加速度力矩项的系数，表示角速度力矩项的系数，I是3×3的矩阵，表示卫星相 对其质心的转动惯量，e_m表示欧拉轴，ω_oI表示轨道角速度，H₀表示初始时刻卫星系统 的角动量，T_wmax是标量，为飞轮力矩包络半径；

步骤三、将步骤二中涉及的角速度力矩项的系数分解为和两部分，并获 得和的模值：

$\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|=\frac{N_{\stackrel{·}{\Phi }}·N_{\stackrel{··}{\Phi }}}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(7\right),$

$\left|N_{\stackrel{·}{\Phi }\perp }\right|=\sqrt{{\left|N_{\stackrel{·}{\Phi }}\right|}^2-{\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|}^2}---\left(8\right);$

其中，是与平行的角速度力矩项系数，是与垂直的角速度力矩项系数；

步骤四、由规划轨迹的最大角速度表示规划轨迹的最大角加速度

${\stackrel{··}{\Phi }}_{max}\approx \frac{T_{wmax}-\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|{\stackrel{·}{\Phi }}_{max}-\frac{{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2}{2T_{w\max }}{{\stackrel{·}{\Phi }}_{\max }}^2}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(9\right);$

步骤五、获得受飞轮最大角动量限制的最大角速度的表达式：

${\stackrel{·}{\Phi }}_{\max \_H}\approx \frac{H_{wmax}}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(10\right),$

式中，H_wmax表示飞轮角动量包络半径，下角标H表示最大角速度受飞轮角动量上限 约束；

步骤六、计算最短机动时间t_m：

$t_m=\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{max}}+\frac{{\stackrel{·}{\Phi }}_{max}}{{\stackrel{··}{\Phi }}_{max}}---\left(11\right),$

并联合式(9)得到以规划轨迹的最大角速度为自变量的机动时间t_m的函数：

$t_m\left({\stackrel{·}{\Phi }}_{\max }\right)=\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{\max }}+\frac{2T_{w\max }\left|N_{\stackrel{··}{\Phi }}\right|{\stackrel{·}{\Phi }}_{\max }}{2{T_{w\max }}^2-2T_{w\max }\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|{\stackrel{·}{\Phi }}_{\max }-{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2{{\stackrel{·}{\Phi }}_{\max }}^2}---\left(12\right),$

结合如下约束：$\left(\begin{array}{c}{\stackrel{·}{\Phi }}_{\max }\le {\stackrel{·}{\Phi }}_{\max \_H}\\ \frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{\max }}\ge \frac{2T_{w\max }{N}_{\stackrel{··}{\Phi }}{\stackrel{·}{\Phi }}_{\max }}{2{T_{w\max }}^2-2T_{w\max }{N}_{\stackrel{·}{\Phi }\left|\right|}{\stackrel{·}{\Phi }}_{\max }-{N_{\stackrel{·}{\Phi }\perp }}^2{{\stackrel{·}{\Phi }}_{\max }}^2}---\left(13\right),\end{array}\right)$

求出使机动时间t_m取最小值时规划轨迹的最大角速度

步骤七、将步骤六获得的规划轨迹的最大角速度代入式(9)，求出规划轨迹的 最大角加速度

本发明的有益效果为：

本发明实现了一种设计卫星规划轨迹的最大角加速度和最大角速度的方法，通过限定 项的取值范围，使得飞轮的力矩得以充分利用，通过获得的规划轨迹的最大角速 度保证所设计的规划轨迹充分利用飞轮的能力，以使得机动时间最短。

本发明方法在分析飞轮力矩空间的利用过程中设计各明确物理含义的物理量，易于理 解；使得本发明方法简单易行，简化计算过程，使得卫星在轨运行时快速获取地面目标的 信息和相应。

附图说明

图1为本发明的流程图；

图2为具体实施方式步骤三、四所述在飞轮力矩包络内的力矩示意图；

图3为具体实施方式步骤六所述的规划轨迹的角加速度曲线和角速度曲线；且图中上 半部分图示为规划轨迹的角加速度曲线示意图，下半部分图示为规划轨迹的角加速度曲线 示意图；

图4为算例的规划轨迹角加速度曲线；

图5为算例的规划轨迹角速度曲线；

图6为算例的卫星姿态角速度曲线；

图7为算例的卫星姿态角曲线；

图8为算例的飞轮力矩曲线；

图9为算例的飞轮力矩模值曲线；

图10为算例的飞轮角动量曲线；

图11为算例的飞轮角动量模值曲线；

具体实施方式

具体实施方式一：

本实施方式的一种基于最大角加速度和最大角速度卫星规划轨迹方法，结合图1至图 3所示内容，所述方法通过以下步骤实现：

步骤一、设定与目标姿态对应的目标坐标系o_bx_ty_tz_t，根据：

$e_m=\frac{1}{sin\frac{{\Phi }_m}{2}}{q}_m---\left(1\right),$

Φ_m＝2arccosq_m0                              (2)，

计算卫星由初始姿态机动至目标姿态的欧拉轴e_m和转角Φ_m；式中，q_m0表示卫星由初始 姿态机动至目标姿态的机动四元数Q_m的标部，q_m表示卫星由初始姿态机动至目标姿态的 机动四元数Q_m的矢部，Q_m为四维矢量Q_m＝[q_m0 q_m]^T，且

式中，Q_i表示初始姿态四元数，Q_t为表示目标姿态四元数；

步骤二、获得规划轨迹的最大角加速度和最大角速度的约束方程：

${\stackrel{··}{\Phi }}_{\max }{N}_{\stackrel{··}{\Phi }}+{\stackrel{·}{\Phi }}_{\max }{N}_{\stackrel{·}{\Phi }}\approx T_{w\max }---\left(4\right),$

其中，$N_{\stackrel{··}{\Phi }}=-{\mathrm{Ie}}_m---\left(5\right),$

$N_{\stackrel{·}{\Phi }}={{\mathrm{Ie}}_m}^{\times }{\omega }_{oI}-{e_m}^{\times }{H}_0---\left(6\right),$

式中，表示规划轨迹的最大角加速度，表示规划轨迹的最大角速度，表 示角加速度力矩项的系数，表示角速度力矩项的系数，I是3×3的矩阵，表示卫星相 对其质心的转动惯量，e_m表示欧拉轴，ω_oI表示轨道角速度，H₀表示初始时刻卫星系统 的角动量，T_wmax是标量，为飞轮力矩包络半径；

步骤三、将步骤二中涉及的角速度力矩项的系数分解为和两部分，并获 得和的模值：

$\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|=\frac{N_{\stackrel{·}{\Phi }}·N_{\stackrel{··}{\Phi }}}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(7\right),$

$\left|N_{\stackrel{·}{\Phi }\perp }\right|=\sqrt{{\left|N_{\stackrel{·}{\Phi }}\right|}^2-{\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|}^2}---\left(8\right);$

其中，是与平行的角速度力矩项系数，是与垂直的角速度力矩项系数；

步骤四、由规划轨迹的最大角速度表示规划轨迹的最大角加速度

${\stackrel{··}{\Phi }}_{max}\approx \frac{T_{wmax}-\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|{\stackrel{·}{\Phi }}_{max}-\frac{{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2}{2T_{w\max }}{{\stackrel{·}{\Phi }}_{\max }}^2}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(9\right);$

步骤五、获得受飞轮最大角动量限制的最大角速度的表达式：

${\stackrel{·}{\Phi }}_{\max \_H}\approx \frac{H_{wmax}}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(10\right),$

式中，H_wmax表示飞轮角动量包络半径，下角标H表示最大角速度受飞轮角动量上限 约束；

步骤六、计算最短机动时间t_m：

$t_m=\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{max}}+\frac{{\stackrel{·}{\Phi }}_{max}}{{\stackrel{··}{\Phi }}_{max}}---\left(11\right),$

并联合式(9)得到以规划轨迹的最大角速度为自变量的机动时间t_m的函数：

$t_m\left({\stackrel{·}{\Phi }}_{max}\right)=\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{max}}+\frac{2T_{wmax}\left|N_{\stackrel{··}{\Phi }}\right|{\stackrel{·}{\Phi }}_{max}}{2{T_{w\max }}^2-2T_{wmax}\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|{\stackrel{·}{\Phi }}_{max}-{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2{{\stackrel{·}{\Phi }}_{\max }}^2}---\left(12\right),$

结合如下约束：$\left(\begin{array}{c}{\stackrel{·}{\Phi }}_{\max }\le {\stackrel{·}{\Phi }}_{\max \_H}\\ \frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{\max }}\ge \frac{2T_{w\max }{N}_{\stackrel{··}{\Phi }}{\stackrel{·}{\Phi }}_{\max }}{2{T_{w\max }}^2-2T_{w\max }{N}_{\stackrel{·}{\Phi }\left|\right|}{\stackrel{·}{\Phi }}_{\max }-{N_{\stackrel{·}{\Phi }\perp }}^2{{\stackrel{·}{\Phi }}_{\max }}^2}---\left(13\right),\end{array}\right)$

求出使机动时间t_m取最小值时规划轨迹的最大角速度

步骤七、将步骤六获得的规划轨迹的最大角速度代入式(9)，求出规划轨迹的 最大角加速度

具体实施方式二：

与具体实施方式一不同的是，本实施方式的一种基于最大角加速度和最大角速度卫星 规划轨迹方法，步骤二所述获得规划轨迹的最大角加速度和最大角速度的约束 方程过程为，忽略外界干扰力矩的影响，卫星的动力学方程为： 通常，卫星的一次机动时间远小于其轨道周期，那么， 在机动过程中，卫星本体系相对轨道系的角速度ω_bo的模值远大于轨道角速度ω_o的模值， 再结合角动量守恒定理：Iω_bI+H_w＝H₀，动力学方程简化为：

$I{\stackrel{·}{\omega }}_{bo}\approx T_w+{{\mathrm{I\omega }}_{bo}}^{\times }{\omega }_{oI}-{{\omega }_{bo}}^{\times }{H}_0---\left(15\right),$

假设卫星能够对规划轨迹进行很好地跟踪，则：

${\omega }_{bo}\approx \stackrel{·}{\Phi }\left(t\right)e_m---\left(16\right),$

${\stackrel{·}{\omega }}_{bo}\approx \stackrel{··}{\Phi }\left(t\right)e_m---\left(17\right),$

将式(16)，式(17)代入式(15)，得到直观体现角加速度和角速度与飞轮力矩T_w关系 的动力学方程：$\stackrel{··}{\Phi }(-{\mathrm{Ie}}_m)+\stackrel{·}{\Phi }({{\mathrm{Ie}}_m}^{\times }{\omega }_{oI}-{e_m}^{\times }{H}_0)\approx T_w---\left(18\right),$令最大角加速度和最大角速度的取值使得式(18)右端的飞轮力矩T_w的模值达到上限， 以充分利用飞轮的力矩空间，即：${\stackrel{··}{\Phi }}_{max}(-{\mathrm{Ie}}_m)+{\stackrel{·}{\Phi }}_{max}({{\mathrm{Ie}}_m}^{\times }{\omega }_{oI}-{e_m}^{\times }{H}_0)\approx T_{wmax}---\left(19\right),$

(-Ie_m)和(Ie_m^×ω_oI-e_m^×H₀)为已知量，令(-Ie_m)由表示，(Ie_m^×ω_oI-e_m^×H₀)由表 示，将式(19)简化为：${\stackrel{··}{\Phi }}_{max}{N}_{\stackrel{··}{\Phi }}+{\stackrel{·}{\Phi }}_{max}{N}_{\stackrel{·}{\Phi }}\approx T_{wmax}---\left(4\right).$

具体实施方式三：

与具体实施方式一或二不同的是，本实施方式的一种基于最大角加速度和最大角速度 卫星规划轨迹方法，步骤三所述获得和的模值的过程为，如图2所示，由于由平行于的和垂直于的两部分组成，从而有：

$\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|=\frac{N_{\stackrel{·}{\Phi }}·N_{\stackrel{··}{\Phi }}}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(7\right)$

$\left|N_{\stackrel{·}{\Phi }\perp }\right|=\sqrt{{\left|N_{\stackrel{·}{\Phi }}\right|}^2-{\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|}^2}---\left(8\right).$

具体实施方式四：

与具体实施方式三不同的是，本实施方式的一种基于最大角加速度和最大角速度卫星 规划轨迹方法，步骤四所述由规划轨迹的最大角速度表示规划轨迹的最大角加速度 的过程为，由图2和式(19)可知，最大角速度的存在使飞轮力矩T_w无法被充 分利用，导致规划轨迹的最大角加速度的可用力矩包络从T_wmax缩减为T_wmax-△T，则 减少量△T表示为：$\Delta T={\stackrel{·}{\Phi }}_{max}\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|+T_{wmax}(1-cos\alpha )---\left(20\right),$

令sinα≈α，以充分利用飞轮力矩T_w的空间，加快机动速度，使角加速度力矩项 尽量大，则式(20)近似的表示为：

$\Delta T\approx {\stackrel{·}{\Phi }}_{max}\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|+\frac{{{\stackrel{·}{\Phi }}_{\max }}^2{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2}{2T_{wmax}}---\left(21\right),$

根据图2可知，规划轨迹的最大角加速度可表示为：

${\stackrel{··}{\Phi }}_{\max }=\frac{T_{w\max }-\Delta T}{\left|N_{\stackrel{··}{\Phi }}\right|}\approx \frac{T_{w\max }-\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|{\stackrel{·}{\Phi }}_{\max }-\frac{{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2}{2T_{w\max }}{{\stackrel{·}{\Phi }}_{\max }}^2}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(9\right).$

具体实施方式五：

与具体实施方式一、二或四不同的是，本实施方式的一种基于最大角加速度和最大角 速度卫星规划轨迹方法，步骤五所述获得受飞轮最大角动量限制的最大角速度的表 达式的过程为，

由于卫星在机动过程中，卫星本体角动量的变化量是由飞轮吸收获得，所以，规划轨 迹的最大角速度不能超过飞轮转速的上限；根据动量守恒定理：Iω_bI+H_w＝H₀，并 假设在机动过程中ω_bo的模值远大于ω_o的模值，得到飞轮角动量上限对最大角速度的限制 方程的表达式：${\stackrel{·}{\Phi }}_{\max \_H}\approx \frac{H_{wmax}}{\left|N_{\stackrel{··}{\Phi }}\right|}---\left(10\right).$

具体实施方式六：

与具体实施方式五不同的是，本实施方式的一种基于最大角加速度和最大角速度 卫星规划轨迹方法，步骤六所述求出使机动时间t_m取最小值的规划轨迹的最大角速度 的过程为，

第一，根据图3得机动时间t_m表示为：$t_m=\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{max}}+\frac{{\stackrel{·}{\Phi }}_{max}}{{\stackrel{··}{\Phi }}_{max}}---\left(11\right),$

将式(9)代入式(11)，得到以规划轨迹的最大角速度为自变量的机动时间t_m的 函数：

$t_m\left({\stackrel{·}{\Phi }}_{\max }\right)=\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{\max }}+\frac{2T_{w\max }\left|N_{\stackrel{··}{\Phi }}\right|{\stackrel{·}{\Phi }}_{\max }}{2{T_{w\max }}^2-2T_{w\max }\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|{\stackrel{·}{\Phi }}_{\max }{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2{{\stackrel{·}{\Phi }}_{\max }}^2}---\left(12\right),$

控制最大角速度满足

第二，由图3所示的角速度和角加速度的变化规律知，规划轨迹的最大角 速度和规划轨迹的最大角加速度的取值满足即： $\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{\max }}\ge \frac{2T_{w\max }\left|N_{\stackrel{··}{\Phi }}\right|{\stackrel{·}{\Phi }}_{\max }}{2{T_{w\max }}^2-2T_{w\max }\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|{\stackrel{·}{\Phi }}_{\max }-{\left|N_{\stackrel{·}{\Phi }\perp }\right|}^2{{\stackrel{·}{\Phi }}_{\max }}^2},$避免在跟踪过程中角速度未增加到最大值 前开始减速；

第三，利用Matlab工具箱中的非线性优化方法，得到带有约束的机动时间t_m的最小 值及其对应的规划轨迹的最大角速度

实施例：

结合图4～11说明本具体实施例。

进行卫星姿态轨迹跟踪的仿真，说明本发明设计的规划轨迹的最大加速度和最大速度 可以充分利用执行机构，使机动时间最短。卫星初始姿态为[0° 10° 10°]^T，目标姿态为 [0° 30° 40°]^T。飞轮采用四斜装构型，单个飞轮可产生的力矩和角动量分别为0.2Nm， 15Nm·s，则飞轮的力矩包络半径和角动量包络半径分别为0.267Nm，20Nm·s。卫星转 动惯量为I＝diag(2000 1000 3000)kg·m²，轨道角速度为ω_oI＝[0 -0.001 0]^Trad/s， 系统初始角动量为H₀＝[0 0 0]^TNm·s。姿态角速度测量误差为三轴10^-4°/s，姿态角测 量误差为三轴5×10^-4°。干扰力矩为$\left(\begin{array}{c}2\times {10}^{-4}sin\left({\omega }_{o}t\right)\\ 1.5\times {10}^{-4}\cos \left({\omega }_{o}t\right)\\ 2\times {10}^{-4}sin\left({\omega }_{o}t\right)\end{array}\right)+\left(\begin{array}{c}1\times {10}^{-4}\\ 1\times {10}^{-4}\\ 1\times {10}^{-4}\end{array}\right)Nm.$PD控制器参数为 K_p＝0.04I₃，K_d＝0.36I₃，I₃为3×3的单位矩阵。

根据步骤一，得到e_m＝[-0.287 0.544 0.789]^T，Φ_m＝0.627rad。

根据步骤二、三，得到$\left|N_{\stackrel{·}{\Phi }\left|\right|}\right|=0.357,\left|N_{\stackrel{·}{\Phi }\perp }\right|=\mathrm{1.563.}$

根据步骤五，得到${\stackrel{·}{\Phi }}_{\max \_H}=0.0105rad/s.$

根据步骤六，得到t_m＝134.702s，

根据步骤七，得到${\stackrel{··}{\Phi }}_{max}=1.384\times {10}^{-4}rad/s^2.$

由图5～7可以得出，卫星的姿态角速度和姿态角能够很好的跟踪规划轨迹，而且， 卫星在加速过程结束后，经过短暂的匀速过程，立马进入减速过程，显然能使机动时间尽 量缩短。由图8～9可以看出，飞轮力矩模值为0.237Nm，而力矩包络为0.267Nm，力矩 空间能够较为充分的利用。由图10～11可以看出，飞轮角动量模值为13.031Nm·s，而角 动量包络为20Nm·s，角动量空间的利用不如力矩空间利用得充分，原因在于，根据角动 量约束求得的${\stackrel{·}{\Phi }}_{\max \_H}=0.0105rad/s,$而根据$\frac{{\Phi }_m}{{\stackrel{·}{\Phi }}_{max}}\ge \frac{{\stackrel{·}{\Phi }}_{max}}{{\stackrel{··}{\Phi }}_{max}}$约束求得的${\stackrel{·}{\Phi }}_{max}=0.009rad/s,$的值减少了15％，所以使得角动量空间不能得到有效利用。