季新宇，李廷豪

(南京海河測繪科技有限公司，江蘇 南京 211101)

0 引 言

長江總長約6 300 km，是僅次于亞馬孫河和尼羅河的世界第三大長河，長江流域總面積約為190萬平方公里。糧食產(chǎn)量約占中國的一半，人口約占中國的三分之一，在中國社會發(fā)展進程中發(fā)揮著重要的作用。隨著長江流域社會經(jīng)濟的快速發(fā)展，水資源問題越來越成為一個十分突出的因素[1]。

GRACE(Gravity Recovery and Climate Experiment)重力衛(wèi)星計劃是由美國宇航局和德國航空航天中心聯(lián)合發(fā)射，提供了一種基于空間測量來監(jiān)測全球時變重力場的手段。GRACE重力衛(wèi)星對陸地水儲量的變化同樣敏感，可探測平均小于1 cm的陸地水儲量變化(TerrestrialWater Storage Change,TWSC)。在流域尺度研究方面，胡小工等利用2002年4月至2003年12月之間的數(shù)據(jù)分析了長江流域水儲量的變化，結果表明長江流域水儲量周年變化幅度可達3.4 cm的等效水高，最大值出現(xiàn)在春季和秋季[2]。Chen等通過研究成功地監(jiān)測了亞馬孫流域2005年的干旱和2009年的洪水[3-4]。馮偉等研究了亞馬孫流域在2010年經(jīng)歷的一次干旱的過程[5]。尼勝楠等利用GRACE時變引力場數(shù)據(jù)并結合水文模型數(shù)據(jù)研究了長江流域和黃河流域的陸地水儲量變化，發(fā)現(xiàn)長江、黃河流域的水儲量變化與降雨量在不同尺度和頻段存在較好的時空變化一致性[6]。

由于GRACE衛(wèi)星軌道的設計和傳感器誤差以及解算時背景模型誤差等一系列因素的影響，利用GRACE數(shù)據(jù)解算的原始時變引力場模型中包含明顯的南北方向條帶誤差，在一定程度上掩蓋了真實的陸地水儲量變化信號。為了研究不同去相關方法對反演長江流域陸地水儲量變化的影響，本文分別利用SWENSON、P3M6、P4M15三種方法對2003年1月至2014年5月之間長江流域的陸地水儲量變化進行反演，與不進行去條帶處理結果及GLDAS水文模型計算的結果進行比較，從而分析不同去條帶誤差方法對GRACE反演長江流域陸地水儲量變化的空間和時間特征的影響。

1 反演陸地水儲量基本原理

利用引力場模型的球諧系數(shù)變化計算地球表面質量變化(以等效水高的形式表示)的模型為：

(1)

(2)

式中，a指地球平均半徑，r表示高斯濾波半徑，本文r取值500 km。

由于衛(wèi)星軌道設計和觀測精度的制約，直接利用位模型反演的陸地水儲量變化網(wǎng)格值存在明顯的南北向條帶誤差[8]，對于某固定次，偶(奇)數(shù)階之間存在明顯的相關性，并提出利用滑動多項式消除相關性。以Cl,m為例，對Cl-2α,m,…Cl-2,m,Cl,m,Cl+2,m,…Cl+2α,m進行4階多項式擬合，然后從原始值中扣除擬合值，從而得到去相關之后的值。關于平滑窗口ω(用于多項式擬合的球諧系數(shù)的個數(shù))和α之間的關系為ω=2α+1。(Swenson & Wahr,2006)對5次以上的球諧系數(shù)進行了去條帶處理，其平滑窗口寬度與球諧系數(shù)的次數(shù)m之間的關系為：

(3)

式中，max()表示兩者之的最大值，文獻[8]里選擇的經(jīng)驗參數(shù)值為A=30和K=10。上述方法(下文簡稱為SWENSON法)是利用階數(shù)來確定滑動窗口的寬度，并進行滑動多項式擬合，而Chen等[3]在計算Sumatra-Andaman地震引起的重力變化時，采用了P3M6方法，即未選擇滑動窗口，而是直接對同一次的所有偶(奇)數(shù)階進行多項式擬合，即，保持球諧系數(shù)的前m階保持不變，然后對m階及m階以上的系數(shù)用n階多項式進行擬合，奇數(shù)項和偶數(shù)項分開擬合，將擬合后的值從原始值扣除之后便得到曲線管制后的值，這種方法簡稱為PnMm方法。同樣的，Chambers等[9]在利用GRACE RL04數(shù)據(jù)和RL05數(shù)據(jù)比較研究海底壓力時使用了P4M15方法。

2 實驗驗證與分析

2.1 數(shù)據(jù)來源

實驗數(shù)據(jù)采用美國德克薩斯大學空間研究中心(CSR)發(fā)布的GRACE GSM RL05引力場模型數(shù)據(jù)，研究時段為2003年1月至2014年5月(缺失9個月的數(shù)據(jù)，后續(xù)計算出的陸地水儲量變化時間序列對缺失月份數(shù)據(jù)進行線性插值)。該產(chǎn)品已經(jīng)扣除了大氣、海洋、海潮、固體潮等等信號，主要反映了由于地表質量重新分布引起的時變重力信號。SLR獲得的C20項和一階項代替GRACE數(shù)據(jù)的C20項和一階項[10-11]。

為了驗證GRACE反演結果的有效性，采用全球陸地數(shù)據(jù)同化系統(tǒng)(Global Land Data Assimilation System，GLDAS)中的Noah陸地表面模型進行驗證。該系統(tǒng)是美國宇航局和美國國家環(huán)境預報中心(National Centers for Environmental Prediction，NCEP)國家海洋大氣局聯(lián)合項目，時間分辨率為一個月，1°×1°的空間分辨率，選取四層土壤濕度數(shù)據(jù)、雪深水當量數(shù)據(jù)模擬陸地水儲量。

2.2 時間序列分析

根據(jù)式(1)，分別采用SWENSON、P3M6、P4M15、NONE(不進行去相關)和500 km的高斯濾波組合的方法計算全球陸地水儲量變化的等效水高格網(wǎng)值，再依據(jù)長江流域輪廓，對流域內的格網(wǎng)值數(shù)據(jù)按照緯度余弦加權的方法獲得長江流域平均陸地水儲量變化的時間序列：

(4)

為了更合理地比較GRACE反演結果與水文模型計算結果，首先對GLDAS-Noah的格網(wǎng)數(shù)據(jù)做球諧展開，并截斷至60階，然后進行半徑為500 km的高斯濾波處理，最后同樣采用緯度余弦加權的方法計算長江流域的TWSC時間序列。結果如圖1所示。

由圖1可知，GRACE反演的長江流域陸地水儲量變化與GLDAS水文模型計算的結果均表現(xiàn)出明顯的季節(jié)性變化，水儲量變化極小值出現(xiàn)在春季，極大值出現(xiàn)在秋季，在季節(jié)尺度上，GRACE反演結果與GLDAS模擬結果具有較好的一致性，SWENSON、P3M6、P4M15、NONE 4種情況下反演的TWSA與GLDAS模擬結果之間的相關系數(shù)分別達到0.886、0.878、0.883、0.869。另外，從上圖明顯看出，利用GRACE位模型反演的陸地水儲量變化值均高于GLDAS模擬結果，這是由于GLDAS-Noah模擬的陸地水儲量變化只包含了雪深水當量變化值和土壤深度為0～200 cm之間水儲量變化，而GRACE監(jiān)測的陸地水儲量變化不僅包含了土壤水和冰雪水變化，還包含了地下水、地表水等一系列變化[12]。

為了進一步比較GRACE不同去條帶誤差方法反演的陸地水儲量變化與GLDAS模擬結果在季節(jié)上的變化規(guī)律[13]，對時間序列進行了擬合分析，結果如表1所示。

表1 GRACE不同去相關方法反演結果與GLDAS模擬長江流域TWSC的周年項和半周年項的振幅、相位

由表1可知，基于GRACE數(shù)據(jù)，利用SWENSON、P3M6、P4M15去相關方法計算的長江流域陸地水儲量變化的周年振幅關系是：P4M15> P3M6>NONE>SWENSON，半周年振幅：P3M6>P4M15>NONE>SWENSON，GRACE反演結果的周年、半周年振幅均大于GLDAS模擬結果。圖2為采用以上幾種方法得到的長江流域陸地水儲量年平均變化空間分布趨勢圖。

由圖2可知，在長江流域中下游，利用SWENSON方法反演的陸地水儲量變化趨勢的空間分布與其余3種情況相差較大，這可能是由于SWENSON方法在使用時對于A與K兩個經(jīng)驗值在選擇上，采用的是文獻[8]里的原始參數(shù)，并未進行改動，容易產(chǎn)生較大誤差。另外，就P3M6與P4M15比較而言，在長江流域上游，P4M15反演出的陸地水儲量變化趨勢幅值明顯大于P3M6的反演結果，而在中下游地區(qū)，P3M6反演出的陸地水儲量變化趨勢幅值明顯大于P4M15的反演結果。

3 結 語

(1)GRACE反演的陸地水儲量變化時間序列的與GLDAS-Noah模擬的結果較為一致，相關系數(shù)均達到0.8以上，另外，GRACE反演結果的幅值均大于GLDAS模擬的幅值，這是由于GRACE反演的陸地水儲量變化不僅包括了土壤水和雪水當量變化，還包含了地下水、地表水等一些列水儲量變化。

(2)在時間序列上，不同去相關方法反演的陸地水儲量變化周年項的半周年項的振幅接近，在空間分布上，不同去相關方法反演的陸地水儲量變趨勢的差異較大，但P3M6和P4M15反演結果的分布十分接近，而與SWENSON相差較大。在長江流域上游，P4M15方法反演的陸地水儲量變化趨勢信號強度最強，而在中下游地區(qū)，SWENSON方法反演的陸地水儲量變化趨勢信號強度最強。