图 2 模型的系数阶误差RMSFig. 2 Degree error RMS of the model coefficients
图选项
3 SGG-UGM-1超高阶重力场模型的求解
EGM2008重力场模型在构建时使用了精度较高且覆盖较全的地面重力数据集[8]，本文使用其计算球面上的全球重力异常格网平均值，并联合该数据和GOCE卫星数据解算了2159阶次重力场模型SGG-UGM-1。
3.1 超高阶重力场模型解算的数据处理流程
本文构建超高阶重力场模型的解算流程如图 3所示，具体描述如下：
图 3 SGG-UGM-1模型的解算流程Fig. 3 Calculation process of SGG-UGM-1 model
图选项
(1) 使用EGM2008模型计算参考球面(球半径为6 378 136.3 m)上5′×5′格网分辨率的重力异常数据ΔgEGM；由该数据使用BDLS方法构建并求解2159阶块对角法方程，得到2~2159阶位系数，这部分系数命名为BDLS，使用BDLS方法构建法方程时认为格网数据的权相等。
(2) 根据步骤(1)求得系数中2~100阶以及251~2159阶系数计算重力异常，然后在数据ΔgEGM中减去这部分重力异常得到残余重力异常数据ΔgR。ΔgR数据中仅有101~250阶系数的贡献，再使用ΔgR数据构建101~250阶系数的块对角法方程。
(3) 使用GOCE卫星的卫星重力梯度(satellite gravity gradiometry，SGG)数据和卫星跟踪卫星(satellite-to-satellite tracking，SST)数据构建220阶卫星观测法方程，并求解得到220阶的卫星重力场模型，其数据处理方法的描述见下节。
(4) 步骤(2)建立了101~250阶系数的法方程，步骤(3)建立了2~220阶系数的法方程，使用1.2节描述的联合平差方法，联合这两个法方程得到2~250阶系数的联合观测法方程。3.3节中有关于法方程联合过程更为详细的描述。
(5) 求解步骤(4)建立的2~250阶系数联合观测法方程，可以得到2~250阶系数，这部分系数与步骤(1)中求解的系数中251~2159阶部分一同组成2~2159阶重力场模型系数。
3.2 GOCE卫星观测法方程的构建
本文使用的GOCE卫星重力观测数据为2009年11月1日—2012年5月31日期间的SGG数据和2009年11月1日—2010年7月5日期间的SST数据，数据的采样间隔为1 s。其中包括：卫星重力梯度观测数据、梯度仪坐标系相对于惯性系之间的姿态数据(四元素)、精密卫星轨道数据、地固系与惯性系之间的转换矩阵、加速度观测数据、运动学轨道的方差-协方差数据等[18]。
基于上面的观测数据，建立GOCE卫星观测法方程的简要描述如下：首先，对原始数据进行粗差探测与剔除、时变重力场信号的分离、数据内插分段等预处理；直接利用沿轨的引力梯度对角线数据(Vxx，Vyy，Vzz)建立观测方程及相应的法方程，再使用方差分量估计方法估计引力梯度张量对角线分量各自的单位权中误差，3个分量的单位权中误差分别为2.5、3.3和4.3 mE。考虑到引力梯度观测值的有色噪声特点，在最小二乘过程中引入通带为0.005~0.041 Hz的ARMA(auto regressive moving-average)滤波器。采用点域加速度方法由PKI(precise kinematic orbit)轨道观测值建立相应的SST法方程，求解SST法方程得到观测值残差，并基于残差估计了单位权方差。基于SGG各分量的单位权方差以及SST的单位权方差，得到了它们之间的加权因子，最后基于加权因子将SGG数据各分量的法方程以及SST法方程相加得到GOCE卫星观测法方程，并求解得到相应的卫星重力场模型，求解时采用Tikhonov正则化方法解决极空白问题引起的法方程不稳定性问题。文献[18]中有更为详细的数据处理方法和过程的描述。
3.3 GOCE卫星观测法方程和格网重力异常法方程的联合解算
联合求解GOCE卫星法方程与重力异常块对角法方程可以得到最优的模型解，但是由于两个法方程在形式、阶次、频谱精度上有较大差异，联合解算时需要采取有效的策略。根据文献[9]，本文制定的法方程联合解算策略如图 4所示。
图 4 SGG-UGM-1的法方程联合解算示意图Fig. 4 The schematic diagram of normal equation combination in SGG-UGM-1
图选项
使用EGM2008模型的系数误差信息，可以计算得到模型在2190阶的重力异常累计误差约为5 mGal。本文使用5 mGal作为重力异常数据的单位权中误差，并基于该单位中误差以及上文中GOCE数据的单位权方差得到重力异常法方程与GOCE观测法方程的加权因子。
联合观测法方程分为两部分，其中251~2159阶系数法方程与重力异常法方程中对应部分相同，而2~250阶系数法方程是使用1.2节中的最小二乘联合平差方法联合101~250阶重力异常数据块对角法方程和GOCE全阶次法方程(220阶)得到的。构建重力异常数据101~250阶块对角法方程时需要在原始重力异常数据中移除2~100阶及251~2159阶系数的贡献。
模型的220~250阶系数是联合GOCE数据和重力异常数据联合求解的，这样可以保证联合模型的这部分系数及其验后方差的频谱有更好的连续性。重力异常数据求解低阶系数的精度较差，重力异常低阶次系数法方程与卫星法方程联合可能会“污染”卫星重力的低阶系数，所以为减弱其影响，2~250阶联合法方程中只采用了重力异常的101~250阶系数块对角法方程。
4 SGG-UGM-1超高阶重力场模型的精度分析4.1 SGG-UGM-1模型的系数误差
EIGEN-6C4模型是当前精度较高的超高阶重力场模型[11]。本文以此为参考模型，认为该模型的系数为“真值”计算SGG-UGM-1、GOSG-EGM、EGM2008和EIGEN-6C2[22]模型的“误差”，并对它们的误差进行对比分析，以在频域内分析SGG-UGM-1的特性，评价模型的精度与合理性。
GOSG-EGM模型是利用本文计算得到的220阶卫星重力场模型与EGM2008重力场模型系数直接拼接得到的，其2~220阶系数来自卫星重力场模型，221~2159阶系数来自EGM2008模型。图 5给出了SGG-UGM-1、GOSG-EGM、EGM2008和EIGEN-6C2模型系数的阶误差RMS。
图 5 重力场模型的系数阶误差RMSFig. 5 Degree error RMS of the gravity field models
图选项
由图 5可知，在2~70阶的频率范围内，EIGEN-6C2的误差最小，EGM2008次之，而SGG-UGM-1和GOSG-EGM较大。这是由于EIGEN-6C2与EIGEN-6C4的2~70阶系数的计算方法和使用数据类似，都是联合LAGEOS、GRACE和GOCE等数据求解的[9, 22]，EGM2008的2~70阶次系数主要是来自于GRACE观测数据的贡献；而SGG-UGM-1和GOSG-EGM的2~70阶次系数主要由GOCE数据求解，因此其与EIGEN-6C4模型的差异大于EGM2008、EIGEN-6C2。
在70~220阶，SGG-UGM-1和EIGEN-6C2与EIGEN-6C4更为接近，主要原因是两个模型在此频段内均有GOCE和地面重力数据的贡献，而GOSG-EGM仅有GOCE卫星的贡献，从170阶开始它的误差迅速增大。同时，在170~220阶部分，因为SGG-UGM-1是联合GOCE和地面重力数据解算的，因此比GOSG-EGM模型的系数误差平滑，在220阶没有出现突变现象。
图 6为SGG-UGM-1等重力场模型的大地水准面阶误差，反映了系数的误差。图 6中模型的大地水准面系数误差与系数相对参考模型的系数误差呈现出了一定的相似性。由图 6可知，SGG-UGM-1的系数误差在220阶附近达到最大，然后开始下降，它的误差谱连续，没有出现突变现象。相比EGM2008、GOCE-EGM、BDLS模型，SGG-UGM-1模型低阶次的系数误差得到了明显改善，这也说明了本文模型联合求解方法的合理性。
