Open Access
Issue
JNWPU
Volume 37, Number 2, April 2019
Page(s) 283 - 290
DOI https://doi.org/10.1051/jnwpu/20193720283
Published online 05 August 2019

© 2019 Journal of Northwestern Polytechnical University. All rights reserved.

Licence Creative Commons
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

目前对于飞卫星的任务设计[4]、系统设计[2]、能耗设计[5]和传感器设计[6]等领域的研究已经逐渐丰富。文献[7]提出了一种使用MEMS控制力矩陀螺(control moment gyroscope, CMG)设计的飞卫星姿态控制执行器, 对MEMS CMG进行建模并详细介绍了一个的飞卫星设计方案。为了降低制造成本, 文献[2]采用商业化元件选型, 对嵌入式卫星的质量、能源、硬件、软件、散热问题做了详细的介绍, 但是对姿态控制系统分析不足, 未能针对芯片的性能提出软件设计方案。文献[8]重点分析了PCBSat和WikiSat 2颗飞卫星在功能、结构、功耗等方面的特点。PCBSat设计有2个太阳敏感器、1个GPS导航接收机和被动气动控制系统, 而WikiSat包含4个光学传感器、1个三轴加速度计、1个三轴陀螺仪以及磁力距器用于姿态控制。但是由于飞卫星的质量、体积和功耗方面的限制, 目前飞卫星很少设计姿态和轨道控制任务[4]。

与纳卫星和皮卫星姿态控制相比, 飞卫星姿态控制系统的MCU计算能力偏弱, 且存储空间更小。例如,现有卫星PCBSAT[8]采用了ATmega 128L微控制器,FLASH容量为128 kB;RyeFemSat[9]飞卫星采用了32 kB FLASH的CC2510控制器;WikiSat[8]飞卫星所采用的控制器ATmega328P也仅有32 kB FLASH。因此考虑工程实现的可行性, 飞卫星姿态控制建模也应该与纳卫星和皮卫星姿态控制建模不尽相同[7]。工程上一般使用历表的方法解决此类复杂模型计算问题[10]。本文根据飞卫星处理器计算性能较弱的特性, 对卫星姿态控制系统中计算复杂度较高、存储空间较大的地磁场模型进行简化, 设计了均匀模型和非均匀模型, 并通过将模型应用到卫星姿态确定中对比分析, 最后总结并选取适用于飞卫星MCU的80点非均匀地磁场历表模型。

1 均匀地磁场历表模型

1.1 地磁场标准模型

通常可将地磁模型分为高空全球模型和区域模型, 全球模型主要以国际地磁参考IGRF为代表, 其模型的建立利用地面、海洋、高空和卫星地磁测量数据, 通过Gauss理论而建立。自1968年起, 国际地磁与高空物理协会(IAGA)每隔5年会发布国际地磁参考模型(IGRF)。

本文的地磁场模型采用IAGA提供的IGRF2010模型。使用球谐波函数表示的地球磁位势为

式中:Re为地球半径; r为格林尼治起算的东经; λθ分别为地理经度和地心余纬; gnmhnm为基本磁场的高斯系数; Pnm(cosθ)为nm阶缔合勒让德函数。

地球磁势位在三轴方向的梯度即为地磁场矢量:

式中,Pnm(θ)和使用递归方式计算。

仿真分析所用到的地磁矢量由IGRF模型模拟得到, 但是受限于星载计算机计算速度和星上存储容量, 实际星上加载地磁模型的级数不能取过高[11]。为了本文的仿真分析, 本文选择了10阶球谐级数的地磁场模型, 模拟了轨道高度为650 km, 轨道倾角为97°, 升交点赤经10°的卫星轨道上的地磁场矢量随时间变化曲线如图 1所示。

thumbnail 图1

卫星轨道上地磁场矢量随时间变化

1.2 地磁场历表模型

地磁场模型只有在精度比较高时, 才能获得满意的精度, 这种实时迭代解算的计算量是星上计算机很难承受的。以IGRF2010模型为例, 该模型包含2003年最新修订的195个高斯系数, 在弱太阳活动时间段内, 计算精度可达100 nT以下[12]。历表模型以损失精度、增加存储量为代价, 保证了地磁场模型的计算效率和双矢量算法的正常运行。

本文设计地磁场历表模型的主要思想为:在存在有效GPS数据的情况下, 通过卫星的轨道位置(GPS)推算轨道6根数中的真近点角, 进而利用真近点角和地磁场强度的对应关系进行建表和查表, 然后对查到的历表数据进行线性插值, 最终解算出卫星所处轨道位置的磁场强度数据。算法设计流程如图 2所示。在没有有效GPS数据的情况下, 通过卫星运行轨道推算卫星轨道位置的真近点角, 然后查表、插值得到卫星所处轨道位置的磁场强度。

thumbnail 图2

地磁场模型算法流程图

1.3 历表设计

本文设计的地磁场历表模型充分考虑了地磁场模型在不同应用场景下可能存在的不同精度需求, 即考虑了同一轨道不同表格长度的历表模型。本文设计的地磁场历表模型的存储格式由三部分组成:起始真近点角θ0、表格长度L和地磁场数据M0图 3表述了不同起始真近点角和表格长度的地磁场历表模型下的表格存储格式。

本文存储实数数据采用单精度浮点格式, 因此存储初始真近点角和表格长度均使用4字节, 而存储轨道上每一点的三轴地磁场需要占用12字节的存储空间, 因此存储L表格长度为的历表模型需要12L+8字节。

thumbnail 图3

不同磁场历表模型的存储示例

1.4 选取表格长度

卫星轨道上的地磁场是连续分布的, 理论上来说可以采集到无穷多点的地磁场数据, 但是实际上表格长度选取的过大无疑会增加存储负担, 同时冗余数据的存储也无助于提高算法精度。相反, 表格数据过小的情况下, 查表法历表模型本身存在较大的误差甚至错误, 导致无法利用这种历表模型对姿态确定的双矢量算法精度进行测试。因此选择一个合适的表格长度显得尤为重要。

thumbnail 图4

表格长度和相应真近点角间隔示意图

1.5 解算真近点角

在卫星飞行的过程中, 认为卫星的运行轨迹为圆形, 轨道高度固定, 其圆点在地心。轨道6根数是经典万有引力定律描述天体按圆锥曲线运动时所必需的6个参数, 如图 5所示。

轨道上任意一点的GPS数据均可由轨道6根数计算确定, 轨道6根数确定卫星经度和纬度的计算如公式(3)所示。

式中:λ0为0时刻的升交点赤经;λs(t)为经度;φs(t)为纬度;i为轨道倾角;θ为真近点角;ωe为地球自转角速度。

当卫星所处的轨道和轨道位置即λs(t)和φs(t)的值已知时, 可以通过(3)式求解出真近点角θ的值, 真近点角的取值范围是0~360°。

在未给出GPS数据的情况下, 卫星的运动可近似看成匀速圆周运动, 即卫星绕地心飞行的角速度ωs为定值, 真近点角θ随时间均匀变化。因此, 在给定初始真近点角θ和飞卫星飞行时间Δt的情况下, 可以推算出当前时刻卫星所处轨道位置的真近点角为

使用推算真近点角的方法, 只需要卫星轨道参数、飞行角速度和至少一个初始真近点角或者GPS数据, 即可推算之后任意时刻卫星所处的轨道位置和该处地磁场强度的大小。

thumbnail 图5

轨道6根数

2 非均匀地磁场历表模型

观察图 1发现磁场数据并不是随着时间均匀变化, 而是在某些时间段磁场曲线弯曲程度大, 另外一些磁场数据弯曲程度小, 因此采用均匀采样的方法势必造成曲率小的时间段内磁场数据冗余采样, 而在曲率较大的时间段内磁场数据不足造成精度损失。为了解决磁场曲线弯曲程度和采样点数之间的矛盾, 本文借用曲率的概念, 在曲线曲率大的地方多采样, 曲率小的地方少采样。地磁场曲线的曲率计算过程中需要对公式(2)求二阶微分, 因此计算复杂度较高, 本文采用离散曲率的方法代替直接对地磁场直接求导的过程, 从而实现计算过程的简化。

2.1 极坐标归一化

考虑到磁场曲线的横坐标表示时间, 纵坐标表示磁场强度, 直接对图 1曲线求曲率是没有意义的, 因此在求曲率之前先把磁场曲线转换成极坐标形式。这种使用一维曲线来表示二维目标轮廓的方法在形状特征提取方面也有应用[13]。

图 1转换成极坐标的形式考虑到三轴磁场曲线是相互独立的, 于是再对每轴进行单独的归一化。需要注意的是, 必须将极坐标下归一化磁场曲线转换成直角坐标系下归一化地磁场曲线的形式之后, 才可以对曲线求导。其原因在于直角坐标转换成极坐标的过程中图形不同区域曲率的相对大小发生了改变, 损失了曲率的相对信息。

2.2 U弦长曲率

U弦长曲率是一种离散曲率计算方法[14], 与现有的离散曲率计算方法相比, U弦长曲率具有更强的抗旋转性和抗噪性, 适用于完成曲线匹配等对曲率计算稳定性要求高的一类任务。方法基本思想是:对于输入的参数U, 按照欧氏距离在曲线当前点处确定一个支持领域, 并应用文献[15]中的曲线精化策略改进计算精度, 由此计算离散曲率。

记:l={pi:(xi, yi, i=1, 2, …, N)}为一条数字曲线, 计算pi处的U弦长曲率时, 应用与支持领域前后臂矢量夹角相关的一个余弦值作为离散曲率, 具体计算公式为

式中,Qipi点处的离散曲率; si=sign[(xi-xib)(yif-yib)-(xif-xib)(yi-yib]为离散曲率值的符号; Di=pifpippifpib这两点间的欧氏距离; (xif, yif), (xib, yib)分别是pifpib的坐标; U为支持邻域的长度, 具体参见文献[11]。

2.3 K曲率加权非均匀采样

K曲率加权的方法是本文提出的一种对每个离散点进行加权的方法, K代表离散点的固定权重, 通过实验的方法选取最优的K值。令地磁场数据为1列N个有序的离散点, 在pi处的离散曲率为Qi, 则所有离散点曲率模的和为

N个有序离散点的离散曲率进行归一化处理, 得到pi点处的离散曲率模为|Qi|/Qsum。在对地磁场数据进行非均匀离散化处理时, 若只考虑每个点的曲率权重, 则容易造成所采样的点过度地聚集到曲率较大的位置处。所以对于每个离散点除了考虑曲率权重外, 还应附加自身固有的权重, 用以平衡这种过度聚集的情况。本文将固有权值K赋值给N个点, 则pi点的权重为K/N。因此, K曲率加权的采样方法需要满足

式中, n为采样点数, j=1, 2, …, n, mj表示采样点位置对应标准模型中的位置。

易知, K的选取决定了曲线的非均匀采样结果, 当K无限接近于正无穷时, K曲率加权的采样方法便蜕化成了均匀采样。K曲率加权采样方法是曲率和均匀采样的融合, 不同的K值会得到不同的非均匀采样结果, 通过实验方法筛选出最优的K, 使得非均匀采样的点能够更好地描述原曲线。

图 6是采用K曲率加权非均匀采样方法得到的80点地磁场模型, 和文章的设计意图相同, 采样密度较高的地方集中在曲率较大的区域。可以看出随着采样点数的减小, 三轴地磁场图形变得越来越简陋, 离散点所表达的信息也越来越少, 图形的轮廓变得越来越模糊。

K曲率加权采样选取了采样点之后, 使用线性内插获得地磁场数据lK={BiK:(1=1, 2, …, N)}, 与标准地磁场数据lS={BiS:(i=1, 2, …, N)}之间的平均误差记作

本文选取的一个轨道的地磁场数据的N=5 875。E用于描述2条磁场曲线的偏离程度, E越小, 表示2条曲线越贴近。

图 7表示采样点数为80时KE的关系, 从图中可以看出, 采样点数相同时, 均匀采样的E为常数, 记作EJ, 而K曲率加权采样的EK变化, 取最小的E记作EminK。总能找到存在EminK < EJ, 表示K曲率加权采样得到的磁场曲线和原磁场曲线更加贴近, 从而证明了算法的有效性。另外, 随着采样点数的减小, EminKEJ都在增加, 表示采样点数减小, 2种方法的拟合曲线和原曲线之间的误差都增加。

thumbnail 图6

80点K曲率加权采样

thumbnail 图7

采样点数为80

3 飞卫星姿态确定实验

3.1 双矢量姿态确定模型

本文中飞卫星的姿态确定系统采用太阳敏感器和磁强计的组合作为敏感器, 并基于矢量算法进行飞卫星姿态确定算法的设计。双矢量定姿的思路如下:

1) 根据轨道坐标系中的单位矢量S0, B0, 在轨道坐标系中建立新的正交坐标系L, 各坐标表轴的单位矢量为:

2) 根据卫星本体坐标系中的单位矢量Sb, Bb在卫星本体系中建立一个正交坐标系S, 各坐标轴的单位矢量为:

3)定义矩阵ML=[l1   l2   l3], MS=[s1    s2   s3], 可得Ms=Ab0ML, 则存在唯一的正交姿态矩阵, 满足

图 8描述了双矢量姿态确定的流程, 其中对地磁场数据的使用包括图中虚线框内的3个流程。首先, 测定轨道位置、加载轨道系地磁场数据; 其次, 初始化卫星的轨道位置并进行轨道推演; 最后使用地磁场模型解算对应的轨道系地磁场矢量数据。图中磁场强度模型和太阳能电池板电压模型引自文献[16]。

thumbnail 图8

双矢量姿态确定流程

3.2 蒙特卡洛模拟

在模拟地磁场数据和太阳矢量数据的过程中, 由于随机误差的引入, 往往导致模拟的量测数据和实际量测的数据之间存在偏差, 为了降低此类偏差以提高模拟精度, 本文采用蒙特卡洛模拟降低此类误差的影响。蒙特卡洛方式是一种概率统计理论为指导的数值计算方法, 指的是使用随机数来解决很多计算问题的方法。

本文使用蒙特卡洛模拟执行M次单轨姿态确定解算, 将IGRF地磁场模型作为轨道系地磁场数据解算出的欧拉角记作Qi, jS, 其中i=1, 2, …, N, j=1, 2, …, M。同样, 使用均匀采样地磁场模型解算的欧拉角记作Oi, jU, 使用K曲率加权地磁场模型解算的欧拉角记作Oi, jK。为描述Oi, jUOi, jK分别与Oi, jS之间的误差, 本文引入了平均均方根误差θARMSE的概念, 如公式(13)~(14)所示。其中, 上标SU表示均匀采样, 上标SK表示非均匀采样。θARMSE越小说明Oi, jUOi, jK分别与Oi, jS之间的差异越小, 即采样结果越精确。

3.3 实验结果和分析

实验中, 选取蒙特卡洛模拟次数M=50, 轨道周期N=5 875 s, 地球半径为6 385 km, 轨道高度为650 km, 升交点赤经为10°, 近地点幅角为10°, 离心率为0.1, 轨道倾角为97°, 轨道常量为3.986×1014, 轨道角速度由轨道常量和轨道半径计算得出; 地磁场模型噪声的均值为0, 标准差为10-8 T; 太阳矢量模型的均值为0, 标准差为0.01。

这里以X轴为例, 从图 9中可以看出θARMSESK < θARMSESU始终成立, 说明使用K曲率加权的地磁场模型得到的姿态确定结果和使用标准姿态确定结果更加接近, 从而证明了K曲率加权的方法比均匀采样的方法更加优越; 另外, 当采样点数较少时, θARMSESKθARMSESU之间差异更明显, 说明采样点数越少, K曲率加权采样的优点突出的更明显; 最后, 当采样点数为80左右时, 三轴θARMSESK均趋于稳定, 因此在不影响姿态确定精度的前提下, 本文选取80点作为地磁历表模型的采样个数较为合理。

图 10为使用80采样点的K曲率加权地磁场模型得到的姿态确定结果,从图中可以看出其姿态确定精度维持在7°以内。

thumbnail 图9

XθARMSE对比

thumbnail 图10

80采样点的K曲率加权姿态确定

3.4 应用前景分析

本文的研究表明,对于一条轨道(5 875 s)三轴地磁场数据,使用80采样点的地磁场历表模型基本可以替代IGRF模型在姿态确定中的作用。文献[10]提出了网格化地磁场模型,即使用经纬度数据去映射存储磁场数据,若以4字节的float格式存储磁场数据,则共需要128.2 kB。由于飞卫星没有变轨功能,因此网格化磁场模型中存储了大量的数据冗余。相比之下,若IGRF模型和80采样点的地磁场数据以4字节的float格式存储,则各需要70.5 kB和0.96 kB。由此可见,本文提出的基于K曲率加权的地磁场历表模型不仅占用存储空间较小,满足在前文提及的3种飞卫星的FLASH中轻量化存储的需求,同时可实现将FLASH中的地磁场数据一次性加载到MCU的RAM中,从而达到对地磁数据快速访问的目的。所以本文设计的磁场历表模型不仅极大程度上节约了存储空间,而且提高了程序的运行速度,从而缩短了算法的控制周期。

4 结论

本文致力于解决地磁场模型简化计算量和存储空间的问题,给出了地磁场历表模型的表格存储格式;并提出了K曲率加权的非均匀采样地磁场模型,在采样点数固定的情况下,给出非均匀点的选取方法;最后,为验证K曲率加权模型的正确性,文章采用蒙特卡洛模拟方法对卫星姿态确定过程进行分析和验证。K曲率加权模型在满足卫星姿态确定精度的前提下,尽量减少选取的采样点数,从而降低所需的存储空间,在计算能力偏弱、存储空间较小的飞卫星姿态确定领域具有非常旷阔的应用前景。

References

  1. Su Ruifeng, Zhang Keke, Song Haiwei. Summarization of Very Small Satellite Development Spacecraft Engineering, 2013, 22(6): 104–111 (in Chinese) [Article] [Google Scholar]
  2. Perez T R Subbarao K. A Survey of Current Femto-Satellite Designs, Technologies, and Mission Concepts[J]. Journal of Small Satellites, 2016, 5(3): 467–482 [NASA ADS] [Google Scholar]
  3. Hadaegh F Y, Chung S J, Manohara H M. On Development of 100-Gram-Class Spacecraft for Swarm Applications[J]. IEEE Systems Journal, 2014, 10(2): 1–12 [Article] [Google Scholar]
  4. Geoffrey R, Mc Vittie Krishna Dev Kumar. Design of a COTS Femto-Satellite and Mission[C]//AIAA SPACE 2007 Conference & Exposition, Long Beach, California, 2007: 1–7 [Article] [Google Scholar]
  5. Barnhart D J, Vladimirova T, Baker A M, et al. A Low-Cost Femtosatellite to Enable Distributed Space Missions[J]. Acta Astronautica, 2009, 64(11/12): 1123–1143 [Article] [NASA ADS] [CrossRef] [Google Scholar]
  6. Izquierdo L, Tristancho J. Next Generation of Sensors for Femto-Satellites Based on Commercial-of-the-Shelf[C]//IEEE Digital Avionics Systems Conference, 2011: 81–87 [Article] [Google Scholar]
  7. Post M, Bauer R, Li J, et al. Study for Femto Satellites Using Micro Control Moment Gyroscope[C]//IEEE Aerospace Conference, Montana, 2016: 1–8 [Article] [Google Scholar]
  8. Tahri N, Hamrouni C, Alimi A M. Study of Current Femto-Satellite Approches[J]. International Journal of Advanced Computer Science & Applications, 2013, 4(5): 148–153 [Article] [CrossRef] [Google Scholar]
  9. Stuurman B, Kumar K. RyefemSat: Ryerson University Femtosatellite Design and Testing[C]//Spaceops 2010 Conference, Huntsviue, Alabama, 2010 [Google Scholar]
  10. Meng Tao. Design and Implementation of Attitude Determination and Control System for Pico Satellite[D]. Hangzhou, Zhejiang University, 2006 (in Chinese) [Google Scholar]
  11. Liu Haiying. Research on Key Technologies of Micro Satellite Attitude Control System[D]. Nanjing, Nanjing University of Aeronautics and Astronautics, 2008 (in Chinese) [Google Scholar]
  12. Li Dong. Research on Attitude Control System Key Technologies for Micro-Satellite[D]. Shanghai, Shanghai Institute of Microsystem and Information Technology Chinese Academy of Sciences, 2005 (in Chinese) [Google Scholar]
  13. Zhou Zhengjie Wang Runsheng. A Contour-Based Method of Feature Extraction Computer Engineering and Applications, 2006, 42(14): 92–94 (in Chinese) [Article] [Google Scholar]
  14. Guo Juanjuan Zhong Baojiang. U-Chord Curvature:A Computational Method of Discrete Curvature Pattern Recognition and Artificial Intelligence, 2014, 27(8): 683–691 (in Chinese) [Article] [Google Scholar]
  15. Zhong Baojiang Liao Wenhe. Corner Detection Based on Accumulative Chord Length of Refined Digital Curves Journal of Computer-Aided Design and Computer Graphics, 2004, 16(7): 939–943 (in Chinese) [Article] [Google Scholar]
  16. Liu Yong, Liu Kunpeng, Li Yilan, et al. A Ground Testing System for Magnetic-Only ADCS of Nano-Satellites[C]//IEEE Navigation and Control Conference, Nanjing, 2016 [Article] [Google Scholar]

All Figures

thumbnail 图1

卫星轨道上地磁场矢量随时间变化

In the text
thumbnail 图2

地磁场模型算法流程图

In the text
thumbnail 图3

不同磁场历表模型的存储示例

In the text
thumbnail 图4

表格长度和相应真近点角间隔示意图

In the text
thumbnail 图5

轨道6根数

In the text
thumbnail 图6

80点K曲率加权采样

In the text
thumbnail 图7

采样点数为80

In the text
thumbnail 图8

双矢量姿态确定流程

In the text
thumbnail 图9

XθARMSE对比

In the text
thumbnail 图10

80采样点的K曲率加权姿态确定

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.