2. 中国科学院大学, 北京 100049
2. University of Chinese Academy of Sciences, Beijing 100049, China
GRACE (Gravity Recovery And Climate Experiment)卫星自2002年发射至今,在恢复地球重力场的精度和空间分辨率方面均有巨大进展,一直是卫星重力领域的研究热点.GRACE卫星对地球重力场的巨大贡献得益于星上同时装载了星载GPS接收机、K波段测距系统(KBR)、加速度计以及星敏感器等主要载荷.其中,K波段测距系统是恢复地球重力场模型的核心载荷,用于测量两颗卫星之间的距离及距离变率;星载GPS接收机主要用于确定GRACE卫星精密轨道;加速度计用于高精度测量卫星的非保守力;星敏感器用于确定卫星的精确姿态.以上载荷提供高低卫卫跟踪和低低卫卫跟踪观测数据,成功应用于反演静态地球重力场和探测全球质量变化信息.基于GRACE观测资料的时变重力场模型在全球、区域水储量变化研究,极区冰盖质量变化研究等方面取得了丰硕的成果(Chambers et al., 2004; Chen et al., 2006, 2009),对地球物理学、大地测量学、海洋学和冰川学等领域的发展具有重要的意义.随着我国低低卫卫跟踪重力卫星计划的提出,自主获取高精度的时变地球重力场模型的相关研究势在必行.
如何从GRACE卫星观测数据中提取出所需的地球质量变化信息,是卫星重力领域的研究重点.基于重力卫星观测数据恢复地球重力场模型主要有以下方法:能量法、加速度法、短弧法、天体力学法和动力学积分法.到目前为止,国际上主要有美国德克萨斯大学空间研究中心CSR、德国地学研究中心GFZ、美国宇航局喷气推进实验室JPL、法国国家空间局CNES、德国波恩大学理论大地测量研究所ITG、荷兰代尔夫特地球观测与空间系统研究所DEOS以及瑞士伯尔尼大学天文研究所AIUB等机构,分别采用不同方法成功反演时变地球重力场模型(Mayer-Gürr, 2006; Bruinsma et al., 2010; Liu et al., 2010; Meyer et al., 2012; Bettadpur et al., 2012; Watkins and Yuan, 2012; Dahle et al., 2013).国内同行也为GRACE模式下利用实测数据反演地球重力场模型做了大量的工作.目前采用短弧法和动力学法得到与国际上相同精度水平的时变重力场模型,如冉将军等(2014)、Chen等(2015)和Shen等(2015)采用短弧法,利用KBR星间观测数据结合精密轨道成功反演了时变重力场模型;王正涛(2005)、周旭华(2005)、肖云(2006)、张兴福(2007)、王庆宾(2009)等人联合GRACE精密轨道和KBR星间观测数据,利用动力学积分法反演地球重力场的工作,取得了显著的研究进展.Zhao等(2010)和Wang等(2015)采用动力学积分法,利用KBR星间测速数据与低轨卫星精密轨道数据成功恢复了时变地球重力场模型.
在恢复地球重力场模型的方法中,能量法、加速度法和短弧法均需要已知重力卫星的精密轨道.基于GRACE卫星观测数据的动力学积分法可分为同解法和两步法,两步法是利用GRACE卫星精密轨道作为伪观测值,同时结合KBR星间观测数据反演地球重力场模型,国内相关机构已采用两步法成功恢复时变地球重力场模型(Zhao et al., 2010; Wang et al., 2015);同解法(Zou et al., 2010)即直接采用原始观测数据同时确定GRACE卫星精密轨道和地球重力场模型.同解法使观测数据和观测模型达到更大程度的同化,理论上更加严密,实现更优的参数估计.目前只有少数国际知名机构(如GFZ、CSR、CNES等)利用该方法成功反演地球重力场模型,因此基于星载GPS和KBR星间观测数据采用同解法恢复时变地球重力场模型的相关研究对我国发展自主卫星重力测量系统以及卫星重力模型的研制理论和实际应用都有积极的意义.
本文利用同解法反演月尺度时变重力场模型,即采用GRACE星载GPS观测数据和KBR星间测速观测数据同时求解GRACE卫星精密轨道和60阶月尺度地球重力场模型,成功提取地球重力场的时变信号,并详细分析了观测数据拟合残差、轨道确定精度、地球重力场时变信号以及特征区域水储量变化等计算结果.
2 基于动力学法反演地球重力场数学模型GRACE卫星装载了星载GPS接收机和K波段测距系统,同时获取星载GPS观测数据和KBR星间测速观测数据.星载GPS观测数据用于确定低轨卫星精密轨道并恢复地球重力场模型,仅采用星载GPS观测数据在恢复中高阶的重力场精度具有一定的局限性,通过结合KBR星间测速观测数据可以提高中高阶部分的精度.本文采用动力学积分法结合两类原始观测数据实现同时确定GRACE卫星的精密轨道和地球重力场模型.
2.1 星载GPS观测数据的观测方程基于动力学积分法反演地球重力场模型和精密定轨,根据卫星的运动方程将卫星的瞬间状态矢量表示为初始状态矢量和重力位系数及其他力学模型参数的函数,具体表示如下:
(1) |
式中:
对(1) 式用泰勒级数展开,则其线性方程表示如下:
(2) |
式中:δXj0, δP, δUj分别表示为各参数的改正量.
星载GPS观测数据包括双频伪距和相位观测值,由于伪距观测相对相位观测精度低,为了减小计算量,本文在恢复重力场时不采用伪距观测值.为了消除电离层影响形成无电离层组合相位观测值,将其作为进行精密轨道确定和反演地球重力场模型的基本观测量.建立基于双频无电离层相位观测值的观测方程,见(3) 式,具体公式推导过程参考郭南男(2015).
(3) |
其中ρoL3为低轨卫星相对于GPS卫星的无电离层组合相位观测值;ρc为低轨卫星相对于GPS卫星的几何距离,
(4) |
由(4) 式,采用最小二乘法即可得基于星载GPS相位观测值的法方程.星载GPS观测数据中涉及到的各类观测误差可通过模型进行改正,最后仅解算模糊度和接收机钟差参数,其中接收机钟差为历元参数.最终法方程中的解算参数包括:初始轨道状态、加速度尺度因子和偏差、重力场模型球谐系数、模糊度参数和接收机钟差.除地球重力场模型参数为全局参数外,其余参数均为局部参数,若用Q表示局部参数,P表示全局参数,可以将法方程简化表示为:
(5) |
KBR星间测距系统的主要观测量是星间距离ρ及其一次变率
采用动力学积分法利用KBR星间测速观测值的具体观测方程如下:
(6) |
星间测速观测值是卫星瞬时状态矢量的函数,表示为:
(7) |
星间测速观测值对于两颗卫星瞬间状态矢量的偏导数如下:
式中
其中
(8) |
联合星载GPS相位观测和KBR星间测速观测法方程,本文根据两类观测数据残差确定两类观测值的权,相位观测设置σPHA=1 cm,星间测速观测设置为σKBRR=0.2 μm·s-1,叠加两类观测的法方程具体表示为:
(9) |
对(9) 式约化局部参数向量,得到仅含地球重力场位系数的法方程:
(10) |
最后的法方程由叠加每一弧段对应的约化法方程形成,求解法方程即可得到重力场位系数,将全局参数回带到每个弧段的法方程中可分别得到两颗低轨卫星的精密轨道.
3 力学模型配置以及参数设置在数据处理过程各类求解参数的具体设置以及力学模型的配置是决定解算精度的关键,不同的研究机构采用的模型配置和参数设置各不相同,本文具体的模型配置和参数设置详见表 1.
加速度尺度因子和偏差参数的设置对GRACE卫星轨道和地球重力场模型的精度有着重要影响.本文采用Bettadpur (2009)计算的加速度尺度因子参考值作为初始值,但由于加速度尺度因子与动力学模型中的CPR (once-per-rev)参数有较高的相关性,因此文中不设置该参数.仅每一小时解算一组加速度偏差参数.
KBR经验参数(也称之为运动学参数),包括常数项、一阶项以及周期项参数,其中周期项参数可避免轨道周期相关的误差,参数设置具体表示为(11) 式,其中a, b每45 min求解一组,c, d每1.5 h求解一组.
(11) |
式中a, b, c, d为求解参数,T为卫星轨道周期,在这里为1.5 h.
4 结果分析为了证明上述算法的可靠性,下面将从以下五个方面对解算结果进行分析,包括相位观测拟合残差、KBRR拟合残差、低轨卫星定轨精度、重力场模型系数阶方差以及反演的陆地水储量变化的空间分布.
4.1 观测数据拟合残差(1) 相位数据拟合残差
基于星载GPS观测数据确定卫星精密轨道和地球重力场模型位系数的整体解算精度,可以通过相位观测值的拟合残差反映.以2007年为例,选取两颗卫星一年的相位观测值以天为单位进行统计.图 1为2007年4月1日一天的GRACE卫星相位残差,相位残差的RMS为5.8 mm,从图中可以看出相位残差分布不存在系统误差,呈现随机分布.图 2为以天为单位统计2007年的相位残差RMS计算结果,从图中可以看出每天的拟合残差的RMS均优于1 cm,计算一年的相位残差RMS的平均值为6.18 mm,与Bruinsma等(2010)计算的相位拟合残差RMS的平均值6.1 mm精度相当.
(2) KBRR观测数据拟合残差
KBRR观测值作为恢复地球重力场模型的关键数据,其拟合残差可以反映解算重力场模型精度.图 3以2007年4月1日为例,给出一天KBRR拟合残差计算结果,拟合残差的RMS为0.18 μm·s-1.
从图 3中可以看出KBRR拟合残差分布不存在系统误差,呈现随机分布.图 4为统计2007年KBRR拟合残差结果,计算一年KBRR拟合残差RMS的平均值为0.201 μm·s-1,与Bruinsma等(2010)计算的KBRR拟合残差的RMS平均值0.199 μm·s-1精度相当,且其变化趋势一致.
轨道误差影响地球重力场模型恢复精度.通过对轨道质量评估,不仅可以验证算法的正确性,同时可以避免较大的轨道误差影响地球重力场模型反演精度.以JPL提供的高精度的GRACE约化动力学精密轨道为标准轨道,将2007年4月1日定轨结果与之比较.图 5表示GRACE-A和GRACE-B在地固坐标系下与JPL发布的约化动力学精密轨道在三个方向的位置之差.GRACE-A在X、Y、Z三个方向位置之差的RMS分别为1.81 cm、2.02 cm、1.40 cm,三维位置之差的RMS为3.05 cm;GRACE-B在X、Y、Z三个方向位置之差的RMS分别为2.02 cm、1.54 cm、1.82 cm,三维位置之差的RMS为3.13 cm,与目前GRACE卫星每个方向的位置精度为2 cm的定轨精度相当(Jaggi et al., 2007;Kang et al., 2009).
通过对比国际上知名机构CSR、JPL、GFZ发布的RL05模型,分析本文解算的60阶月尺度时变重力场模型精度.为了方便表述,本文解算的结果在文中均命名为SHA,而其他三个机构的RL05模型均用机构代号表示,分别表示为CSR、JPL、GFZ.分别计算包括本文在内的4个机构月尺度时变重力场模型的时变信号的阶方差.由于GRACE卫星本身的限制,不能准确地测定C20项,不同机构解算的C20项存在一定的差异.为方便比较,图 6给出2007年4月和10月从3~60阶的时变信号的阶方差结果.已有研究结果表明地球重力场模型的时变信号主要集中在前30阶.从图 6中可以看出,在前30阶部分,包括本文在内的四组模型的计算结果基本一致.随着阶次的增加,噪声会逐步增大.由于不同机构采用的具体数据处理策略不同导致噪声引入程度不同,30阶以后四组模型的时变信号阶方差均出现微小差异.本文结果在高阶部分阶方差略小于其他三组模型的结果,但从总体上本文模型时变信号阶方差与其他三组模型具有很好的一致性.
为了进一步分析本文解算的时变重力场模型对全球质量变化信号的恢复程度,将地球时变重力信号用等效水高表示,并对比CSR、JPL、GFZ发布的RL05模型.具体处理过程如下:首先用卫星激光测距解算的C20项替换时变重力场模型结果,采用P4M6去相关滤波去除同一次的奇数阶和偶数阶位系数的相关性,减少南北条带误差影响(Swenson and Wahr, 2006),并采用500 km滤波半径的Gauss滤波(Jekeli, 1981)对时变重力场模型进行处理,计算月尺度全球水储量信号.图 7和图 8分别表示2007年4月和10月四组模型的全球时变重力场信号空间分布情况.对比CSR、GFZ、JPL模型计算结果,本文模型与其他三个机构反演的陆地水时变信号在全球范围内的空间分布十分接近,具有较高的空间一致性.
为了量化比较本文恢复的月尺度时变地球重力场模型与CSR、GFZ和JPL三个机构发布的RL05的差异,本文针对亚马逊流域和长江流域,分别对包括本文在内的四个机构的时变重力场在上述两个地区质量变化信号进行对比.同样对四个机构的月尺度时变重力场模型进行相同的高斯滤波和去相关处理.亚马逊流域选取[21°S, 5°N],[45°W, 80°W]区域范围,长江流域选取[24°N, 35°N],[90°E, 122°E],分别针对该区域范围的纬度进行加权平均得到该区域的质量变化对应的等效水高.
如图 9所示,在亚马逊流域,四个机构结果都表现出一致的周期性变化,其主要反映了陆地水储量的变化.计算SHA、CSR、GFZ和JPL四个机构2005—2007年亚马逊流域平均水储量变化的周年振幅分别为26.2±3.3 cm、27.2±3.9 cm、27.3±4.0 cm和27.0±3.9 cm等效水柱高,SHA与其他三组产品的相关系数均大于0.98,呈现出较强的相关性.
如图 10所示,在长江流域,2005—2007年SHA得到的等效水柱高与其他三个机构的周年信号基本一致.时变信号对应的等效水高的周期和振幅与CSR、GFZ和JPL这三个机构符合得较好,计算的相关系数分别为0.96、0.94和0.95.证明本文模型与其他三个机构的陆地水时变信号在区域范围的等效水高十分接近,具有较高的一致性.
反演高精度地球重力场模型是卫星重力测量中的核心内容.本文基于卫星动力学方法,采用2005—2007年GRACE卫星的Level 1B资料,在建立卫星重力场和卫星轨道参数联合解算基础上,给出一组完整的60阶月平均地球重力场模型和卫星精密轨道.从观测数据拟合残差、轨道确定精度、重力场模型精度三个方面详细分析解算结果,得到如下结论:
(1) 统计2007年卫星重力场和卫星轨道参数联合解算后的星载GPS相位数据和KBR星间测速数据残差,分别为5~8 mm、0.18~0.30 μm·s-1,与国际上公布的一些高精度结果接近.
(2) 与JPL发布的约化动力学精密轨道相比,本文计算的GRACE卫星轨道的X、Y、Z三个方向精度约为2 cm,三维位置精度优于5 cm,表明本文确定的GRACE卫星轨道与JPL公布的精密轨道精度相当.
(3) 与CSR、GFZ和JPL发布的RL05模型相比,本文解算的60阶月平均地球重力场模型与RL05模型精度相当,且在反演陆地水储量变化的区域等效水高差异很小,并具有很好的空间分布一致性.如亚马逊和长江区域,本文计算2005—2007年水储量变化与CSR、GFZ、JPL的RL05模型结果差异很小.对于亚马逊流域,其相关系数均大于0.98;对于长江流域,其相关系数分别为0.96、0.94和0.95,且年变化周期明显.
综上所述,本文联合GRACE卫星的星载GPS和KBRR星间测速数据,反演的时变重力场模型达到国外最新的时变重力场模型精度,且具有同时获得高精度卫星精密轨道的能力,为我国未来卫星重力反演和轨道确定提供另一重要技术途径.
致谢感谢CSR、GFZ、JPL以及IGS提供的相关数据.
Bettadpur S. 2009. Recommendation for a-priori bias & scale parameters for Level-1B ACC data (version 2). GRACE Technical Note-02 Version 2. | |
Bettadpur S. 2012. Gravity Recovery and Climate Experiment UTCSR Level-2 processing standards document for Level-2 product release 0005. Austin:Center for Space Research, University of Texas at Austin. | |
Bruinsma S, Lemoine J, Biancale R, et al. 2010. CNES/GRGS 10-day gravity field models (release 2) and their evaluation. Adv. Space. Res., 45(4): 587-601. DOI:10.1016/j.asr.2009.10.012 | |
Chambers D P, Wahr J, Nerem R S. 2004. Preliminary observations of global ocean mass variations with GRACE. Geophys. Res. Lett., 31(3). DOI:10.1029/2004GL020461 | |
Chen J L, Wilson C R, Blankenship D D, et al. 2006. Antarctic mass rates from GRACE. Geophys. Res. Lett., 33(11). DOI:10.1029/2006GL026369 | |
Chen J L, Wilson C R, Blankenship D, et al. 2009. Accelerated Antarctic ice loss from satellite gravity measurements. Nat. Geosci., 2(12): 859-862. DOI:10.1038/NGEO694 | |
Chen Q J, Shen Y Z, Zhang X F, et al. 2015. Monthly gravity field models derived from GRACE Level 1B data using a modified short arc approach. J. Geophys. Res. Solid Earth, 120(3): 1804-1819. DOI:10.1002/2014JB011470 | |
Dahle C, Flechtner F, Gruber C, et al. 2013. GFZ GRACE level-2 processing standards document for level-2 product release 0005:revised edition. Potsdam:Deutsches GeoForschungsZentrum. | |
Guo N N, Zhou X H, Wu B. 2015. Research on rapid precise orbit determination of Haiyang-2A using on-board GPS data. Journal of Astronautics, 36(7): 797-803. DOI:10.3873/j.issn.1000-1328.2015.07.008 | |
Jaggi A, Hugentobler U, Bock H, et al. 2007. Precise orbit determination for GRACE using undifferenced or doubly differenced GPS data. Advances in Space Research, 39(10): 1612-1619. DOI:10.1016/j.asr.2007.03.012 | |
Jekeli C. 1981. Alternative methods to smooth the Earth's gravity field. Columbus:The Ohio State University. | |
Kang Z G, Tapley B D, Bettadpur S V, et al. 2009. Quality of GRACE orbits using the reprocessed IGS products.//Proceedings of 2009 AGU Fall Meeting. Moscone Convention Center. | |
Liu X, Ditmar P, Siemes C, et al. 2010. DEOS Mass Transport model (DTM-1) based on GRACE satellite data:methodology and validation. Geophys. J. Int., 181(2): 769-788. DOI:10.1111/j.1365-246X.2010.04533.x | |
Mayer-Gürr T. 2006. Gravitationsfeldbestimmung aus der Analyse kurzer Bahnböegen am Beispiel der Satellitenmissionen CHAMP und GRACE[Ph. D. thesis]. Bonn:Universitat Bonn. | |
Meyer U, Jäggi A, Beutler U. 2012. Monthly gravity field solutions based on GRACE observations generated with the celestial mechanics approach. Earth Planet. Sci. Lett., 345: 72-80. DOI:10.1016/j.epsl.2012.06.026 | |
Ran J J, Xu H Z, Zhong M, et al. 2014. Global temporal gravity field recovery using GRACE data. Chinese Journal of Geophysics, 57(4): 1032-1040. DOI:10.6038/cjg20140402 | |
Shen Y Z, Chen Q J, Xu H Z. 2015. Monthly gravity field solution from GRACE range measurements using modified short arc approach. Geodesy and Geodynamics, 6(4): 261-266. DOI:10.1016/j.geog.2015.05.009 | |
Swenson S, Wahr J. 2006. Post-processing removal of correlated errors in GRACE data. Geophys. Res. Lett., 33(8): L08402. DOI:10.1029/2005GL025285 | |
Wang C Q, Xu H Z, Zhong M, et al. 2015. An investigation on GRACE temporal gravity field recovery using the dynamic approach. Chinese Journal of Geophysics, 58(3): 756-766. DOI:10.6038/cjg20150306 | |
Wang Q B. 2009. Study on the recovery of gravitational potential model with dynamical method[Ph. D. thesis] (in Chinese). Zhengzhou:The PLA Information Engineering University. | |
Wang Z T. 2005. Theory and methodology of earth gravity field recovery by satellite-to-satellite tracking data[Ph. D. thesis] (in Chinese). Wuhan:Wuhan University. | |
Watkins M, Yuan D N. 2012. GRACE JPL level-2 processing standards document for level-2 product release 05. GRACE document 327-744. | |
Xiao Y. 2006. Research on the earth gravity field recovery from satellite-to-satellite tracking data[Ph. D. thesis] (in Chinese). Zhengzhou:School of Geodesy and Geomatics, Information Engineering University. | |
Zhang X F. 2007. The earth's field model recovery on the basis of satellite-to-satellite tracking missions[Ph. D. thesis] (in Chinese). Shanghai:Tongji University. | |
Zhao Q, Guo J, Hu Z G, et al. 2011. GRACE gravity field modeling with an investigation on correlation between nuisance parameters and gravity field coefficients. Adv. Space Res., 47(10): 1833-1850. DOI:10.1016/j.asr.2010.11.041 | |
Zhou X H. 2005. Study on satellite gravity and its application[Ph. D. thesis] (in Chinese). Wuhan:Institute of Geodesy and Geophysics, Chinese Academy of Sciences. | |
Zou X C, Li J C, Jiang W P, et al. 2010. Research on the simultaneous solution method for satellite gravity data analysis and its simulation. Acta Geodaetica et Cartographica Sinica, 39(4): 344-348. | |
郭南男, 周旭华, 吴斌. 2015. 利用星载GPS数据进行海洋2A卫星快速精密定轨. 宇航学报, 36(7): 797–803. DOI:10.3873/j.issn.1000-1328.2015.07.008 | |
冉将军, 许厚泽, 钟敏, 等. 2014. 利用GRACE重力卫星观测数据反演全球时变地球重力场模型. 地球物理学报, 57(4): 1032–1040. DOI:10.6038/cjg20140402 | |
王长青, 许厚泽, 钟敏, 等. 2015. 利用动力学方法解算GRACE时变重力场研究. 地球物理学报, 58(3): 756–766. DOI:10.6038/cjg20150306 | |
王庆宾. 2009. 动力法反演地球重力场模型研究[博士论文]. 郑州: 解放军信息工程大学. | |
王正涛. 2005. 卫星跟踪卫星测量确定地球重力场的理论与方法[博士论文]. 武汉: 武汉大学. | |
肖云. 2006. 基于卫星跟踪卫星数据恢复地球重力场的研究[博士论文]. 郑州: 信息工程大学测绘学院. | |
张兴福. 2007. 应用低轨卫星跟踪数据反演地球重力场模型[博士论文]. 上海: 同济大学. | |
周旭华. 2005. 卫星重力及其应用研究[博士论文]. 武汉: 中国科学院测量与地球物理研究所. | |
邹贤才, 李建成, 姜卫平, 等. 2010. 卫星重力资料分析的同解法研究及其仿真. 测绘学报, 39(4): 344–348. | |