基于罗德里格斯参数的星敏感器几何线性姿态确定方法

摘要

一种基于罗德里格斯参数的星敏感器几何线性姿态确定方法，该方法基于姿态旋转几何分析，建立参考矢量、方位矢量、姿态参数之间的线性投影关系，通过优化的误差评估准则获得最小估计误差下的姿态估计结果，并进一步采用预旋转方法消除计算奇异，从而建立基于罗德里格斯参数的几何线性姿态确定方法。本发明可以很好的快速估计载体姿态，并具有形式简单、计算量小、便于递归扩展等优点，可以满足星敏感器姿态测量数据处理需要。

说明书

技术领域

本发明涉及的是一种卫星定位领域的技术，具体是一种基于罗德里格斯参数的星敏感器 几何线性姿态确定方法。

背景技术

面向轻质化、小型化敏捷卫星对地观测多种扫描模式，提出适应快速机动的几何线性姿 态确定方法，为星体姿态敏捷机动提供精确的测量数据，具有重要的应用价值。姿态确定是星 敏感器数据处理的重要环节，其输出结果直接利用在姿态控制、遥感应用计算中。目前国际上 普遍采用上世纪八十年代提出的四元数估计方法，需要求解迭代求解一元四次方程的根，计算 量较大。因此提出新的快速线性姿态确定方法来解决敏捷机动下的姿态测量是十分必要的。

经过对现有技术的检索发现，中国专利文献号CN104020671A公开(公告)日2014.09.03， 公开了一种量测干扰下用于飞行器姿态估计的鲁棒递推滤波方法。该技术包括：采集飞行器运 动过程中陀螺与星敏感器的输出数据；建立量测干扰下基于姿态四元数的飞行器非线性状态空 间模型；进行时间更新，求得一步状态预测和预测方差的上界；进行鲁棒递推滤波量测更新， 求得最优的滤波增益，进而求出k+1时刻的状态估计值和方差的上界，将k+1时刻的状态估计 中四元数部分进行强制的归一化约束；输出姿态四元数及陀螺漂移的结果，完成姿态估计。该 技术采用基于最小方差的鲁棒滤波设计，能够实现在最小方差意义下最优滤波增益设计，有利 用提高系统的姿态估计精度。但该技术无法解决的技术问题包括：1、该方法基于飞行器运动 状态方法来解决姿态估计问题，但并没有从考虑测量矢量和参考矢量之间的几何关系，以及利 用这种关系来建立新的方法以提高姿态估计精度；2、鲁棒递推滤波方法计算量大，且存在飞 行器运动建模不准确、传感器噪声特性难以获取等问题。

发明内容

本发明针对现有技术存在的上述不足，提出一种基于罗德里格斯参数的星敏感器几何线 性姿态确定方法，可以很好的快速估计载体姿态，并具有形式简单、计算量小、便于递归扩展 等优点，可以满足星敏感器姿态测量数据处理需要。

本发明是通过以下技术方案实现的：

本发明基于姿态旋转几何分析，建立参考矢量、方位矢量、姿态参数之间的线性投影关 系，通过优化的误差评估准则获得最小估计误差下的姿态估计结果，并进一步采用预旋转方法 消除计算奇异，获得非奇异的辅助姿态并最终得到星敏感器的最优估计姿态。

所述的方位矢量通过采集星敏感器的图像处理信息获得。

所述的姿态旋转几何分析是指：根据欧拉旋转理论，姿态余弦矩阵R与旋转轴a、旋转 角θ之间满足：R(a,θ)＝cosθI₃+(1-cosθ)aa^T-sinθ[a^×]，其中：[a^×]为旋转轴a生成的 斜对称矩阵，$\left[a^{\times }\right]=\left(\begin{array}{ccc}0& -a_z& a_y\\ a_z& 0& -a_x\\ -a_y& a_x& 0\end{array}\right).$

所述的线性投影关系是指：在观测坐标系中，参考矢量和测量矢量在旋转轴上的投影相 等、姿态罗德里格斯参数g与y轴垂直以及x轴与姿态罗德里格斯参数g的叉乘等于y轴。

所述的观测坐标系为右手正交坐标系，具体为：对于惯性系下的参考矢量r以及在星敏 感器本体系中的测量矢量b，构建观测坐标系$\left(\begin{array}{c}x=(b+r)/2\\ y=(b-r)/2\\ z=x^{\times }y\end{array}\right),$坐标轴的长度为$\left(\begin{array}{c}x=\left|\right|x\left|\right|\\ y=\left|\right|y\left|\right|\\ z=\left|\right|z\left|\right|\end{array}\right),$对应 所述的线性投影关系具体是指：$\left(\begin{array}{c}{a}^{T}b=a^{T}R(a,\theta )r\\ =a^T[\cos \mathrm{\theta I}+(1-\cos \theta ){\mathrm{aa}}^T-\sin \theta [a^{\times }]r\\ =a^{T}r\end{array}\right),$即如图1所示。

所述的优化的误差评估准则是指：对于每一组测量{r_i,b_i}及其测量坐标系{x_i,y_i,z_i}， 对应的罗德里格斯参数的误差向量为由于测量误差的存在，使得误差向量e_i不 为零。以加权累计误差为性能函数，建立误差评估准则：其中：ξ_i为归一化加权值，满足

所述的姿态估计结果通过以下方式得到：将优化的误差评估准则整理成关于向量的二次 型函数J＝g^TMg-2g^Tz+c，其中：$\left(\begin{array}{c}M=-\underset{i=1}{\overset{n}{\Sigma }}{\xi }_i{\left[x_i^{\times }\right]}^2\\ \upsilon =\underset{i=1}{\overset{n}{\Sigma }}{\xi }_i\left[x_i^{\times }\right]y_i\end{array}\right);$通过求解优化问题，获得线性姿态 确定方程，即求解二次型优化问题min J＝g^TMg-2g^Tz+c，由于矩阵M是一系列具有对称 性的3×3矩阵之和，具有半正定性M≥0，因此姿态估计结果为

所述的预旋转方法消除计算奇异，其步骤包括：

i)改进的罗德里格斯参数是一种无奇异姿态描述方式，用来表示姿态估计值，则有

$\hat{\sigma }=\frac{1}{1+\sqrt{1+{\left|\right|\hat{g}\left|\right|}^2}}\hat{g}=\frac{1}{1+\sqrt{1+{\upsilon }^T{M}^{-1}{M}^{-1}\upsilon }}{M}^{-1}\upsilon$

为了避免矩阵求逆容易导致计算奇异，上式分子和分母同时乘以矩阵M的行列式 γ＝det M，则有其中：M^*为矩阵M的伴随矩阵。

ii)采用预旋转方法，维持γ一直在较高的值附近，即选择辅助旋转轴：

$a^{\oplus }=\left(\begin{array}{cc}\underset{i=1}{\overset{n}{\Sigma }}{\xi }_i\frac{x_i}{{\tilde{x}}_i},& x_1\ge y_1\\ \underset{i=1}{\overset{n}{\Sigma }}{\xi }_i\frac{y_i}{y_i},& \mathrm{otherwise}\end{array}\right),$并构建辅助参考向量r_i^⊕＝R(a^⊕,π)r_i＝(r_i^Ta^⊕)a^⊕-r_i， 并用r_i^⊕代替r_i，就可以计算非奇异的辅助姿态为

iii)根据辅助姿态和辅助旋转轴a^⊕，计算最终的最优估计姿态为

$\hat{\sigma }=\frac{1}{1-2{\hat{\sigma }}^{\oplus T}{a}^{\oplus }+{\hat{\sigma }}^{\oplus T}{\hat{\sigma }}^{\oplus }}[(1-{\hat{\sigma }}^{\oplus T}{\hat{\sigma }}^{\oplus })I-2[{\hat{\sigma }}^{\oplus \times }\left]\right]a^{\oplus }.$

附图说明

图1是本发明的姿态旋转几何关系图。

图2是本发明方法流程示意图。

图3是本发明实施例1的最优估计姿态示意图。

具体实施方式

下面对本发明的实施例作详细说明，本实施例在以本发明技术方案为前提下进行实施， 给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。

实施例1

如图2所示，本实施例包括以下步骤：

1.1)通过采集星敏感器的图像处理信息获得多颗高星等恒星方位矢量b_i，i＝1,2,...,N； 通过查询星历表获得这些恒星在J2000惯性坐标系的参考矢量r_i，i＝1,2,...,N。

本实施例以低轨卫星对地观测为例模拟产生1000个采样时刻的方位矢量和参考矢量， 每组方位矢量和参考矢量各有4个，也就是N＝4，采样步长为0.1s，卫星轨道高度为800km。

本实施例中参考矢量固定为r₁＝[1,0,0]^T、r₂＝[0,1,0]^T、r₃＝[0,0,1]^T和t_k时刻的方位矢量为b_k,1＝A_kr₁+δb_k,1、b_k,2＝A_kr₂+δb_k,2、b_k,3＝A_kr₃+δb_k,3和 b_k,4＝A_kr₄+δb_k,4，其中：A_k＝cos(nt_k)I₃+(1-cos(nt_k))aa^T-sin(nt_k)[a^×]为t_k时刻从参考坐 标系到卫星本体坐标系的方向余弦矩阵，n＝0.001s^-1为卫星轨道角速度，a可取为[0,-1,0]^T； δb_k,1、δb_k,2、δb_k,3和δb_k,4为星敏感器测量噪声，不同噪声之间具有统计无关性，且各分量方 差为σ²＝1×10^-4。

1.2)对于t_k时刻的每一组测量{r_i,b_k,i}，i＝1,2,3,4，计算坐标轴$\left(\begin{array}{c}{x}_{k,i}=(b_{k,i}+r_i)/2\\ y_{k,i}=(b_{k,i}-r_i)/2\end{array}\right),$以及$\left(\begin{array}{c}{M}_k=-\underset{i=1}{\overset{n}{\Sigma }}{\xi }_{k,i}{\left[x_{k,i}^{\times }\right]}^2\\ {\upsilon }_k=\underset{i=1}{\overset{n}{\Sigma }}{\xi }_{k,i}\left[x_{k,i}^{\times }\right]y_{k,i}\end{array}\right).$

2.1)设置奇异避免临界角为θ_c＝π/2，则判别条件为当满足判别条件则直接解算姿态参数否则执行步骤2.2)；

2.2)计算t_k时刻的辅助旋转轴为$a_k^{\oplus }=\left(\begin{array}{cc}\underset{i=1}{\overset{n}{\Sigma }}{\xi }_{k,i}\frac{x_{k,i}}{\left|\right|x_{k,i}\left|\right|},& \left|\right|x_{k,1}\left|\right|\ge \left|\right|y_{k,1}\left|\right|\\ \underset{i=1}{\overset{n}{\Sigma }}{\xi }_{k,i}\frac{y_{k,i}}{\left|\right|y_{k,i}\left|\right|},& \mathrm{otherwise}\end{array}\right);$构建辅助参考向 量为r_i^⊕＝(r_i^Ta^⊕)a^⊕-r_i，将r_i^⊕代替r_i，计算非奇异的辅助姿态为 ${\hat{\sigma }}_k^{\oplus }=\frac{1}{{\gamma }_k^{\oplus }+\sqrt{{\gamma }_k^{\oplus 2}+{\upsilon }_k^{\oplus }{M}_k^{\oplus *}{M}_k^{\oplus *}{\upsilon }_k^{\oplus }}}{M}_k^{\oplus *}{\upsilon }_k^{\oplus };$根据辅助姿态和辅助旋转轴a^⊕，计算最优估计 姿态为${\hat{\sigma }}_k=\frac{1}{1-2{\hat{\sigma }}_k^{\oplus T}{a}_k^{\oplus }+{\hat{\sigma }}_k^{\oplus T}{\hat{\sigma }}_k^{\oplus }}[(1-{\hat{\sigma }}_k^{\oplus T}{\hat{\sigma }}_k^{\oplus })I-2[{\hat{\sigma }}_k^{\oplus \times }\left]\right]a_k^{\oplus }.$在100s的仿真时间内，采用本 方法得到的姿态估计精度如图3所示。