陳家樑，張寶劍

(福建省地震局,福建福州，350003)

0 引言

福建地震臺網(wǎng)中心從1971年開始建設(shè)，從最初3個測震臺站的模擬觀測及人工震中定位，發(fā)展到2020年底88個測震臺的數(shù)字觀測臺網(wǎng)，得益于國家的發(fā)展壯大。福建地震臺網(wǎng)包含了測震臺、強震臺、烈度臺及GNSS臺網(wǎng)，臺站密度走在了全國前列。其中測震臺網(wǎng)中部署有寬頻帶地震計的臺站占據(jù)84%，為組建虛擬臺網(wǎng)我中心接入了鄰省臺站及16個臺灣臺，福建本省測震臺站的平均臺間距約37千米，震級的觀測范圍為0-6級。臺站數(shù)量的增加提高了地震定位的精度，而地震數(shù)據(jù)的可靠性決定了地位精度的準確性，因此有必要對地震數(shù)據(jù)質(zhì)量進行多方面的監(jiān)測。方位角應(yīng)用于包括SKS剪切波分裂，接收函數(shù)研究等多項科學(xué)研究，如果方位角存在偏差將會影響相關(guān)研究結(jié)果的可靠性。以往都是在相關(guān)研究出現(xiàn)異常時才會發(fā)現(xiàn)臺站有故障，這樣的故障維護就會出現(xiàn)延遲，且會造成該時間段內(nèi)的數(shù)據(jù)不可用。因此需要有專門軟件對臺站方位角進行實時監(jiān)測，發(fā)現(xiàn)問題及時處理，這樣不僅可以提高故障修復(fù)的時效性，且可以讓儀器維護人員的工作更具有針對性。

圖1 福建臺網(wǎng)及臺灣交換臺站分布

1 程序設(shè)計

方位角的檢測需要用到震中距30度以上的遠震，所以該程序為準實時系統(tǒng)。根據(jù)具體的操作步驟，程序分為地震事件選取模塊、數(shù)據(jù)預(yù)處理模塊、臺站屬性讀取模塊、最小能量計算模塊、相關(guān)系數(shù)計算模塊及虛擬臺站構(gòu)建模塊。為簡化代碼數(shù)量，減少程序運行中所占用內(nèi)存，本程序使用Python中的Obspy這個開源的地震數(shù)據(jù)處理框架。

1.1 地震事件選取模塊

本程序中計算方位角偏差的方法為利用P波質(zhì)點的運動特性法，需要震中距在30度以上的遠震，因此在準備地震事件時需要計算震中距是否達到30度。福建本省面積較小，對于遠震而言省內(nèi)不同臺站對于同一地震事件的震中不會有太大差別，因此我們選取位于福建較為中心的永安燕西作為所有臺站的參照臺，只要計算該臺與地震事件的震中距即可，這樣可以減少運算量又不影響數(shù)據(jù)的可靠性。計算震中距我們引入Obspy.geodetics中的gps2dist_azimuth及kilometers2degrees兩個方法，前者用于計算震中距，其返回值是一個數(shù)組，震中距為第一個元素且單位為千米，所以需要使用kilometers2degrees將單位從千米轉(zhuǎn)化為度。程序?qū)逊险鹬芯嘁蟮牡卣鹗录丶罢鹬芯噍敵鲋林付窂较碌腡XT文件中。

圖2 程序架構(gòu)

1.2 數(shù)據(jù)預(yù)處理

因?qū)掝l帶地震計中有60秒與120秒甚至更高的頻帶，為保證帶寬一致所以需要進行帶通濾波。臺站的采樣率為100HZ，為減少計算量程序采取降采樣處理。利用Obspy.Stream類循環(huán)讀取地震離線波形數(shù)據(jù)，全部讀取后利用該類中的resample（sampling_rate）方法即可進行降采樣處理。Obspy.Stream.filter可以對數(shù)據(jù)進行帶通濾波處理，用該方法離線波形數(shù)據(jù)進行帶通濾波處理的具體參數(shù)設(shè)置如下所示 filter(‘bandpass’,freqmin=0.02,freqmax=0.2),參數(shù)表示使用帶通濾波處理數(shù)據(jù)且設(shè)置了低頻與高頻的具體數(shù)值。

1.3 臺站屬性讀取

臺站屬性指的是臺站名、經(jīng)緯度及高程，這些數(shù)據(jù)均存儲于MySQL數(shù)據(jù)庫中。對于省級地震臺網(wǎng)中心，有固定的數(shù)據(jù)庫用于配置各臺站經(jīng)緯度、高程及靈敏度等參數(shù)，所以我們用xml文件配置好數(shù)據(jù)庫的地址等，程序就可以通過xml文件進入數(shù)據(jù)庫進行臺站屬性的讀取。Xml文件的內(nèi)容形式如下所示：

數(shù)據(jù)庫密碼

數(shù)據(jù)庫名

圖片存儲路徑

該模塊中我們創(chuàng)建一個專門讀取xml文件的類，在Python中引用xml.dom.minidom類包，就可以按照該類讀取xml文件的步驟方法，將xml文件中的各配置內(nèi)容讀取到指定的數(shù)組中并返回該數(shù)組。

獲取上述數(shù)據(jù)庫的配置信息后連接MySQL數(shù)據(jù)庫，寬頻帶地震計的通道命名有顯著要求，即開頭字母為B，其三通道的命名為BHN、BHE、BHZ。因此在選取寬頻帶臺站的相關(guān)屬性時對數(shù)據(jù)庫的操作語句可寫作“select Sta_code from Channel_info where Chn_code like ‘%BHE%’”,其 中Sta_code為臺站名，Channel_info為表名稱，Chn_code為通道存儲字段名。通過以上操作可以獲取全部寬頻帶地震計臺站名，而后運用循環(huán)方法獲取所有的臺站屬性后封裝到dictionary中返回。

1.4 最小能量計算

通過以上以上步驟計算每個臺站對應(yīng)每個地震事件的P波到時，選取到時前后五秒內(nèi)經(jīng)過旋轉(zhuǎn)的能量，再以1度為間隔進行180度的旋轉(zhuǎn)計算每一度對應(yīng)的能量，再進行對比選出最小能量對應(yīng)的旋轉(zhuǎn)度數(shù)，認為此時對應(yīng)的角度為與正北方向的角度偏差。再將此時的偏差角與經(jīng)過經(jīng)緯度經(jīng)過計算出的理論偏差值進行對比，即可得出此臺對應(yīng)該地震事件時間點上的角度偏差。將每個地震事件的偏差值輸出至String類型的數(shù)組，計算完所有事件的偏差值后，再利用“SNR”法計算所有事件中每一度的能量值后選出最小值對應(yīng)的偏差角，即認為是最終的地震計與正北方向的偏差值。

圖3 時間序列圖

為更清晰的顯示每個事件對應(yīng)的方位角偏差，在此步驟中將所得到的數(shù)據(jù)進行繪圖顯示。因利用了GMT繪圖軟件所以需要在主機上安裝GMT5.0以上版本，否則會出現(xiàn)繪圖命令錯誤或不兼容的情況。時間序列圖的縱軸單位為的形式，所以需要用到GMT中的轉(zhuǎn)義字符“@～\152\244\050\260\051@～”。將GMT繪圖命令以String形式寫入bat文件后，執(zhí)行該文件即可進行繪圖操作，繪制后的圖片將按xml文件中配置的路徑存儲。

1.5 虛擬臺陣

為更加精準的分析異常方位角偏差臺站，該程序還設(shè)計了異常臺站與正常臺站所構(gòu)建的虛擬臺陣對應(yīng)通道的相關(guān)系數(shù)分析模塊。此模塊以故障臺站為中心，計算所有臺站與此臺的臺間距，并選取臺間距最近的15個正常臺站構(gòu)成虛擬臺陣。線性疊加臺陣內(nèi)所有臺站相同通道的數(shù)據(jù)后，與故障臺相同通道進行相關(guān)系數(shù)計算，如果故障臺不出現(xiàn)儀器方向擺放錯誤或極性反轉(zhuǎn)等重大問題，則其相關(guān)系數(shù)會接近1相反則會接近0或負數(shù)。其相關(guān)系數(shù)的計算引入scipy.stats類中的pearsonr方法，可以迅速計算出兩個相同通道的相關(guān)系數(shù)。為了便于觀察，此模塊同樣具備輸出波形形狀及相關(guān)系數(shù)值的功能。

圖4 上圖為故障臺，下圖為正常臺與虛擬臺陣的相關(guān)系數(shù)波形圖

2 程序在福建臺網(wǎng)中的應(yīng)用

通過該程序?qū)Ω＝ǖ卣鹋_網(wǎng)所屬寬頻帶地震計進行計算，地震事件為2015年1月到2019年8月之間的遠震，得到74個寬頻帶地震計臺站中有72個臺的方位角偏差值在-7°～7°之間，這樣的偏差值在該方法的理論中可以忽略不計，因此福建臺網(wǎng)寬頻帶地震計臺站的方位角總體合格率在97%如圖5所示。其中兩個存在較大偏差值的臺站均為BBVS-60地震計,此類型地震計沒有明顯指北標志，所以導(dǎo)致維護人員在更換過程中無法準確的進行指北操作，經(jīng)過實地校正后已使方位角的偏差值減小到要求值1°以內(nèi)。

圖5 方位角偏差臺站分布

3 結(jié)論

通過對程序模塊化開發(fā)，可以理順每個模塊的具體功能，使程序整體研發(fā)更清晰，更易調(diào)試出現(xiàn)問題的部分代碼進行小范圍的修復(fù)，而不至于出現(xiàn)牽一發(fā)動全身式的更改。因程序為處理地震數(shù)據(jù)而開發(fā)，所以使用到Obspy庫，因其幾乎支持地震學(xué)界內(nèi)通常使用的所有波形格式的讀寫，可簡化大量的數(shù)據(jù)格式之間轉(zhuǎn)換的代碼，同時他還集成了大量的地震學(xué)及數(shù)學(xué)所用的專有庫，對于開發(fā)者只需要對其進行引入就可使用。程序應(yīng)用在福建臺網(wǎng)中心后，經(jīng)過實地校核發(fā)現(xiàn)該程序的計算結(jié)果可靠。該程序的應(yīng)用彌補了福建臺網(wǎng)中心在方位角檢測方面的不足，且提高了方位角偏差的檢測效率及修復(fù)故障的時效性。