考慮橢圓模型MW≥5.5地震震源參數時域計算方法

王江, 馬強*, 陶冬旺, 章旭斌, 解全才, 林德昕

1 中國地震局工程力學研究所地震工程與工程振動重點實驗室, 哈爾濱 150080

2 地震災害防治應急管理部重點實驗室, 哈爾濱 150080

0 引言

隨著密集臺網設施布設和地震監測技術系統的迅速發展,以地震預警(Earthquake Early Warning,簡稱EEW)為代表的實時防震減災技術得到了長足的發展,將由地震預報、地震危險性分析、結構抗震和應急救援組成的防震減災救災鏈條擴展到了震時階段.震時地震信息(發震時刻,地震位置,震級和地震動場)的快速獲取是地震預警系統的核心任務.現行EEW系統采用的策略是通過有限臺站、有限波形信息快速定位震源位置并估算地震震級(Kanamori,2005;Olson and Allen,2005;Wu and Zhao,2006;Zollo et al.,2006;馬強,2008;金星等,2012);當獲取發震位置和震級后,采用“點源”震源模型,將其代入地震動預測方程(Ground Motion Prediction Equation,簡稱GMPE)預測震中附近目標場地的地震動強度(Wu and Kanamori,2005,2008;Wu and Zhao, 2006;馬強,2008;B?se et al.,2012;Peng et al.,2014).預測地震動強度和空間分布的可靠性關乎EEW系統震時防震減災的成敗,直接決定地震警報的發布范圍和預警等級.而破壞性地震往往震源規模較大,破裂由一點開始破裂,沿著主破裂方向傳播形成較大尺度斷層,并向外輻射地震波.因此,破壞性地震的震源能量(震級)、震源區應力降和破裂尺度及方向等震源參數共同控制著致災地震動的強度和空間分布.EEW系統為提高計算效率采取大量簡化處理,將震源視為點源,不考慮震源區應力降和破裂尺度等震源參數,難以保證中、大地震事件中地震動預測結果的可靠性.例如,EEW系統對日本2011年太平洋沿岸的東北地震(MW9.1)的地震動預測出現了嚴重低估(如Hoshiba and Iwakiri,2011;Hoshiba and Ozaki,2014),導致預警信息和警報發布范圍出現較大偏差,其一方面原因是采用有限波形估計的震級偏低,另一方面原因是未考慮震源尺度和方向等關鍵震源參數.此外,EEW系統估算震級采用的經驗統計關系對具有區域特征地震個性化和地震動傳播衰減差異考慮不足,這也會對地震動預測造成影響.日本緊急地震速報系統和我國川滇地震預警先行先試區處理結果表明,目前地震預警系統地震動幅值的預測精度多在0.5個數量級,對地震烈度的預測偏差多數在±1度,地震動預測精度明顯偏低(Fujinawa and Noda,2013;Kodera et al.,2016;Peng et al.,2020).

確定震源參數是地震學研究中一項基礎工作,震后常規獲取地震矩、震源尺度和應力降等震源參數的處理方法可分為非參數化方法和參數化方法.非參數化方法通過經驗衰減關系估算地震震級等震源參數(Castro et al.,1990;Picozzi et al.,2017).參數化方法采用圓形或矩形斷層破裂輻射的體波,通過擬合震源位移譜獲取地震矩、震源尺度和應力降(Abercrombie,1995;Zollo et al.,2014;Supino et al.,2019).擬合使用的震源譜是從觀測體波譜中消除路徑和場地效應的無偏譜,最常用的校正方法是經驗格林函數法(Shearer et al.,2019).目前已發展多種校正方法,例如譜分解和全局經驗格林函數擬合法(Shearer et al.,2006;Trugman and Shearer,2017),譜比擬合方法(Abercrombie et al.,2017;Ruhl et al.,2017),但是采用不同方法對相同數據計算的應力降結果并不一致,表現出一定程度的不確定性.由于非參數化方法和參數化方法均為震后收集記錄,需將時域信號轉至頻域處理,難以實現震源參數實時計算,所以無法用于EEW系統.近幾年,時域快速計算(估算)震源參數的方法也有報道,一種方法是針對中、小地震,根據P波位移增長特征估計震源破裂尺度、應力降和地震矩(Zollo et al.,2022);另一種是根據震源譜ω-2衰減模型,由“圓盤”震源(Brune,1970)中震源參數與地震動均方根參數的約束關系實現震源參數快速計算(Lior and Ziv,2017,2018,2020),該時域方法為震源參數實時計算和考慮震源參數的地震動預測建立了橋梁.以上時域方法研究主要面向中、小震事件,采用“圓盤”震源模型實現震源參數計算.破壞性大震震源斷層一般沿著破裂方向的傳播距離大于斷層的寬度,因此,圓形震源模型不適用描述中、大震震源的破裂面.此外,Kaneko和Shearer(2015)在分析地震應力降和拐角頻率等參數估計結果不確定性的原因時,也表明了震源幾何形狀對拐角頻率和應力降的結果存在一定影響.

由此本文引入了符合中、大震震源幾何特征的橢圓震源模型,推導了橢圓模型的地震矩、拐角頻率、應力降和破裂面積等震源參數與地震動加速度、速度和位移均方根的關系,發展了用于中、大震的震源參數時域計算方法(簡稱“地震動均方根法”).本文采用日本Kik-net臺網82次MW≥5.5地震井下和地表臺強震動觀測數據驗證該方法的可行性,并討論了場地效應,地震動傳播衰減,震源效應和濾波頻帶等因素對地震動均方根法的影響,分析地震動均方根法的特點和適用性.

1 地震動均方根法

(1)

(2)

其中,W是斷層窄邊的一半,常數CG是與斷層幾何形狀和泊松比有關的常數(Madariaga,1977).由式(1)和(2)可得應力降：

(3)

中、小地震震源可采用圓形模型(Brune,1970),其震源位移譜的零頻極限水平Ω0和拐角頻率f0與地震矩和圓形半徑的關系如式(4)和(5)所示(Eshelby,1957; Sato and Hirasawa,1973; Madariaga,1976),

(4)

(5)

其中,Uφθ表示輻射花樣因子,FS是自由表面因子,ρ是介質密度,C是體波(P波或S波)速度,CS是S波速度,R是震源距,a是圓形震源半徑,k是與特定震源理論模型的破裂速度有關的常數.在圓形震源模型中,CG=16/(7π),應力降Δτ和拐角頻率f0為

(6)

(7)

大震的震源破裂受地殼厚度的限制,沿著斷層走向的尺度一般大于斷層的寬度,因此圓形模型不適用描述大震的斷層破裂面幾何特征.大震震源模型常視為矩形,但觀測研究表明,通過測量震源譜的拐角頻率僅能近似估算出震源斷層長度,很難同時給出矩形震源的長度和寬度(Madariaga,1977).而在橢圓模型中,震源譜拐角頻率f0和應力降Δτ與震源尺度關系如式(8)和式(9)所示(Eshelby,1957;Dong and Papageorgiou,2003;Kaneko and Shearer,2015).

(8)

(9)

其中,a,b分別為橢圓形震源的長半軸和短半軸.橢圓面積S=πab,因此當拐角頻率已知時,便可近似計算橢圓形震源面積：

(10)

從圓形模型(a=b)到矩形模型(a>b),其常數CG從最小值16/(7π)增至最大值π/2,數值變化不大,約為2倍(Madariaga,1977).橢圓介于圓形和矩形之間,矩形斷層中CG值大于任何一個橢圓形斷層,因此為簡化計算,本文令CG=π/2.將式(9)代入到式(8)可推導出橢圓模型震源譜的拐角頻率和應力降關系：

(11)

震源破裂過程控制地震動的強度和空間分布(Hanks,1979;McGuire and Hanks,1980;Luco,1985;Lior et al.,2016),基于圓盤震源模型(Brune,1970)的地震動位移,速度和加速度均方根可由震源譜零頻極限水平Ω0和拐角頻率f0近似表示為(Lior and Ziv,2017,2018)：

(12a)

(12b)

(12c)

其中T是記錄的時長,κ0是場地效應因子.由式(12a)、(12b)和(12c),可導出震源譜參數Ω0和f0的表達式為

(13)

(14)

拐角頻率f0的另一種表達為

(15)

由此,圓盤震源模型的震源譜參數Ω0和f0可由位移均方根和速度均方根計算.

參考Kaneko和Shearer(2015)對圓形和橢圓形震源的位移譜的擬合方法,將式(13)和式(15)代入到橢圓模型震源參數表達式(4)和式(11),修正橢圓模型的地震矩M0和應力降Δτ的表達式：

(16)

(17)

其中常數

(18)

與Lior和Ziv(2020)圓盤震源模型震源參數計算公式對比可知,本文針對中、大地震引入橢圓模型修正后,僅改變了應力降計算公式(式17)中的系數,而地震矩和拐角頻率計算公式并無變化.因此,仍可在時域內采用地震動均方根參數計算震源參數,計算過程不涉及具有區域特征的先驗參數,適用于高時效性的地震預警系統.

2 數據

2.1 數據及來源

鑒于日本Kik-net臺網同時配置了地表和井下三分量強震儀,同一臺站記錄可分析地震動的場地放大效應,本文采用Kik-net井下和地表強震動記錄驗證地震動均方根法的可靠性.淺源的中、大震破壞力強,影響范圍大,破裂過程復雜,因此篩選地震事件和強震動記錄的原則如下：地震事件矩震級MW≥5.5,震源深度小于60 km;5.5～6.5級地震記錄的震源距設置為100 km,6.6～6.9級地震記錄的震源距為150 km,7.0及以上地震記錄的震源距為200 km;每個地震事件至少包含4組記錄.數據集共選取82次地震,剔除異常記錄(如P波頭丟失記錄和噪聲引起P初至難以識別記錄等)后,井下記錄共計2035組,地表記錄共計2031組,每組記錄包括兩個水平向分量和一個豎向分量.地震的目錄矩震級(MWcata)和破裂半持時(Half-duration)由全球矩張量項目(www.globalcmt.org)提供.地震震中(彩色圓形)和臺站(灰色三角形)位置分布圖和震源深度及記錄震源距分布見圖1.

圖1 地震震中及臺站分布圖(a)和震源深度(b)及記錄震源距統計圖(c)

2.2 數據處理

各分向原始加速度記錄減掉記錄前5 s均值完成基線校正.每組記錄的P波初至通過豎向分量拾取,采用長短時比值粗撿拾和AIC精撿拾結合方法執行自動撿拾(馬強等,2013),人工檢查P波拾取結果并修正因噪聲干擾引起的誤撿拾.S波到時使用理論到時差確定,其中P波和S波波速見表1,為避免地震定位和走時的誤差影響,將理論到時差的90%定為最終P-S到時差.基線校正后加速度記錄分別經過一次和兩次積分得到速度時程和位移時程.為消除因長周期地脈動等噪聲引起的漂移,采用雙向巴特沃斯高通濾波器對速度和位移時程濾波,低頻截止頻率fl為0.02 Hz.加速度、速度和位移時程的三個分量分別代入下式計算地震動均方根：

(19)

其中,參數Prms分別為加速度均方根Arms,速度均方根Vrms和位移均方根Drms.UD,NS和EW分別為豎直向,南北向和東西向地震動.n為輸入記錄的采樣點數.近場P波記錄是地震預警系統最常用的數據,所以將P波的地震動均方根參數代入式(15)、(16)、(17)和(10),即可計算震源譜拐角頻率、地震矩、應力降和破裂面積等震源參數,其中每次地震破裂面積使用其拐角頻率平均值計算.

(20)

(21)

其中,fc是缺失信號的截止頻率,本文令fc=max(1/T,fl).T是記錄的時長,本文計算震源參數均采用下式計算,

(22)

其中,等號右側第一項視為震源破裂持時(Hanks and McGuire,1981),f0由式(15)計算;第二項是在數據長度中加入震源距相關項,考慮隨地震波傳播距離增大引起波包的擴展(Boore and Thompson,2014).

(23)

震源破裂機制復雜,不同震級地震的破裂規模和破裂持時差異巨大,本文主要采用在P波初至后逐步增加的波形獲取隨破裂演化的震源參數變化過程,進而分析震源破裂的演化特征.地震動均方根法涉及常數按表1取值,每組強震動記錄的震源參數計算過程如圖2所示.

表1 計算震源參數采用的常數取值表Table 1 Constant values for source parameters

圖2 MW6.6地震FKOH05臺記錄的震源參數計算過程示例(a)—(c) 分別為逐步增加P波的加速度均方根,速度均方根和位移均方根曲線; (d)—(f) 分別為矩震級,應力降和拐角頻率計算結果; (d) 中水平虛線為目錄矩震級.

3 結果

3.1 震源參數計算

Kik-net臺網井下儀器位于孔底,記錄噪聲低,默認為基巖場地,記錄只包含震源和傳播路徑信息.本節分析井下記錄中P波獲取的矩震級,應力降和拐角頻率等震源參數結果.為消除震源輻射花樣的影響,每次地震震源參數取為該地震所有記錄計算結果的平均值.圖3為MW5.6,MW6.4,MW7.0和MW8.3地震通過地震動均方根法獲取的矩震級,應力降和拐角頻率的示例.由圖可知,地震動均方根法可提供矩震級、應力降Δτ和拐角頻率f0隨破裂傳播的變化過程.隨著P波波形的增加,四個震例的矩震級均逐漸趨近于目錄矩震級(MWcata),其中前三次地震的矩震級,應力降Δτ和拐角頻率f0在破裂半持時前均達到穩定值.而對MW8.3巨震,其震源破裂半持時約為33.5 s,在20 s P波窗內震源參數仍在上升,未完全穩定,這表明對破裂過程異常復雜的大震,通過破裂早期部分波形難以約束整個震源的破裂特征.

圖3 四個地震震例的震源參數計算結果灰線為各臺記錄計算結果,紅線為所有記錄的平均值;水平虛線表示目錄矩震級,豎向虛線表示破裂半持時.

圖4匯總了由地震動均方根法計算所有地震的矩震級、應力降、拐角頻率和破裂面積隨破裂發展的變化曲線.顯然,隨目錄矩震級的增大,震源參數曲線越復雜,達到穩定值的所需P波時長也增加.不同震級地震初始破裂的震源參數存在一定波動,但隨破裂發展逐漸穩定,穩定后的矩震級、拐角頻率和破裂面積參數隨震級的分布體現出明顯的相關性,其中最終矩震級幾乎完全趨近于目錄震級.在圖4a中的矩震級增長曲線體現出震源能量的釋放過程,下節將對比其與震源時間函數的關系.一般震源區應力降與地震規模無關,本文地震動均方根法獲取的早期破裂的應力降曲線波動較大,穩定后的數值分布在0.3～20 MPa范圍,82次淺源地震的應力降與震級未表現出相關性.值得注意的是不同震級的震例中初始破裂的應力降多高于后續穩定值(圖4b),這可能是由于采用P波時窗長(tW)遠小于震源持時,對震源過程采樣不足,計算結果可能出現不穩定的情況.Lior和Ziv(2020)認為在這種情況下,應力降表達式可能返回對應于拐角頻率f0小于1/tW的非物理結果.但不能排除地震成核后的初始破裂或早期破裂區域的應力降大于整個破裂區的平均值的情況,因為相同均方根參數獲取得矩震級增長過程與震源時間函數計算的結果基本一致(圖7).

拐角頻率變化曲線表明從地震早期破裂輻射P波獲取拐角頻率均遠大于整個破裂的均值;拐角頻率的穩定值隨震級關系體現出地震的基本規律,地震越大,輻射地震動的低頻成分越豐富,拐角頻率越低(圖4c).由拐角頻率均值計算的破裂面積隨震級的增加而增大,因地震破裂過程復雜,震時破裂面幾何形狀并不規則,因此通過拐角頻率近似計算的震源破裂面積精度可能較低.

對比2 s、4 s和8 s P波的矩震級結果發現,通過擴大時窗,地震動均方根法獲取的矩震級均值與目錄震級的分布逐漸接近于1∶1線,提高矩震級結果的可靠性(圖5).采用全P波段記錄使得大震震級的低估問題得到有效改善.因此,地震動均方根法能夠利用P波記錄快速可靠地獲取表征地震絕對能量的矩震級.

圖6進一步展示了所選地震在全P波記錄中獲取的拐角頻率、應力降和破裂面積參數隨震級的分布.作為描述地震震源譜的關鍵參數,拐角頻率隨震級增大而降低,目錄矩震級(MWcata)換算的地震矩和拐角頻率的擬合關系式為

(24)

拐角頻率的冪次在-4與-3之間,符合地震自相似比例關系(Archuleta et al.,1982;Abercrombie, 1995;Hiramatsu et al.,2002).82次地震的震源區應力降主要分布在0.3～15 MPa之間,接近日本地區地震的統計范圍(Oth et al.,2010).按震源深度分類,深度在30 km以內視為地殼事件,大于30 km的為地殼下事件,那么由全P波段計算的震源參數中,相同或相近震級的地殼地震的拐角頻率和應力降普遍大于深部的地殼下事件,其中地殼地震的應力降均值約為1.45 MPa,地殼下地震的應力降均值約為3.96 MPa.此外,震級接近的淺源殼震的破裂面積幾乎都大于地殼下地震,這也符合日本地區地震基本規律(Oth et al.,2010).

圖4 井下P波記錄計算的震源參數結果

圖5 2 s、4 s、8 s時窗和全P波段的矩震級結果灰點為各組記錄計算的矩震級,菱形為每次地震所有記錄計算的矩震級平均值,灰色虛線為1∶1參考線.

3.2 地震能量釋放過程

中、大地震震源破裂物理過程復雜,震后常用遠場測震記錄反演獲取地震震源能量釋放過程,即震源時間函數(STF).目前全球中、大震的震源時間函數SCARDEC數據庫已經建立(Vallée and Douet,2016).為驗證地震動均方根法獲取矩震級曲線是否反映地震能量的釋放過程,本節將其與SCARDEC數據庫提供的STF進行對比,其中由STF計算的矩震級定義為MW-STF.圖7展示了觀測記錄較豐富的四個地震事件(MW6.6、MW7.0、MW7.4和MW8.3)的震源能量釋放過程.由圖可知,MW6.6、MW7.4和MW8.3地震的計算矩震級曲線與MW-STF曲線基本一致,地震動均方根法能通過P波記錄可靠量化地震能量釋放過程,但MW7.0地震的矩震級曲線在局部時間段內(1～6 s)的偏離MW-STF曲線.四個地震震中位置分別對應圖1中的地震序號1、2、3和4,對比觀測臺站位置發現,MW7.0地震的強震臺集中分布在震中方位角的一個象限(270°～360°),而其他三次地震的記錄至少分布在地震震中方位角的兩個象限.地震動均方根法的矩震級計算偏差很可能受輻射花樣和(或)破裂方向性效應的影響.由此可知,地震動均方根法對觀測臺站震中方位角分布比較敏感,當地震的臺站方位角分布均勻時,能夠通過P波記錄快速可靠地量化地震能量釋放過程.

圖6 由全P波記錄獲取的拐角頻率、應力降和破裂面積結果灰色圓點為各臺記錄獲取的震源參數結果,彩色圓形為各地震所有記錄獲取參數的均值.

圖7 地震能量釋放過程灰線為地震動均方根法由每組記錄獲取的矩震級曲線,黑實線為所有記錄的矩震級均值曲線,粉虛線為STF計算的矩震級曲線.

4 影響因素分析

震級表征了地震破裂能量,是地震預警發布的關鍵地震信息,便于從地震目錄中獲取,本節以地震動均方根法的矩震級計算結果為例,分析場地效應,地震動傳播衰減,震源效應和濾波頻帶等因素對該方法的影響.

場地效應對地震動影響十分復雜,土層厚度、密實程度、地形地貌都可能會對地震動不同程度的放大,對不同頻率成分地震動的放大也有差異.為明確地震動均方根法對場地效應敏感程度,本節簡單對比采用井下和地表記錄8 s P波獲取矩震級的誤差.經統計,MW7.0及以下地震震源破裂半持時均小于8s,對該震級段的地震,地震動均方根法計算的矩震級接近于目錄矩震級(圖8).其中井下記錄矩震級結果同目錄矩震級的誤差介于±0.5之間,平均值為0.07個震級單位(圖8a);地表記錄矩震級計算誤差介于±0.5之間,平均值為0.02個震級單位(圖8b).然而大震(MW>7.0)破裂半持時大于8 s,震源破裂規模大,導致從井下和地表8 s P波波形計算的矩震級都偏低.圖8c是由井下和地表記錄獲取的矩震級之差,井下記錄計算的矩震級略大于地表記錄,矩震級之差介于±0.2之間,平均值為0.09個震級單位.由此可知,地震動均方根法在計算矩震級時對場地效應不敏感,采用井上和地表記錄都可獲取可靠的地震震級.因此,在EEW系統中,綜合考慮時效性和精度要求,地震動均方根法可忽略場地效應影響直接使用地面記錄計算矩震級.

圖8 井下和地表記錄計算矩震級的誤差統計(a) 井下記錄計算誤差; (b) 地面記錄計算誤差; (c) 井下記錄和地面記錄計算誤差.

地震動傳播衰減是EEW系統Pd方法估算震級的關鍵修正項,而且往往存在區域特征.震動記錄匱乏和地震風險水平較低的地區難以獲取Pd經驗關系,可能會限制EEW系統應用.若對這類地區采用來自其他區域的先驗關系,則需開展專項驗證測試工作.地震動均方根法考慮了地震動傳播衰減項,即震源距R(式(16)—(17)).圖9結果表明地震動均方根法采用井下記錄4 s和8 s P波計算的矩震級隨震源距分布未出現明顯的衰減趨勢.因此,地震動均方根法無需再進行地震動傳播衰減的校正,放寬了在缺少強震記錄地區的限制,有效提升地震預警系統的適用范圍.

圖9 井下記錄獲取的矩震級隨震源距的分布

地震破裂在時間和空間上擴展傳播的不對稱行為,即地震破裂方向性,常出現于大震中(McGuire et al.,2002;Wen et al.,2015),部分中、小地震的破裂方向性效應也有報道(Kane et al.,2013;Wang et al.,2019).破裂方向性效應往往會引起在破裂傳播方向上的地震動持時變短、幅值升高、視拐角頻率升高,破裂后方地震動持時延長、幅值降低、視拐角頻率降低.地震震源對地震幅值調制的另一種常見的機制為震源輻射花樣,常用于量化輻射花樣的模型是四葉模型,對震中附近地震動幅值調制多為四象限花樣分布.本節選取了4次陸域地震(MW5.5,MW5.8,MW6.8和MW6.9),其臺站方位角分布相對均勻,震中位置分別對應圖1地震序號5、6、7和8.4次地震井下記錄的4 s和8 s P波計算的矩震級隨方位角的分布均體現出正(余)弦形式的規律性波動(圖10).鑒于傳播距離對地震矩影響較小,因此推測這種調制結果可能反映震源的一種或綜合效應.其中兩次中強震(MW6.8和MW6.9)的矩震級計算結果隨方位角波動較大,震源效應更強.因此,如果臺站集中分布在較窄方位角范圍內,很可能因震源效應造成較大的矩震級計算偏差,如圖7b.若計算矩震級隨臺站方位角的律動是破裂方向性的效應引起的,那么當震中周圍臺站方位角分布均勻時,可通過矩震級隨方位角分布特征,如方向性效應系數Cd,估計出震源破裂方向和破裂方式(單側破裂還是雙側破裂),這可為實時地震動預測提供至關重要的震源斷層參數.

圖10 井下記錄獲取的矩震級隨臺站方位角的分布及其擬合曲線

強震動記錄積分處理時采用了高通濾波以消除因長周期地脈動等噪聲積分引起的漂移,但由此也可能濾掉地震動中部分有效低頻成分,尤其是破壞性大震.將高通濾波器的低頻截止頻率設為0.02 Hz和0.05 Hz,分析頻帶差異對地震動均方根法獲取矩震級結果的影響.圖11表明井下8 s P波經過fl=0.02 Hz的高通濾波處理能夠獲取更可靠的矩震級效果.因此,建議選用低頻截止頻率為0.02 Hz的高通濾波處理記錄以保留地震動中有效的低頻成分.

圖11 濾波頻帶對地震動均方根法矩震級結果的影響

5 實時震級計算效果

地震動均方根法可通過P波波形快速獲取震源參數,為分析其矩震級計算結果的可靠性,本節選用地表記錄分別通過估算震級經驗方法和地震動均方根法計算矩震級.EEW系統一般采用P波特征參數(峰值參數Pd和周期參數τc等)的經驗關系快速估算地震震級.在離線數據中,Pd參數估計震級的精度表現更好(金星等,2012).本節對比全球經驗Pd方法(式(25))和地震動均方根法的實時震級估算結果,

MPd=1.23log(Pd)+1.38log(R)+5.39,

(25)

+1.105(±0.126),

(26)

其中MPd是Pd震級,對應日本地震震級Mj,R是震源距(Kuyuk and Allen,2013).為統一震級標度,Pd方法估算的震級MPd通過式(26)轉換為矩震級MW(Oth et al.,2010).Pd方法采用垂直向記錄,地震動均方根法仍采用三分向記錄,用于速度和位移時程的濾波頻帶與2.2節一致.本節選擇4 s和8 s P波對比兩種方法的矩震級結果(圖12).結果表明,對統一的P波時窗,地震動均方根法計算矩震級的離散性和準確度均優于Pd方法.其中8 sPd參數估算的矩震級普遍高于目錄矩震級(圖12c),說明估算震級的經驗方法需要考慮不同區域地震動衰減規律的差異.因此,在地震風險低和強震動記錄匱乏的地區,若采用統一先驗關系估算震級,其適用性需要進一步驗證.

圖12 Pd方法和地震動均方根法計算矩震級結果對比(a)和(c)分別是4 s和8 s時窗的Pd法計算結果;(b)和(d)分別是4 s和8 s時窗的地震動均方根法計算結果,誤差棒表示95%置信區間.

6 討論

鑒于震源幾何形狀一定程度上影響著拐角頻率和應力降等震源參數(Kaneko and Shearer,2015),本文在Lior和Ziv(2020)提出的震源參數算法的基礎上,針對中、大震引入了橢圓震源模型,發展了橢圓模型下震源參數的時域算法,其中主要修正項是應力降計算公式,矩震級和拐角頻率算法并未改變.

Kaneko和Shearer(2015)開展了震源幾何特征對震源譜和拐角頻率及應力降等震源參數影響的研究,其對圓形和橢圓形震源輻射體波的位移譜均采用了Brune(1970)提出的ω2震源譜擬合.借鑒這種方法,本文從Brune震源譜出發,在圓形震源模型震源參數與地動參數之間的關系基礎上,嘗試通過橢圓模型震源參數間的關系,修正中、大震的震源參數計算公式(式(10)、(15)、(16)和(17)),進而實現由時域地震動參數計算震源參數.也可理解為直接利用Brune震源譜擬合中、大震震源譜獲取Ω0和f0,然后根據橢圓模型修正公式計算震源參數.

事實上,利用Brune震源譜擬合大震震源譜可能對量化低頻幅值水平(Ω0)存在一定誤差(Joyner,1984).由于Brune震源譜模型在中低頻段高估觀測譜幅值,因此在Brune模型的基礎上發展了適用于大震的雙拐角震源譜模型(如陶夏新和王國新,2003),以修正Brune模型Ω0的擬合誤差.與低頻幅值水平Ω0相對應的震源參數是地震矩(式(4)),本文推導的地震矩計算公式(式(16))在圓形和橢圓形震源模型下是一致的,也就是引入橢圓模型并未改變該公式.為明確由地震矩計算公式(式(16))對大震適用性,本節補充了日本82次中、大地震由全波列記錄計算的矩震級結果(圖13),具體計算方法參考Lior和Ziv(2020),相關常數按表1取值.計算結果表明當采用足夠長度記錄計算的矩震級和目錄矩震級基本一致,未出現明顯偏差,證實了本文方法的可靠性.但是我們也注意到不論是采用P波還是全波列計算82次地震的矩震級結果確實存在略微偏離1∶1線的趨勢(圖5和圖13),且隨震級增大,計算矩震級向下偏移趨勢越明顯.這可能是采用Brune震源譜模型擬合中、大震震源譜參數Ω0的誤差引起的.

地震動均方根法推廣最大的限制性因素是場地效應.第4節已明確了場地效應對矩震級結果的影響較小,在此分析場地效應對拐角頻率和應力降和破裂面積的影響.圖14為地表和井下全P波記錄計算震源參數的比值.對比結果表明,場地效應使拐角頻率放大約1～3倍,震級越大,影響越弱;應力降的變化范圍較大,放大約為1～12倍,同樣震級越大,場地效應的影響降低.而由地表記錄獲取的破裂面積普遍小于井下記錄結果.當震級相近時,場地效應對深度更大地震的拐角頻率和應力降影響略高于淺源地震.根據式(11)可知,當場地效應使得拐角頻率增加2倍時,計算應力降將增加8倍,這是應力降放大倍數遠大于拐角頻率的主要原因.傳統采用頻域參數法擬合震源譜獲取的應力降結果差異也比較大,表現出一定的不確定性,這主要是因為應力降是高度依賴于模型的(Shearer et al.,2019).值得注意的是,對震級大于6.5的地震,場地效應對震源參數的影響變弱,這可能是由于大震輻射地震波頻率偏低,與場地的卓越頻率相差較大,地震記錄的視拐角頻率受其影響弱化引起的.因此對于大震,本文發展方法可以通過地面記錄快速獲取相對可靠的震源參數.

圖13 井上臺10 s、20 s、40 s和全波列的矩震級計算結果灰色圓點為臺站記錄計算矩震級結果,菱形為各地震事件的平均矩震級.灰色虛線為矩震級1∶1線.

圖14 82次地震地表和井下全P波記錄計算拐角頻率、應力降和破裂面積的比值

以上僅對地震動均方根法主要影響因素展開討論,地震動均方根法推導過程還涉及一些簡化假定,如殼幔介質中P波和S波波速隨著深度增加而增大,而本文采用固定數值,勢必引入計算誤差.P波常數k的取值是來自破裂傳播速度Cr=0.9CS的假設,實際地震破裂傳播速度存在差異,因此,在震源參數計算中常數k也會引入誤差.除此之外,實際破裂面的幾何形狀要遠比橢圓模型復雜,大震破裂面的長度可能遠大于寬度,通過橢圓模型估計應力降和破裂面積可能存在一定偏差.盡管如此,本文結果與日本地區地震震源參數的統計結果對比可知,地震動均方根法能夠通過時域近場波形快速獲取可靠的震源參數,不僅為EEW系統提供了新方法,還能夠量化震源參數隨破裂發展的演化過程,為分析地震破裂特征提供了新的途徑.

7 結論

考慮中、大震震源幾何長短軸差異,在圓盤模型震源參數與地震動均方根的約束關系的基礎上,本文推導了橢圓震源模型中震源地震矩、拐角頻率和應力降同地震動加速度,速度和位移均方根的約束關系,提出了一種適用于中、大型地震震源參數時域快速估計方法.該方法繞過了傳統震源參數的反演或者在地震信號頻率域測量的過程,不涉及地震區域特征先驗參數,計算方法簡便,無需人工參與計算,非常適用于高時效性的地震預警系統.為驗證地震動均方根法的可靠性,本文采用日本Kik-net臺網MW≥5.5中、大地震的井下和地表強震動記錄進行驗證.通過對比分析,得到以下結論：

(1)考慮橢圓模型的地震動均方根法能夠采用P波波形快速可靠地計算中、大震矩震級、應力降和拐角頻率及破裂面積等震源參數;

(2)由全P波段和全波列記錄計算的矩震級與目錄矩震級基本一致,大震的震級低估得以有效緩解;

(3)由地震動均方根法計算82次地震的拐角頻率隨地震矩的分布滿足地震自相似規律,震源區應力降數值符合日本地區的統計范圍(0.3～15 MPa),驗證了該方法的可靠性;

(4)地震動均方根法考慮地震動傳播距離,其震源參數計算結果可不再做路徑衰減項校正;

(5)地震動均方根法對震源效應比較敏感.當震中附近臺站方位角分布均勻時,能夠通過P波記錄快速可靠地量化地震能量釋放過程;

(6)場地效應對地震動均方根法影響復雜,其中矩震級計算結果對場地效應不敏感;場地效應對拐角頻率和應力降存在一定影響,但隨震級的增大其影響變弱.因此對于大震,地震動均方根法可以通過地面記錄快速獲取相對可靠的震源參數.

致謝感謝日本Kik-net數據中心(https：∥www.kyoshin.bosai.go.jp/)為本文研究提供的強震動數據和全球矩張量項目(https：∥www.globalcmt.org/CMTsearch.html)提供的震源參數.感謝匿名審稿專家為本文提出的寶貴意見與建議.文中部分圖片使用通用制圖工具GMT(Genetic Mapping Tools)繪制.