利用噪聲互相關(guān)技術(shù)計(jì)算福建臺(tái)網(wǎng)臺(tái)站絕對(duì)鐘差

胡夢(mèng)迪 ，金 星 ，李 軍

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

0 引言

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

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

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

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

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

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

圖1 臺(tái)站空間分布圖Fig.1 Station spatial distribution map

圖2 臺(tái)間距分布Fig.2 Station spacing distribution

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

1.1 重采樣

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

1.2 去均值去趨勢(shì)

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

1.3 帶通濾波

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

1.4 時(shí)域歸一化

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

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

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

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

2 方法及其原理

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

一個(gè)彌散場(chǎng)中臺(tái)站A，B的背景噪聲記錄在時(shí)域下的互相關(guān)函數(shù)C（τ）計(jì)算式如下式：

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

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

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

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

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

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

2.2 求解單臺(tái)鐘差

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

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

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

3 單臺(tái)鐘差結(jié)果

3.1 每日鐘差結(jié)果

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

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

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

表1 地震事件與走時(shí)殘差

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

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

3.2 方法改進(jìn)

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

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

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

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

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

4 討論與結(jié)論

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

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

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