深空探测器三程多普勒建模与算法实现

曹建峰，黄 勇，刘 磊，张 宇，郑爱武，胡松杰

(1. 航天飞行动力学技术重点实验室，北京100094；2. 北京航天飞行控制中心，北京100094；3. 中科院上海天文台，上海 200030)

深空探测器三程多普勒建模与算法实现

曹建峰1,2,3，黄 勇3，刘 磊1,2，张 宇1,2，郑爱武1,2，胡松杰1,2

(1. 航天飞行动力学技术重点实验室，北京100094；2. 北京航天飞行控制中心，北京100094；3. 中科院上海天文台，上海 200030)

针对高精度深空探测多普勒应用，提出了一种间接使用中心天体、航天器的速度与加速度进行多普勒建模的方法。避免了计算机字长截断误差、星历表插值误差对建模精度的影响。使用该算法对环火星探测器进行了连续1周的多普勒计算，通过常规的距离差分方法检验了算法的正确性，对计算的1小时多普勒数据进行10阶多项式拟合，拟合后数据残差小于0.002mm/s. 分析结果表明该算法计算多普勒结果正确，计算精度较常规距离差分算法在精度上得到1个量级的提升。

深空探测；多普勒；观测建模

0 引 言

深空网(Deep space network, DSN)在深空探测任务中发挥着至关重要作用，测量类型主要包括多普勒测速与测距[1-2]。测距数据因为不可避免的存在系统性偏差，环绕飞行阶段的航天器测轨工作更依赖于高精度的DSN多普勒测量。

美国推进喷气实验室(Jet propulsion laboratory, JPL)于1958年负责筹建了首个深空网，完成了早期的深空探测任务。近60年深空探测任务的不断发展使得深空网测量能力与精度得到了大幅提升。早期Mariner 9、Viking 1-2的测轨数据主要为S波段双程多普勒，60s积分双程多普勒测量的噪声水平约为1mm/s[3]。欧空局(European space agency, ESA)发射火星探测器——火星快车(Mars express, MEX)的测轨主要依赖于S和X波段的双程、三程多普勒测量，其1s积分周期的S和X波段多普勒测量噪声水平分别约为1.2mm/s和0.2mm/s[4]。日本臼田站在试验跟踪中的多普勒精度优于0.1mm/s[5]。嫦娥二号(Chang’E-2, CE-2)小行星探测试验后期的测轨完全依赖于DSN和中国VLBI网，受限于测控频段以及星载转发设备的能力，在小行星探测试验期间，CE-2在S频段的测量精度约0.7mm/s(1 s积分周期)[6]，郑为民在对火星快车跟踪试验中获取了0.17mm/s的多普勒测量精度[7]，这一状况在未来的自主火星探测任务必将得到改善，测量精度有望达到ESA水准。

使用多普勒测量数据对航天器开展精密定轨计算工作中，对多普勒测量进行观测建模是轨道计算的基础，通常建模精度需高于测量数据精度1个量级以上时可以忽略建模引入的误差。目前对于多普勒的建模主要使用距离差分建模，对于深空跟踪数据在算法实现时必须采用四精度浮点数。哥达德分析中心指出双精度浮点数引入的截断误差，大行星位置计算误差都是制约多普勒建模精度的重要因素。对于积分周期较大的多普勒建模，如60积分周期，可以采用距离差分建模，文中称为常规算法。对于较小的积分周期，必须在算法实现上进行优化处理。GEODYN软件针对深空多普勒建模问题，在行星历表位置计算中重新对插值引数进行了计算[8]。欧空局在对多普勒观测建模进行了算法上的改进，基于泰勒展开，但是仍然必须依赖四精度浮点数，其建模精度可以优于0.01mm/s。论文对多普勒建模算法进行了进一步改进，使用双精度浮点数定义可以实现优于0.01mm/s的建模精度。

1 三程多普勒测量模型

根据无线电信号上下行链路构型的差异，地基测站对航天器的跟踪模式可以分为单程、双程和三程模式[9]。图1给出了不同测量模式下的跟踪示意图。单程多普勒测量模式通过星载晶振产生下行信号，由地面台站直接接收跟踪信号，该模式只有下行链路，其测量精度受限于星载晶振的稳定性，因此目前很少在实际深空探测任务中使用，我国“萤火一号”火星探测任务原计划采用该测量模式[10]。双程模式由地面台站向探测器发射上行信号，星上转发器接收上行信号，并产生相干的下行信号，再由地面跟踪站接收。三程模式与双程模式相类似，但是上行发射站与下行接收站不同，属于开源测量。在观测建模算法中，三程测量与双程测量的原理相同，仅仅是上下行站不同，因此本文仅针对三程多普勒测量进行建模分析，将上下行站设置为相同则得到双程多普勒建模。

Moyer给出了三程多普勒计算公式[11]，

(1)

(2)

式中：fR为接收频率，fS为发射频率，M为星上转发比，c表示光速，ΔT为接收端的积分周期，UTC为协调世界时。下标s,e分别表示积分周期的起始与结束，下标1表示上行站发射信号，2表示卫星转发信号，3表示下行站接收星上转发信号，12表示信号上行路径，23表示信号下行路径(下文中的下标意义相同)。

为便于计算，对式(1)～(2)进行变换，三程多普勒可以表示为，

[(TDB1-UTC1)e-(TDB1-UTC1)s]-

(3)

式中：TDB表示质心动力学时，大括号内第1项表示1个积分周期内信号传输的光行时差异，第2、3项为时间系统的差异。

2 三程多普勒测量模型算法实现

受限于计算机字长、行星历表精度，当积分周期较小时(如积分周期为1s)，通过式(3)直接计算深空探测器多普勒精度较低，无法满足优于0.01mm/s的测量需求。通常做法是使用四精度浮点数定义来解决计算机字长截断问题，但仍不足以克服星历表插值误差，另外四精度浮点数较双精度浮点数在数值计算效率上会大大降低。JPL在建立行星历表DE414时发现使用双精度与四精度混合定义会比使用双精度定义的效率降低30倍[12]。

根据式(3)，将多普勒计算步骤分2部分完成：首先计算第一部分的光行时，然后计算第二部分的时间系统差异，并通过算法改进提高计算精度。

2.1 光行时的计算

式(3)中大括号内第1项根据其物理意义可以表示为，

(TDB3-TDB1)e-(TDB3-TDB1)s=

[(l23)e-(l23)s]+[(l12)e-(l12)s]+

[(T23)e-(T23)s]+[(T12)e-(T12)s]

(4)

式中：l表示信号传输的直线距离，T表示引力延迟。积分开始时刻的直线距离可以表示为，

(5)

式中：rE表示信号接收时刻(发射时刻)地球在质心天球参考系中的位置矢量，rsta表示该时刻测站在天球参考系位置矢量，rP表示信号转发时刻中心天体在质心天球参考系中的位置矢量，rsat表示卫星相对于中心天体的质心天球参考系位置矢量，式(5)的计算需要进行光行时迭代。

积分结束时刻的距离为，

(6)

式中：ΔrE,ΔrP表示表示地球、中心天体从积分开始至积分结束在质心天球参考系中的位置矢量变化，Δrsta表示测站在一个积分周期内相对于地球的位置矢量变化，Δrsat表示航天器在1个积分周期内相对于中心天体的位置矢量变化。

标记，

Δr= (ΔrE+Δrsta)-(Δrsat+ΔrP)=

(7)

式中：x,y,z为r的坐标分量，Δx,Δy,Δz为Δr的坐标分量。则le-ls在r处可做泰勒展开，丢弃高阶项有，

(8)

为确保计算精度，式(7)中ΔrE与ΔrP不宜直接使用行星历表计算，而使用质心天球参考系中地球与中心天体的速度计算，

(9)

引力延迟部分可以表示为，

(10)

2.2 时间系统差异的计算

式(3)中大括号内第2、3项表示一个积分周期内TDB与UTC差的变化，可以分解为,

[(TDB-UTC)e-(TDB-UTC)s]=

[(TDB-TT)e-(TDB-TT)s]+

[(TT-TAI)e-(TT-TAI)s]+

[(TAI-UTC)e-(TAI-UTC)s]

(11)

TDB与TT的转换关系一般表示为下式，

TDB-TT=wf+wj+wt

(12)

式中：wf为Fairhead模型，wt为与地理位置相关的项，wj为基于行星历表给出的改正项。

经分析，行星历表改正项对多普勒建模的影响小于10-7m/s，可以不予考虑。Fairhead模型项展开表达式包含787项，

(13)

ΔT表示从J2000.0起算的时间，单位为1000年，

(14)

与地理位置相关的项，

wt= 2.9×10-14usin(UT1+λ+L-LSa)+

1.0×10-13usin(UT1+λ-2M)+

1.33×10-13usin(UT1+λ-D)+

专名屈折变化之非本真性的意义在于殊相相对于语言的独立性：虽然我们对语言的不同使用方式创造了各种不同的对象，然而在某种意义上，殊相可以被视为语言不得不原封不动地接纳的前语言“自在客体”，而作为其语言上的对应物，本真专名在充当句子之专名时在很多语言中都不经历变化。

1.33×10-13usin(UT1+λ+L-LJ)-

2.29×10-13usin(UT1+λ+2L+M)-

2.2×10-12vcos(L+M)+

5.312×10-12usin(UT1+λ-M)

-1.3677×10-11usin(UT1+λ+2L)-

1.31840×10-10vcos(L)+

3.17679×10-10usin(UT1+λ)

(15)

式中：M表示地月系质心轨道运动的平近点角，UT1为世界时，λ表示测站地理经度，D为日月平角距，L为太阳平黄经，LJ为木星的日心平黄经，LSA为土星的日心平黄经。u是测站距地球自转轴的距离(单位km)，v是测站距赤道平面北向距离。

时间系统差异对积分周期的导数则可以表示为

(16)

在1950-2050时间范围内，考虑0.003ps/s精度需要(对应于多普勒约0.001mm/s)，式(16)需要考虑wf的部分系数(A系数170项，B系数4项，C系数1项)，wt的全部项，而wj完全不需要考虑，A,B,C系数参见SOFA文档。

2.3TDB至UTC的迭代计算

式(11)的计算中涉及到TDB至UTC的转换，而SOFA标准算法仅提供了TT计算TDB的算法，其说明是可以使用TDB替代TT计算TDB-TT，因而严格来讲，TDB至UTC需要进行迭代计算，具体迭代步骤如图2。

图3给出了不迭代计算与迭代1次计算的差异，可以看到两者的差异小于0.5ps，且主要为周年变化，属于长周期项，对于短时间的多普勒计算，该差异可以忽略，因此在实际算中，可以直接使用TDB代替TT计算，该简化引入的误差对计算精度的影响可以忽略。

至此，利用式(3)、 (4)、 (11)可以实现多普勒的理论值计算。

3 精度验证

使用常规的距离差分的方法已经对观测建模的正确性进行了检验[13-14]，本节的计算主要验证本文方法计算的精度。这里的常规算法是通过距离差分实现，其优点是原理与算法简单，但由于航天器与测站距离的直接使用，计算中截断误差直接影响建模的精度；另外，采用JPL星历表插值计算大行星位置时，其精度仅能达到0.05mm，该误差并不能通过距离差分而完全消除。本文提供的算法相对复杂，避开了直接使用大行星位置，间接使用了航天器相对于中心天体的速度与加速度。

以火星快车轨道为例，基于双精度浮点数，分别使用常规数值算法与本文提供的算法进行多普勒测速的计算，计算的弧长为7天，考虑积分周期为1s，为了便于比较，不考虑测站与目标之间的遮挡关系，假设卫星全弧段可见，然后比较两者之间的差异(图4)。两者的一致性非常高，绝对差异小于0.2mm/s，且差异表现为噪声。

理论上每个观测弧段内多普勒测速值应该是光滑连续的，因此使用多项式拟合后的残差可以在一定程度上反映多普勒测速建模的误差，分别对两种算法计算的结果选取1小时观测弧段数据进行拟合。图5(a)给出了拟合示意图，使用的是10阶多项式，(B)为常规算法拟合后的残差，最大偏差达±0.1mm/s，图5(c)为使用本文算法计算多普勒数据拟合后的残差，最大偏差为±0.002mm/s。比较图5(b)与图5(c)可知,图中算法的差异主要由常规算法计算的误差引起。

4 结 论

论文提出了一种适用于深空探测任务的高精度多普勒建模算法，该算法使用环绕中心天体与航天器的位置、速度、加速度进行多普勒测速的计算，避免了直接使用测距建模引入的截断误差。使用常规距离差分建模算法验证了算法的正确性，通过多项式拟合的方法检验了算法计算的精度。

[1] 李平，张纪生．NASA深空网(DSN)的现状及发展趋势[J]．飞行器测控学报，2003, 22(4): 10-17．[LiPing,ZhangJi-sheng.ThecurrentstatusandthefutureofNASAdeepspacenetwork[J].JournalofSpacecraftTT&CTechnology, 2003, 22(4): 10-17.]

[2] 于志坚，于益农，董光亮，等．欧空局深空网的现状和发展趋势[J]．飞行器测控学报，2004, 12(3)：1-6．[YuZhi-jian,YuYi-nong,DongGuang-liang,etal.StatusandtrendofdevelopmentofESAdeepspacenetwork[J].JournalofSpacecraftTT&CTechnology, 2004, 12(3): 1-6.]

[3]BornGH,ChristensenEJ,FerrariAJ,etal.ThedeterminationofthesatelliteorbitofMariner9[J].CelestialMechanics, 1974, 9: 395-414.

[4]RosenblattP,LaineyV,MaistreSL,etal.AccurateMarsexpressorbitstoimprovethedeterminationofthemassandephemerisoftheMartianmoons[J].PlanetaryandSpaceScience, 2008, 56(7): 1043-1053.

[5]KikuchiF,LiuQ,HanadaH,etal.PicosecondaccuracyVLBIofthetwosubsatellitesofSELENE(KAGUYA)usingmultifrequencyandsamebeammethods[J].RadioScience, 2009, 44(2): 413-421.

[6] 曹建峰，胡松杰，刘磊，等．“嫦娥二号”小行星探测试验定轨计算与分析[J]．北京航天航天大学学报，2014, 40(8): 1095-1101．[CaoJian-feng,HuSong-jie,LiuLei,etal.OrbitdeterminationandanalysisforChang'E-2asteroidexploration[J].JournalofBeijingUniversityofAeronauticsandAstronautics, 2014, 40(8): 1095-1101.]

[7] 郑为民，马茂莉，王文彬．深空探测器被动式高精度多普勒测量方法与应用[J]. 宇航学报, 2013, 34(11): 1462-1467. [ZhengWei-min,MaMao-li,WangWen-bin.High-precisionpassivedopplermeasurementmethodanditsapplicationsindeepspaceexplorater[J].JournalofAstronautics, 2013, 34(11): 1462-1467.]

[8]RowlandsDD,MarshallJA,MccarthyJ,etal.GeodynIIsystemdescription.volume1 [R].HughesSTXCorp.Greenbelt,Maryland, 1993.

[9]JianN,ShangK,ZhangS,etal.Adigitalopen-loopdopplerprocessingprototypefordeep-spacenavigation[J].ScienceinChinaSeriesG:Physics,Mechanics&Astronomy, 2009, 52(12): 1849-1857.

[10] 吴季，朱光武，赵华，等．萤火一号火星探测计划的科学目标[J]．空间科学学报，2009, 29(5): 449-455．[WuJi,ZhuGuang-wu,ZhaoHua,etal.OverviewofscientificobjectivesofChina-RussiajointMarsexplorationprogramYH-1 [J].Chin.J.SpaceSCI, 2009, 29(5): 449-455.]

[11]MoyerTD.Formulationforobservedandcomputedvaluesofdeepspacenetworkdatatypesfornavigation[M].California:JetPropulsionLaboratory, 2000.

[12]StandishEM.JPLplanetaryephemerisDE414，IOM343R-06-002[R].California:JetPropulsionLaboratory, 2006: 1-8.

[13]CaoJF,HuangY,HuXG,etal.MarsexpresstrackingandorbitdeterminationtrialswithChineseVLBInetwork[J].ChineseScienceBulletin, 2010, 55(32): 3654-3660.

[14] 曹建峰，黄勇，胡小工，等．深空探测中多普勒的建模与应用[J]．宇航学报，2011, 32(7): 1583-1589． [CaoJian-feng,HuangYong,HuXiao-gong,etal.Modelingandapplicationofdopplerdataindeepspaceexploration[J].JournalofAstronautics, 2011, 32(7): 1583-1589.]

通信地址：北京市5130信箱120分箱(100094)

电话：13811206835

E-mail:jfcao@foxmail.com

胡松杰(1973-)，男，博士，研究员，主要从事轨道动力学研究。本文通信作者。

通信地址：北京市5130信箱120分箱(100094)

E-mail:husongjie@aliyun.com

(编辑：张宇平)

Modeling and Algorithm Realization of Three-Way Doppler for Deep Space Exploration

CAO Jian-feng1,2,3, HUANG Yong3, LIU Lei1,2, ZHANG Yu1,2, ZHENG Ai-wu1,2, HU Song-jie1,2

(1. Science and Technology on Aerospace Flight Dynamics Laboratory, Beijing 100094, China; 2.Beijing Aerospace Control Center, Beijing 100094, China; 3. Shanghai Astronomical Observatory, Shanghai 200030, China)

Aiming at the high precision Doppler application for deep space exploration, a new algorithm for Doppler modeling using the velocity and acceleration of the central object and the spacecraft is proposed. The algorithm avoids the truncation error induced by the direct use of ranging and the positional uncertainty of ephemerides interpolation. Continuous Doppler computation for a Martian orbiter lasting one week is performed, and the general computation algorithm using ranging difference is used to test the correctness of the new algorithm. In addition, one-hour length of the Doppler data is fitted with a polynomial of order 10, and the residual error is less than 0.002mm/s. The analysis results show that the new algorithm is correct and is approximately one order of magnitude improvement in accuracy compared with the conventional modeling algorithm.

Deep space exploration; Doppler; Observation modeling

2016-08-01;

2016-11-21

国家自然科学基金(11203003, 61573049, 11373013, 61304233, 11473056)

V412.4+1

A

1000-1328(2017)03-0304-06

10.3873/j.issn.1000-1328.2017.03.011

曹建峰(1982-)，男，博士，工程师，主要从事轨道动力学研究。