GPS垂向年周期運動的精確獲取與方法評估*

李經緯，占 偉，梁洪寶，王 勇，郭南男

(1.中國地震局第一監測中心，天津 300180；2.天津城建大學 地質與測繪學院，天津 300384)

0 引言

近年來，我國許多區域(如川滇地區)GPS垂向存在顯著的周期性運動(其中年周期運動最為突出)(張飛鵬等，2002；劉任莉等，2013；梁洪寶等，2015；王巖等，2018)，盛傳貞等(2014)認為滇西地區GPS坐標變化時間序列的垂向分量中，普遍包含有明顯的年周期非構造形變波動，高值可達12 mm；Zhan等(2017)基于2010—2015年云南27個GPS連續站垂向時間序列得出該時段云南地區GPS連續站垂向年周期振幅平均值為9.7 mm，認為云南地區季節性的降雨與GPS垂向年周期運動相關性較好。GPS垂向周期性運動研究有重要的研究意義和實際應用價值，目前我國GPS流動觀測的觀測數據量太少、無法消除周期性運動的影響，導致得到的垂向線性運動速率結果精度不高，GPS垂向速度場應用較少(方穎等，2014；朱爽等，2017；李長軍等，2018)，因此精確獲取垂向周期性運動信息能夠提高GPS流動觀測的垂向線性運動速率精度，極大地促進GPS垂向速度場的應用。

一些學者使用了不同的方法獲取GPS流動站的垂向周期性運動參數(Liangetal，2013；占偉等，2016)，但是對這些方法的效果評估、分析較少，尤其是中國大陸幅員廣闊，不同區域影響垂向周期性運動的主要因素不一，不同方法的適用性也不同，因此對GPS垂向周期性運動參數獲取方法進行全面的評估十分必要。本文以云南為研究區域，選用4種垂向年周期運動參數獲取方法，以研究區域內2010—2015年陸態網絡26個GPS連續站觀測數據得到的垂向年周期運動結果，評估這4種方法的效果，綜合評估各方法的優缺點，對不同的區域如何獲取垂向年周期運動參數給出建議。

1 方法介紹

常用的垂向年周期運動獲取方法分為2類，一類是根據一些物理模型(如負荷模型或者GRACE)直接求解某一空間位置的垂向年周期運動參數；另一類是對區域內GPS連續站實測的垂向年周期運動參數進行空間內插。本文選用4種GPS垂向年周期運動獲取方法，分別為：

(1)方法1：負荷模型方法。地球表面各種負荷變化是GPS垂向周期性運動的主要影響因素，而在云南地區，陸地水的質量變化是該區域GPS垂向周期性運動的主要影響因素(姜衛平等，2013；Jiangetal，2017)，因此筆者從EOST Loading Service(1)http：//loading.-ustrasbg.fr/ditrf_form.php.下載了2010—2015年全球陸地水文負荷模型(GLDAS)解算的垂向位移格網數據(CM)，數據空間分辨率為0.5°×0.5°，時間分辨率為3 h。基于云南地區各格網點垂向位移時間序列，采用線性+年周期+半年周期方法擬合得到各格網點年周期運動正弦波振幅A和余弦波振幅B，為了與GPS連續站結果比對，采用多核函數法空間插值獲取云南地區26個GPS連續站所在空間位置的年周期運動正弦波振幅A1和余弦波振幅B1。

(2)方法2：GRACE模型方法。利用GRACE數據給出的地球重力場球諧系數(Bruinsmaetal，2010)和負荷勒夫數(Farrell，1972)，可以根據彈性負荷形變理論獲取負荷引起的形變信息(Wahretal，1998；Kusche，2005；Fuetal，2012)。筆者按照Hao等(2016)的方法處理了CSR提供的GRACE-RL05數據：使用GCM05S模型去掉平均重力場；用SLR觀測的C20項系數項和一階項系數替換模型中相應值；采用高斯濾波方法去除高階模型噪聲，最終得到了26個GPS連續站所在空間位置的垂向位移結果(采樣率為每月)，采用線性+年周期+半年周期方法擬合得到各GPS連續站的年周期運動正弦波振幅A2和余弦波振幅B2。

(3)方法3：空間加權法。基于Zhan等(2017)給出的2010—2015年云南地區26個GPS連續站垂向位移時間序列(進行了預處理，包括去除了儀器更換等原因產生的突跳)，采用線性+年周期+半年周期方法擬合得到各連續站年周期運動正弦波振幅A0和余弦波振幅B0。對某一個GPS連續站(記為插值點)，根據其余25個連續站的年周期運動正弦波振幅A0和余弦波振幅B0，采用空間加權法插值得到插值點的年周期運動正弦波振幅A3和余弦波振幅B3：

(1)

式中：dj為插值點與第j個GPS連續站的水平距離；A0j和B0j為第j個GPS連續站的年周期運動正弦波和余弦波振幅。

(4)方法4：多核函數法。該方法計算過程與方法3類似，對于每一個GPS連續站，都是其它25個GPS連續站垂向年周期運動參數(A0和B0)空間插值得到插值點的年周期運動參數(A4和B4)，唯一的區別就是空間插值采用多核函數法。多核函數法由美國Hardy(1971)提出，其空間插值計算過程如下：

(2)

式中：f(λ，φ)為多核函數；sj(λ，φ，λj，φj)為多核函數；dj1.1為球面上2點間的大地線長度(單位為km)；cj為待定系數(楊國華等，2012)。將所有可用的數據點作為核函數點，可以獲得數學解析式，計算任一空間位置的結果。

2 效果評估

本文依據4種方法與GPS連續觀測得到的垂向年周期運動振幅差異和一致性評估各方法的效果。

2.1 振幅差異

圖1給出了GPS連續觀測得到的垂向年周期運動與4種方法結果之間的差異，從中可以看出：①方法1和方法2與GPS連續觀測的垂向年周期運動振幅差異較大，且方向上大體一致，說明這2種方法得到的垂向年周期運動與GPS連續觀測結果存在系統的時間差；②方法3和方法4與GPS連續觀測的垂向年周期運動振幅差異比前2種方法偏小，且差異的方向較為雜亂、沒有表現出較為明顯的規律，說明后2種方法更為接近GPS連續站得到的垂向年周期運動。表1給出了4種方法與GPS實測結果的振幅差異統計結果，4種方法與GPS連續站實測的振幅差異平均值分別為4.97，3.31，1.38，1.43 mm，而26個GPS連續站實測的垂向年周期運動振幅平均值為9.87 mm，可見方法1和方法2與GPS實測結果的差異較為顯著(分別達到實測振幅的50.35%和34.95%)。

圖1 4種方法與GPS連續觀測的年周期運動差異Fig.1 Difference of crustal annual motions between continuous GPS and four methods

表1 4種方法與GPS連續觀測的年周期運動振幅差異Tab.1 Difference of crustal annual motion amplitude between continuous GPS and four methods (單位：mm)

2.2 一致性評估

本文用WRMS減速比(Fu，Freymueller，2012)來評估各方法得到的垂向年周期參數與GPS連續站時間序列的一致性，定義為：

(3)

式中：RatioW表示WRMS減速比；WRMSG是2010—2015年26個GPS連續站垂向時間序列去除線性后WRMS值；WRMSm是GPS垂向時間序列去除線性后再去除年周期運動后的殘差WRMS值；WRMSf是GPS垂向時間序列去除線性和年周期運動后的殘差WRMS值。

RatioW能夠反映GPS連續站的垂向時間序列和各方法求得的年周期運動相位和振幅的一致性。圖2給出了4種方法的RatioW結果，從中可以看出后2種方法的絕大多數測站RatioW大于0.9，遠遠優于前2種方法，但方法3和4中YNZD站RatioW相對較小(圖2c，d)。表2給出了4種方法的RatioW統計結果，其平均值分別為0.70，0.84，0.95，0.94，表明方法3和4效果大體相當，優于方法1和2，相對而言，方法2的效果優于方法1。

表2 4種方法的RatioW統計結果Tab.2 Static results of WRMS reduction ratio of the four methods

圖2 4種方法的RatioW
Fig.2 WRMS reduction ratio of four methods

3 討論

從4種方法的整體評估結果可以看出，方法1和方法2得到的垂向年周期運動無論是振幅大小還是運動方向，均與GPS連續觀測的結果存在較大的差異，而方法3和方法4利用GPS連續觀測進行空間內插，得到的垂向年周期運動與GPS連續觀測結果一致性較好。

由于全球負荷模型是基于全球觀測數據建立的模型，模型建立時使用的中國大陸范圍內的實測數據較少，因此在中國大陸尤其是云南這樣的局部區域與實測數據一致性較差，使得全球負荷模型得到的垂向年周期運動與GPS實測結果差異較大。

方法2的GRACE模型得到的垂向周期性信號反映的是區域空間尺度達300～400 km的整體形變，GPS垂向周期性信號不僅反映了區域整體形變，還可能疊加一些局部效應(姜衛平等，2013；Yanetal，2016)，因此GRACE得到垂向年周期運動振幅與GPS實測結果存在一定差異。

需要說明的是：不同的全球負荷模型(方法)在中國大陸的符合程度存在差異，同時GRACE數據處理方法(方法2)也存在一定的差異，這2種方法在中國大陸的應用優化值得更為深入的研究。

云南地區GPS連續站較多，且區域垂向周期性運動空間一致性較好，各點的垂向年周期運動無論是振幅大小還是初相均較為一致(Zhanetal，2017)，因此2種空間插值的方法效果要優于前2種方法，而且方法3和方法4評估結果大體相當，說明只要區域內垂向周期性運動空間一致性較好，即使選用不同的空間內插方法，對結果影響不大。值得注意的是，在26個GPS連續站中，YNZD站空間位置最為突出(最北端)，周邊的GPS連續站與之距離較遠，因此在方法3和方法4中，該站利用其余25個測站空間內插得到的垂向年周期運動參數與實測結果差異最大，這說明2種空間插值方法要求GPS連續站的分布盡量均勻。

上述4種方法各有優劣，雖然在本文的實驗中，方法1和方法2的效果不如方法3和4，但前2種方法可以直接利用已有的全球負荷模型和GRACE給出的任意一點垂向位求得任意一點的垂向年周期運動。方法3和方法4雖然效果較好，但是對區域內GPS連續站的數量和分布有一定的要求，在實際應用中可能受到限制。事實上，目前中國大陸仍有一些區域GPS連續站較少，以陸態網絡為例，目前中國大陸范圍內只有260個GPS連續站，平均站間距為130 km，僅少數區域如首都圈、南北地震帶GPS連續站密度相對較高(云南地區GPS連續站平均站間距79 km)，在西藏、東北、華南等區域GPS連續站較少，因此方法1和2在這些GPS連續站較少的區域仍然有一定的實際使用價值。

4 結論

本文選用4種獲取GPS垂向年周期運動參數的方法，以GPS連續觀測得到的垂向年周期運動結果評估了這4種方法的效果，結果表明：由于云南地區GPS連續站較多且垂向年周期運動空間一致性較好，利用GPS連續觀測的垂向年周期運動參數進行空間插值(方法3和方法4)都能夠取得較好的效果，GRACE模型方法(方法2)效果稍差，負荷模型方法(方法1)效果最差，因此對于獲取區域內(如GPS流動站)的垂向年周期運動參數，在GPS連續觀測站較多的區域，建議采用空間加權法(方法3)或者多核函數法(方法4)；在GPS連續站較少的區域，建議采用GRACE模型方法。

由于本文的主要目的是驗證垂向年周期運動信息的構建模型有效性，因此只考慮了季節性運動中最主要的年周期運動，但本方法也適用于其它周期性運動(如半年周期運動)。

感謝NASA提供GRACE數據、EOST提供全球陸地水文負荷模型垂向位移數據。