基于GRACE KBRR數(shù)據(jù)的動(dòng)力積分法反演時(shí)變重力場(chǎng)模型

羅志才， 周浩， 李瓊， 鐘波

1 武漢大學(xué)測(cè)繪學(xué)院， 武漢　430079　2 華中科技大學(xué)物理學(xué)院， 武漢　430074　3 華中科技大學(xué)地球物理研究所， 武漢　430074　4 武漢大學(xué)地球空間環(huán)境與大地測(cè)量教育部重點(diǎn)實(shí)驗(yàn)室， 武漢　430079

?

基于GRACE KBRR數(shù)據(jù)的動(dòng)力積分法反演時(shí)變重力場(chǎng)模型

羅志才1,4， 周浩2,3*， 李瓊2,3， 鐘波1,4

1 武漢大學(xué)測(cè)繪學(xué)院， 武漢4300792 華中科技大學(xué)物理學(xué)院， 武漢4300743 華中科技大學(xué)地球物理研究所， 武漢4300744 武漢大學(xué)地球空間環(huán)境與大地測(cè)量教育部重點(diǎn)實(shí)驗(yàn)室， 武漢430079

摘要基于動(dòng)力積分法恢復(fù)了一組60階的時(shí)變重力場(chǎng)模型WHU-Grace01s，且在位系數(shù)解算過(guò)程中僅使用KBRR數(shù)據(jù).通過(guò)與CSR、GFZ和JPL發(fā)布的Release 05模型的階方差和位系數(shù)誤差譜對(duì)比可知，WHU-Grace01s模型在高階次部分的階方差較小，且對(duì)軌道共振現(xiàn)象不敏感.將WHU-Grace01s時(shí)變重力場(chǎng)模型與CSR、GFZ、JPL、DEOS、Tongji、ITG、AIUB和GRGS等8家機(jī)構(gòu)發(fā)布模型通過(guò)相同的濾波處理，獲得了全球地表質(zhì)量變化的時(shí)空分布，從結(jié)果可以看出：各個(gè)模型計(jì)算的時(shí)變信號(hào)在空域上分布十分接近，且WHU-Grace01s模型計(jì)算的太平洋中心和撒哈拉沙漠區(qū)域的質(zhì)量變化較?。粚?duì)比幾個(gè)典型質(zhì)量變化區(qū)域，WHU-Grace01s模型和JPL模型計(jì)算的長(zhǎng)江流域和珠江流域時(shí)變信號(hào)呈強(qiáng)相關(guān)，其相關(guān)系數(shù)分別為0.948和0.976，且與上述8個(gè)模型計(jì)算的兩個(gè)流域時(shí)變信號(hào)的相關(guān)系數(shù)均達(dá)到0.9以上；在南極區(qū)域和格陵蘭島，WHU-Grace01s模型和其他各個(gè)模型均能反映區(qū)域冰川質(zhì)量的積累或消融，且各模型計(jì)算獲得的長(zhǎng)期趨勢(shì)變化結(jié)果相當(dāng).研究結(jié)果表明，WHU-Grace01s模型和國(guó)內(nèi)外已發(fā)布機(jī)構(gòu)模型具有很好的一致性，且受到軌道共振影響較小.

關(guān)鍵詞時(shí)變重力場(chǎng)模型； GRACE； 動(dòng)力積分法； 衛(wèi)星重力

1引言

GRACE衛(wèi)星任務(wù)的成功實(shí)施為人類了解地球內(nèi)部構(gòu)造及淺層物質(zhì)運(yùn)移開啟了新紀(jì)元.自Wahr等(2004)發(fā)布第一組月重力場(chǎng)模型以來(lái)，基于GRACE衛(wèi)星時(shí)變重力場(chǎng)模型的全球或局域陸地水儲(chǔ)量變化研究(Chen et al.，2005a；Hu et al.，2006；Zhong et al.，2009；李瓊等，2013；Feng and Zhong，2015)、極區(qū)冰蓋質(zhì)量變化研究(Luthcke et al.，2006a；Chen et al.，2008；Luo et al.，2012；Ju et al.，2014)等均取得了豐碩的研究成果，對(duì)地球物理學(xué)、大地測(cè)量學(xué)、海洋學(xué)和冰川學(xué)等領(lǐng)域的發(fā)展具有重要意義.

基于GRACE L1b數(shù)據(jù)獲取時(shí)變重力場(chǎng)模型需要充分考慮各個(gè)載荷的誤差特性和綜合影響，且數(shù)據(jù)觀測(cè)量大，時(shí)間跨度長(zhǎng)，解算復(fù)雜.迄今為止，國(guó)際上僅有少數(shù)幾家知名機(jī)構(gòu)完成該項(xiàng)復(fù)雜工作，分別是德克薩斯大學(xué)空間研究中心(CSR，Center for Space Research)、美國(guó)宇航局噴氣推進(jìn)實(shí)驗(yàn)室(JPL，Jet Propulsion Laboratory)、德國(guó)地學(xué)研究中心(GFZ，GeoForschungsZentrum Potsdam)、波恩大學(xué)理論大地測(cè)量研究所(ITG，Institute of Theoretical Geodesy)、代爾夫特地球觀測(cè)與空間系統(tǒng)研究所(DEOS，Delft institute of Earth Observation and Space system)、伯爾尼大學(xué)天文研究所(AIUB，Astronomical Institute, University of Bern)和法國(guó)空間研究中心空間大地測(cè)量研究所(GRGS，Groupe de Recherche de Geodesie Spatiale)等.國(guó)內(nèi)的學(xué)者也長(zhǎng)期致力于GRACE衛(wèi)星任務(wù)恢復(fù)時(shí)變重力場(chǎng)模型，具有代表意義的有Shen等(2013)、冉將軍等(2014)分別基于短弧法解算了時(shí)變重力場(chǎng)模型Tongji-GRACE01和IGG-CAS，且與國(guó)際知名機(jī)構(gòu)發(fā)布模型的精度相當(dāng).表 1給出了國(guó)內(nèi)外已發(fā)布時(shí)變重力場(chǎng)模型的統(tǒng)計(jì)信息，大部分機(jī)構(gòu)選用動(dòng)力積分法進(jìn)行解算.在國(guó)內(nèi)，周旭華(2005)、肖云(2006)、鄒賢才(2007)、張興福(2007)和王慶賓(2009)等在動(dòng)力積分法恢復(fù)地球重力場(chǎng)靜態(tài)信息方面做了大量研究.Zhao等(2011)聯(lián)合軌道和KBRR觀測(cè)數(shù)據(jù)，比較了經(jīng)驗(yàn)參數(shù)的處理方式對(duì)時(shí)變重力場(chǎng)模型解算精度的影響.本文將在上述研究基礎(chǔ)上編寫一套動(dòng)力積分法軟件，以實(shí)現(xiàn)時(shí)變重力場(chǎng)模型反演.此外，Luthcke等(2006b)關(guān)于基線法的研究結(jié)果表明，由于軌道幾何特性和共振周期等綜合影響，僅使用KBRR(K-Band Range Rate)觀測(cè)能夠使某些特殊階次的時(shí)變重力場(chǎng)模型精度更高.因此，本文僅采用KBRR觀測(cè)值解算時(shí)變重力場(chǎng)模型.

表1　已發(fā)布的時(shí)變重力場(chǎng)模型

2利用動(dòng)力積分法反演時(shí)變重力場(chǎng)模型的基本原理

動(dòng)力積分法是將衛(wèi)星星歷擾動(dòng)作為地球重力場(chǎng)信息的泛函，通過(guò)建立軌道攝動(dòng)與地球重力場(chǎng)位系數(shù)之間的關(guān)系，精密獲取地球重力場(chǎng)模型的方法.

2.1基于HL-SST技術(shù)的動(dòng)力積分法

(1)

考慮到力模型誤差的累積影響，動(dòng)力積分法通常分弧段進(jìn)行.待解算的未知數(shù)中，與弧段相關(guān)的變量稱為局部變量，如整周模糊度和鐘差等與GPS觀測(cè)相關(guān)的變量、初始狀態(tài)和加速度計(jì)校準(zhǔn)參數(shù)等與軌道積分相關(guān)的變量；影響整個(gè)觀測(cè)時(shí)段內(nèi)的變量稱為全局變量，這里主要指重力場(chǎng)位系數(shù).根據(jù)對(duì)局部參數(shù)處理方式的不同，動(dòng)力積分法可分為一步法、兩步法和三步法.其中，一步法是基于公式(1)分弧段積分并構(gòu)建法方程，然后通過(guò)法方程約化和疊加同時(shí)求解所有局部變量和全局變量；二步法首先利用載波和偽距完成定軌，然后將與軌道積分相關(guān)的局部變量和全局變量一并求解；三步法是依次求得局部變量，最后疊加弧段法方程求解全局變量.Luthcke等(2006b)采用三步法實(shí)現(xiàn)了基于KBRR觀測(cè)解算時(shí)變重力場(chǎng)模型的基線法，取得了很好的計(jì)算結(jié)果.鑒于此，本文也將采用三步法實(shí)現(xiàn)動(dòng)力積分法反演時(shí)變重力場(chǎng)模型.

精密校準(zhǔn)星載加速度計(jì)是基于動(dòng)力積分法獲取高精度地球重力場(chǎng)模型的關(guān)鍵.本文引入加速度計(jì)三個(gè)觀測(cè)方向的尺度因子k、偏差因子b.

(2)

其中，fobs為加速度計(jì)觀測(cè)的非保守力，R為星固系轉(zhuǎn)換至慣性系的旋轉(zhuǎn)矩陣.在計(jì)算參數(shù)敏感矩陣時(shí)需要非保守力相對(duì)于校準(zhǔn)因子的偏導(dǎo)數(shù)，依次表示為：

(3)

考慮到力模型誤差對(duì)軌道積分的累積影響，通常在加速度計(jì)校準(zhǔn)過(guò)程中引入尺度和偏差因子的先驗(yàn)參數(shù)，該參數(shù)可由地面數(shù)據(jù)處理中心獲得(Bettadpur，2009).鄒賢才(2007)的研究結(jié)果表明，靜態(tài)重力場(chǎng)模型對(duì)基于參考模型的加速度計(jì)校準(zhǔn)影響較大.因此，本文選用最新發(fā)布的純GRACE重力場(chǎng)模型GGM05s為參考重力場(chǎng)模型.

2.2基于LL-SST技術(shù)的動(dòng)力積分法

(4)

(5)

(6)

(7)

距s為：

(8)

考慮到地球重力場(chǎng)同時(shí)作用于兩顆衛(wèi)星，根據(jù)偏導(dǎo)數(shù)計(jì)算的鏈?zhǔn)椒▌t，星間距離變率對(duì)球諧系數(shù)c的偏導(dǎo)數(shù)可表示為：

(9)

特別地，對(duì)于加速度計(jì)校準(zhǔn)參數(shù)等與單個(gè)衛(wèi)星飛行狀態(tài)有關(guān)的力模型參數(shù)，距離變率對(duì)該類參數(shù)偏導(dǎo)數(shù)的計(jì)算方式與距離變率對(duì)初始狀態(tài)的偏導(dǎo)數(shù)相一致；對(duì)于潮汐模型參數(shù)、三體引力參數(shù)等同時(shí)作用于兩顆衛(wèi)星的力模型參數(shù)，距離變率對(duì)該類參數(shù)偏導(dǎo)數(shù)的計(jì)算方式與距離變率對(duì)地球重力場(chǎng)位系數(shù)的偏導(dǎo)數(shù)相一致.

基于上述討論，即可根據(jù)GRACE衛(wèi)星LL-SST技術(shù)獲取初始狀態(tài)、加速度校準(zhǔn)因子以及位系數(shù)等參數(shù)的最佳估值，獲取精密動(dòng)力學(xué)軌道，并更新地球重力場(chǎng)模型.需要注意的是，本文在更新地球重力場(chǎng)模型時(shí)，僅采用KBRR觀測(cè)信息，未考慮軌道的影響.最后，利用各個(gè)月份觀測(cè)數(shù)據(jù)解算重力場(chǎng)模型，減去平均重力場(chǎng)模型，即可獲得各個(gè)月份的時(shí)變信息.

3基于GRACE L1b數(shù)據(jù)的時(shí)變重力場(chǎng)模型反演

本文利用動(dòng)力積分法反演軟件確定了GRACE時(shí)變重力場(chǎng)模型(周浩, 2015；Zhou et al.,2015)，所有解算過(guò)程均采用聯(lián)合OpenMP和MPI的并行函數(shù)庫(kù)(周浩等，2011，2015)，軌道積分和變分方程解算均采用Gauss-Jackson數(shù)值積分器(羅志才等，2013)，力模型包括地球中心引力、非球形引力、三體攝動(dòng)力、固體潮、海潮、大氣潮、極潮、相對(duì)論效應(yīng)、低階次位系數(shù)的長(zhǎng)周期變化和非保守力，非保守力由加速度計(jì)觀測(cè)數(shù)據(jù)提供，地固系與慣性系的坐標(biāo)旋轉(zhuǎn)采用IERS2003協(xié)議，具體描述見表2.

表2　各項(xiàng)攝動(dòng)力模型及其主要參數(shù)

本文計(jì)算均采用JPL提供的GRACE軌道數(shù)據(jù)、加速度計(jì)數(shù)據(jù)、星象儀數(shù)據(jù)和星間距離變率數(shù)據(jù).考慮到觀測(cè)數(shù)據(jù)存在空白和粗差等情況(鄭偉等，2009)，采用GRACE用戶手冊(cè)中提供的內(nèi)插、重采樣和剔除粗差的方法進(jìn)行觀測(cè)數(shù)據(jù)的預(yù)處理(Wu et al.，2006).以每個(gè)月的觀測(cè)數(shù)據(jù)為一個(gè)單元，每24小時(shí)積分一個(gè)弧段，步長(zhǎng)選為5 s，解算一組截?cái)嚯A次為60的月重力場(chǎng)模型WHU-Grace01s.

4時(shí)變重力場(chǎng)模型分析

4.1模型解算精度

為了分析WHU模型(為表述方便，后文中均采用各個(gè)模型的發(fā)布機(jī)構(gòu)代替相應(yīng)的時(shí)變重力場(chǎng)模型名稱)的解算精度，首先從2005年時(shí)變模型中減去全年的平均重力場(chǎng)信號(hào)，然后計(jì)算剩余時(shí)變信號(hào)的階方差，與CSR、JPL、GFZ的Release05模型比較結(jié)果如圖1所示.計(jì)算結(jié)果表明：第一，僅僅采用KBRR觀測(cè)值解算的WHU模型在前20階與其他三種模型的階方差一致.第二，各個(gè)模型在2階的階誤差RMS(Root Mean Square Error)差異較為明顯，且WHU模型最小.Luthcke等(2006b)采用基線法的反演結(jié)果表明，僅使用KBRR觀測(cè)值計(jì)算時(shí)變重力場(chǎng)模型的低階次信號(hào)小于聯(lián)合軌道和KBRR觀測(cè)值的情況.第三，WHU模型在高階次的階方差逐步優(yōu)于CSR、JPL和GFZ模型，全年的階方差在60階時(shí)均優(yōu)于上述三個(gè)模型.Zhao等(2011)的研究結(jié)果表明，GRACE時(shí)變重力場(chǎng)模型的信號(hào)主要集中在前30階次左右，隨著階次的增加，條帶誤差會(huì)淹沒(méi)時(shí)變信號(hào).由于沒(méi)有引入軌道觀測(cè)誤差，僅采用KBRR觀測(cè)值解算的WHU模型能夠更好地抑制高頻誤差.

圖1　WHU和CSR、JPL、GFZ的階方差(2005年)Fig.1　Degree variance of CSR, JPL, GFZ and WHU (2005)

圖2　WHU和CSR、JPL、GFZ的位系數(shù)誤差譜(2005年9月)Fig.2　Coefficient error spectral of WHU, CSR, JPL and GFZ (September 2005)

圖3　各種時(shí)變重力場(chǎng)模型計(jì)算的全球等效水高(2005年4月)Fig.3　Global EWH computed from different time-variable gravity field models (April 2005)

為了進(jìn)一步分析模型解算精度，給出了2005年9月WHU和JPL、CSR、GFZ時(shí)變重力場(chǎng)模型的位系數(shù)誤差譜，如圖2所示.整體而言：由于GRACE衛(wèi)星采用的南北飛行模式對(duì)帶諧系數(shù)敏感，四種模型在帶諧項(xiàng)誤差小，扇諧項(xiàng)精度較差，且GFZ模型的誤差最大.具體而言：在低階次部分，WHU模型的位系數(shù)誤差小于CSR、JPL和GFZ模型；在40階次后，CSR、JPL和GFZ模型的相關(guān)誤差較大，而WHU模型的相關(guān)誤差也不如其他三個(gè)模型明顯；特別地，由于僅采用KBRR觀測(cè)數(shù)據(jù)，WHU模型在15、30、45階等與軌道共振關(guān)系密切的階次在位系數(shù)譜圖中沒(méi)有明顯偏大的現(xiàn)象，而CSR、JPL和GFZ均有出現(xiàn)該現(xiàn)象.

4.2時(shí)變信號(hào)的空間分布

時(shí)變重力場(chǎng)模型可反映某個(gè)特定時(shí)間分辨率內(nèi)全球或局部地表質(zhì)量的變化特征，通常將時(shí)變重力場(chǎng)模型轉(zhuǎn)換到空間域，利用等效水高(EWH，Equivalent Water Height)反映地表質(zhì)量遷移的分布狀況.由于觀測(cè)誤差以及先驗(yàn)地球物理場(chǎng)模型誤差等綜合影響，直接采用時(shí)變重力場(chǎng)模型計(jì)算的全球等效水高將會(huì)存在明顯的南北向條帶誤差和高頻誤差(Swenson and Wahr，2006)，而Gauss濾波、Fan濾波可有效削弱高頻誤差的影響(Wahr et al.，2004；Zhang et al.，2009)，去相關(guān)濾波可以有效削弱南北向條帶誤差(Swenson and Wahr，2006).為了更加有效地削弱這部分影響，本文選用Fan濾波和去相關(guān)濾波的組合濾波方式(Luo et al.，2012)，其中，F(xiàn)an濾波半徑取為500 km，去相關(guān)濾波選用P3M6.此外，由于對(duì)C2,0項(xiàng)不敏感，本文采用激光測(cè)衛(wèi)數(shù)據(jù)解算的二階項(xiàng)替代GRACE時(shí)變重力場(chǎng)模型中的C2,0項(xiàng)(Chen et al.，2005b).

將國(guó)內(nèi)外9個(gè)機(jī)構(gòu)發(fā)布的時(shí)變重力場(chǎng)模型轉(zhuǎn)換到空間域，采用上述處理方式計(jì)算了2005年4月的全球等效水高變化，如圖 3所示.與其他機(jī)構(gòu)解算的時(shí)變重力場(chǎng)模型相比，WHU模型在全球大部分區(qū)域的時(shí)變信號(hào)分布趨于一致.結(jié)合Liu(2008)對(duì)流域的劃分，各個(gè)模型計(jì)算的亞馬遜流域的等效水高最大，處在雨季；非洲的剛果河和贊比西河等效水高均為正值；北美洲的密西西比河、馬更些河流域，歐洲的多瑙河流域，俄羅斯境內(nèi)的鄂畢河和葉尼塞河流域均為正值；北美洲的育空河流域，印度北部的恒河流域以及中國(guó)境內(nèi)的珠江流域均為負(fù)值.

此外，選取了太平洋中心(146°W—134°W，17°S—23°S)和撒哈拉沙漠(10°E—20°E，20°N—30°N)等全年質(zhì)量變化較小的區(qū)域，以評(píng)估時(shí)變重力場(chǎng)模型的反演精度.以2005年1月至2006年12月為例，計(jì)算上述兩個(gè)區(qū)域在不加入任何濾波和加入組合濾波的EWH，截?cái)嚯A次分別為20階和60階，最后給出了各個(gè)模型在該時(shí)間段內(nèi)EWH的RMS(如表3所示).結(jié)果表明：濾波前，DEOS、WHU和GRGS模型的精度較高；濾波后，DEOS、Tongji、WHU和GRGS模型的精度較高.

4.3時(shí)變信號(hào)在特征流域的檢核

表3　各種時(shí)變重力場(chǎng)模型計(jì)算的EWH RMS(單位：cm)

為了更詳細(xì)地分析WHU時(shí)變重力場(chǎng)模型與其他模型的一致性，計(jì)算了2003年1月至2008年12月的月重力場(chǎng)模型.以我國(guó)的長(zhǎng)江流域和珠江流域?yàn)槔酶鱾€(gè)時(shí)變重力場(chǎng)模型分別計(jì)算了其水儲(chǔ)量變化的時(shí)間序列(如圖4和圖 5所示)，基于各時(shí)間序列得到了各個(gè)流域水儲(chǔ)量變化的周年振幅A、周年相位φ和變化趨勢(shì)P(如表 4所示)，并比較了WHU模型與其他模型計(jì)算各個(gè)流域水儲(chǔ)量變化的相關(guān)系數(shù)(如表5所示).其中，長(zhǎng)江流域和珠江流域的相關(guān)系數(shù)平均值分別為0.900和0.949，WHU模型和JPL模型的相關(guān)系數(shù)最高，分別達(dá)到了0.948和0.976，且與CSR、GFZ、JPL、ITG、Tongji和GRGS在兩個(gè)流域的相關(guān)系數(shù)均超過(guò)了0.90；以兩個(gè)流域相關(guān)系數(shù)平均值為基準(zhǔn)，WHU模型與各個(gè)機(jī)構(gòu)模型的相關(guān)系數(shù)從大到小依次為JPL、Tongji、GRGS、CSR、ITG、GFZ、AIUB和DEOS.整體而言，WHU模型計(jì)算的流域水儲(chǔ)量變化與其他模型的結(jié)果較為一致.

表4　各種時(shí)變重力場(chǎng)模型計(jì)算各個(gè)流域的統(tǒng)計(jì)信息

表5　WHU模型與其他模型計(jì)算各個(gè)流域水儲(chǔ)量變化的相關(guān)系數(shù)

圖4　各種時(shí)變重力場(chǎng)模型計(jì)算的EWH(長(zhǎng)江流域)Fig.4　EWH of Yangtze river basin computed from different time-variable gravity field models

圖5　各種時(shí)變重力場(chǎng)模型計(jì)算的EWH(珠江流域)Fig.5　EWH of Pearl river basin computed from different time-variable gravity field models

4.4時(shí)變信號(hào)在冰蓋區(qū)域的檢核

南極和格陵蘭島作為南北兩極具有代表性意義的冰蓋區(qū)域，是地球物理學(xué)、大地測(cè)量學(xué)、海洋學(xué)和冰川學(xué)等領(lǐng)域的熱點(diǎn)研究區(qū)域，也是基于GRACE衛(wèi)星探測(cè)全球質(zhì)量變化的重要組成部分，因此，本文分別利用9個(gè)時(shí)變重力場(chǎng)模型計(jì)算了南極和格陵蘭島區(qū)域的地表質(zhì)量變化率(如圖6和圖7所示).由于僅用于評(píng)價(jià)各個(gè)時(shí)變重力場(chǎng)模型，冰后回彈和泄露誤差等對(duì)解算結(jié)果的影響是一致的，計(jì)算中并未考慮這些改正項(xiàng).

對(duì)于南極區(qū)域，西南極Amundsen海灣(圖6中A點(diǎn))均表現(xiàn)為明顯的質(zhì)量減少趨勢(shì)，是整個(gè)南極地區(qū)質(zhì)量變化率最大的區(qū)域；B點(diǎn)Ronne冰架均表現(xiàn)為質(zhì)量累積；C點(diǎn)北方為橫貫?zāi)蠘O的山脈，各個(gè)模型結(jié)果均表明該點(diǎn)附近的質(zhì)量呈現(xiàn)增加趨勢(shì)；各個(gè)模型在D點(diǎn)Wilkes Land的冰川呈現(xiàn)消融趨勢(shì)，其中CSR和WHU模型的質(zhì)量減少趨勢(shì)最明顯.

圖6　各種重力場(chǎng)模型計(jì)算的南極地表質(zhì)量變化率(2003年至2005年)Fig.6　The rate of Antarctic mass change recovered from released time-variable gravity field models (2003 to 2005)

圖7　各種時(shí)變重力場(chǎng)模型計(jì)算的格陵蘭島地表質(zhì)量變化率(2003年至2005年)Fig.7　The rate of Greenland mass change recovered from released time-variable gravity field models (2003 to 2005)

對(duì)于格陵蘭島區(qū)域，各個(gè)模型計(jì)算的地表質(zhì)量變化率表明，冰蓋整體呈現(xiàn)消融的趨勢(shì)，且中部和南部消融速度快，東北部的消融速度相對(duì)較??；格陵蘭島區(qū)域附近島嶼的地表質(zhì)量變化率也趨于一致.特別的，WHU模型和JPL模型計(jì)算的格陵蘭島東北部有小部分出現(xiàn)質(zhì)量增加現(xiàn)象，具體原因有待進(jìn)一步查明.此外，利用CSR、GFZ、JPL、DEOS、WHU、Tongji、ITG、AIUB和GRGS等模型計(jì)算了2003年至2005年格陵蘭島冰蓋質(zhì)量年變化率，依次為-101 Gt/a、-98 Gt/a、-93 Gt/a、-81 Gt/a、-95 Gt/a、-91 Gt/a、-96 Gt/a、-95 Gt/a、-102 Gt/a，均與Luthcke 等(2006a)中采用Mascon解算結(jié)果-99 Gt/a相當(dāng).

5結(jié)論

本文實(shí)現(xiàn)了基于GRACE L1b數(shù)據(jù)反演時(shí)變重力場(chǎng)模型的動(dòng)力積分法，僅采用KBRR觀測(cè)數(shù)據(jù)解算了一組時(shí)間分辨率為1個(gè)月、截?cái)嚯A次為60的時(shí)變重力場(chǎng)模型WHU-Grace01s，并通過(guò)多組數(shù)值計(jì)算比較了該模型與國(guó)內(nèi)外知名機(jī)構(gòu)發(fā)布模型的精度，得到如下結(jié)論：

(1) 通過(guò)與GFZ、JPL和CSR的Release 05模型對(duì)比可知，WHU模型在5至20階次左右與上述三個(gè)模型的階方差一致，而在高階次部分的階方差較?。幌啾扔贕FZ、JPL和CSR模型，WHU模型的位系數(shù)誤差譜在15及其整數(shù)倍的階次誤差沒(méi)有出現(xiàn)明顯偏大的現(xiàn)象.

(2) 通過(guò)Fan濾波和去相關(guān)濾波的組合濾波處理后的全球等效水高表明，僅采用KBRR觀測(cè)值解算的WHU模型在時(shí)變信號(hào)的空間分布上與其他模型相似，可明顯探測(cè)大尺度范圍的質(zhì)量變化，且太平洋中心和撒哈拉沙漠等質(zhì)量變化較小區(qū)域的精度較高.

(3) WHU模型與CSR、GFZ、JPL、DEOS、Tongji、ITG、AIUB和GRGS等8個(gè)模型計(jì)算的長(zhǎng)江流域和珠江流域質(zhì)量變化一致性較高，在上述兩個(gè)流域的相關(guān)系數(shù)平均值分別可達(dá)0.900和0.949；WHU模型和JPL模型的相關(guān)系數(shù)最高，分別達(dá)到了0.948和0.976；與CSR、GFZ、JPL、ITG、Tongji和GRGS在兩個(gè)流域的相關(guān)系數(shù)均超過(guò)了0.90.

(4) WHU模型與CSR、GFZ、JPL、DEOS、Tongji、ITG、AIUB和GRGS等8個(gè)模型計(jì)算的南極和格陵蘭島冰蓋區(qū)域的質(zhì)量變化率一致性較高.在南極，上述9個(gè)模型均能反映出Amundsen海灣、Ronne冰架、橫貫?zāi)蠘O山脈西南部部分區(qū)域、Wilkes Land等區(qū)域的質(zhì)量增加或者減少；在格陵蘭島，上述9個(gè)模型計(jì)算的質(zhì)量變化率依次為-91 Gt/a、-101 Gt/a、-98 Gt/a、-93 Gt/a、-81 Gt/a、-91 Gt/a、-96 Gt/a、-95 Gt/a、-102 Gt/a，均與Luthcke 等(2006a)僅使用KBRR觀測(cè)數(shù)據(jù)解算的-99 Gt/a相當(dāng).

綜上所述，本文解算的WHU-Grace01s模型與目前國(guó)內(nèi)外主要機(jī)構(gòu)發(fā)布的時(shí)變重力場(chǎng)模型精度相當(dāng).致謝感謝JPL提供解算所需的GRACE衛(wèi)星觀測(cè)數(shù)據(jù)，感謝兩位匿名審稿專家對(duì)本文提出的寶貴修改意見.

References

Bettadpur S. 2009. Recommendation for a-priori bias and scale parameters for Level-1B ACC Data (Version 2). GRACE TN-02.

Chen J L, Wilson C R, Tapley B D, et al. 2005a. Seasonal global mean sea level change from satellite altimeter, GRACE, and geophysical models.J.Geodesy, 79(9): 532-539.

Chen J L, Rodell M, Wilson C R, et al. 2005b. Low degree spherical harmonic influences on Gravity Recovery and Climate Experiment (GRACE) water storage estimates.GeophysicalResearchLetters, 32(14): L14405.

Chen J L, Wilson C R, Tapley B D, et al. 2008. Antarctic regional ice loss rates from GRACE.EarthPlanet.Sci.Lett., 266(1-2): 140-148. Feng W, Zhong M. 2015. Global sea level variations from altimetry, GRACE and Argo data over 2005—2014.GeodesyandGeodynamics, 6(4): 274-279.

Hu X G, Chen J L, Zhou Y H, et al. 2006. Seasonal water storage change of the Yangtze River basin detected by GRACE.ScienceinChinaSeriesD, 49(5): 483-491.

Ju X L, Shen Y Z, Zhang Z Z. 2014. GRACE RL05-based ice mass changes in the typical regions of Antarctica from 2004 to 2012.GeodesyandGeodynamics, 5(4): 57-67.

Li Q, Luo Z C, Zhong B, et al. 2013. Terrestrial water storage changes of the 2010 southwest China drought detected by GRACE temporal gravity field.ChineseJ.Geophys. (in Chinese), 56(6): 1843-1849, doi: 10.6038/cjg20130606.Liu X L. 2008. Global gravity field recovery from satellite-to-satellite tracking data with the acceleration approach [Ph. D. thesis]. Delft: Delft University of Technology.

Luo Z C, Li Q, Zhang K, et al. 2012. Trend of mass change in the Antarctic ice sheet recovered from the GRACE temporal gravity field.ScienceChinaEarthSciences, 55(1): 76-82.

Luo Z C, Zhou H, Zhong B, et al. 2013. Analysis and validation of Gauss-Jackson integral algorithm.GeomaticsandInformationScienceofWuhanUniversity(in Chinese), 38(11): 1364-1368.

Luthcke S B, Zwally H J, Abdalati W, et al. 2006a. Recent Greenland ice mass loss by drainage system from satellite gravity observations.Science, 314(5803): 1286-1289.

Luthcke S B, Rowlands D D, Lemoine F G, et al. 2006b. Monthly spherical harmonic gravity field solutions determined from GRACE inter-satellite range-rate data alone.GeophysicalResearchLetters, 33(2): L02402.Ran J J, Xu H Z, Zhong M, et al. 2014. Global temporal gravity filed recovery using GRACE data.ChineseJ.Geophys. (in Chinese), 57(4): 1032-1040, doi: 10.6038/cjg20140402.Shen Y Z, Chen Q J, Xu H Z, et al. 2013. A modified short arc approach for recovering gravity field model. ∥GRACE Science Team Meeting. Austin, TX: Center for Space Research, University of Texas.Swenson S, Wahr J. 2006. Post-processing removal of correlated errors in GRACE data.GeophysicalResearchLetters, 33(8): L08402.

Wahr J, Swenson S, Zlotnicki V, et al. 2004. Time-variable gravity from GRACE: First results.GeophysicalResearchLetters, 31(11): L11501. Wang Q B. 2009. Study on the recovery of gravitational potential model with dynamical method [Ph. D. thesis](in Chinese). Zhengzhou: PLA Information Engineering University.

Wu S C, Kruizinga G, Bertiger W. 2006. Algorithm theoretical basis document for GRACE level-1B data processing V1. 2. Jet Propulsion Laboratory, California Institute of Technology.

Xiao Y. 2006. Research on the earth gravity field recovery from satellite-to-satellite tracking data [Ph. D. thesis] (in Chinese). Zhengzhou: PLA Information Engineering University.

Zhang X F. 2007. The earth′s field model recovery on the basis of satellite-to-satellite tracking missions [Ph. D. thesis] (in Chinese). Shanghai: Tongji University. Zhang Z Z, Chao B F, Lu Y, et al. 2009. An effective filtering for GRACE time-variable gravity: Fan filter.GeophysicalResearchLetters, 36(17): L17311.

Zhao Q L, Guo J, Hu Z G, et al. 2011. GRACE gravity field modeling with an investigation on correlation between nuisance parameters and gravity field coefficients.AdvancesinSpaceResearch, 47(10): 1833-1850.

Zheng W, Xu H Z, Zhong M, et al. 2009. Effective processing of measured data from GRACE key payloads and accurate determination of Earth′s gravitational field.ChineseJ.Geophys. (in Chinese), 52(8): 1966-1975, doi: 10.3969/j.issn.0001-5733.2009.08.003.

Zhong M, Duan J B, Xu H Z, et al. 2009. Trend of China land water storage redistribution at medi- and large-spatial scales in recent five years by satellite gravity observations.ChineseScienceBulletin, 54(5): 816-821.

Zhou H, Zhong B, Luo Z C, et al. 2011. Application of parallel algorithms based on OpenMP to satellite gravity field recovery.JournalofGeodesyandGeodynamics(in Chinese), 31(5): 123-127.

Zhou H. 2015. Study on the determination of the earth′s gravity field from the combination of multi-type satellite gravimetry data [Ph. D. thesis] (in Chinese). Wuhan: Wuhan University.

Zhou H, Luo Z C, Zhong B, et al. 2015. MPI parallel algorithm in satellite gravity field model inversion on the basis of least square method.ActaGeodaeticaetCartographicaSinica(in Chinese), 44(8): 833-839. Zhou H, Luo Z C, Zhong B. 2015. WHU-Grace01s: a new temporal gravity field model recovered from GRACE KBRR data alone.GeodesyandGeodynamics, 6(5): 316-323.

Zhou X H. 2005. Study on satellite gravity and its application [Ph. D. thesis] (in Chinese). Wuhan: Institute of Geodesy and Geophysics, Chinese Academy of Sciences.

Zou X C. 2007. Theory of satellite orbit and earth gravity field determination [Ph. D. thesis] (in Chinese). Wuhan: Wuhan University.

附中文參考文獻(xiàn)

李瓊, 羅志才, 鐘波等. 2013. 利用GRACE時(shí)變重力場(chǎng)探測(cè)2010

羅志才, 周浩, 鐘波等. 2013. Gauss-Jackson積分器算法分析與驗(yàn)證. 武漢大學(xué)學(xué)報(bào)(信息科學(xué)版), 38(11): 1364-1368.

冉將軍, 許厚澤, 鐘敏等. 2014. 利用GRACE重力衛(wèi)星觀測(cè)數(shù)據(jù)反演全球時(shí)變地球重力場(chǎng)模型. 地球物理學(xué)報(bào), 57(4): 1032-1040, doi: 10.6038/cjg20140402.

王慶賓. 2009. 動(dòng)力法反演地球重力場(chǎng)模型研究[博士論文]. 鄭州: 中國(guó)人民解放軍信息工程大學(xué).

肖云. 2006. 基于衛(wèi)星跟蹤衛(wèi)星數(shù)據(jù)恢復(fù)地球重力場(chǎng)的研究[博士論文]. 鄭州: 中國(guó)人民解放軍信息工程大學(xué).

張興福. 2007. 應(yīng)用低軌衛(wèi)星跟蹤數(shù)據(jù)反演地球重力場(chǎng)模型[博士論文]. 上海: 同濟(jì)大學(xué).

鄭偉, 許厚澤, 鐘敏等. 2009. GRACE衛(wèi)星關(guān)鍵載荷實(shí)測(cè)數(shù)據(jù)的有效處理和地球重力場(chǎng)的精確解算. 地球物理學(xué)報(bào), 52(8): 1966-1975, doi: 10.3969/j.issn.0001-5733.2009.08.003.

周浩, 鐘波, 羅志才等. 2011. OpenMP并行算法在衛(wèi)星重力場(chǎng)模型反演中的應(yīng)用. 大地測(cè)量與地球動(dòng)力學(xué), 31(5): 123-127.

周浩. 2015. 聯(lián)合多類衛(wèi)星重力數(shù)據(jù)反演地球重力場(chǎng)的研究[博士論文]. 武漢: 武漢大學(xué).

周浩, 羅志才, 鐘波等. 2015. 利用最小二乘直接法反演衛(wèi)星重力場(chǎng)模型的MPI并行算法. 測(cè)繪學(xué)報(bào), 44(8): 833-839.

周旭華. 2005. 衛(wèi)星重力及其應(yīng)用研究[博士論文]. 武漢: 中國(guó)科學(xué)院測(cè)量與地球物理研究所.

鄒賢才. 2007. 衛(wèi)星軌道理論與地球重力場(chǎng)模型的確定[博士論文]. 武漢: 武漢大學(xué).

(本文編輯何燕)

基金項(xiàng)目國(guó)家自然科學(xué)基金項(xiàng)目(41131067，41374023，41474019，41504014，41504015)；大地測(cè)量與地球動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室開放基金項(xiàng)目(SKLGED2015-1-3-E)；中國(guó)博士后科學(xué)基金項(xiàng)目(2016M592337)資助.

作者簡(jiǎn)介羅志才，男，1966年生，工學(xué)博士，教授，博士生導(dǎo)師，現(xiàn)主要從事物理大地測(cè)量學(xué)和衛(wèi)星重力學(xué)研究.E-mail:zhcluo@sgg.whu.edu.cn*通訊作者周浩，男，1987年生，工學(xué)博士，主要從事衛(wèi)星重力數(shù)據(jù)處理及應(yīng)用研究.E-mail:lengyeshanren@126.com

doi:10.6038/cjg20160606 中圖分類號(hào)P223 年中國(guó)西南干旱陸地水儲(chǔ)量變化. 地球物理學(xué)報(bào), 56(6):1843-1849, 10.6038/cjg20130606.

收稿日期2015-01-27，2015-06-15收修定稿

A new time-variable gravity field model recovered by dynamicintegral approach on the basis of GRACE KBRR data alone

LUO Zhi-Cai1,4， ZHOU Hao2,3*， LI Qiong2,3， ZHONG Bo1,4

1SchoolofGeodesyandGeomatics,WuhanUniversity,Wuhan430079,China2SchoolofPhysics,HuazhongUniversityofScienceandTechnology,Wuhan430074,China3InstituteofGeophysics,HuazhongUniversityofScienceandTechnology,Wuhan430074,China4KeyLaboratoryofGeospaceEnvironmentandGeodesy,MinistryofEducation,WuhanUniversity,Wuhan430079,China

AbstractA new monthly time-variable gravity field model WHU-Grace01s truncated to 60 degrees and orders is determined solely from GRACE KBRR data using dynamic integral approach in this paper. Compared with the Release 05 models from CSR, GFZ and JPL, the degree variance of WHU-Grace01s are smaller in high degrees, and it is not very sensitive to orbit resonance phenomenon. After applying the same filter to the temporal signals derived from the time-variable gravity field models published by WHU, CSR, GFZ, JPL, DEOS, Tongji, ITG, AIUB and GRGS, we can get several useful conclusions as follows, the distribution of temporal signal computed by WHU-Grace01s is similar to the other models, and its signal in the Pacific Ocean center and Sahara desert is relative smaller. The correlation coefficients of the temporal signal in Yangtze river and Pearl river basin between WHU-Grace01s and JPL can reach to 0.948 and 0.976 respectively, the average correlations coefficients between WHU-Grace01s and the aforementioned 8 time-variable gravity field models is bigger than 0.90. All of the aforementioned time-variable gravity field models can detect the mass loss or accumulation in the Antarctica and Greenland, and the mass loss rate is similar to each other. In short, the accuracy of WHU-Grace01s has a good consistency with the previously published GRACE solutions, while the impact of orbit resonance on WHU-Grace01s is unconspicuous.

KeywordsTime-variable gravity field model; GRACE; Dynamic integral approach; Satellite gravity

羅志才， 周浩， 李瓊等. 2016. 基于GRACE KBRR數(shù)據(jù)的動(dòng)力積分法反演時(shí)變重力場(chǎng)模型.地球物理學(xué)報(bào)，59(6):1994-2005,doi:10.6038/cjg20160606.

Luo Z C, Zhou H, Li Q, et al. 2016. A new time-variable gravity field model recovered by dynamic integral approach on the basis of GRACE KBRR data alone.ChineseJ.Geophys. (in Chinese),59(6):1994-2005,doi:10.6038/cjg20160606.