環渤海地區GPS臺站時間序列分析

邱亞輝，薄志毅，郎 博，岳彩亞(. 北京工業職業技術學院，北京 0004； . 中國礦業大學(北京)，北京 00037)

環渤海區域主要包括燕山、中朝、魯東海和華北等4個二級板塊的部分區域，其西鄰鄂爾多斯塊體，南鄰寬闊穩定的華南塊體，北鄰中蒙塊體，東鄰太平洋。由于該區域復雜多變的環境和氣候，以及不同地方的地質、地貌差異較大，并且長時間受到周圍板塊的相互作用，因而使得該地區板塊較為活躍。而GPS基準站坐標時間序列為大地測量學及地球動力學研究提供了寶貴的基礎數據[1-2]。隨著GPS定位技術的不斷更新與發展，其定位精度目前可以達到10-8～10-9量級。

高精度GPS站點坐標時間序列是對站點連續觀測數據經過精密解算后得到的高精度位置信息序列，它既是對站點位置連續性的描述，也是對站點穩定的反映。研究測站坐標時間序列可以獲得臺站噪聲模型并將其剔除，有利于分析板塊水平方向運動特性和板塊高程方向運動規律，以及對地震的預測。前人已經對該領域做了許多研究，Johnson等分析了隨機漫步噪聲給GPS站點速率估計可能帶來的影響[3]；符養等通過對GPS坐標時間序列計算分析，獲取了全球高程季節變化和整體振蕩運動的特性[4]；袁林果等對香港GPS基準站時間序列特征進行分析，證實了香港GPS連續運行參考站網存在明顯的本地季節性周期變化[5]；張鵬等研究分析了中國地殼運動觀測網絡基準站的坐標位置變化規律，發現位于中國境內的基準站的時間序列具有一定的周期性，在高程方向分量的周期性最為明顯，并且時間序列的時頻特性表現出明顯的區域性[6]。本文基于中國陸態網絡中環渤海區域連續觀測得到的臺站坐標時間序列數據，分析了臺站主要受到的噪聲類型和噪聲頻率特性，研究了環渤海區域地殼運動特性。

1 GPS坐標時間序列噪聲譜類型分析原理

1.1 Lomb周期圖法

通常采用快速離散傅里葉變換、譜分析、Lomb周期圖法及最大似然估計對臺站坐標序列進行噪聲類型分析。其中快速傅里葉變換和譜分析要求觀測數據無間隔，在漫長的觀測時間段中，很難滿足這種要求，必須對缺失天數進行插值，因而會帶來一些誤差。因此本文使用Lomb周期圖法[7-8]，它可以計算缺失天數的坐標序列功率譜。

典型的Lomb周期圖法公式基于離散傅里葉變換

(1)

式中，fn=n/T為采樣頻率，T為采樣總長度；xi為GPS站點坐標時間序列。典型Lomb周期圖法的使用和譜分析的方法類似，采樣必須是無間隔的。

簡化Lomb周期圖法公式為

(2)

式中，τ可通過式(3)求出。

(3)

式中，ω為采樣頻率；ti為所取觀測天數，時間單位為1 d。簡化Lomb周期圖法可用來計算缺失天數的GPS測站坐標時間序列噪聲譜，除此法之外，近年來發展起來的小波分析也可以做到。正是由于這個原因，Lomb周期圖法和小波分析方法受到了該方向研究者的一致好評，避免了插值帶來的誤差影響。

1.2 GPS坐標時間序列譜指數反演

地球物理現象功率譜函數通常表示為1/fα,而GPS坐標時間序列功率譜可以表示為如下形式[8-9]

P(f)=P0f-α

(4)

式中，P0和α為待求量，P0為常量；α表示譜指數。自然界中現象多表現為-3<α<-1,此時為非平穩過程，即分形布朗運動，其中包括隨機漫步噪聲(RWN)。對于平穩過程現象，當α為-1時，稱為分形白噪聲，當α為0，則為白噪聲(WN)；當α為1時，則為閃爍噪聲(FN)；除此之外的冪指數噪聲系統稱為有色噪聲[9-11]。

在計算GPS站點坐標時間序列譜時，通常的計算方法是將式(4)兩邊取對數，形式如下

lnP(f)=lnP0-αlnf

(5)

式中，P(f)可通過式(2)計算得到；P0和α可通過曲線擬合或最小二乘擬合得到；f為對應時間序列頻率。

在GPS坐標時間序列的功率譜分析中，還可以用到基于小波的多分辨率分析方法，這種方法不僅可以分析出序列的譜性質，還可以分離出序列的頻率特性，可以一目了然地看出組成信號噪聲的頻率。除此之外，多分辨率分析方法可以提取坐標時間序列高程方向的周期項，據以前學者研究，大部分臺站高程方向時間序列存在半周年項和周年項，少部分臺站存在兩周年項。

2 實例分析與討論

數據來源：由中國地震局與總參測繪局、中科院、原國家測繪地理信息局、中國氣象局、教育部合作共建的中國大陸構造環境監測網絡(以下簡稱陸態網絡)已基本完成，其由260個連續GNSS基準站、30個連續相對重力觀測站和100個絕對重力觀測站組成[12]。自2007年以來積累了將近10年的觀測數據，為監測中國大陸地殼運動、重力場形態及變化，以及現代大地測量基準系統的建立和維持提供了充足的數據資源。取中國陸態網絡中位于環渤海區域自2010年至2014年的GPS臺站時間序列數據。

陸態網絡中的環渤海地區測站布設情況為：主要行政區域包括北京、天津、遼寧、河北和山東共三省兩市，共同構成一個“C”形區域。在這個“C”形區域內布設了密密麻麻的GPS測站點。陸態網絡的布設特點是緊緊圍繞著地震頻發地區和地殼相對不穩定地區，這樣能夠精確地監測中國大陸板塊整體運動和區域地殼板塊運動，為區域地塊的進一步劃分提供了數據基礎，從而可更加精確地預測地震。

2.1 臺站時間序列分析

本文使用的是環渤海地區的GPS連續觀測站的觀測數據，使用GMT功能繪出臺站的N、E、U這3個方向的時間序列差值圖。限于篇幅限制，選擇了位于河北省的HETS臺站和遼寧省的LNYK臺站進行分析，如圖1所示。

圖1 臺站HETS和臺站LNYK時間序列差值

圖1中N代表南北方向，北方向為正；E代表東西方向，東方向為正；U代表高程方向，向上的方向為正。橫軸是觀測的天數，以年積日標記，并且年與年之間采用相疊加的方式，縱軸是時間段內的單天解與第一天解的差值，圖1(a)為臺站HETS時間序列差值圖，(b)為臺站LNYK時間序列差值圖。從時序圖中可以看出臺站的垂直方向波動比水平方向波動要大，并且具有明顯的周期特性，類似于正弦函數，波峰多出現在每年的6、7月，波谷多出現在每年的12月與次年的1月，這說明環渤海區域的GPS臺站周年運動較為規則，并且季節性是影響該區域臺站高程方向運動的主要因素。據以前學者研究，臺站高程方向波動的周期性常表現出周年性和半周年性及季節性，有些臺站會表現出2年的周期性。這種波動特性是由地球周期性的膨脹和收縮，還是由周圍地塊間的相互作用，還是由該區域復雜環境和氣候造成，需要進一步的研究。

臺站的水平方向同樣在做微弱的運動，雖然海洋潮汐負荷與海洋非潮汐負荷、大氣潮汐負荷與大氣非潮汐負荷及固體潮負荷等對測站水平方向的運動有影響，但影響較小，因此臺站坐標的水平方向的變化仍主要以板塊運動為主。從時序圖中可以看出臺站有向北和西的微弱運動趨勢，運動速度分別為1.10和3.15 mm/a，從而可得在2010—2014年，該地區板塊有往西北向的微弱運動，但運動量級較小，為3.34 mm/a，并且該區域內所細分的次級塊體間相對運動不是很活躍。總而言之，位于該區域的GPS臺站在水平方向比垂直方向穩定。

2.2 臺站噪聲特性分析

對臺站HETS和臺站LNYK做lomb周期圖法計算并畫出頻譜和解算天數的對應關系，如圖2所示。功率譜密度的計算結果一般會隨所使用的方法和計算中選定的參數不同而有差別，但總體上都可以反映臺站的頻譜關系。

圖3是使用Lomb周期圖法經計算后繪出的頻譜圖。N、E是水平面在南北、東西兩個方向的分量，U是高程方向。從圖中可看出N、E兩個方向的坐標序列表現出了微弱的周期特性，而U方向的坐標序列則表現出了明顯的周期特性。正如上所述，臺站水平方向和高程方向的時間序列周期都存在一年和半年兩種形式，除此之外高程方向還呈現出了季節特性，這說明臺站垂直方向比水平方向容易受到季節性影響，這與張鵬等用小波分析法對IGS站坐標序列頻譜分析所得出的結論一致[6]。

測站在水平方向的周期性弱于測站高程方向的周期性，是因為板塊的運動基本可分解為兩個方向上的運動：一是在水平方向的運動，二是在垂直方向的運動。板塊在水平方向的運動主要表現為整體的移動，除此之外由于受到區域內次級板塊的相互作用以及其他因素的影響還夾雜著微弱的周期振蕩運動。板塊在垂直方向主要表現周期性的上升或者下沉，這種運動與重力激發、熱效應、水文動力學等有關，由于該區域測站位于渤海附近，容易受到海洋潮汐負荷的影響，因此，該區域臺站的周年和半周年運動特性與海洋潮汐密切相關。

圖2 HETS頻譜分析圖

圖3 LNYK頻譜分析圖

從圖3中各測站振幅分析可得，臺站HETS水平方向分量在2010—2013年間振幅幾乎呈現正弦形式的遞增和遞減，并且周期約為1年，在2013年后的時間段內臺站HETS在南北向上一直在做急劇的向北運動，而東西向上在做向西的運動，則反映出該臺站所在地塊可能受到了某些外力作用發生了明顯的運動趨勢變化。對臺站LNYK在水平方向分析可得，2010—2011年時間段內，該站所在的地塊在做不明顯的周期運動，且整體表現出向西南方向的運動趨勢，但在2011年后，該地區板塊呈現了明顯的東北向運動。由此可得兩個測站雖都位于環渤海地區，但位于兩個不同的次級塊體上，充分反映了環渤海區域地質、地貌復雜，斷裂帶較多，次級塊體的運動差異較大的特性。據不完全統計，僅在北京、天津、遼寧、河北、山東等地出現了多達十幾處斷裂。

利用Lomb周期圖法計算出功率譜后，可通過式(5)對常數P0和噪聲類型參數α采取線性擬合，見表1。表1中列出了HETS和LNYK這2個測站擬合后的值。由表1可得出α的值在-1和1之間，則為平穩過程，且α值越大，所代表的隨機過程越平滑。由上述可知環渤海地區的臺站主要受分形白噪聲、白噪聲、閃爍噪聲的影響，而隨機漫步噪聲(α=-2)的成分可能性不大，這與黃立人對中國地殼運動觀測網絡研究得出的結論基本一致[13]。但是α值在0附近波動，因此可從三者的權重分析得出臺站受到的白噪聲影響較大，而分形白噪聲和閃爍噪聲所占的成分相對較小。因此可得，臺站主要受到白噪聲的影響。

3 結 論

文中分析了環渤海地區GPS臺站速度場及各臺站的坐標時間序列信息，得出以下結論：

表1 坐標時間序列擬合結果

(1) 該區域臺站高程方向呈現明顯的周期振蕩特性，周期約為1年，振幅最大可達2 cm，波峰和波谷分別出現在每年的6、7月和12月及下一年的1月；從U方向可以看出，在所觀測的時間段內，該區域臺站高程方向的波動較為規則，即此處的板塊相對而言較為穩定。

(2) 該處的地殼有向西北方向運動的趨勢，但在2010—2014年間運動量級較小，運動值僅在3 mm左右。因此可判定這段時間內該區域塊體較為穩定。

(3) 使用Lomb周期圖法對環渤海地區的臺站坐標時間序列進行了分析，得出臺站水平方向和高程方向都表現出了一年、半年的周期性，但高程方向的周期性相對更加明顯，并且高程方向同時表現出了季節性變化，這是水平方向不具有的，進一步證實了結論(1)。

(4) 通過對α和P0線性擬合，得出α主要在0附近波動，并對影響臺站的噪聲做了類型分析，結果表明臺站主要受到白噪聲、分形白噪聲、閃爍噪聲的影響，在這3種噪聲中，最主要的是白噪聲的影響。但此處所有的臺站幾乎不會受到隨機漫步噪聲的影響?？偟膩碚f，在求高精度測站坐標和速度場信息時，為了消除季節性對測站造成的影響，需要選擇相同季節進行GPS重復觀測和解算；考慮影響臺站高程方向的因素表現出周年性和半周年性，因此在對環渤海地區高精度地殼運動監測時，最佳方案是通過長期跟蹤觀測，從而可以消除周期性誤差。

參考文獻：

[1] 姜衛平，夏傳義，李昭，等.環境負載對區域GPS基準站時間序列的影響分析[J].測繪學報，2014,43(12)：1217-1223.

[2] 陳俊勇，張鵬，武軍酈，等.關于在中國構建全球導航衛星國家級連續運行站系統的思考[J].測繪學報，2007,36(4)：366-369.

[3] 符養.中國大陸現今地殼形變與GPS坐標時間序列分析[D].上海：中國科學院上海天文臺，2002.

[4] JOHNSON H O，AGNEW D C.Monument Motion and Measurement of Crustal Velocities[J].Geophysical Research Letters，1995,22(21)：2905-2908.

[5] 袁林果，丁曉利，陳武，等.香港GPS基準站坐標序列特征分析[J].地球物理學報，2008,51(5)：1372-1384.

[6] 范朋飛.高精度GPS站點坐標時間序列分析與應用[D].西安：長安大學，2013：22-23.

[7] SCARGLE J D.Studies in Astronomical Time Series Analysis.Ⅱ-statistical Aspects of Spectral Analysis of Unevenly Spaced Data[J].The Astrophysical Journal,1982,263:835-853.

[8] 張鵬，蔣志浩，秘金鐘，等.我國GPS跟蹤站數據處理與時間序列特征分析[J].武漢大學學報，2007,32(3)：251-254.

[9] 朱文耀，符養，李彥.GPS高程導出的全球高程振蕩運動及季節變化[J].中國科學，2003,33(5)：470-481.

[10] AGNEW D C.The Time-domain Behavior of Power-law Noises[J].Geophys Res Lett，1992,19:333-336.

[11] 田云鋒.GPS位置時間序列中的中長期誤差研究[D].武漢：中國地震局地質研究所，2011.

[12] 趙斌，聶兆生，黃勇.大規模GPS揭示的華北地區現今垂直運動[J].大地測量與地球動力學，2014，34(5)：35-39.

[13] 黃立人.GPS基準站坐標分量時間序列的噪聲特性分析[J].大地測量與地球動力學，2006,26(2)：31-38.