胡夢迪 ，金 星 ，李 軍

（1.中國地震局工程力學(xué)研究所，黑龍江 哈爾濱 150080；2.中國地震局地震工程與工程振動重點實驗室，黑龍江 哈爾濱 150080；3.福建省地震局，福建 福州 350003）

0 引言

地震記錄的時間準(zhǔn)確性對一些時間精度要求較高的地震研究至關(guān)重要，特別是對于地震預(yù)警、地震定位以及地殼結(jié)構(gòu)層析成像等研究。地震臺站授時主要由兩部分完成，一是通過GPS（全球定位系統(tǒng)）進行授時，二是由地震觀測設(shè)備自帶的晶振計時。一般而言，臺站定期接收GPS授時信號，時鐘穩(wěn)定，然而GPS信號會受到雷雨等不良天氣或是GPS模塊故障的影響導(dǎo)致無法通過外部授時信號校正時間，此時只能依靠內(nèi)部晶振計時[1]，而晶振易受溫度、濕度等因素影響導(dǎo)致計時偏差。因此，鐘差是實際工作中經(jīng)常遇到的問題，如何判別鐘差是一項重要工作。

目前，鐘差判別通常是基于地震事件定位得出臺站走時殘差。此方法需人工標(biāo)注震相，依賴于工作人員的經(jīng)驗和地震事件的發(fā)生以及地震事件的大小。由于遠震震相往往不夠清晰導(dǎo)致無法準(zhǔn)確標(biāo)注，因此一般只適用于近震。同時，該方法只能判別出地震事件發(fā)生時有無鐘差，無法得知鐘差起止時間。可見利用地震時間判別鐘差的方法所受限制因素較多。

Derode等[2-3]發(fā)現(xiàn)臺站間背景噪聲的互相關(guān)格林函數(shù)與真實格林函數(shù)只有幅值上的差異，而在相位上是一致的。而在鐘差研究中，我們關(guān)注的核心也是相位變化，由此可利用背景噪聲互相關(guān)研究鐘差。Stehly等[4]用11年的連續(xù)記錄提取面波格林函數(shù)，用噪聲互相關(guān)方法分析了位于加利福尼亞南部的三個寬頻臺站對（PAS-PFO、GSC-PAS、PFO-GSC） 的互相關(guān)函數(shù)走時信息，計算出臺站絕對鐘差；Sch?nfelder等[5]在對Merapi火山臺陣的研究中，通過重建格林函數(shù)，構(gòu)建相對鐘差超定方程組求得臺陣中四個臺站的絕對鐘差。P.Gouédard等[6]分別利用基于背景噪聲互相關(guān)技術(shù)的時間對稱分析法和虛擬重復(fù)地震法校正OBS（海底地震儀）時鐘，結(jié)果表明OBS時鐘不論是線性還是非線性漂移都可以通過背景噪聲互相關(guān)計算得到。鄭宏等[7]基于背景噪聲互相關(guān)技術(shù)校正雅浦俯沖帶OBS時鐘，分別進行了線性校正和考慮鐘差閉合性的計算。

本文利用背景噪聲互相關(guān)技術(shù)就2018年福建省地震臺網(wǎng)單臺鐘差進行計算，選取每天的格林函數(shù)作為短期格林函數(shù)。

1 數(shù)據(jù)及預(yù)處理

本文利用福建臺網(wǎng)88個測震臺2018年的記錄進行研究，選取三分量中的豎向分量。88個臺站在全省分布大體均勻（圖1），臺間距主要集中在 50～300km，最大超過 500km（圖 2）。

圖1 臺站空間分布圖Fig.1 Station spatial distribution map

圖2 臺間距分布Fig.2 Station spacing distribution

臺站原始地震記錄中包含一些對于背景噪聲來說是干擾的信號，如地震信號、儀器響應(yīng)等，為得到純凈的噪聲信號，需要對原始數(shù)據(jù)進行預(yù)處理，預(yù)處理步驟如下。

1.1 重采樣

本文研究的臺站其采樣率為100Hz，采樣率大不利于數(shù)據(jù)的處理及讀取存儲，以減少數(shù)據(jù)點，因此本文進行采樣率為5Hz的重采樣。

1.2 去均值去趨勢

臺站記錄會受到基線漂移、突跳以及儀器本身的長周期振動的影響，需要去除這類影響，通過去均值處理可以達到目的。去趨勢可以消除地震儀器獲取數(shù)據(jù)時產(chǎn)生的偏移對后期計算產(chǎn)生的影響。

1.3 帶通濾波

地震信號中不同頻率成分包含不同信息。噪聲大體上可分為高頻（1Hz以上）、低頻（1～10s） 及長周期（10～50s） 三部分[8]，其中高頻噪聲主要受人類活動影響[9]，低頻是由于海洋波與海底或海岸線的非線性作用引起的[10]，長周期主要是風(fēng)、急流等自然因素引起[11]。本文利用來源于海洋性質(zhì)穩(wěn)定的低頻成分，選取1～10s帶通濾波。

1.4 時域歸一化

為剔除記錄中夾雜的地震信號，需做時域歸一化，本文選取滑動窗絕對均值歸一化法。給定一個時間序列dj，對于計算時間點n處的歸一化權(quán)重，計算如下式（1）：

圖3波形數(shù)據(jù)預(yù)處理Fig.3 Waveform data preprocessing

圖3 中，a為SHLX臺2018年11月26日臺灣海峽地震豎向分量原始記錄波形；b為經(jīng)重采樣、去均值去趨勢后的波形；c為1～10s帶通濾波后的波形；d為經(jīng)滑動窗絕對均值歸一化后的噪聲波形圖。

至此，經(jīng)過以上四個步驟處理后，數(shù)據(jù)預(yù)處理即完成，得到的噪聲記錄可以用于互相關(guān)計算。

2 方法及其原理

2.1 通過互相關(guān)函數(shù)獲取走時信息

一個彌散場中臺站A，B的背景噪聲記錄在時域下的互相關(guān)函數(shù)C（τ）計算式如下式：

式中u、u*分別為A，B兩臺記錄到的噪聲。

噪聲源的空間分布及噪聲水平只會影響互相關(guān)格林函數(shù)的振幅[13]。理論上，對于臺站A與B，不論其所在場地噪聲源是否為均勻分布，AB臺站對得到的噪聲互相關(guān)函數(shù)峰值到時總是對稱的，背景噪聲互相關(guān)函數(shù)與格林函數(shù)在相位上一致，只有幅值上的差異。因此，可利用背景噪聲互相關(guān)函數(shù)提取格林函數(shù)。本文選用臺站垂直分量計算背景噪聲格林函數(shù)，垂直分量記錄的是瑞利面波，即提取的是瑞利面波格林函數(shù)。

圖4不同鐘差情況時格林函數(shù)圖Fig.4 Green’s function graph for different clock conditions

圖4 中，a為噪聲源均勻分布，不存在相對鐘差時的互相關(guān)格林函數(shù)；b為存在正相對鐘差時的互相關(guān)格林函數(shù)；c為存在負相對鐘差時的互相關(guān)格林函數(shù)。虛線代表走時對稱軸，具有相對鐘差時，對稱軸發(fā)生偏移。

臺站間噪聲做互相關(guān)，其結(jié)果包含臺站間相對鐘差信息（圖5），為提取鐘差信息，需要無鐘差的格林函數(shù)作為參考。鐘差是偶爾存在的，長期的噪聲互相關(guān)結(jié)果疊加可以忽略鐘差影響，且隨疊加時長的增加，互相關(guān)格林函數(shù)的信噪比會增大，因此本文選取一年疊加（本文將噪聲以不重疊的半小時時窗劃分，一年疊加即是17520個窗疊加）結(jié)果作為參考格林函數(shù)。將參考格林函數(shù)與短期互相關(guān)格林函數(shù)再做互相關(guān)，得到的互相關(guān)函數(shù)峰值的位置反映臺站間的相對鐘差值。左右分支疊加能提高信噪比，本文以兩分支疊加的結(jié)果計算相對鐘差。圖5中峰值位置與0s的偏移值即為鐘差值，不同顏色代表不同時段的結(jié)果。

圖5 互相關(guān)函數(shù)左右分支疊加后結(jié)果Fig.5 Cross-correlation function left and right branches superimposed results

2.2 求解單臺鐘差

由2.1節(jié)得到各臺站間的相對鐘差，臺站i，j間的相對鐘差與絕對鐘差可表示為：

其中ΔSi、ΔSj分別為臺站i，j的絕對鐘差，ΔSi-j為臺站i，j間的相對鐘差。所有臺站間的相對鐘差關(guān)系可整合成矩陣形式，如式（4），本文方法的核心即是構(gòu)建此矩陣。

此矩陣簡記為Gm=d，其中m是由ΔS1、ΔS2……組成的列向量，即臺站1、臺站2……的絕對鐘差。d是由ΔDS1-S2、ΔDS1-S3……組成的列向量，即臺站1與2的相對鐘差、臺站1與3的相對鐘差......。88個臺站共組成3828個臺站對，因此m的行數(shù)為3828行，列數(shù)為88列，即臺站個數(shù)。此方程組為超定方程組，系數(shù)矩陣非方陣，所以m不可逆，求解此方程只能借助系數(shù)矩陣的廣義逆，即m=G+d，尋找其滿足最小二乘的解。

3 單臺鐘差結(jié)果

3.1 每日鐘差結(jié)果

本文以不重疊的，窗長為半小時的時窗將噪聲記錄進行劃分，半小時窗做互相關(guān)得到的格林函數(shù)信噪比約為1～3，考慮到其信噪比不夠高會影響計算結(jié)果準(zhǔn)確性。信噪比與疊加時長有關(guān)，疊加時長越長信噪比越大。為提到信噪比，采用一天內(nèi)48個半小時窗互相關(guān)格林函數(shù)疊加作為每天的互相關(guān)格林函數(shù)，信噪比最大達到21左右。下面以全年無鐘差臺站AXCK及CTCX、PTDT、ZPLA三個存在明顯鐘差的臺站為例進行分析，其鐘差計算結(jié)果如圖6所示，其中紅、藍、綠、紫紅色虛線分別表示1，2，3，4號地震發(fā)生時間，灰色部分為誤差棒，代表標(biāo)準(zhǔn)差。

圖6 每日鐘差計算結(jié)果Fig.6 Calculated daily clock difference result

驗證鐘差計算結(jié)果的準(zhǔn)確性，需要與地震事件得到的走時殘差對比。選取四個地震事件，震中位于島陸和近海，對于此類地震，單純型法定位優(yōu)于HYOSAT、LOCSAT法[14]。因此本文用單純型法定位獲得臺站走時殘差（表1）。

表1 地震事件與走時殘差

圖6中，四個臺都在7月存在斷記，可能是由于供電設(shè)備故障、數(shù)采故障、地震計無輸出、通信鏈路中斷等造成，但短時間斷記不影響本文研究結(jié)果。

對于不存在鐘差的臺站AXCK，鐘差變化幅度不大，年平均鐘差在0.0092s±0.1157s范圍內(nèi)波動，無較大突跳，時鐘良好。對于存在明顯鐘差的臺站，每日鐘差計算結(jié)果與地震事件得到的走時殘差存在一定偏差，但鐘差的線性漂移及突跳都比較清晰地顯現(xiàn)出來。如PTDT臺在11月23日，每日鐘差計算結(jié)果為+1.3s，與走時殘差+2.0s相差0.7s；在4月份，PTDT臺存在約-0.015s/天的線性漂移。

3.2 方法改進

用本文方法求解一個臺站數(shù)眾多的網(wǎng)絡(luò)中的單臺絕對鐘差，結(jié)果的準(zhǔn)確性會受限于互相關(guān)格林函數(shù)的信噪比，信噪比越大，計算結(jié)果越準(zhǔn)確。信噪比除了與疊加時長有關(guān)外，還與臺間距相關(guān)，臺間距越小信噪比越大。短期互相關(guān)格林函數(shù)取疊加時長為一天時，臺間距50km以下信噪比達到21，臺間距300km以上信噪比不及9。臺網(wǎng)中有超過500對臺站間距大于300km，其互相關(guān)格林函數(shù)信噪比不夠高，求得的相對鐘差值不夠準(zhǔn)確，會對超定方程產(chǎn)生擾動，進而使單臺鐘差結(jié)果存在偏差。

因此本文做出如下改進，為提高臺站對信噪比整體水平，選取CTCX臺站方圓50km范圍內(nèi)7個臺站組成小型網(wǎng)絡(luò)進行計算，最大臺間距85km，最小24km，平均間距48km，取每天的疊加結(jié)果作為短期互相關(guān)格林函數(shù)，此臺站網(wǎng)絡(luò)平均信噪比約21。求解超定方程組得單臺鐘差結(jié)果，CTCX鐘差如圖7所示。

圖7 CTCX臺鐘差結(jié)果Fig.7 Station CTCX clock difference result

2018年5月25日福建漳平地震CTCX臺走時殘差為-1.0s，CTCX臺鐘差約-0.9s，兩者對應(yīng)較好；同樣6月6日臺灣海峽南部地震該臺站走時殘差為+2.7s，計算鐘差約+2.5s，相對誤差不及8%，在9月初—12月末是時鐘正常的時段，平均鐘差為-0.1187s±0.3157s，是可以接受的誤差。由此可見，選用間距較近的臺站后，得到的結(jié)果更加準(zhǔn)確。

本文求得福建臺網(wǎng)88個測震臺2018年全年的鐘差變化?？傮w來看，將時間服務(wù)存在故障的臺站數(shù)按月統(tǒng)計，發(fā)現(xiàn)8月與10月兩個月存在鐘差的臺站最多，均為15個，鐘差臺站數(shù)占比不足20%，多數(shù)臺站時鐘正常。時鐘偏差有線性漂移與鐘差突跳兩種形式，多數(shù)存在鐘差的臺站，在鐘差突跳或線性漂移結(jié)束后都會恢復(fù)到正常狀態(tài)，原因是GPS授時恢復(fù)正?；蚴桥_站受人工維護。

4 討論與結(jié)論

圖6中AXCK臺是不存在鐘差的臺站，但是其計算結(jié)果依然存在小幅度波動。王俊等[15]在對江蘇省臺站鐘差計算時同樣發(fā)現(xiàn)計算結(jié)果的小幅波動，進行分析后發(fā)現(xiàn)江蘇地區(qū)夏季噪聲源能量優(yōu)勢方向來自太平洋，而冬季則來自于印度洋?？紤]到福建省地理位置與江蘇省類似，同為東南沿海，同樣會受到太平洋與印度洋活動的影響，推測這可能是引起計算結(jié)果小幅波動的原因。AXCK臺年平均鐘差在±0.15s內(nèi)，在可以接受的范圍，不影響對臺站時間偏差的評估。在驗證鐘差計算結(jié)果時，最好的辦法是與GPS對鐘日志進行對比。遺憾的是，GPS對鐘日志未被長期保存，在緩存中被不斷地刷新，所以無法得到2018年的GPS對鐘記錄，只能將單純型法定位得到的走時殘差用以驗證，雖然精度無法與GPS對鐘記錄相比，但由于地震事件的發(fā)生具有隨機性，所以由地震事件定位得到的走時殘差還是可信的。

在計算鐘差的過程中發(fā)現(xiàn)互相關(guān)格林函數(shù)信噪比對于鐘差計算結(jié)果的準(zhǔn)確性至關(guān)重要，如何增大信噪比是提高鐘差計算準(zhǔn)確性的關(guān)鍵。本文嘗試將大的臺站網(wǎng)絡(luò)分為若干個小的網(wǎng)絡(luò)，如3.2節(jié)以CTCX臺為中心，方圓50km內(nèi)的臺站組成一個小型網(wǎng)絡(luò)，使得臺網(wǎng)內(nèi)互相關(guān)格林函數(shù)信噪比整體水平更高，結(jié)果更好。除了選用臺間距較近的臺站對外，在今后的研究中還可以通過尋找良好的加權(quán)方法來提高信噪比。

目前運用此方法能實現(xiàn)精度相對較高的鐘差監(jiān)測，將其整合到地震監(jiān)測系統(tǒng)上，能提高臺網(wǎng)時間服務(wù)質(zhì)量，實現(xiàn)近實時的鐘差監(jiān)測，而且通過歷史地震數(shù)據(jù)可以計算任何時間段任何臺站的鐘差，不受傳統(tǒng)方法的局限性制約。對位于高緯度、山區(qū)等GPS信號接收條件不良地區(qū)的臺站具有重要的現(xiàn)實意義。