田 平,蔣 策,陳 瑛,梁 明
(1.廣東省地震局,廣州 510070;2.中國(guó)地震局地震監(jiān)測(cè)與減災(zāi)技術(shù)重點(diǎn)實(shí)驗(yàn)室,廣州 510070)
廣東省新豐江水庫(kù)自1959 年截流蓄水以來(lái),水庫(kù)區(qū)內(nèi)的地震活動(dòng)明顯增加,并在1962 年觸發(fā)了6.1 級(jí)地震,成為華南地區(qū)地震活動(dòng)最為活躍的地區(qū)之一。為研究庫(kù)區(qū)的地震活動(dòng)特征及發(fā)震機(jī)制,廣東省地震局在新豐江地區(qū)開(kāi)展了大量的研究工作,不斷加密測(cè)震臺(tái)網(wǎng),目前該地區(qū)的水平定位誤差小于1千米,達(dá)到百米級(jí)別的定位能力[1]。水庫(kù)區(qū)密集的測(cè)震臺(tái)站可以高質(zhì)量記錄到區(qū)域內(nèi)發(fā)生的大量小震微震事件,對(duì)這些地震活動(dòng)及其時(shí)空分布特征進(jìn)行精確分析,完善地震目錄,有利于對(duì)新豐江地區(qū)的發(fā)震斷層結(jié)構(gòu)的精細(xì)探測(cè)、水庫(kù)蓄水對(duì)地震的觸發(fā)作用以及地震預(yù)測(cè)等方面研究工作的開(kāi)展。
波形模板匹配(matched filter)方法是一種有效的微震檢測(cè)方法,該方法將模板地震的各分量波形與連續(xù)波形進(jìn)行互相關(guān)計(jì)算,通過(guò)對(duì)互相關(guān)波形進(jìn)行疊加來(lái)判斷地震事件[2-4]。波形模板匹配方法對(duì)微弱的地震信號(hào)敏感,即使在低信噪比的情況下也可以檢測(cè)到震級(jí)極小的微震事件,廣泛應(yīng)用于地震前震、余震檢測(cè)[5-8]、觸發(fā)地震[9-10]以及非天然地震檢測(cè)[11-12]等方面。Zhang 等[13]基于波形模板匹配方法發(fā)展了匹配定位方法(Match and Locate,簡(jiǎn)稱(chēng)M&L 方法),在相關(guān)波形疊加前考慮待檢測(cè)事件與模板事件之間的相對(duì)位置,在檢測(cè)微震的同時(shí),獲取高精度的發(fā)震位置信息。Zhang等[14]利用該方法檢測(cè)了日本Ontake 火山2007 和2014 年兩次噴發(fā)前的密集微震活動(dòng),得到30 余倍日本氣象廳的地震目錄結(jié)果;馬曉靜等[15]將該方法應(yīng)用于2014 年6 月河源地區(qū)的小震群事件研究中,結(jié)果顯示模板檢測(cè)到的地震事件數(shù)目相較于地震目錄結(jié)果增加了一倍;王志偉等[16]利用匹配定位方法對(duì)2008 年10 月至2011 年7 月間重慶榮昌地區(qū)的微震進(jìn)行檢測(cè)與定位,將該地區(qū)地震目錄的震級(jí)下限從ML1.0 降低至ML0.3;Liu 等[17]在匹配定位方法的基礎(chǔ)上,開(kāi)發(fā)了GPU 版本的匹配定位計(jì)算程序,顯著地提高了模板匹配計(jì)算的速度;杜瑤等[18]在水庫(kù)區(qū)爆破事件的檢測(cè)研究中,發(fā)現(xiàn)匹配定位方法可以有效區(qū)分不同類(lèi)型的地震活動(dòng),進(jìn)而提高對(duì)疑爆事件的識(shí)別率。
在本研究中,我們利用M&L 方法檢測(cè)新豐江河源地區(qū)的微震事件,通過(guò)數(shù)字地震臺(tái)網(wǎng)中心數(shù)據(jù)處理系統(tǒng)建立適合新豐江地區(qū)的微震模板事件庫(kù),實(shí)現(xiàn)對(duì)地震臺(tái)網(wǎng)記錄的連續(xù)波形的自動(dòng)檢測(cè),將其應(yīng)用于臺(tái)網(wǎng)的地震目錄完善工作。
M&L 方法不需要假定待檢測(cè)地震與模板地震發(fā)生在同一位置。該方法將各臺(tái)站分量上記錄的模板波形與包含潛在地震信號(hào)的連續(xù)波形進(jìn)行互相關(guān)計(jì)算,在對(duì)相關(guān)波形進(jìn)行疊加前,考慮待檢測(cè)地震與模板地震發(fā)生位置之間的差異,在模板地震周?chē)娜S空間進(jìn)行網(wǎng)格搜索,計(jì)算模板地震與每個(gè)假定待檢測(cè)地震位置之間的走時(shí)差,在此基礎(chǔ)上對(duì)互相關(guān)波形進(jìn)行走時(shí)矯正和疊加,計(jì)算疊加后的互相關(guān)波形中的平均互相關(guān)(Cross-Correlation,CC)系數(shù)和信噪比(Signal-Noise ratio,SNR)。
兩個(gè)距離很近的地震波形之間的歸一化互相關(guān)可表示為:
互相關(guān)系數(shù)與地震波形的相似度成正比,可以作為事件檢測(cè)的閾值標(biāo)準(zhǔn)。但僅使用較低互相關(guān)系數(shù)作為閾值時(shí),盡管可以增加檢測(cè)事件的數(shù)量,同時(shí)也會(huì)極大提高了事件的誤檢概率;信噪比的引用可以對(duì)其進(jìn)行約束,在較低閾值標(biāo)準(zhǔn)設(shè)定的情況下,檢測(cè)到更多地震,同時(shí)避免增加誤檢概率。
M&L 方法將平均互相關(guān)系數(shù)和信噪比二者聯(lián)立作為地震事件的判斷標(biāo)準(zhǔn)。當(dāng)其值超過(guò)設(shè)定的相關(guān)系數(shù)閾值和信噪比閾值時(shí),即判斷為檢測(cè)到地震事件。之后,根據(jù)所有臺(tái)站三分量波形與參考事件波形的振幅比的中位數(shù)確定地震事件的震級(jí)。
2018—2019 年期間,廣東測(cè)震臺(tái)網(wǎng)記錄到新豐江地區(qū)發(fā)生的天然地震共2967件(圖1),其中震級(jí)ML≤2.0 的地震事件共2935 件,震級(jí)ML≤1.0 的地震事件共2625 件,微震小震事件占新豐江地區(qū)地震事件總數(shù)的88%以上。
圖1 2018—2019年間新豐江地區(qū)地震目錄Fig.1 Earthquake catalog in Xinfengjiang area from 2018 to 2019
上述地震目錄全部為人工瀏覽連續(xù)地震波形識(shí)別、編目得到,由于人工識(shí)別誤差、震級(jí)未達(dá)到編目要求(ML<0.0)等原因,還存在大量微震事件未被記錄在內(nèi)。本研究將地震目錄中的部分地震事件作為地震模板,檢測(cè)2018—2019 年間新豐江河源地區(qū)的地震事件,將其與人工地震目錄結(jié)果比較,測(cè)試M&L方法在新豐江河源地區(qū)的適用性。
M&L 方法的檢測(cè)結(jié)果依賴(lài)于地震模板的豐富程度,無(wú)法利用有限的模板識(shí)別出所有的微震事件。地震模板的種類(lèi)和數(shù)量越豐富,檢測(cè)的結(jié)果越豐富越接近真實(shí)情況。相應(yīng)地,模板的大量使用也極大增加計(jì)算時(shí)間,進(jìn)而降低檢測(cè)效率。在使用模板匹配檢測(cè)方法時(shí),一般將研究時(shí)間范圍內(nèi)的所有事件進(jìn)行篩選作為地震模板,對(duì)連續(xù)波形進(jìn)行掃描。新豐江庫(kù)區(qū)地震數(shù)量較多,為測(cè)試在該地區(qū)模板檢測(cè)的能力及效率,對(duì)使用的地震模板進(jìn)行了人為震級(jí)篩選。最終確定選定2019 年間新豐江河源地區(qū)(23.700°~23.766°N,114.57°~114.75°E),震級(jí)范圍0.5≤ML≤1.0,共125條地震事件作為地震模板(圖2中紅點(diǎn)所示)。新豐江水庫(kù)區(qū)測(cè)震臺(tái)站較為密集,經(jīng)過(guò)對(duì)臺(tái)站數(shù)據(jù)的記錄質(zhì)量以及相對(duì)研究區(qū)域的位置分布,選取其中的10 個(gè)臺(tái)站(圖4)記錄的數(shù)據(jù)進(jìn)行分析處理。M&L 方法受地殼速度模型影響較小,本文采用PREM 速度模型[20],利用TauP 工具包[23]計(jì)算理論震相到時(shí)。微震事件的判斷閾值參考采用9 倍背景平均相關(guān)系數(shù)。另外,新豐江水庫(kù)區(qū)的地震深度普遍較淺,深度分布比較集中,為提高檢測(cè)效率,在地震檢測(cè)中沒(méi)有對(duì)震源深度進(jìn)行搜索,固定深度為模板地震的深度。
圖2 2018—2019年間利用模板檢測(cè)得到的地震事件M-T分布結(jié)果Fig.2 M-T distribution results of earthquake events obtained by template detection from 2018 to 2019
通過(guò)對(duì)連續(xù)波形的模板匹配掃描,共檢測(cè)到2018—2019 年間新豐江河源地區(qū)的地震事件共2160 條(圖2),最大震級(jí)ML3.08(圖3),最小震級(jí)ML0.8。對(duì)檢測(cè)結(jié)果的震級(jí)分布進(jìn)行統(tǒng)計(jì)分析,震級(jí)0.0≤ML≤1.0 的微震事件共625 條,占地震事件總數(shù)的28.9 %;震級(jí)ML<0.0 的地震事件共1438 條,占地震事件的66.5%。檢測(cè)到的地震事件在空間上分布在地震模板的發(fā)震位置附近。2018—2019 年地震目錄結(jié)果(圖4)顯示新豐江河源區(qū)域記錄地震共1764 件,震級(jí)0.0≤ML<0.5 的微震事件共1163條,震級(jí)0.5≤ML≤1.0 的微震事件共425 條,震級(jí)1.0≤ML≤2.0 的地震事件共159 條,ML>2.0 的地震事件17條。
圖3 模板檢測(cè)結(jié)果示例(最大震級(jí)為ML 3.08,發(fā)震時(shí)刻2018-07-25 08:53:45)Fig.3 Example of template detection results(the maximum magnitude is ML 3.08,which occurred at 08:53:45 2018-07-25)
上述僅利用少量震級(jí)0.5≤ML≤1.0 的地震模板,對(duì)2018-2019年間新豐江河源地區(qū)的地震事件的檢測(cè)結(jié)果表明,地震模板事件全部被自檢出,相關(guān)系數(shù)接近于1.0。將這些事件進(jìn)行比較,發(fā)現(xiàn)檢測(cè)結(jié)果的發(fā)震時(shí)刻和震級(jí)與編目結(jié)果具有較小的差異,具有較高的準(zhǔn)確度。但是,由于缺少2018 年0.5≤ML≤1.0 的地震模板,對(duì)2018 年期間檢測(cè)出的震級(jí)0.5≤ML≤1.0 的地震事件僅92 件,少于廣東臺(tái)網(wǎng)編目結(jié)果300條,M&L方法的檢測(cè)結(jié)果受地震模板的豐富程度影響較大,僅利用部分模板無(wú)法對(duì)所有地震事件進(jìn)行檢測(cè)。
上述的模板檢測(cè)結(jié)果表明,M&L 方法適用于新豐江水庫(kù)區(qū)的微震檢測(cè),而且可以對(duì)檢測(cè)到的微震事件進(jìn)行定位以及震相標(biāo)記,可將其應(yīng)用于測(cè)震臺(tái)網(wǎng)的地震編目工作中。在數(shù)字臺(tái)網(wǎng)中心的服務(wù)器上進(jìn)行軟件部署,利用數(shù)字臺(tái)網(wǎng)中心數(shù)據(jù)處理軟件系統(tǒng)JOPENS[21-22]獲取地震事件以及連續(xù)波形,實(shí)時(shí)擴(kuò)充地震模板庫(kù),并且對(duì)其調(diào)用在連續(xù)波形上自動(dòng)進(jìn)行微震檢測(cè),達(dá)到完善地震目錄的目的。
根據(jù)M&L 模板匹配計(jì)算程序的處理流程(圖5),為達(dá)到對(duì)地震波形的自動(dòng)檢測(cè)的要求,需要改進(jìn)的主要部分為流程中最主要的輸入部分:即“連續(xù)波形”和“地震模板”。
圖5 模板檢測(cè)程序的流程Fig.5 Process of template detection program
對(duì)于“連續(xù)波形”部分,通過(guò)數(shù)據(jù)處理系統(tǒng)獲取到連續(xù)波形數(shù)據(jù)后,需要對(duì)其進(jìn)行格式轉(zhuǎn)換,完成波形的拼接、濾波以及降采樣等數(shù)據(jù)處理工作。我們利用JOPENS 的數(shù)據(jù)導(dǎo)出功能,將其設(shè)置為定時(shí)將新入庫(kù)的連續(xù)波形數(shù)據(jù)導(dǎo)出為seed 文件,保存在臺(tái)網(wǎng)中心服務(wù)器端。檢測(cè)程序部署在信息中心機(jī)群,利用腳本從臺(tái)網(wǎng)中心服務(wù)器下載波形數(shù)據(jù)seed 文件,并自動(dòng)將其轉(zhuǎn)換為sac 格式數(shù)據(jù)文件;對(duì)連續(xù)波形進(jìn)行1~20 Hz 的濾波以及降采樣處理工作,統(tǒng)一命名放置于M&L 檢測(cè)程序內(nèi)部的Trace目錄,每日對(duì)其進(jìn)行更新。
對(duì)于“檢測(cè)模板”部分,臺(tái)網(wǎng)編目報(bào)告的地震事件seed 文件以及地震目錄catalog 文件經(jīng)數(shù)據(jù)處理系統(tǒng)導(dǎo)出上傳至中心機(jī)群,放置于目標(biāo)目錄下。系統(tǒng)定時(shí)在該目錄下進(jìn)行搜索,發(fā)現(xiàn)地震事件seed 文件并將其轉(zhuǎn)換為sac 格式數(shù)據(jù)文件,對(duì)事件波形進(jìn)行1-20 Hz 的濾波以及降采樣處理,保證與連續(xù)波形數(shù)據(jù)格式相同。經(jīng)上得到的新的地震模板文件,統(tǒng)一命名放置于M&L 檢測(cè)程序內(nèi)部的Template目錄,實(shí)時(shí)對(duì)地震模板庫(kù)進(jìn)行更新。
上述利用震級(jí)范圍0.5≤ML≤1.0 的地震模板進(jìn)行檢測(cè)的結(jié)果表明,僅利用部分地震模板不能檢測(cè)出潛在的所有地震,地震模板的豐富程度對(duì)檢測(cè)結(jié)果產(chǎn)生非常大的影響。在保證檢測(cè)效率的同時(shí),為了檢測(cè)盡量多的微震事件,在部署自動(dòng)檢測(cè)軟件中,將待檢測(cè)長(zhǎng)度為24 小時(shí)的連續(xù)波形包括其前24小時(shí)以及后24小時(shí),共72小時(shí)內(nèi)臺(tái)網(wǎng)編目的所有地震事件作為模板,對(duì)連續(xù)波形進(jìn)行檢測(cè)。若此時(shí)間段內(nèi)不存在地震事件,自動(dòng)將調(diào)取時(shí)間范圍進(jìn)行延長(zhǎng);另一方面,系統(tǒng)也可隨時(shí)根據(jù)需要調(diào)用數(shù)據(jù)庫(kù)中的其他地震模板,以保證充足的地震模板進(jìn)行微震檢測(cè)。自相關(guān)系數(shù)閾值等參數(shù)設(shè)置采用與上文相同設(shè)置。
基于對(duì)程序的上述配置,重新檢驗(yàn)其在新豐江地區(qū)地震檢測(cè)中的適用性。2023 年2 月11 日10時(shí)41分廣東省河源市源城區(qū)發(fā)生了M4.3地震,震后發(fā)生一系列余震,廣東數(shù)字地震臺(tái)網(wǎng)中心的編目結(jié)果顯示,截至2月12日0時(shí),共記錄到該地區(qū)余震171 件,其中,0≤ML<1.0,143 件;1.0≤ML<2.0,23件;ML≥2.0,5件;最大余震ML3.4。
利用上述改進(jìn)的程序配置,對(duì)其在該時(shí)間段內(nèi)對(duì)連續(xù)波形的自動(dòng)檢測(cè)結(jié)果進(jìn)行分析。結(jié)果顯示,截至2023 年2 月12 日0 時(shí),自動(dòng)檢測(cè)配置共檢測(cè)到河源M4.3 地震的余震265 件,與人工余震目錄相比,新檢測(cè)到的地震事件中,有15 件事件的震級(jí)在ML0.0 左右,最大的檢測(cè)震級(jí)為ML0.64;其余震級(jí)均在ML<0.0(圖6),最小的震級(jí)為ML-0.65。
圖6 檢測(cè)出的地震事件(ML 0.37)Fig.6 One detected earthquake event with magnitude of ML0.37
地震目錄記錄的171件地震事件,共檢測(cè)出其中的167件,經(jīng)人工對(duì)波形校對(duì)發(fā)現(xiàn),未檢測(cè)出的4 件地震事件均為雙震事件,震級(jí)在0≤ML<1.0 范圍。通過(guò)對(duì)模板的檢查發(fā)現(xiàn),雙震中的兩起事件發(fā)生間隔時(shí)間過(guò)短(圖7),生成的模板同時(shí)包含了兩起事件,在連續(xù)波形計(jì)算過(guò)程中對(duì)另一事件無(wú)法檢測(cè),導(dǎo)致對(duì)雙震事件的檢測(cè)失效。
圖7 未檢測(cè)出的雙震事件Fig.7 Undetected double-earthquake event
將地震模板中被重新檢測(cè)出中被地震目錄收錄的167件地震事件的發(fā)震時(shí)刻、震級(jí)與人工結(jié)果進(jìn)行比較(圖8),發(fā)現(xiàn)大部分震級(jí)與人工編目的震級(jí)結(jié)果相近,差值一般小于ML0.5;發(fā)震時(shí)刻與人工編目的發(fā)震結(jié)果相近,差值小于1.0 s。
圖8 編目報(bào)告與檢測(cè)結(jié)果中同事件的震級(jí)(a)與發(fā)震時(shí)刻(b)的差值統(tǒng)計(jì)比較Fig.8 Statistical comparison of the difference between the magnitude(a)and the occurrence time(b)of the same event in the catalog report and the detection results
利用上述配置對(duì)河源M4.3 的余震自動(dòng)檢測(cè)結(jié)果顯示,多個(gè)雙震事件沒(méi)能被掃描識(shí)別出來(lái)。經(jīng)過(guò)對(duì)地震目錄的分析,發(fā)現(xiàn)這幾次雙震事件發(fā)生的時(shí)間間隔很短,普遍在2 s 以?xún)?nèi)。由于系統(tǒng)的配置未對(duì)地震事件進(jìn)行篩選,全部處理為地震模板并應(yīng)用于模板匹配,導(dǎo)致部分單個(gè)地震模板內(nèi)包含了多次地震事件。在對(duì)波形進(jìn)行掃描時(shí),直接將其與連續(xù)波形進(jìn)行滑動(dòng)互相關(guān)并輸出匹配結(jié)果,當(dāng)雙震的震級(jí)差比較大時(shí),輸出為掃描時(shí)間窗內(nèi)震級(jí)較大的地震事件。
在利用M&L 方法進(jìn)行匹配定位時(shí),一般需要對(duì)地震模板進(jìn)行篩選,尤其是多個(gè)地震事件的發(fā)震時(shí)間間隔較近,前一事件的s 震相淹沒(méi)在后一事件更大的p 震相中的情況。盡管在處理地震模板時(shí),可以控制地震模板的波形長(zhǎng)度,在一定程度上消除另一事件的影響,但這也可能導(dǎo)致其他地震模板的長(zhǎng)度過(guò)短,不能包含所有必需的震相信息;另外,雙震事件作為模板也將導(dǎo)致檢測(cè)出的事件的定位結(jié)果和震級(jí)產(chǎn)生較大的誤差,雙震以及多地震事件不適合作為地震模板進(jìn)行匹配定位檢測(cè)。在對(duì)地震模板庫(kù)進(jìn)行自動(dòng)更新,需結(jié)合地震目錄,剔除使用發(fā)震間隔較近、且發(fā)震位置也較近的地震事件。
本文與前人[15]在新豐江河源地區(qū)的微震檢測(cè)中,對(duì)于平均互相關(guān)系數(shù)的參數(shù)設(shè)定存在差異,這種差異對(duì)微震的檢測(cè)結(jié)果數(shù)量產(chǎn)生影響。在不同的研究區(qū)域,需要對(duì)研究區(qū)域內(nèi)擬使用的地震模板進(jìn)行自檢來(lái)估計(jì)平均背景相關(guān)系數(shù),以此確定檢測(cè)標(biāo)準(zhǔn)。數(shù)字臺(tái)網(wǎng)中心利用M&L 方法對(duì)不同區(qū)域的微震進(jìn)行檢測(cè)時(shí),需要重新確定檢測(cè)標(biāo)準(zhǔn)以保證檢測(cè)質(zhì)量。
M&L 方法在地震模板處理中,利用TauP 工具包[23]根據(jù)給定地殼速度模型標(biāo)記理論震相信息截取地震模板。通過(guò)計(jì)算模板事件與待檢測(cè)事件的可能位置在同一臺(tái)站上的走時(shí)差來(lái)進(jìn)行定位計(jì)算,對(duì)模型的依賴(lài)較小。本文利用PREM 模型生成的地震模板,發(fā)現(xiàn)其震相標(biāo)注與實(shí)際震相的到時(shí)差異較大,這對(duì)于檢測(cè)微震以及對(duì)其進(jìn)行定位的結(jié)果沒(méi)有影響。但是,地震目錄包含精確的震相標(biāo)記信息,可以利用這些震相信息替代理論到時(shí)標(biāo)記,對(duì)地震模板進(jìn)行截取、標(biāo)記處理,在匹配定位檢測(cè)的同時(shí),輸出比較準(zhǔn)確的震相信息,更加適合地震編目的完善工作。
另外,M&L 程序已經(jīng)發(fā)布GPU(Graphics Processing Unit-Based)版本,該程序的計(jì)算速度提高到M&L 程序的4.5倍,同時(shí)在地震模板的各分量波形上引入權(quán)重系數(shù)顯著提高檢測(cè)結(jié)果的穩(wěn)定性[13]。這可以在一定程度上緩解由于地震模板過(guò)多導(dǎo)致的檢測(cè)效率過(guò)低的問(wèn)題,適合臺(tái)網(wǎng)處理大量微震檢測(cè)工作。
本文利用匹配定位M&L 方法,通過(guò)對(duì)匹配定位程序與數(shù)據(jù)處理系統(tǒng)軟件進(jìn)行配置和參數(shù)設(shè)置,在新豐江河源地區(qū)建立地震模板數(shù)據(jù)庫(kù),進(jìn)行微震事件自動(dòng)檢測(cè),并對(duì)河源M4.3 地震震后余震的自動(dòng)檢測(cè)結(jié)果進(jìn)行分析,得到以下結(jié)論:
(1)在新豐江地區(qū),存在有大量未被記錄的微震事件,匹配定位方法可以準(zhǔn)確地對(duì)這些事件進(jìn)行檢測(cè)以及定位;
(2)利用數(shù)據(jù)處理系統(tǒng)軟件將數(shù)據(jù)庫(kù)內(nèi)的地震事件自動(dòng)轉(zhuǎn)換為地震模板,建立新豐江地區(qū)的地震模板庫(kù)并對(duì)其進(jìn)行實(shí)時(shí)擴(kuò)充,方便在微震檢測(cè)中進(jìn)行調(diào)用;
(3)根據(jù)對(duì)匹配定位程序的配置和參數(shù)設(shè)置,檢測(cè)出河源M4.3 地震震后大量ML<0.0 的地震事件,可以達(dá)到完善地震目錄的要求。
致謝:感謝張淼博士提供地震模板匹配程序(Match and Locate),本文部分圖件使用GMT繪制。