《测绘学报》
构建与学术的桥梁 拉近与权威的距离
梁伟1, 徐新禹1,2
1. 武汉大学测绘学院, 湖北 武汉 430079;
2. 武汉大学地球空间环境与大地测量教育部重点实验室, 湖北 武汉 430079;
3. 国家测绘地理信息局卫星测绘应用中心, 北京 100048
收稿日期：2017-05-22；修回日期：2017-10-28
基金项目：国家自然科学基金（41774020；41210006；41404020）；国家高分专项高分遥感测绘应用示范系统项目（AH1701-4）
第一作者简介：梁伟(1991—)，男，博士生，研究方向为地球重力场模型的构建。E-mail: wliang@whu.edu.cn
通信作者：徐新禹, E-mail: xyxu@sgg.whu.edu.cn
摘要：本文研究了联合卫星观测数据和重力异常数据确定超高阶重力场模型的理论方法，并使用EGM2008模型重力异常和GOCE（gravity field and ocean circulation explorer）观测数据构建了重力场模型SGG-UGM-1。重点研究了由球面格网重力异常快速构建超高阶重力场模型的块对角最小二乘方法，将OpenMP技术引入到块对角最小二乘中以提高计算效率，并基于模拟数据验证了方法及算法和软件模块的正确性。采用本文制定的联合解算策略，利用GOCE重力卫星观测数据构建的220阶次法方程和EGM2008模型重力异常构建的2159阶次块对角法方程，联合求解了2159阶次的重力场模型SGG-UGM-1。将SGG-UGM-1与EGM2008、EIGEN-6C2、EIGEN-6C4等超高阶模型在频谱域内进行了比较分析，结果表明SGG-UGM-1相对参考模型的系数误差较小，且在220阶次内的系数精度相比EGM2008模型有了提高。采用中国与美国的GPS/水准数据和毛乌素测区的航空重力观测数据对这些模型进行了外符合精度的检验。检核结果表明，在中国区域，SGG-UGM-1模型大地水准面的精度在EIGEN-6C2和EIGEN-6C4两个模型之间，优于GOSG-EGM模型和EGM2008模型，与美国区域几个模型的精度相当。利用毛乌素测区的航空重力数据对几个模型进行了检核，结果表明SGG-UGM-1模型计算的重力扰动精度与EGM2008、EIGEN-6C4模型相当，优于GOSG-EGM模型和EIGEN-6C2模型。
The Determination of an Ultra-high Gravity Field Model SGG-UGM-1 by Combining EGM2008 Gravity Anomaly and GOCE Observation Data
LIANG Wei1, XU Xinyu1,2, LI Jiancheng1,2, ZHU Guangbin3
Abstract: The theory and methods of the determination of an ultra-high gravity field model by combination of satellite observation data and gravity anomaly data are studied.And an ultra-high gravity field model named SGG-UGM-1 is computed using EGM2008 derived gravity anomaly and GOCE observation data.The block-diagonal least squares (BDLS) method for quickly estimating an ultra-high gravity field model is researched and the corresponding software module is validated by numerical experiments.OpenMP technique is introduced into the BDLS program, which improves computing efficiency dramatically.An ultra-high gravity model SGG-UGM-1 complete to degree and order 2159 is derived using the proposed calculation strategies.The fully occupied normal equation system up to degree and order 220 formed by GOCE satellite data and the block-diagonal normal equation system up to degree and order 2159 formed by EGM2008 gravity anomaly data are used for the combination.Comparison among the models SGG-UGM-1 and EGM2008, EIGEN-6C2, EIGEN-6C4, GOSG-EGM in frequency domain has been done, which shows that SGG-UGM-1 is close to the reference models and that the coefficients lower than degree 220 of SGG-UGM-1 are more accurate than that of EGM2008.The models are also validated by GPS-leveling data in China and America and airborne gravity data in Maowusu surveying area.The results show that in China the accuracy level of SGG-UGM-1 derived geoid is between EIGEN-6C2 and EIGEN-6C4, and better than GOSG-EGM and EGM2008, while in America they are almost the same.In Maowusu area the accuracy level of SGG-UGM-1 derived gravity disturbance is almost at the same accuracy level with EGM2008 and EIGEN-6C4 and better than GOSG-EGM and EIGEN-6C2.
Key words: SGG-UGM-1 ultra-high gravity field model block diagonal least squares method OpenMP parallel computing technique
构建超高阶地球重力场模型一直是物理大地测量学领域的热点研究问题。采用重力梯度测量技术的GOCE重力卫星于2009年成功发射，仅利用其观测数据恢复的全球重力场模型阶次达到了280阶[1]，模型中长波分量的精度更是有了2个数量级的提高，尤其在地面重力数据空白区，相比由GRACE和地面观测数据联合反演的EIGEN-5C等模型，在100 km空间分辨率上，GOCE任务反演的模型精度有了明显提升[2]。这一进展使得联合GOCE卫星观测数据和地面观测数据获得全球更高精度的超高阶重力场模型成为可能。反演的超高阶重力场模型可用于区域高精度大地水准面的确定[3]、全球高程基准的统一[4]、地球内部结构反演[5]等科学研究及应用。
基于球面重力异常数据反演超高阶全球重力场模型常用的方法主要有数值积分(numerical quadrature，NQ)方法、最小二乘配置(least squares collocation，LSC)方法和最小二乘(least squares，LS)方法[6]。NQ方法是基于面球谐函数的正交性，采用积分方法独立求解模型的每个系数，其计算简单，但无法直接给出模型系数的验后方差。相比NQ方法，LS方法不仅可评估模型系数的精度，而且可以联合解算多类观测数据。同时，在LS方法中，如果位系数按照特定的顺序排列，法方程矩阵会具有块对角特点，因此可以分块求解位系数，这种方法称为块对角最小二乘(block diagonal least squares，BDLS)方法[6]。该方法可避免超高维矩阵的构建及求逆，计算量大大减少，已成为目前求解超高阶重力场模型的主要方法，国际上精度得到认可的EGM96[7]、EGM2008[8]、EIGEN-6C4[9]等模型在构建时均使用了该方法。文献[10]分析了BDLS方法在重力场模型确定中的应用，并通过模拟试验比较了BDLS方法与NQ方法求解360阶重力场模型的精度。
考虑到卫星重力观测数据和重力异常数据在反演全球重力场模型时存在频谱互补性，前者对重力场的中长波信号更为敏感，后者对短波及甚短波信号更为敏感，因此联合卫星重力数据和重力异常数据可以反演得到高精度高分辨率的重力场模型，这也是当前求解超高阶重力场模型的主要方法。EGM2008模型构建时联合了GRACE卫星观测法方程和重力异常数据[8]，EIGEN-6C4模型不仅使用了这些数据，还使用了GOCE观测数据，相比EGM2008模型在一些区域的精度有明显提高[9]。文献[11]联合EGM2008模型和GOCE-TIM-R5模型解算了2160阶次重力场模型GECO；文献[12]采用全矩阵求逆的最小二乘方法解算了720阶次的GOCO05c模型，模型精度在15′的空间分辨率上(720阶次)与EIGEN-6C4相当；文献[13]基于卫星、陆地、海洋等观测数据，使用球谐分析方法确定了2159阶次的UGM08模型；文献[14]采用谱权组合方法联合EIGEN-6S卫星重力场模型和格网重力异常构建了2160阶次的重力场模型。
本文将研究BDLS方法的基本原理，并将OpenMP并行计算技术引入到BDLS方法的求解中，研制相应的软件模块，并基于数值模拟试验，分析BDLS方法求解2160阶重力场模型的精度及软件模块的可靠性。然后联合由GOCE任务独立建立的卫星观测法方程和EGM2008模型格网重力异常法方程，采用最小二乘方法估计一个2159阶次的重力场模型，并利用独立的GPS/水准数据对模型精度进行分析与检验。
1 基本原理1.1 由重力异常数据快速求解模型系数的BDLS方法
格网平均重力异常观测值的观测方程为[6]
(1)
式中，Δσi=Δλ(cos θi-cosθi 1)；
