一种非球形非匀质SLR（人卫激光测距）卫星太阳光压建模方法

本发明涉及一种基于box-wing物理分析的卫星各部件信息太阳光压建模方法，尤其涉及一种slr非球形非匀质材料卫星的太阳光压模型建模方法，属于slr精密定轨和数据处理领域。

背景技术：

1、slr精密定轨正在向毫米级精度迈进，需要更为精细化的摄动力，以此来满足高精度slr精密定轨和时空基准应用的需求。太阳光压摄动是slr卫星的主要摄动力，太阳光压摄动精准建模直接影响了slr卫星精密定轨精度。slr观测卫星除球形动力学卫星外，还有其他用途的非球形观测卫星，如grace-fo、jason-1/2/3、海洋-2卫星等，而在现有阶段slr数据处理过程中都是把slr观测卫星的太阳光压摄动力当作简单的球模型进行处理，这并不适合复杂结构或者非均匀材料的卫星。因此，研究slr复杂结构卫星或者非匀质材料卫星的太阳光压摄动力精准建模十分重要。项目组已研制出了一套适用于gnss/bds卫星系统的纯物理和半经验半物理的高精度太阳光压模型，具有一定的科研基础，本实用发明拟将建立的适用于gnss的纯物理和半经验半物理的高精度太阳光压模型引入精密slr数据处理中，进行精密定轨和时空基准扩展应用研究，提高其精度和应用水平。目前、国内外还没有针对非球形非匀质材料slr卫星，建立高精度的纯物理和半经验半物理分析太阳光压摄动力模型。

技术实现思路

1、本发明的技术解决问题是：克服现有slr卫星精密定轨中太阳光压摄动球模型的不足，提供一种适合高、中、低轨非球形非匀质slr卫星太阳光压建模方法，就可给出slr卫星受到的太阳光压摄动和各轴向的加速度及其偏导数，可利用计算出的摄动修正卫星的状态从而提高slr非球形非匀质材料卫星精密定轨的精度和轨道预报精度，从而提高其应用水平。

2、本发明的技术解决方案是：一种非球形非匀质slr卫星太阳光压建模方法，其特征是适合任何形状和结构材料卫星进行太阳光压建模，步骤如下：

3、(1)记录待分析卫星id号，获取待分析卫星的三维卫星结构图，对不同结构或者材料的部件进行划分，并对待分析卫星所有部件进行编号，记为i，i＝1,…n；

4、(2)获取待分析卫星的各表面详细参数信息和光学特性信息，包括：表面积ai，法向量发射系数ρi，散射系数δi；

5、步骤(1)和步骤(2)设计为pre_boxwing函数，输入卫星id号，输出为卫星整体的各面信息，包括质量、受照信息，具体如图2所示。

6、(3)对slr形状复杂的卫星可以根据部件组成及物理特性将整体分解为多个形状面积不同的部件并分别计算太阳辐射压加速度，然后，进行矢量相加获得所有部件各个表面的太阳辐射压之和，即为该颗卫星的太阳光压摄动加速度，可写为：

7、

8、当且仅当cosθi≥0时

9、式中，p为卫星处的太阳辐射流量，取p＝4.56×10-6n/m2；

10、m为卫星的质量；

11、λ为卫星的蚀因子；

12、ai为平面i的面积；

13、αi为平面i的方向因子，cosθi<0时为0，cosθi>0时为1；

14、分别为平面i的法向矢量和卫星到太阳的方向矢量；

15、θi为平面i的法向与卫星到太阳方向之间的夹角；

16、ρi,δi分别为平面i的反射系数和散射系数。

17、其中，卫星的蚀因子λ计算方法为：

18、

19、

20、此公式能获得待分析卫星受到的太阳光压加速度，代表了待分析卫星在轨道运行过程中整颗星所受太阳光压加速度大小。

21、(4)输入步骤(2)获取的卫星各部件信息，根据步骤(3)的公式可输出整星太阳光压加速度，即为考虑各个部件表面材料物理特性信息建立的太阳光压摄动模型，将其作为卫星动力学模型的一个摄动项，引入slr数据处理，可用于slr非球形非匀质卫星精密定轨、eop和地球参考架的精密确定等。

22、其中slr精密定轨详细步骤为：

23、卫星在地心天球参考坐标系(gcrs)中的受摄运动方程为：

24、

25、式中，gme为地球引力常数；

26、r为卫星到地球距离；

27、(或)、分别为卫星在gcrs中的位矢、速度与加速度；

28、为卫星所受的摄动加速度，包括n体摄动、固体潮、地球非球形引力、太阳辐射压摄动(即为步骤(3)中asolar)、地球辐射压摄动(地球反照辐射和红外热辐射)、相对论摄动力、大气阻力摄动、经验力摄动、其他摄动；

29、t为时刻；

30、上标t为矩阵转置，下标0为初始时刻；

31、为待估计的动力学模型参数矢量，即可取它可以包括引力场系数、大气阻力系数、太阳光压系数等，如对光压系数cr、大气阻力系数cd进行估计。

32、不出现在卫星运动方程中的待估参数，称为几何参数，记为可以包括台站坐标、地球自转参数等，如考虑估计地球自转参数极移xp0，yp0，和日长变化dr，可取为：

33、

34、其中xp0，yp0，满足：

35、

36、对言，显然有以下微分方程：

37、

38、定义状态矢量，令：

39、

40、状态方程定义为：

41、

42、式中，状态函数为

43、如果状态真值x(t)与参考状态x*(t)在一定的时间间隔内(t0≤t≤t1)足够接近，可将其在参考状态处展开，即每个点处的非线性问题可线性化展开为：

44、

45、省略二阶及以上小量，设x(t)为参考状态x*(t)的偏差，即满足那么

46、

47、未知状态量x0被状态偏差xk代替，同时非线性的参数估计转换为线性状态参数估计。设l为观测历元数，n为参数个数，m为观测总数，初始历元的轨道偏差x0作为xk，那么，为了将线性估计系统中的l×m+n个未知数转化为m+n个，需要将所有的状态偏差xi表达成在某一特定历元的状态偏差xk的函数，即可通过状态转移矩阵将状态偏差转换到某一特定历元：

48、xi＝φ(ti,tk)xk

49、其中，φ(ti,tk)为状态转移矩阵，具有如下性质：

50、(1)φ(ti,ti)＝φ(tk,tk)＝i

51、(2)φ(ti,tk)＝φ(ti,tj)·φ(tj,tk)

52、(3)φ(ti,tj)＝φ-1(tj,ti)

53、这样，状态转移矩阵的微分方程可写为：

54、

55、通过数值积分方法可得到φ(ti,tk)的解，这样方程xi＝φ(ti,tk)xk可写为：xi＝φ(ti,t0)x0。利用观测数据对卫星初始轨道和有关参数进行迭代求解，达到轨道改进的目的。

56、slr观测方程中，令:y＝ρ0,g(x,t)＝ρc，ρ0,ρc分别为卫星到slr测站的观测距离，理论距离，则观测与状态的关系式可表示为：

57、yi＝g(xi,ti)+εi,i＝1,...,l

58、式中，xi、yi、εi分别为ti时刻的状态矢量、观测量和以观测噪声。

59、对观测量g(xi,ti)进行线性化：

60、

61、其中，g(xi,ti)*是先验值，省略二阶及以上小量即为：

62、

63、线性化后的观测方程可写为：

64、令则最终的观测方程被表达为：

65、yi＝hix0+εi

66、最终，通过上述一系列转换，将一个动态系统的状态参数确定转换成对初始历元未知参数的估计问题，即通过不断的观测数据对初始轨道以及未知状态量进行改进，通过不断迭代不断逼近真实轨道。

67、对于观测方程yi＝hix+εi：

68、令

69、通过最小二乘法进行求解，可得到状态偏差的最佳估值

70、

71、对每个观测值进行定权后，令权阵为w，最终的解可写为：

72、

73、其中：

74、

75、若给定一组初始条件x*(t0)，给定先验的和协方差阵以及观测值权阵w，使用批处理算法进行计算的解为：

76、

77、令则协方差矩阵p为:

78、

79、这样就获得了最终状态矢量的改正数及其协方差矩阵。若求得了某一时刻tk弧段的对于等待处理的下一时刻的tk+1弧段，可将利用状态转移矩阵传播到tk+1，作为先验值参与下一弧段的轨道改进，公式为：

80、

81、这样通过迭代就可获得最优轨道。继而完成轨道确定。

82、本发明专利与现有技术相比的有益效果是：本专利通过获取待分析卫星的各表面详细参数信息，根据box-wing太阳光压分体建模的基本原理，对各表面太阳光压辐射力的矢量进行合成，获得整颗星的太阳光压辐射加速度。该方法建立的slr复杂卫星或者非匀质材料卫星高精度太阳光压模型，能够弥补传统slr精密定轨中太阳光压摄动模型当做球模型的不足，进一步提高slr非球形非匀质卫星精密定轨的精度，也可将非球形非匀质卫星纳入slr时空基准确定中，打破目前仅由5颗匀质球形卫星进行slr时空基准确定的局面，拓展了高精度的slr时空基准源，可以提高多种空间大地测量技术进行时空基准确定的精度和稳定性。