GRACE时变重力场滤波方法比较

赵元元，苏宗跃，刘 元，穆大鹏，3

(1.华北地质勘查局五一九大队，河北保定071051;2.山东科技大学测绘科学与工程学院，山东青岛266590;3.中国科学院遥感与数字地球研究所，北京100094)

一、引 言

GRACE采取卫星跟踪卫星的观测模式，能够确定地球重力场的时变特征，其时间分辨率为10 d到30 d，空间分辨率为 300 km 到 400 km［1］。GRACE卫星极大的提高了人类对地球表层质量迁移和再分布的认识，尤其是两极地区冰川消融［2］和大尺度流域水文变化［3］。

由于受到卫星轨道误差、观测误差、模型误差以及数据处理造成的误差等影响，使用GRACE数据研制的时变地球重力场模型含有较多噪声，尤其是位系数的高阶部分。为了抑制GRACE含有的噪声，一般通过设计平滑核函数，对位系数采取某种限制措施，比如经典的高斯平滑核函数［4］，通过降低高阶位系数的权重，来降低GRACE的随机噪声;Fan滤波是在高斯滤波基础上［5］，不仅降低高阶位系数的权重，而且也降低高次位系数的权重，从而进一步降低GRACE的噪声;DDK2滤波是通过贝叶斯估计［6］，利用GRACE两颗卫星的轨道误差来设计误差矩阵和模型信息设计信号矩阵，来提高位系数解的精度。GRACE滤波的方法还有很多，比如将卫星观测误差和泄露误差最小化来设计适合研究区域的平滑核函数［7］，基于高斯滤波的其他非各向同性滤波［8］，基于时间序列的统计方法［9］和经验正交函数方法［10］，以及被广泛应用的经验去相关滤波等［11］。

本文将使用高斯滤波、Fan滤波和DDK2滤波这3种方法来反演地表质量变化，并分析比较它们的结果异同。GRACE数据使用的是由CSR(Center for Space Research)发布的Level-2 RL05版本，用户可以从http:∥isdc.gfz-potsdam.de/grace进行下载，DDK2滤波的结果由Kushce提供，下载网址是http:∥icgem.gfz-potsdam.de/ICGEM/TimeSeries.html。

二、滤波方法

GRACE确定的时变重力场可以反演地球表层的质量变化，由于其位系数含有较大噪声，直接反演的结果很难分辨出需要的信号，需要对位系数进行某种限制。经过平滑之后，地表的质量异常可以由下式计算［12］

式中，θ和λ分别为地心纬度和地心经度;Δσ表示质量异常;a为地球平均半径;ρave为地球平均密度(5517 kg/m3)，˜Plm是规格化缔合勒让德函数，kl表示负荷勒夫数，Wlm为平滑核函数，ΔClm和ΔSlm为GRACE位系数与平均值的差值。Δσ除以水密度ρw(=1000 kg/m3)就得到了以等效水高表示的质量异常。

对于平滑核函数Wlm而言，其不同的设计构成了不同的滤波器。Wahr最早引入高斯平滑核函数来抑制高阶位系数的噪声，该平滑核函数仅对位系数的阶部分起作用，即Wlm退化为

由于高斯滤波的实质是对位系数的不同阶赋予不同的权重，同一阶下的不同次的权重一样，文献［5］在高斯滤波基础上，构造了Fan滤波

式中，Wl和Wm均为高斯平滑核函数，这样Fan滤波对位系数的阶和次同时起作用。

文献［6］通过贝叶斯估计，使用GRACE两颗卫星的几何轨道误差来近似估计误差矩阵E，使用模型确定的信息来估计GRACE的信号协方差矩阵S，这样，平滑核函数可以由下式得到

式中，a表示正则化因子。

三、结果分析

本文使用CSR发布的Level-2 RL05数据，时间跨度为2003年1月到2012年12月，共120个月数据，缺失月份的数据通过线性内插得到，利用这120个月的数据计算GRACE位系数的平均值，再计算每个月位系数的残差。利用式(1)，使用高斯滤波、Fan滤波和DDK2滤波计算了2010年4月和10月的质量异常，并转化成等效水高值，其中高斯滤波和Fan滤波的半径均为400 km。选择4月份和10月份的原因是，在赤道南北两侧附近，这两个月份的质量异常变化较其他月份强烈，并且趋势相反。

从等效水高图看，在4月份，3种滤波方法均可以在亚马逊流域、刚果河流域、澳洲北部以及南极和格陵兰岛等地区观察到强烈的质量变化信号。对于残余的南北条带状噪声，高斯滤波要明显大于Fan滤波和DDK2滤波。Fan滤波不仅抑制高阶的位系数，而且抑制高次的位系数，相对于高斯滤波，Fan滤波能更有效的降低噪声，但同时也会造成信号的损失;在10月份，3种滤波方法的南北条带状噪声均要比各自4月份的结果要大，但Fan滤波和DDK2滤波的噪声残余仍然比高斯滤波小(见表1)。此外，在亚马逊地区，Fan滤波结果的振幅明显要小于高斯滤波和DDK2滤波，这也说明Fan滤波虽然能有效抑制高阶次的噪声，但也会造成该部分信号的损失。

表1 3种滤波方法结果统计

3种滤波的统计结果也有较大差异，如图1所示，高斯滤波和Fan滤波结果较为接近，后者的最大值、最小值以及均方根都要略小于前者，这也是前面指出的Fan滤波对于高阶次位系数的抑制。对于DDK2滤波，在4月份，其最大值要高于高斯滤波48 mm，高于Fan滤波71 mm，而最小值差距则超过了260 mm，均方根也较前两者有9 mm和12 mm的差距。在10月份，DDK2滤波与高斯滤波和Fan滤波的结果差异进一步拉大，其最小值是Fan滤波的2倍，而均方根的差距则达到了13 mm和24 mm。造成这种较大差异的原因在于，DDK2滤波使用的是贝叶斯估计方法，依赖于使用的先验信息，也就是GRACE卫星的几何轨道误差矩阵和模型确定的信号矩阵，它在有效降低南北条带噪声的同时，也可能会造成信号的失真。

图1 2010年等效水高图

四、结 论

本文使用高斯滤波、Fan滤波和DDK2滤波反演了2010年4月和10月GRACE时变重力场模型确定的地表质量异常，并作了比较分析。高斯滤波只对GRACE位系数的不同阶起降权作用，也就是同一阶里面的不同次的权重一样;Fan滤波的构造是在高斯滤波基础上，它不仅对位系数的阶起作用，而且对次也同时起作用，进一步压缩GRACE高阶次的噪声;DDK2滤波是基于贝叶斯估计，使用先验信息构造平滑核函数，其滤波的结果依赖于先验信息的近似程度，文献［6］使用GRACE两颗卫星的轨道误差来设计误差矩阵，使用模型来计算信号矩阵，有效地降低了南北条带噪声。

对于4月份和10份的反演结果，从等效水高图和统计结果看，不同方法之间有着较大差异:Fan滤波和DDK2滤波之后的南北噪声残余要小于高斯滤波，3种方法的最大值、最小值和均方根也有着较大区别，尤其是DDK2与前两者之间;此外，对于4月份和10月份这两个不同的时间，3种滤波方法的对于南北条带噪声的抑制作用也很不同，4月份要明显优于10月份，可能原因是GRACE轨道运行以及地面跟踪、GPS跟踪在不同时间的误差水平不同。

［1］TAPLEY B D，BETTADPUR S，RIES J C，et al.GRACE Measurements of Mass Variability in the Earth System［J］.Science，2004，305(5683):503-505.

［2］LUTHCKE SB，ZWALLY H J，ABDALATI W，et al.Recent Greenland Ice Mass Loss by Drainage System from Gravity Observations［J］.Science，2006，314(5803):1286-1289.

［3］冯伟，莱莫尼 JM，钟敏，等.利用重力卫星GRACE监测亚马逊流域2002-2010年的陆地水变化［J］.地球物理学报，2012，55(3):814-821.

［4］詹金刚，王勇，郝晓光.GRACE时变重力位系数误差的改进去相关算法［J］.测绘学报，2011，40(4):442-446.

［5］ZHANG Z Z，CHAO B F，LU Y，et al.An Effective Filtering for GRACE Time-variable Gravity:Fan Filter［J］.Geophysical Research Letters，2009(36):L17311.

［6］KUSCHE J.Approximate De-correlation and Non-isotropic Smoothing of Time-variable GRACE-type Gravity Field Models［J］.Journal of Geodesy，2007(81):733-749.

［7］SWENSON S，WAHR J.Methods for Inferring Regional Surface-mass Anomalies form Gravity Recovery and Climate Experiment(GRACE)Measurements of Time-variable Gravity［J］.Journal of Geophysical Research，2002，107(B9):109-117.

［8］HAN SC，SHUM CK，JEKELI C，et al.Non-isotropic Filtering of GRACE Temporal Gravity for Geophysical Signal Enhancement［J］.Geophysical Journal International，2005，163(1):18-25.

［9］DAVISJM，TAMISIEA E，ELOSEGUI P，et al.A Statistical Filtering Approach for Gravity Recovery and Climate Experiment(GRACE)Gravity Data［J］.Journal of Geophysical Research，2008，113(01):117-202.

［10］WOUTERS B，SCHRAMA E O.Improved Accuracy of GRACE Gravity Solutions through Empirical Orthogonal Function Filter of Spherical Harmonics［J］.Geophysical Research Letters，2007(34):711-715.

［11］SWENSON S，WAHR J.Post-processing Removal of Correlated Errors in GRACE Data［J］.Geophysical Research Letters，2006(33):553-561.

［12］WAHR J，MOLENAAR M，BRYAN F.Time Variability of the Earth’s Gravity Filed:Hydrological and Oceanic Effects and Their Possible Detection Using GRACE［J］.Journal Geophysical Research，1998，103(B12):30205-30229.