CN103674023A - 一种基于陀螺精确角度关联的星敏感器动态测姿方法 - Google Patents

Abstract

本发明公开了一种基于陀螺精确角度关联的星敏感器动态测姿方法，属于导航定位领域。在现有技术对星敏感器每个测量曝光帧进行动态补偿处理并得到含有动态误差和噪声影响的恒星匹配矢量矩阵的基础上，通过与星敏感器固联的三个陀螺的组合体精确测量出星敏感器各相邻测量帧之间的变换矩阵，由这种变换矩阵将星敏感器相邻测量帧的匹配矢量矩阵关联起来，最后用一系列的关联测量帧建立一个关联的测量方程，相当于将一系列的测量帧等同于一个单测量帧处理，用最小二乘法求解出姿态矩阵，可有效抑制动态误差和噪声影响，从而得到高精度的运动载体的动态姿态及其变化信息。

Description

一种基于陀螺精确角度关联的星敏感器动态测姿方法

技术领域

本发明涉及导航定位领域,特别是惯性测量与星敏感器高精度复合导航定位领域，具体是一种基于陀螺精确角度关联的星敏感器动态测姿方法。

背景技术

陀螺是一种惯性测量器件，可以测量出运动载体相对于惯性空间的转动角速度，角速度的测量精度是表征陀螺性能高低的关键参量。角度是对角速度的积分，随着积分时间的增长，角度的测量误差将越来越大，这是包括陀螺在内所有惯性测量器件存在的共性问题；星敏感器（SS：Star Sensor)是一种依赖于恒星实现运动载体导航定位的设备，特别适合于运动载体的姿态测量，其主要特点是测量精度不随时间而变化，但受运动载体的动态影响很大，特别是运动载体的转动运动，当运动载体的转动角速度越大时，测得的恒星在星敏感器上的成像模糊就越严重，测量精度越差，与星敏感器在静止时得到的测量精度相差甚远。所以，在高动态环境下如何保持星敏感器姿态测量的高精度是目前导航定位亟待解决的难题。

2012年8月发表在《Mathematical Problems in Engineering》（《工程数学问题》）上的论文“A Novel Approach Based on MEMS-Gyro’s Data Deep Coupling for Determining theCentroid of Star Spot”（一种基于MEMS陀螺数据深耦合进行星点提取的新方法）介绍了一种采用三个MEMS陀螺实时测量星敏感器曝光时间内的角速度，根据角速度计算出星点在星敏感器像平面上的长度和方向，确定星点质心提取算法的窗口大小并用质心算法提取星点中心坐标的方法，该方法基于扩展卡尔曼滤波（EKF）算法，建立了陀螺实测角速度与星点质心提取算法融合的EKF方程，得到了动态下星敏感器单次测量星点坐标的最优估计，对5等星8゜/s的动态情况，提取精度优于1个像素。

中国专利“一种基于星敏感器和陀螺的高精度卫星姿态确定方法”（CN201010194288.2）主要针对卫星姿态确定系统中的状态方程的非线性问题，提出了非线性预测滤波方法估计模型误差，对状态方程进行补偿后利用插值滤波进行最优状态估计，得到卫星的姿态。该专利并未涉及星敏感器动态图像的处理问题。

1999年发表在《Proceedings of the Fourth ESA International Conference on SpacecraftGuidance,Navigation and Control Systems,Netherlands,1999》（第四届ESA国际航天器引导、导航与控制系统大会，荷兰，1999）上的论文“Dynamical Binning for High Angular Rate StarTracking”（高角速度下实施星跟踪的动态同步像素移位技术）介绍了一种与主动像素传感（APS：Active Pixel Sensor）技术结合的动态同步像素移位（Dynamical Binning）的方法，通过陀螺实时测量星敏感器的角速度，同步控制APS-CCD中像素的移动，消除动态的影响。该法虽在卫星上使用达到了较大动态下的高精度，但无法消除绕星敏感器光轴转动角速度的影响。

发表在2012年第12期《Sensors》（《传感器》）上的论文“Blurred Star Image Processingfor Star Sensors under Dynamic Conditions”（星敏感器动态模糊星图处理）介绍了一种根据相邻帧得到的角速度对动态拖尾图像进行恢复的方法，该方法先将星敏感器的星图帧进行傅里叶变换到频域，然后在频域根据角速度大小利用维纳滤波反解卷积，最后再进行傅里叶逆变换到空域，可以较好地解决动态拖尾和星像动态模糊的影响，恢复的星点可达亚像素的精度。该方法利用的是相邻帧的信息来预测运动参数，在高动态情况下误差较大，反解卷积得到的原星图与静态星图相差很大。另外，姿态解算时处理的依然是单帧星图，单帧星图视场内有限的有效恒星数量在统计上使得其精度有限。

中国专利“一种星敏感器高动态下的星跟踪方法”（CN200810209622.X）提出根据星敏感器的上一帧预测的理想星像坐标信息来从当前帧星图中提取相应星像坐标，然后利用当前帧提取的星像坐标以及这些星像对应的天球坐标来计算星敏感器当前帧的姿态信息，提供了一种根据上一帧预测的理想星像中心为参考来提取实际星像位置的跟踪算法。中国专利“一种高动态下恒星星像恢复方法”（CN201310053071.3）提出根据前两帧的姿态信息预测当前帧的姿态，预测的依据是前两帧星图的姿态机动四元数信息，该机动四元数被用于当前帧恒星星像的恢复，进而提取星像坐标进行姿态解算。二者的本质均是在没有陀螺的情况下，根据星敏感器测量的前两帧姿态得到角速度来预测和处理当前图像帧，只是具体方法不同而已，因动态对星敏感器的姿态测量精度有很大影响，由此根据这两种方法得到的角速度误差很大，依此角速度进行动态补偿效果较差。同时，二者在进行姿态解算时实际利用的还是单帧星图的星像信息，并没有扩大视场以及增加视场内有效星点的数量，姿态精度受限于单帧星图恒星数量是有限的。

因此，基于陀螺测量得到的角速度对星敏感器动态图像帧进行补偿处理，可有效降低角运动对星点坐标提取的影响，提高星敏感器在动态下的测量精度。但由于星敏感器的噪声、曝光时间、反卷积等因素的共同影响，其星点坐标提取精度与静态相比误差还是很大，还不能满足动态下高精度测量的要求。

发明内容

本发明所要解决的技术问题是，针对现有技术不足，提供一种基于陀螺精确角度关联的星敏感器动态测姿方法，以满足动态下高精度测量的要求。

本发明的思路是：对星敏感器每个测量曝光帧用现有技术进行处理，得到每个测量帧经角速度补偿处理后的星点坐标，与标准星图进行搜索与匹配，得到各星点坐标的匹配矢量对。在此基础上，再利用陀螺精确测量出星敏感器各相邻测量帧之间的角度关系，通过这种精确的角度关系将星敏感器前后测量帧匹配矢量对关联起来。最后用一系列的关联测量帧建立一个关联的测量方程，相当于将一系列的测量帧等同于一个单测量帧处理，用最小二乘法求解出姿态矩阵，从而得到高精度的运动载体的姿态及其变化信息。

定义惯性坐标系O-X_iniY_iniZ_ini，星敏感器坐标系O-X_ssY_ssZ_ss，陀螺组合体（GU：GyroUnit）O-X_guY_guZ_gu，如图1、图2所示。本说明书中，通过给相关物理量添加上标和/或下标“gu”、“ss”、“ini”分别表示该变量为陀螺组合体坐标系、星敏感器坐标系、惯性坐标系和/或从一个坐标系变换到另一个坐标系的变量，以下同。

设星敏感器坐标系到惯性坐标系的姿态矩阵为

星敏感器光轴的单位矢量为I，该矢量对应的惯性坐标系方向矢量L为：

$L=C_{\mathrm{ss}}^{\mathrm{ini}}I---\left(1\right)$

由于动态的影响，单帧中星点位置提取有较大的测量误差，可用其方差σ_ε表征。不失一般性，σ_ε以像元为单位，则单星点角度测量误差E_single可写为：

$E_{\sin \mathrm{gle}}=\frac{A_{\mathrm{FOV}}}{N_P}{\mathrm{\σ}}_{\mathrm{\ϵ}}---\left(2\right)$

其中：N_P为星敏感器图像探测器一维像元数，A_FOV为星敏感器的一维视场角。

提取的星点在星敏感器坐标系中的方向矢量用s表示，误差矢量用e表示。星点所对应的恒星在惯性坐标系中的方向矢量为p，则：

$P=C_{\mathrm{ss}}^{\mathrm{ini}}(s+e)---\left(3\right)$

SS的每一帧星图都有多颗星可供匹配使用，设数量为n，该n颗星在SS坐标系的方向矢量分别为s₁、s₂、…、s_n,误差矢量分别为e₁、e₂、…、e_n；n的平均值

可用下式表示：

$\stackrel{\‾}{n}=6.57e^{1.08M_{\mathrm{ss}}}\frac{1-\cos \left(\frac{A_{\mathrm{FPV}}}{2}\right)}{2}---\left(4\right)$

其中：M_ss为SS星等探测灵敏度。

n颗星在SS坐标系的方向矢量用矩阵S表示，对应的误差矢量用矩阵E表示，则S和E用下式表示：

S=[s₁,s₂,...,s_n]       (5)

E=[e₁,e₂,...,e_n]

经过与标准星图的搜索与匹配，SS坐标系中每一颗星的方向矢量s与惯性坐标系中的唯一方向矢量p对应。n颗恒星在惯性坐标系下的方向矢量分别为p₁、p₂、...、p_n表示,则n颗恒星在惯性坐标系下的方向矢量矩阵P可表示为：

P=[p₁,p₂,...,p_n]           (6)

根据(3)式，P与S满足：

$P=C_{\mathrm{ss}}^{\mathrm{ini}}(S+E)---\left(7\right)$

根据n颗恒星在SS坐标系下的方向矢量S与在惯性坐标系下的方向矢量P，由(7)式，就可以用最小二乘法求出单帧姿态矩阵

其最优解即为单帧最优姿态矩阵。这是一个定姿过程，从而确定了SS坐标系到惯性坐标系之间的关系。可以用三个欧拉角φ_x、φ_y、φ_z表示，即星敏感器坐标系可由惯性坐标系绕Z轴、X轴、Y轴依次旋转角度φ_z、φ_x、φ_y得到，如图1所示，此时

的具体形式如(8)式：

$C_{\mathrm{ss}}^{\mathrm{ini}}=\left[\begin{array}{ccc}\cos {\mathrm{\φ}}_y\cos {\mathrm{\φ}}_z-\sin {\mathrm{\φ}}_y\sin {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_z& -\cos {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_z& \sin {\mathrm{\φ}}_y\cos {\mathrm{\φ}}_z+\cos {\mathrm{\φ}}_y\sin {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_z\\ \cos {\mathrm{\φ}}_y\sin {\mathrm{\φ}}_z+\sin {\mathrm{\φ}}_y\sin {\mathrm{\φ}}_x\cos {\mathrm{\φ}}_z& \cos {\mathrm{\φ}}_x\cos {\mathrm{\φ}}_z& \sin {\mathrm{\φ}}_y\sin {\mathrm{\φ}}_z-\cos {\mathrm{\φ}}_y\sin {\mathrm{\φ}}_x\cos {\mathrm{\φ}}_z\\ -\sin {\mathrm{\φ}}_y\cos {\mathrm{\φ}}_x& \sin {\mathrm{\φ}}_x& \cos {\mathrm{\φ}}_y\cos {\mathrm{\φ}}_x\end{array}\right]---\left(8\right)$

由于误差矩阵E的存在，单帧最优姿态矩阵

与之间存在误差，其值可用SS的光轴指向误差表示。设一帧星图中不同星点的误差是相互独立的，则由(2)式，n颗星定姿后SS的光轴指向误差E为：

$E=\frac{A_{\mathrm{FOV}}}{N_P}\frac{{\mathrm{\σ}}_{\mathrm{\ϵ}}}{\sqrt{n}}---\left(9\right)$

显然，对单帧测量而言，星点位置提取误差σ_ε越大，则光轴指向误差E越大。这是现有星敏感器测量技术普遍存在的问题。

本发明在单帧测量的基础上，将视野拓宽到帧与帧之间的关联上。本发明提出，将三个单轴陀螺组成的相互正交安装的陀螺组合体1与星敏感器2捷联刚性固联，如图2所示。

设SS坐标系O-X_ssY_ssZ_ss与GU坐标系O-X_guY_guZ_gu之间的三个欧拉角为ψ_x、ψ_y、ψ_z，从GU坐标系到SS坐标系的变换矩阵为

的形式同(8)式类似。

和三个欧拉角ψ_x、ψ_y、ψ_z可通过现有技术（具体见王岩于2009年在《战术导弹控制技术》第26卷第3期上发表的“一种星敏感器与捷联惯导高精度安装误差标定方法”）精确测量标定。

因GU与SS刚性固联，故SS在惯性空间转动的角度可同时由GU精确测量出来。动态条件下SS相对于惯性空间的角速度矢量为Ω^ss，通过GU测量得到的三个角速度分量为Ω_x,Ω_y,Ω_z，Ω^gu=[Ω_x,Ω_y,Ω_z]^T。Ω^gu在SS坐标系上的投影Ω^ss为：

设SS的帧与帧之间的曝光时间间隔为τ，则其曝光频率f_ss=1/τ，陀螺组合体测量角速度的采样周期为δτ。为了保证陀螺测量角速度的实时性，δτ要远小于τ，且与曝光时间间隔同步，τ=mδτ，m为远大于1的整数，如图4所示。SS在τ时间内转动的角度，可由Ω^ss按四元数的运动方程求出。姿态四元数q的运动方程可表示为：

$\stackrel{\·}{q}=\frac{1}{2}q\⊗{[0,{\mathrm{\Ω}}^{\mathrm{ss}}]}^T---\left(11\right)$

其中，

为姿态四元数q的微分，

为四元数乘法符号。

对(11)式进行积分并作离散化处理，并记陀螺第λ次采样前后（间隔δτ）的姿态四元数分别为q(λ-1)、q(λ)，则q(λ)可由q(λ-1)表示为：

$\begin{array}{c}q\left(\mathrm{\λ}\right)=q(\mathrm{\λ}-1)\⊗\mathrm{\δq}\\ =\left[\begin{array}{c}{q}_0(\mathrm{\λ}-1),-q_1(\mathrm{\λ}-1),-q_2(\mathrm{\λ}-1),-q_3(\mathrm{\λ}-1)\\ q_1(\mathrm{\λ}-1),-q_0(\mathrm{\λ}-1),-q_3(\mathrm{\λ}-1),-q_2(\mathrm{\λ}-1)\\ q_2(\mathrm{\λ}-1),-q_3(\mathrm{\λ}-1),-q_0(\mathrm{\λ}-1),-q_1(\mathrm{\λ}-1)\\ q_3(\mathrm{\λ}-1),-q_2(\mathrm{\λ}-1),-q_1(\mathrm{\λ}-1),-q_0(\mathrm{\λ}-1)\end{array}\right]\left[\begin{array}{c}\cos \left(\right|\mathrm{\δ\Φ}|/2)\\ (\mathrm{\δ}{\mathrm{\Φ}}_x/|\mathrm{\δ\Φ}\left|\right)\sin \left(\right|\mathrm{\δ\Φ}|/2)\\ (\mathrm{\δ}{\mathrm{\Φ}}_y/|\mathrm{\δ\Φ}\left|\right)\sin \left(\right|\mathrm{\δ\Φ}|/2)\\ (\mathrm{\δ}{\mathrm{\Φ}}_z/|\mathrm{\δ\Φ}\left|\right)\sin \left(\right|\mathrm{\δ\Φ}|/2)\end{array}\right]\end{array}---\left(12\right)$

其中：δΦ为δτ时间内的角增量，δΦ=[δΦ_x,δΦ_y,δΦ_z]^T=Ω^ss(λ)δτ，δΦ_x、δΦ_y、δΦ_z为δΦ沿X轴、Y轴、Z轴的分量，Ω^ss(λ)为陀螺第λ次采样在SS坐标系上的投影值，可通过(10)式计算得到；δq为δτ时间内累计更新的姿态四元数。

在星敏感器第(k-1)帧到第k帧的曝光时间τ内，有m次δτ的采样。q的初始条件为q(0)=[1,0,0,0]^T，利用Ω^ss数据通过(12)式进行m次姿态更新，就可以得到τ时刻的姿态四元数q(τ)。由q(τ)=[q₀,q₁,q₂,q₃]^T就可以得到第k帧到第(k-1)帧的角度关联矩阵：

$C_k^{k-1}=\left\{\begin{array}{ccc}{q}_0^2+q_1^2-q_2^2-q_3^2& 2(q_1{q}_2-q_0{q}_3)& 2(q_1{q}_3+q_0{q}_2)\\ 2(q_1{q}_2+q_0{q}_3)& q_0^2-q_1^2+q_2^2-q_3^2& 2(q_2{q}_3-q_0{q}_1)\\ 2(q_1{q}_3-q_0{q}_2)& 2(q_2{q}_3+q_0{q}_1)& q_0^2-q_1^2-q_2^2+q_3^2\end{array}\right)---\left(13\right)$

SS每隔τ时间间隔获取一帧星图。因为载体的运动，SS测量的一系列星图帧F₁,F₂、…、F_n相对于标准星图（以第一帧F₁为原点）发生姿态变化，如图3所示。

动态条件下，设第k帧坐标系可由第(k-1)帧坐标系按欧拉角θ_z、θ_x、θ_y的旋转顺序得到，则：

$\left\{\begin{array}{c}{\mathrm{\θ}}_x=\arcsin (2q_2{q}_3+2q_0{q}_1)\\ {\mathrm{\θ}}_y=-\arctan \frac{2(q_1{q}_3-q_0{q}_2)}{q_0^2-q_1^2-q_2^2+q_3^2}\\ {\mathrm{\θ}}_z=-\arctan \frac{2(q_1{q}_2-q_0{q}_3)}{q_0^2-q_1^2+q_2^2-q_3^2}\end{array}\right.---\left(14\right)$

因此，帧与帧之间的精确角度关系θ_x、θ_y、θ_z可以由陀螺组合体测量得到。

设SS第(k-1)帧坐标系到惯性坐标系的变换矩阵为

则SS第k帧坐标系到惯性坐标系的变换矩阵为第(k-1)帧和第k帧星图经搜索和匹配的方向矢量矩阵P_k-1、P_k和S_k-1、S_k分别满足下式：

$\left\{\begin{array}{c}{P}_{k-1}=C_{\mathrm{ss}}^{\mathrm{ini}}(k-1)(S_{k-1}+E_{k-1})\\ P_k=C_{\mathrm{ss}}^{\mathrm{ini}}\left(k\right)(S_k+E_k)\\ =C_{\mathrm{ss}}^{\mathrm{ini}}(k-1)C_{\mathrm{ak}}^{k-1}(S_k+E_k)\\ =C_{\mathrm{ss}}^{\mathrm{ini}}(k-1)(C_k^{k-1}{S}_k+C_k^{k-1}{E}_k)\\ =C_{\mathrm{ss}}^{\mathrm{ini}}(k-1)(S_k^{k-1}+E_k^{k-1})\end{array}\right.---\left(15\right)$

由(15)式可以看出，相邻两帧共用相同的姿态变换矩阵

两帧可以合成一帧处理：

$[P_{k-1},P_k]=C_{\mathrm{ss}}^{\mathrm{ini}}(k-1)[(S_{k-1}+E_{k-1}),(S_k^{k-1}+E_k^{k-1})]---\left(16\right)$

这种星图帧称角度关联帧（ACF：Angle Correlated Frame）。与单帧的方程相比，合成帧方程中的恒星矢量数多了一倍，用最小二乘方法求解姿态矩阵

时，误差要降低

倍，等同于两帧星图通过ACF变成了一帧大的星图进行处理。

根据ACF，同样可以将SS测量的一个星图帧序列1、2、…、k关联到第1帧：

$\left\{\begin{array}{c}[P_1,P_2,...,P_k]=C_{\mathrm{ss}}^{\mathrm{ini}}\left(1\right)[(S_1+E_1),(S_2^1+E_2^1),...,(S_k^1+E_k^1)]\\ S_{l_0}^1=\left(\underset{r=2}{\overset{r=l_0}{\mathrm{\Π}}}{C}_r^{r-1}\right)S_{l_0}=C_{l_0}^1{S}_{l_0}\\ E_{l_0}^1=\left(\underset{r=2}{\overset{r=l_0}{\mathrm{\Π}}}{C}_r^{r-1}\right)E_{l_0}=C_{l_0}^1{E}_{l_0}\end{array}\right.---\left(17\right)$

其中：(P₁,S₁)、…、(P_k,S_k)分别为星敏感器测量的第1帧、…、第k帧星图经搜索匹配得到的星图帧中各星点相对于惯性坐标系和星敏感器坐标系中的方向矢量矩阵；E₁、…、E_k为星敏感器测量的第1帧、…、第k帧每帧星图中各星点的误差矢量矩阵；

为星敏感器第1帧坐标系到惯性坐标系的变换矩阵；

为陀螺组合体测量得到的星敏感器第r帧到第(r-1)帧的角度关联矩阵；

为陀螺组合体测量得到的星敏感器第l₀帧直接关联到第1帧的角度关联矩阵；

的最小二乘解因k帧星图的关联，误差要降低倍。只要k足够大，动态误差的抑制效果是非常明显的，使动态星敏感器的测量精度达到与静态相当的程度。

但由于陀螺测量角度的误差会随着时间的积累而增加，对应的星敏感器第k帧直接变换到第1帧的角度变换矩阵

误差也会增大。当

的误差大到与单帧处理得到的姿态矩阵的误差接近时，再增加帧数k已没有意义。因此，帧数k根据星敏感器测量期望精度值E_exp会有一个合适的最大值K_max（帧序列中关联帧数量的最大值），它与陀螺的测量精度ω、星敏感器曝光频率f_ss以及星敏感器测量的期望精度E_exp有关。

SS帧序列最长关联时间：

$T_{\max }=\frac{E_{\mathrm{esp}}}{\mathrm{\ω}}---\left(18\right)$

SS最大关联帧数，即帧数k的最大值为：

$K_{\max }=f_{\mathrm{ss}}{T}_{\max }=\frac{f_{\mathrm{ss}}{E}_{\exp }}{\mathrm{\ω}}---\left(19\right)$

由(19)式可以看出，星敏感器测量的期望精度E_exp越高，星敏感器曝光频率f_ss越低，陀螺精度ω越差，则K_max就越小，反之K_max越大。

K_max的典型值如表1所示：

表1：最大关联帧数表

由(9)式和(19)式可知，基于ACF技术的关联帧序列星敏感器的光轴指向误差为：

$E_{\mathrm{ss}}=\frac{A_{\mathrm{FOV}}}{N_P}\frac{{\mathrm{\σ}}_{\mathrm{\ϵ}}}{\sqrt{n{K}_{\max }}}---\left(20\right)$

从表1可知，在单测量帧动态处理的基础上，为了使SS达到与静态基本相当的精度，0.01°/hr的陀螺是必要且合适的选择，在此精度上可选用的陀螺还是比较多的，如中精度的激光陀螺、高精度的光纤陀螺、高精度的挠性陀螺等。另外，SS一般的曝光频率在10Hz以下，在动态条件下，提高曝光频率既可以降低单帧动态的影响，又可以通过增加关联帧数来提高抑制噪声的能力，从而有利于提高SS的精度。

(17)式给出了一个星图帧序列关联到第1帧的公式，利用最小二乘法求解(17)式得出星敏感器第1帧的最优姿态矩阵

将第1帧的最优姿态矩阵

与各相应帧l₀∈[2,k]的角度关联矩阵相乘，即得到l₀=2、3、…、k各帧的最优姿态矩阵

为第l₀帧到第1帧的角度关联矩阵。根据星敏感器第l_0A帧的最优姿态矩阵

l_0A∈[1,k]，利用(8)式计算出星敏感器第l_0A帧的姿态角φ_x、φ_y、φ_z：

$\left\{\begin{array}{c}{\mathrm{\φ}}_x=\arcsin \left[{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{3,2}\right)\right]\\ {\mathrm{\φ}}_y=-\arctan [{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{3,1}\right)/{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}(3,3)]\\ {\mathrm{\φ}}_z=-\arctan [{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}(1,2)/{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{2,2}\right)]\end{array}\right.---\left(21\right)$

其中，均为星敏感器第l_0A帧的最优姿态矩阵

中的矩阵元素。如

对应中第三行第二列的元素。

基于ACF的星图帧序列关联方法既可以关联到第1帧，也可以关联到其他帧进行滑动处理。不失一般性的，对于星图帧序列1、…、(j-1)、j、(j+1)、…、(j+k)，可通过陀螺测量的角度关联矩阵将第(j+1)帧、…、第(j+k)帧向前关联到第j帧，得到向前关联公式：

$\left\{\begin{array}{c}[P_j,P_{j+1},...,P_{j+k}]=C_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)[(S_j+E_j),(S_{j+1}^j+E_{j+1}^j),...,(S_{j+k}^j+E_{j+k}^j)]\\ S_{j+l}^j=\left(\underset{r=1}{\overset{r=l}{\mathrm{\Π}}}{C}_{j+r}^{j+r-1}\right)S_{j+l}=C_{j+l}^j{S}_{j+l}\\ E_{j+l}^j=\left(\underset{r=1}{\overset{r=l}{\mathrm{\Π}}}{C}_{j+r}^{j+r-1}\right)E_{j+l}=C_{j+l}^j{E}_{j+l},l\∈[1,k]\end{array}\right.---\left(22\right)$

也可以通过陀螺测量的角度关联矩阵将第j帧、…、第(j+k-1)帧向后关联到当前帧(j+k)，得到向后关联公式：

$\left\{\begin{array}{c}[P_{j+k},P_{j+k-1},...,P_j]=C_{\mathrm{ss}}^{\mathrm{ini}}(j+k)[(S_{j+k}+E_{j+k}),(S_{j+k-1}^{j+k}+E_{j+k-1}^{j+k}),...,(S_j^{j+k}+E_j^{j+k})]\\ S_j^{j+l}=\left(\underset{r=l}{\overset{r=1}{\mathrm{\Π}}}{C}_{j+r-1}^{j+r}\right)S_j=C_j^{j+l}{S}_j\\ E_j^{j+l}=\left(\underset{r=l}{\overset{r=1}{\mathrm{\Π}}}{C}_{j+r-1}^{j+r}\right)E_j=C_j^{j+l}{E}_j,l\∈[1,k]\end{array}\right.---\left(23\right)$

其中：(P_j,S_j)、…、(P_j+k,S_j+k)为星敏感器测量的第j帧、…、第(j+k)帧每帧星图经搜索匹配得到的星图帧中各星点相对于惯性坐标系和星敏感器坐标系中的方向矢量矩阵；E_j、…、E_j+k为星敏感器测量的第j帧、…、第(j+k)帧每帧星图中各星点的误差矢量矩阵；

为星敏感器第j帧坐标系到惯性坐标系的变换矩阵，

为星敏感器当前帧(j+k)坐标系到惯性坐标系的变换矩阵；

为陀螺组合体测量得到的星敏感器第(j+r)帧到第(j+r-1)帧的角度关联矩阵，

为星敏感器第(j+r-1)帧到第(j+r)帧的角度关联矩阵，且

为星敏感器第(j+l)帧直接关联到第j帧的角度关联矩阵，

为第j帧直接关联到第(j+l)帧的角度关联矩阵，且

根据(22)式所示的向前关联公式，利用最小二乘法求解出星敏感器第j帧的最优姿态矩阵

将第j帧的最优姿态矩阵

帧到第j帧的角度关联矩阵

相乘，即得第(j+l)帧的最优姿态矩阵

l∈[1,k]；或者，根据(23)式所示的向后关联公式，利用最小二乘法可直接求解出星敏感器当前帧(j+k)的最优姿态矩阵将当前帧(j+k)的最优姿态矩阵

与第(j+l-1)帧到第(j+k)帧的角度关联矩阵

相乘，即得第(j+l-1)帧的最优姿态矩阵

l∈[1,k]。利用最优姿态矩阵

根据(21)式就可以计算出星敏感器第(j+l_A)帧的姿态角φ_x、φ_y、φ_z，l_A∈[0,k]。

上述关联帧算法因大幅度增加了恒星矢量的数量，显然增强了对噪声的抑制能力，与此同时计算量也指数增加。为此，本发明进一步提出了一种基于单帧最优姿态矩阵的简化求解算法—平均最优解：根据(7)式的单帧最优姿态矩阵解

可以得到每帧独立的最优解序列

根据(22)所示的向前关联公式，可得到第j帧的平均最优姿态矩阵

${\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)=\frac{1}{k+1}\{{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)+{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+1)C_j^{j+1}+...+[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)C_{j+k-1}^{j+k}{C}_{j+k-2}^{j+k-1}...C_j^{j+1}\left]\right\}---\left(24\right)$

根据(23)所示的向后关联公式，可得到当前帧(j+k)的平均最优姿态矩阵

${\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)=\frac{1}{k+1}\{{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)+{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k-1)C_{j+k}^{j+k-1}+...+[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)C_{j+1}^j{C}_{j+2}^{j+1}...C_{j+k}^{j+k-1}\left]\right\}---\left(25\right)$

这种平均最优姿态矩阵求解方法的精度和噪声抑制能力与(22)、(23)式基本相当，而计算量大大降低。

另外，根据(24)式，向前关联的平均公式还可进一步简化为递推公式：

$\begin{array}{c}{\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)=\frac{1}{k+1}{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)+\frac{1}{k+1}\{{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+1)+{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+2)C_{j+1}^{j+2}+...\\ +[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k+1)C_{j+k}^{j+k+1}...C_{j+1}^{j+2}]\}{C}_j^{j+1}-\frac{1}{k+1}[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k+1)C_{j+k}^{j+k+1}...C_j^{j+1}]\\ =\frac{1}{k+1}{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)+{\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+1)C_j^{j+1}-\frac{1}{k+1}[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k+1)C_{j+k}^{j+k+1}...C_j^{j+1}]\end{array}---\left(26\right)$

其中：

为星敏感器第(j+1)帧的平均最优姿态矩阵；

为第(j)帧的单帧最优姿态矩阵；

为当前帧(j+k)的下一帧(j+k+1)的单帧最优姿态矩阵。

同样，根据(25)式，向后关联的平均公式也可进一步简化为递推公式：

$\begin{array}{c}{\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)=\frac{1}{k+1}{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)+\frac{1}{k+1}\{{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k-1)+{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k-2)C_{j+k-1}^{j+k-2}+...\\ +[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j-1)C_j^{j-1}...C_{j+k-1}^{j+k-2}]\}{C}_{j+k}^{j+k-1}-\frac{1}{k+1}[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j-1)C_j^{j-1}...C_{j+k}^{j+k-1}]\\ =\frac{1}{k+1}{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)+{\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k-1)C_{j+k}^{j+k-1}-\frac{1}{k+1}[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j-1)C_j^{j-1}...C_{j+k}^{j+k-1}]\end{array}---\left(27\right)$

其中：为星敏感器第(j+k-1)帧的平均最优姿态矩阵；

为第(j+k)帧的单帧最优姿态矩阵；为第(j-1)帧的单帧最优姿态矩阵。

与现有技术相比，本发明所具有的有益效果为：本发明可有效抑制动态误差和噪声影响，从而得到高精度的运动载体的动态姿态及其变化信息，能满足动态下高精度测量的要求。

附图说明

图1为星敏感器坐标系与惯性坐标系的变换图；

图2为星敏感器与陀螺组合体刚性固联示意图；

图3为星敏感器关联帧序列示意图；

图4（a）是陀螺组合体与星敏感器同步采样时序的同步时序图；图4（b）是陀螺组合体与星敏感器同步采样时序的电路结构示意图；

图5为陀螺测量相邻帧变换矩阵流程图；

图6为星敏感器获取恒星匹配矢量序列流程图；

图7为序列关联帧姿态解算流程图；

图8为姿态更新流程图。

具体实施方式

本发明具体实施方式例述如下：

第一步：建立同步采样时序。

为保证陀螺组合体测量角速度和角度解算的精度，GU的采样频率f_gu要远高于SS的曝光频率f_ss。f_gu一般为K Hz量级，而f_ss为10Hz量级。同时，为了确保精确的角度关联，GU的采样频率必须与SS的曝光频率严格同步，二者必须是严格整数关系m，使f_gu=mf_ss，则τ=mδτ，如图4所示。

图4(a)为f_gu与f_ss的同步时序关系示意图，图4(b)为同步电路图。晶体振荡器产生标准时频f₀，第一同步分频器的分频系数为M，第二同步分频器的分频系数为m。则：

$f_{\mathrm{gu}}=\frac{f_0}{M},f_{\mathrm{ss}}=\frac{f_{\mathrm{gu}}}{m}=\frac{f_0}{\mathrm{Mm}}---\left(28\right)$

通过(28)式确保SS与GU的同步关系。

第二步：GU测量相邻帧的角度关联矩阵和角度关系。

以f_gu同步采样GU的三个角速度分量Ω_x，Ω_y，Ω_z，即得GU的角速度矢量Ω^gu=[Ω_x,Ω_y,Ω_z]^T，按(10)式计算星敏感器的角速度

根据(12)式对姿态四元数q进行更新。以f_gu重复采样GU的三个角速度分量，并不断更新四元数q，直至一个测量帧结束。按(13)式求解出相邻帧的角度关联矩阵

按(14)式求解出相邻帧的角度关系。对应的流程图如图5所示。

第三步：按f_ss同步采样处理SS的星像帧。

以f_ss同步采样SS的星像帧，对其按现有方法（刘朝山等编著、2010年国防工业出版社出版的教科书《弹载星敏感器原理及系统应用》）进行运动动态补偿处理，然后为星点提取进行去噪处理，提取星点坐标并构造搜索匹配三角形序列，与标准星图进行搜索匹配，得到与本帧星图相匹配的惯性坐标系下的恒星矢量矩阵P和星敏感器坐标系下的方向矢量矩阵S。对应的流程图如图6所示。

第四步：序列关联帧姿态解算。

重复第二步到第三步，以f_ss同步不断采样处理SS的星像帧序列并解算相邻帧的角度变换矩阵，直至关联帧数等于K_max；按(17)式建立关联帧序列的测量方程，用最小二乘法求解出最优姿态矩阵

将

与各相应帧的变换矩阵相乘，即得第2,…,K_max各帧的最优姿态矩阵

为星敏感器第l帧到第1帧的角度关联矩阵，

按(21)式计算出星敏感器第1帧、…、第K_max帧各帧相对于惯性坐标系的姿态角φ_x、φ_y、φ_z。对应的流程图如图7所示。

第五步：更新SS的星像帧和GU的变换矩阵。

重复第二、三步，得到新帧的匹配星矢量P、S和相邻帧变换矩阵；将新帧与之前的(K_max-1)帧序列按(23)式进行关联，得到向后关联公式，用最小二乘求解姿态矩阵，得到新帧的最优姿态矩阵；按(21)式求解姿态角φ_x、φ_y、φ_z，实时输出姿态。对应的流程图如图8所示。

第六步：重复第五步，不断更新处理测量帧，并实时输出姿态。

Claims (4)

1.一种基于陀螺精确角度关联的星敏感器动态测姿方法，采用由三个相互正交安装的单轴陀螺组成的陀螺组合体以及与所述陀螺组合体固联的星敏感器组成的测姿装置，其特征在于，该方法为：

1）对于星敏感器测量的一个星图帧序列1、…、(j-1)、j、(j+1)、…、(j+k)，通过陀螺组合体测量的角度关联矩阵将第(j+1)帧到当前帧(j+k)向前关联到第j帧，即向前关联公式：

$\left\{\begin{array}{c}[P_j,P_{j+1},...,P_{j+k}]=C_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)[(S_j+E_j),(S_{j+1}^j+E_{j+1}^j),...,(S_{j+k}^j+E_{j+k}^j)]\\ S_{j+l}^j=\left(\underset{r=1}{\overset{r=l}{\mathrm{\Π}}}{C}_{j+r}^{j+r-1}\right)S_{j+l}=C_{j+l}^j{S}_{j+l}\\ E_{j+l}^j=\left(\underset{r=1}{\overset{r=l}{\mathrm{\Π}}}{C}_{j+r}^{j+r-1}\right)E_{j+l}=C_{j+l}^j{E}_{j+l},l\∈[1,k]\end{array}\right.$

或者通过陀螺组合体测量的角度关联矩阵将第j帧到(j+k-1)帧向后关联到当前帧(j+k)，即向后关联公式：

$\left\{\begin{array}{c}[P_{j+k},P_{j+k-1},...,P_j]=C_{\mathrm{ss}}^{\mathrm{ini}}(j+k)[(S_{j+k}+E_{j+k}),(S_{j+k-1}^{j+k}+E_{j+k-1}^{j+k}),...,(S_j^{j+k}+E_j^{j+k})]\\ S_j^{j+l}=\left(\underset{r=l}{\overset{r=1}{\mathrm{\Π}}}{C}_{j+r-1}^{j+r}\right)S_j=C_j^{j+l}{S}_j\\ E_j^{j+l}=\left(\underset{r=l}{\overset{r=1}{\mathrm{\Π}}}{C}_{j+r-1}^{j+r}\right)E_j=C_j^{j+l}{E}_j,l\∈[1,k]\end{array}\right.$

其中：(P_j,S_j)、…、(P_j+k,S_j+k)为星敏感器测量的第j帧、…、第(j+k)帧每帧星图经搜索匹配得到的星图帧中各星点相对于惯性坐标系和星敏感器坐标系中的方向矢量矩阵；E_j、…、E_j+k为星敏感器测量的第j帧、…、第(j+k)帧每帧星图中各星点的误差矢量矩阵；

为星敏感器第j帧坐标系到惯性坐标系的变换矩阵，

为星敏感器当前帧(j+k)坐标系到惯性坐标系的变换矩阵；

为陀螺组合体测量得到的星敏感器第(j+r)帧到第(j+r-1)帧的角度关联矩阵，

为星敏感器第(j+r-1)帧到第(j+r)帧的角度关联矩阵，且

为星敏感器第(j+l)帧直接关联到第j帧的角度关联矩阵，

为第j帧直接关联到第(j+l)帧的角度关联矩阵，且

2）根据上述步骤1）中的向前关联公式，利用最小二乘法求解出星敏感器第j帧的最优姿态矩阵

将第j帧的最优姿态矩阵

与第(j+l)帧到第j帧的角度关联矩阵

相乘，即得第(j+l)帧的最优姿态矩阵

l∈[1,k]；或者，根据上述步骤1）中的向后关联公式，利用最小二乘法直接求解出星敏感器当前帧(j+k)的最优姿态矩阵

将当前帧(j+k)的最优姿态矩阵

与第(j+l-1)帧到第(j+k)帧的角度关联矩阵

相乘，即得第(j+l-1)帧的最优姿态矩阵 ${\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+l-1)={\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)C_{j+l-1}^{j+k},$ l∈[1,k]；

3）利用上述步骤2）得到的星敏感器第(j+l_A)帧的最优姿态矩阵l_A∈[0,k]，计算星敏感器第(j+l_A)帧的姿态角φ_x、φ_y、φ_z：

$\left\{\begin{array}{c}{\mathrm{\φ}}_x=\arcsin \left[{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{3,2}\right)\right]\\ {\mathrm{\φ}}_y=-\arctan [{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{3,1}\right)/{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{3,3}\right)]\\ {\mathrm{\φ}}_z=-\arctan [{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{1,2}\right)/{\stackrel{\≈}{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(\mathrm{2,2}\right)]\end{array}\right.;$

其中，均为星敏感器第(j+l_A)帧最优姿态矩阵

中的矩阵元素。

2.根据权利要求1所述的基于陀螺精确角度关联的星敏感器动态测姿方法，其特征在于，星图帧序列j、(j+1)、…、(j+k)的关联帧数(k+1)的最大值K_max为：

$K_{\max }=\frac{f_{\mathrm{ss}}{E}_{\exp }}{\mathrm{\ω}};$

其中，ω为陀螺测量精度，f_ss为星敏感器的曝光频率，E_exp为星敏感器期望的测量精度。

3.根据权利要求1所述的基于陀螺精确角度关联的星敏感器动态测姿方法，其特征在于，根据所述向前关联公式，利用单帧最优姿态矩阵简化求解第j帧的最优姿态矩阵

得到向前关联的平均公式：利用角度关联矩阵将第j帧、…、第(j+k)帧的单帧最优姿态矩阵

向前变换到第j帧，然后取平均，得到第j帧的平均最优姿态矩阵

${\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)=\frac{1}{k+1}\{{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)+{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+1)C_j^{j+1}+...+[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)C_{j+k-1}^{j+k}{C}_{j+k-2}^{j+k-1}...C_j^{j+1}\left]\right\};$

根据所述向后关联公式，利用单帧最优姿态矩阵简化求解当前帧(j+k)的最优姿态矩阵

得到向后关联的平均公式：利用角度关联矩阵将第j帧、…、第(j+k)帧的单帧最优姿态矩阵

向后变换到第(j+k)帧，然后取平均，得到当前帧(j+k)的平均最优姿态矩阵

${\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)=\frac{1}{k+1}\{{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)+{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k-1)C_{j+k}^{j+k-1}+...+[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)C_{j+1}^j{C}_{j+2}^{j+1}...C_{j+k}^{j+k-1}\left]\right\};$

4.根据权利要求3所述的基于陀螺精确角度关联的星敏感器动态测姿方法，其特征在于，将向前关联的平均公式进一步简化为递推公式：

${\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)=\frac{1}{k+1}{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}\left(j\right)+{\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+1)C_j^{j+1}-\frac{1}{k+1}[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k+1)C_{j+k}^{j+k+1}...C_j^{j+1}];$

其中：

为星敏感器第(j+1)帧的平均最优姿态矩阵；为第j帧的单帧最优姿态矩阵；

为当前帧(j+k)的下一帧(j+k+1)的单帧最优姿态矩阵；

将向后关联的平均公式进一步简化为递推公式：

${\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)=\frac{1}{k+1}{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k)+{\stackrel{\‾}{\tilde{C}}}_{\mathrm{ss}}^{\mathrm{ini}}(j+k-1)C_{j+k}^{j+k-1}-\frac{1}{k+1}[{\tilde{C}}_{\mathrm{ss}}^{\mathrm{ini}}(j-1)C_j^{j-1}...C_{j+k}^{j+k-1}];$

其中：为星敏感器第(j+k-1)帧的平均最优姿态矩阵；

为第(j+k)帧的单帧最优姿态矩阵；为第(j-1)帧的单帧最优姿态矩阵。

Images