李偉偉 鞠曉蕾 沈云中 張子占

（1同濟(jì)大學(xué)測繪與地理信息學(xué)院，上海200092；2同濟(jì)大學(xué)空間信息科學(xué)與及可持續(xù)發(fā)展應(yīng)用中心，上海200092；3中國科學(xué)院測量與地球物理研究所大地測量與地球動力學(xué)國家重點(diǎn)實(shí)驗(yàn)室，湖北武漢430077）

0 引言

利用GNSS跟蹤站長期觀測數(shù)據(jù)確定的地殼運(yùn)動信息（如速度場）對地球動力學(xué)的研究發(fā)揮了重要的作用。影響GNSS跟蹤站速度場估計(jì)精度的主要因素包括跟蹤站坐標(biāo)序列的長度，周期信號的合理建模，有色噪聲正確確定。通常要求坐標(biāo)序列的時(shí)間跨度大于2.5年，考慮周年與半年周期項(xiàng)［1］，顧及閃爍噪聲和隨機(jī)游走噪聲兩類有色噪聲才能準(zhǔn)確估計(jì)GNSS跟蹤站的速度及其精度［2-6］。近年來，衛(wèi)星重力測量（如GRACE）從質(zhì)量變化的角度分析徑向形變也得到了應(yīng)用。Davis等［7］首次分析比較了亞馬遜流域GRACE時(shí)變重力資料與GNSS坐標(biāo)序列，兩者具有較好的一致性。Van Dam等［8］比較分析了歐洲區(qū)域GNSS與GRACE 34個(gè)月估計(jì)得到的地表垂直位移，并將兩者呈現(xiàn)較低相關(guān)性的原因歸結(jié)于GNSS數(shù)據(jù)處理模型的不準(zhǔn)確，且在沿海區(qū)域的站點(diǎn)表現(xiàn)更明顯。隨后，Tregoning等［9］在 Van Dam研究基礎(chǔ)上，重新處理GNSS與GRACE數(shù)據(jù)，兩者的相關(guān)性有了明顯提高，接近25%的站點(diǎn)相關(guān)系數(shù)都大于 0.5。Tesmer等［10］分析全球區(qū)域 115個(gè)站GNSS站觀測坐標(biāo)和GRACE估算得到的徑向位移序列，60%站相關(guān)系數(shù)達(dá)到0.5以上，34%站相關(guān)系數(shù)在0和0.5之間，還有6%站呈現(xiàn)負(fù)相關(guān)。可見，不同區(qū)域GNSS與GRACE相關(guān)程度不同。

基準(zhǔn)站速度主要包括大尺度板塊運(yùn)動，局部構(gòu)造形變和冰后回彈三種地殼形變因素的綜合影響?；鶞?zhǔn)站位于基巖上且在板塊內(nèi)部相對穩(wěn)定，就可消除局部構(gòu)造形變的影響。而板塊運(yùn)動一般沿水平方向，此時(shí)的徑向變化則主要是冰后回彈引起的地殼形變［11］。Bevis等［12］利用西南極 GNSS觀測網(wǎng)坐標(biāo)序列獲得徑向形變值，并比較了 ICE-5G（VM2），IJ05（6A），RF3S20（β＝0.2）以及 HUY09m四種GIA模型的預(yù)測徑向形變值。觀測徑向形變值與預(yù)測徑向形變值之差的經(jīng)驗(yàn)分布函數(shù)表示四種GIA模型均存在低估形變或高估形變的現(xiàn)象。埃爾斯沃斯山脈和東南極沿海區(qū)域的GNSS站徑向形變結(jié)果與前兩種GIA模型預(yù)測值比較同樣有一定的偏差［12］，并按站點(diǎn)分析了其偏差的可能原因。

本文主要利用GNSS數(shù)據(jù)分析南極區(qū)域的噪聲特性和形變速度，并與GRACE時(shí)變重力場結(jié)果進(jìn)行比較分析，研究兩者在徑向位移的相關(guān)性，且比較了GNSS估計(jì)形變與GIA模型（W12a）預(yù)測站點(diǎn)形變。

1 數(shù)學(xué)模型

1.1 GNSS估計(jì)形變

GNSS跟蹤站坐標(biāo)序列通?？梢杂镁€性趨勢項(xiàng)、周年和半年項(xiàng)、階躍變化項(xiàng)和噪聲項(xiàng)表示［14］，即

式中ti是以年為單位的GNSS站點(diǎn)單日歷元，y0為常數(shù)項(xiàng)，r為線性速度，cf和df表示周期性運(yùn)動（f＝1表示年周期項(xiàng)，f＝2表示半年周期項(xiàng)），表示發(fā)生在歷元Tgj大小為gj的階躍，H表示階梯函數(shù)。按照最小二乘參數(shù)估計(jì)的函數(shù)模型，式（1）可以改寫為

其中y為觀測值向量，A為設(shè)計(jì)矩陣，x＝［y0，r，c1，d1，c2，d2，g1，…，gng］為待估的模型參數(shù)，e為觀測噪聲向量。若Σy為觀測值的協(xié)方差矩陣，則基于最小二乘準(zhǔn)則求得的參數(shù)估值及其協(xié)方差陣為

考慮到白噪聲、閃爍噪聲和隨機(jī)游走噪聲的影響，觀測值的協(xié)方差矩陣Σy可表示為

其中U1，U2，U3分別為上述三種噪聲的協(xié)因數(shù)陣，θ1，θ2，θ3為對應(yīng)的噪聲分量［3］，采用方差-協(xié)方差估計(jì)公式可計(jì)算各種噪聲分量的估值［15］。

1.2 GRACE估計(jì)形變

由GRACE估計(jì)的徑向位移可用（6）式表示，

其中，Δr為地球表面徑向位移變化，R為地球的平均半徑（6 378.136 km），θ、λ分別為余緯和經(jīng)度，Nmax為最大階數(shù)為完全規(guī)格化的勒讓德函數(shù)，ΔClm和ΔSlm為每月Stokes系數(shù)與系數(shù)平均值之差，h′l及 k′l為 l階負(fù)荷勒夫數(shù)，Wl為高斯濾波函數(shù)，主要用來減小高階系數(shù)本身誤差導(dǎo)致的影響。本文采用Han等［16］提供的負(fù)荷勒夫數(shù)，取Nmax為60，高斯濾波半徑為 500 km［17］。

2 數(shù)據(jù)處理與結(jié)果分析

2.1 數(shù)據(jù)介紹

2.1.1 GNSS數(shù)據(jù)

由于南極區(qū)域各跟蹤站建站時(shí)間不同，因此各站采集的數(shù)據(jù)時(shí)間段也不同。為了與GRACE數(shù)據(jù)時(shí)間統(tǒng)一，本文選取了南極大陸12個(gè)時(shí)間跨度為2004—2012年的GNSS站點(diǎn)作為分析對象，每站數(shù)據(jù)長度不少于3年。各站經(jīng)緯度和時(shí)間跨度如表1所示。圖1為選取的GNSS站點(diǎn)分布，除CRAR站和MCM4站處于火山地帶外，其他站均建立在基巖上。

表1 南極GNSS站數(shù)據(jù)概況Table 1.Data overview of GNSS stations in Antarctica

圖1 南極GNSS站點(diǎn)分布Fig.1.Distribution of GNSS stations in Antarctica

2.1.2 GRACE數(shù)據(jù)

本文采用的GRACE數(shù)據(jù)為2012年3月CSR（Centre for Space Research）公布的GSM RL05數(shù)據(jù)，其時(shí)間跨度為2004年1月到2012年6月，共99個(gè)月的數(shù)據(jù)。該數(shù)據(jù)采用了新的地球物理模型（海洋、大氣、潮汐等），其空間分辨率、精度和周期性變化特性等都優(yōu)于RL04數(shù)據(jù)［18］。由于GRACE GSM模型已經(jīng)扣除海洋大氣等的影響，因此在與GNSS數(shù)據(jù)進(jìn)行對比時(shí)，需要加回level 2數(shù)據(jù)中反映海洋和大氣質(zhì)量變化引起的時(shí)變重力場信息的GAC產(chǎn)品系數(shù)；此外，為了使GRACE數(shù)據(jù)能夠與GNSS數(shù)據(jù)的參考框架統(tǒng)一［19］，需要考慮GRACE數(shù)據(jù)的一階項(xiàng)［20］。

2.2 結(jié)果分析

2.2.1 噪聲分量

表2列出了各站估計(jì)的白噪聲、閃爍噪聲和隨機(jī)游走噪聲分量，由此可明顯看出南極GNSS跟蹤站的觀測噪聲中白噪聲并不是主要的噪聲，因此對其形變速度進(jìn)行估計(jì)時(shí)只顧及白噪聲的影響是不符合實(shí)際的。另外，除MCM4站以外，大部分站點(diǎn)均不存在隨機(jī)游走噪聲或其分量較小，這可能是由于MCM4站附近的火山運(yùn)動引起基巖不穩(wěn)定造成的［21］。

表2 噪聲分量Table 2.Noise component

2.2.2 相關(guān)性比較

如2.1節(jié)所述，GNSS與GRACE數(shù)據(jù)的時(shí)間分辨率不同，為了比較GNSS和GRACE去除線性趨勢項(xiàng)的形變數(shù)據(jù)，本文對GNSS數(shù)據(jù)進(jìn)行了月平均。圖2為12個(gè)站點(diǎn)的徑向月變化比較圖（藍(lán)色為GNSS，紅色為GRACE），可見GNSS和GRACE確定的形變總體趨勢一致，但前者的變化振幅比后者大。為了研究GNSS和GRACE的相關(guān)性，計(jì)算了兩者相關(guān)系數(shù)（如表3，顯著性水平為0.05）。結(jié)果顯示兩者的相關(guān)系數(shù)最大可達(dá)0.663，總體上呈現(xiàn)較強(qiáng)的相關(guān)性。

2.2.3 徑向形變速度比較

表4列出12個(gè)站點(diǎn)GNSS估計(jì)的徑向形變速度和W12a模型預(yù)測的形變速度。GNSS結(jié)果顯示位于南極半島的OHI3和PALM站的徑向形變速度最大，12個(gè)站中4個(gè)站下降，但變化趨勢均比較緩慢。由于CRAR站和MCM4站位于火山活動地區(qū)，其形變速度不能準(zhǔn)確反映形變趨勢。其他10個(gè)站形變方向基本一致，且6個(gè)站與W12a模型吻合較好，W12a模型在DEVI站和DUM1站較高估計(jì)形變，而在南極半島較低估計(jì)形變。

圖2 GNSS與GRACE數(shù)據(jù)求得的徑向坐標(biāo)月變化Fig.2.Vertical displacements derived from GNSS and GRACE data

表3 GNSS與GRACE數(shù)據(jù)求得的徑向坐標(biāo)相關(guān)系數(shù)Table 3.Correlation coefficient of GNSS and GRACE

表4 站點(diǎn)速度（mm·a－1）Table 4.Velocities of GNSS stations（mm·a－1）

3 結(jié)論

本文計(jì)算了GNSS與GRACE月平均坐標(biāo)的相關(guān)性，顧及有色噪聲估計(jì)南極GNSS站點(diǎn)的徑向形變，并比較分析了GNSS站估計(jì)形變和GIA模型預(yù)測形變。主要得出以下結(jié)論：

（1）南極GNSS坐標(biāo)序列的觀測噪聲主要是白噪聲和閃爍噪聲，然而部分站（如MCM4站）受當(dāng)?shù)丨h(huán)境影響隨機(jī)游走噪聲更明顯，因此在估計(jì)形變速度時(shí)，只考慮白噪聲是不準(zhǔn)確的。

（2）南極區(qū)域GNSS站序列與GRACE序列的相關(guān)系數(shù)最大為0.663，總體呈現(xiàn)較強(qiáng)的相關(guān)性。

（3）雖然研究站點(diǎn)形變趨勢總體與W12a模型相吻合，但是需要更多的站點(diǎn)來分析兩者的吻合度。

1 Blewitt G，Lavallée D.Effect of annual signals on geodetic velocity.Journal of Geophysical Research，2002，107（B7）：2145，doi：10.1029／2001JB000570.

2 Zhang J，Bock Y，Johnson H，et al.Southern California Permanent GPSGeodetic Array：Error analysis of daily position estimates and site velocities.Journal of Geophysical Research，1997，102（B8）：18035—18055.

3 Mao A L，Harrision C G A，Dixon T H.Noise in GPS coordinate time series.Journal of Geophysical Research，1999，104（B2）：2797—2816.

4 Teferle F N.Strategies for Long-term Monitoring of Tide Gauges Using GPS.Nottingham：University of Nottingham，2003.

5 Santamaría-Gómez A，Bouin M N，Collilieux X，et al.Correlated errors in GPS position time series：Implications for velocity estimates.Journal of Geophysical Research，2011，116（B1），doi：10.1029／2010JB007701.

6 LiW W，Shen Y Z.Velocity estimation of GPS base stations considering the coloured noises.Journal of Applied Geodesy，2012，6（3-4）：149—157.

7 Davis JL，Elósegui P，Mitrovica JX，et al.Climate-driven deformation of the solid Earth from GRACE and GPS.Geophysical Research Letters，2004，31（24），doi：10.1029／2004GL021435.

8 Van Dam T，Wahr J，Lavallée D.A comparison of annual vertical crustal displacements from GPS and Gravity Recovery and Climate Experiment（GRACE）over Europe.Journal of Geophysical Research，2007，112（B3），doi：10.1029／2006JB004335.

9 Tregoning P，Watson C，Ramillien G，et al.Detecting hydrologic deformation using GRACE and GPS.Geophysical Research Letters，2009，36（L15401），doi：10.1029／2009GL038718.

10 Tesmer V，Steigenberger P，Van Dam T，etal.Vertical deformations from homogeneously processed GRACE and global GPS long-term series.Journal of Geodesy，2011，85（5）：291—310.

11 朱新慧，孫付平.用甚長基線干涉測量數(shù)據(jù)檢測冰期后地殼回彈.地球物理學(xué)報(bào)，2005，48（2）：308—313.

12 Bevis M，Kendrick E，Smalley R，etal.Geodeticmeasurements of vertical crustal velocity in West Antarctica and the implications for icemass balance.Geochemistry Geophysics Geosystems，2009，10（10），doi：10.1029／2009GC002642.

13 Argus D F，BlewittG，PeltierW R，etal.Rise of the Ellsworthmountainsand partsof the EastAntarctic coastobserved with GPS.Geophysical Research Letters，2011，38（16），doi：10.1029／2011GL048025.

14 Rosanne N.Observation of Geodetic and Seismic Deformation with the Global Positioning System.San Diego：University of California，2002.

15 Li B F，Shen Y Z，Lou L Z.Efficient estimation of variance and covariance components：A case study for GPS stochasticmodel evaluation.IEEE Transactions on Geoscience and Remote Sensing，2011，49（1）：203—210.

16 Han D Z，Wahr J.The viscoelastic relaxation ofa realistically stratified Earth，and a further analysisof post-glacial rebound.Geophysical Journal International，1995，120（2）：287—311.

17 Jekeli C.Alternative methods to smooth the Earth’s gravity field.Department of Geodetic Science and Surveying.Ohio：Ohio State University，1981.

18 Bettadpur S.Product Specification Document.Center for Space Research，The University of Texas at Austin，2007-b.

19 Swenson S，Chambers D，Wahr J.Estimating geocenter variations from a combination of GRACE and ocean model output.Journal of Geophysical Research，2008，113（B8），doi：10.1029／2007 JB0 05338.

20 Crétaux JF，Soudarin L，Davidson F JM，et al.Seasonal and interannual geocentermotion from SLR and DORISmeasurements：comparison with surface loading data.Journal of Geophysical Research，2002，107（B12）：2374，doi：10.1029／2002JB001820.

21 Raymond C A，Ivins E R，Heflin M B，et al.Quasi-continuous global positioning system measurements of glacial isostatic deformation in the Northern Transantarctic Mountains.Global and Planetary Change，2004，42（1-4）：295—303.