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.