李 斌，李自紅，宋美琴，關(guān)鵬虎

(1.太原理工大學(xué) a.礦業(yè)工程學(xué)院，b.地震與地質(zhì)災(zāi)害防治研究所，太原 030024；2.山西省地震局，太原 030021)

北京時(shí)間2016年4月7日04時(shí)49分，在山西省原平市蘇龍口鎮(zhèn)(38.87°N，112.91°E)發(fā)生了ML4.7地震，震源深度11 km，山西原平、代縣、定襄、太原、晉中、陽(yáng)泉以及朔州等地震感明顯。這是自2016年1月1日以來(lái)，繼2016年3月12日運(yùn)城M4.4地震、3月15日運(yùn)城M3.0地震和3月27日運(yùn)城M4.0地震后在山西省境內(nèi)發(fā)生的第4次顯著地震事件，也是2009年3月28日原平ML4.5地震打破山西地震帶近5年的M4地震平靜后在該地區(qū)發(fā)生的又一次顯著地震。以往震例研究表明，山西地震帶中強(qiáng)地震“南北遷移”活動(dòng)特征明顯[1-2]。此次原平ML4.7地震的發(fā)生，是山西地震帶顯著地震活動(dòng)由南向北遷移的表現(xiàn)，對(duì)山西地震帶中強(qiáng)地震具有重要的預(yù)警作用，預(yù)示著山西地震帶可能進(jìn)入了新一輪的4級(jí)地震活躍時(shí)段[3]。因此，深入理解此次顯著地震事件的震源機(jī)制，有助于理解地震發(fā)生時(shí)震源區(qū)的應(yīng)力狀態(tài)與發(fā)震構(gòu)造，以及山西地震帶中強(qiáng)地震的孕育規(guī)律與活動(dòng)特征，也是跟蹤震情發(fā)展趨勢(shì)的重要依據(jù)。

目前顯著地震震源機(jī)制的求解方法較多，如P波初動(dòng)法[4]、振幅比法[5]、P波初動(dòng)聯(lián)合振幅比法[4,6]、CAP全波形反演法[7-8]與矩張量反演法[9-11]等，且各種方法的研究結(jié)果都對(duì)世界上很多地震作出過(guò)較為合理的解釋[12-16]。地震矩張量相比于其他方法，對(duì)于顯著地震的震源機(jī)制解能有較為準(zhǔn)確的描述[11]，且具有一定的優(yōu)勢(shì)。相對(duì)于P波初動(dòng)與振幅比法，矩張量反演對(duì)臺(tái)站布局要求低，且只需數(shù)字波形資料，無(wú)需進(jìn)行震相分析，避免了人工量取震相誤差帶來(lái)的不確定性；反演結(jié)果反映的是整個(gè)破裂過(guò)程的信息，而非僅僅是初始破裂信息[17]。相比于CAP方法，矩張量反演方法的計(jì)算結(jié)果中含補(bǔ)償線性矢量偶極 (CLVD)和各向同性分量(ISO)，因此地震矩張量可以分解為表示斷層面剪切錯(cuò)動(dòng)的雙力偶分量(DC)、震源體膨脹或收縮的各向同性分量(ISO)、震源體優(yōu)勢(shì)方向的張裂或擠壓變形(CLVD)等[9-10]。通常構(gòu)造地震是由斷層兩盤(pán)的剪切位錯(cuò)觸發(fā)，因此矩張量能更好地反映斷層破裂信息[18]。

本文利用山西區(qū)域地震臺(tái)網(wǎng)寬頻帶數(shù)字波形資料，在采用Hypo2000對(duì)地震重新定位的基礎(chǔ)上，采用TDMT_INV時(shí)間域矩張量反演方法對(duì)2016年4月7日山西原平ML4.7地震進(jìn)行矩張量反演計(jì)算，通過(guò)逐個(gè)增減參與反演的臺(tái)站數(shù)量進(jìn)行多次反演，結(jié)合震源機(jī)制隨反演深度的變化特征進(jìn)行參數(shù)的穩(wěn)定性評(píng)價(jià)。將獲得的震源機(jī)制結(jié)果與不同臺(tái)站垂向分量的P波極性在震源球上的分布位置，以及與利用Snoke方法和CAP方法得到的震源機(jī)制解進(jìn)行對(duì)比，以驗(yàn)證參數(shù)的合理性與可靠性。研究結(jié)果將有助于理解山西地震帶中強(qiáng)地震的孕育規(guī)律與活動(dòng)特征。

1 方法與原理

國(guó)內(nèi)外學(xué)者對(duì)地震矩張量的求解方法進(jìn)行了大量的研究。PATTON et al[19]首先將遠(yuǎn)震矩張量反演算法做了修改，進(jìn)行了區(qū)域地震矩張量反演。DREGER et al[9]提出了在時(shí)間域利用區(qū)域地震波形(Pnl)進(jìn)行地震矩張量反演的方法(time domain moment tensor inversion，TDMT_INV).本文采用的是DREGER于2003年改進(jìn)理論波形計(jì)算方法后的TDMT_INV方法[10]，即在計(jì)算理論格林函數(shù)時(shí)加入了SAIKIA[20]改進(jìn)的離散波數(shù)積分法。目前該方法在國(guó)際上已經(jīng)得到廣泛應(yīng)用[21-24]。隨著寬頻帶數(shù)字臺(tái)網(wǎng)的建立，國(guó)內(nèi)許多研究人員也利用該方法反演震源機(jī)制解：趙翠萍等[25]、唐蘭蘭等[26]、屠泓為等[27]將此方法應(yīng)用于伽師震源區(qū)；王勤彩等[28]反演了汶川大地震余震序列88個(gè)地震的矩張量解。

TDMT_INV方法基于雙力偶點(diǎn)源模型，在自由表面距離震源x處t時(shí)刻的實(shí)際觀測(cè)記錄dn(x,t)可表示為[29-31]：

(1)

(2)

式中：oi和si分別為觀測(cè)數(shù)據(jù)與基于格林函數(shù)計(jì)算的時(shí)間域理論數(shù)據(jù)。VR的計(jì)算不僅考慮了波形的相似性，還考慮了絕對(duì)振幅的大?。籚R值越大，表明波形擬合程度越高，求解出來(lái)的地震矩張量越可靠。在不同深度上搜索VR的最大值，以獲得最佳震源深度與震源機(jī)制的最優(yōu)解。

2 資料選取與數(shù)據(jù)處理

2.1 地震波形數(shù)據(jù)

目前，山西區(qū)域測(cè)震臺(tái)網(wǎng)由72個(gè)實(shí)時(shí)記錄與傳輸?shù)臄?shù)字化寬頻帶地震臺(tái)站組網(wǎng)而成(山西省內(nèi)臺(tái)站57個(gè)，鄰省臺(tái)站15個(gè))，基本上均勻分布于山西省及周邊地區(qū)(圖1(a)).現(xiàn)有臺(tái)網(wǎng)條件下，ML≈2.0級(jí)及以上的地震，通常能同時(shí)被周?chē)鷰资畟€(gè)臺(tái)站記錄到[32]。對(duì)于此次原平ML4.7地震，山西區(qū)域測(cè)震臺(tái)網(wǎng)72個(gè)地震臺(tái)站全部記錄到此次地震，38個(gè)臺(tái)站的垂向分量有清晰的P波初動(dòng)。本文共選取了震中距120 km范圍內(nèi)8個(gè)臺(tái)站的波形數(shù)據(jù)參與反演，臺(tái)站方位及分布如圖1(b)所示。選取數(shù)據(jù)的基本標(biāo)準(zhǔn)為：記錄波形質(zhì)量好，信噪比高，二階Butterworth濾波器濾波后波形清晰，滿足矩張量反演數(shù)據(jù)可靠性的要求。

HAVSKOV et al[33]研究表明，對(duì)于質(zhì)量較好的波形數(shù)據(jù)，利用速度波形數(shù)據(jù)或積分至位移數(shù)據(jù)進(jìn)行矩張量反演,其結(jié)果具有較好的一致性。本文在數(shù)據(jù)預(yù)處理過(guò)程中對(duì)比發(fā)現(xiàn)，速度波形數(shù)據(jù)在三個(gè)分量旋轉(zhuǎn)及濾波后波形表現(xiàn)得更加穩(wěn)定。鑒于此，本文選用了速度波形數(shù)據(jù)進(jìn)行后續(xù)的反演計(jì)算。數(shù)據(jù)處理與計(jì)算過(guò)程主要包括：

1) 對(duì)原始波形數(shù)據(jù)進(jìn)行去均值、去傾斜、去儀器響應(yīng)的預(yù)處理。

2) 將3個(gè)分量分別旋轉(zhuǎn)至垂向、切向和徑向分量。

3) 選取二階Butterworth濾波器進(jìn)行濾波。頻帶范圍的選取與震級(jí)有關(guān)；對(duì)于本次ML4.7地震，選取的頻帶范圍為0.02～0.05 Hz.

4) 設(shè)定計(jì)算起始深度1.1 km(避免計(jì)算深度位于速度結(jié)構(gòu)分層界面上)及步長(zhǎng)1.0 km，基于區(qū)域速度結(jié)構(gòu)等參數(shù)計(jì)算各臺(tái)站不同深度范圍的理論格林函數(shù)。

5) 進(jìn)行矩張量反演及反演結(jié)果置信度評(píng)估。

圖1 山西地震臺(tái)網(wǎng)臺(tái)站(a)及參與計(jì)算的8個(gè)地震臺(tái)站(b)分布圖Fig.1 Distribution of stations of Shanxi Seismic Network (a) and the 8 selected stations(b)

2.2 地殼速度結(jié)構(gòu)

以往研究表明，無(wú)論哪種震源機(jī)制求解方法，地殼速度結(jié)構(gòu)都是影響震源機(jī)制解的重要因素之一[17,24]。為了減少地殼速度結(jié)構(gòu)對(duì)矩張量反演結(jié)果的影響，本文在Hpyo2000重新定位與矩張量反演過(guò)程中均采用了山西地震帶最新的1D速度結(jié)構(gòu)模型[34]，地殼分層及數(shù)據(jù)詳見(jiàn)表1.

表1 本文計(jì)算采用的山西地震帶1D速度結(jié)構(gòu)模型Table 1 1D velocity model of Shanxi seismic belt applied in this study

3 計(jì)算結(jié)果與分析

3.1 重新定位結(jié)果

矩張量反演需要首先計(jì)算各參與反演臺(tái)站的理論格林函數(shù)。除地殼速度結(jié)構(gòu)外，震源位置特別是震源深度誤差會(huì)對(duì)理論格林函數(shù)的計(jì)算結(jié)果產(chǎn)生影響，進(jìn)而可能影響矩張量反演的最終結(jié)果，因此準(zhǔn)確的定位結(jié)果是可靠矩張量反演的基礎(chǔ)數(shù)據(jù)之一[6,17]。對(duì)于本次原平ML4.7地震，山西地震臺(tái)網(wǎng)與中國(guó)地震臺(tái)網(wǎng)均給出了定位結(jié)果，見(jiàn)表2.兩個(gè)結(jié)果的震中位置基本一致，但山西臺(tái)網(wǎng)給出的震源深度為11 km，中國(guó)臺(tái)網(wǎng)為16 km，偏差較大。鑒于此，本文基于山西地震帶最新的1D地殼速度結(jié)構(gòu)，采用Hpyo2000[35]定位方法對(duì)本次地震進(jìn)行了重新定位，得到主震震中位置為：38.86°N，112.92°E，震源深度14.7 km. Hpyo2000定位方法由Klein開(kāi)發(fā)，采用傳統(tǒng)的Geiger法基本思路，更適合于觀測(cè)臺(tái)站分布均勻的網(wǎng)內(nèi)近震及地方震。重新定位結(jié)果將為后續(xù)矩張量反演提供保障。

表2 2016年4月7日山西原平ML4.7地震定位結(jié)果Table 2 Location of the Yuanping ML4.7 earthquake on April 7, 2016

3.2 震源機(jī)制結(jié)果

本文共選取了震中距120 km范圍內(nèi)8個(gè)地震臺(tái)站記錄的三分量速度波形數(shù)據(jù)進(jìn)行地震矩張量反演。采用的濾波頻帶為0.02～0.05 Hz，起算深度1.1 km、步長(zhǎng)1.0 km，共反演了1.1～25.1 km范圍內(nèi)不同深度的矩張量結(jié)果。根據(jù)反映觀測(cè)數(shù)據(jù)與理論波形擬合情況的約化方差VR和Var/Pdc值隨深度變化的情況，搜索出擬合程度最好、雙力偶成分相對(duì)較高的解，即為最佳震源機(jī)制解。圖2為矩張量反演過(guò)程中觀測(cè)數(shù)據(jù)與理論波形最佳擬合情況下的反演結(jié)果。參與反演的8個(gè)地震臺(tái)站中，除SHC臺(tái)垂直分量和SHY臺(tái)南北分量因波形有一定畸變?cè)诜囱葸^(guò)程中被刪除外，其余各分量的觀測(cè)波形與理論波形的擬合程度尚好，平均VR=56.3%，結(jié)果符合信度要求[6,24,33]。

圖2 山西原平ML4.7地震震源機(jī)制的波形擬合情況及矩張量反演結(jié)果Fig.2 Waveform fitting and results of moment tensor inversion of Yuanping ML4.7 earthquake

反演得到的最佳雙力偶解節(jié)面I參數(shù)為走向143°，傾向63°，滑動(dòng)角-66°；節(jié)面II走向278°，傾向36°，滑動(dòng)角-128°；震源機(jī)制類型屬于正斷兼有右旋走滑分量。雙力偶分量Percent DC=93，線性補(bǔ)償分量Percent CLVD=7，震源機(jī)制符合雙力偶點(diǎn)源的假設(shè)條件。主壓應(yīng)力P軸方位角為94°，傾角為64°；主張應(yīng)力T軸方位角為216°，傾角為15°.標(biāo)量地震矩M0= 2.492 64×1015N·m，矩震級(jí)結(jié)果為MW=4.2，與MS震級(jí)基本一致[3]。

3.3 震源深度的確定

圖3是約化方差VR與震源機(jī)制解隨不同反演深度的變化圖?？梢钥闯觯S著反演深度的不斷增加，震源機(jī)制結(jié)果逐漸從正斷向正斷兼有走滑分量、再向逆斷趨勢(shì)變化；但在5～16 km深度范圍內(nèi)震源機(jī)制性質(zhì)及結(jié)果基本穩(wěn)定，表現(xiàn)為正斷兼有右旋走滑分量，且在11～14 km范圍內(nèi)VR相對(duì)穩(wěn)定且取得最大值，為擬合程度最好、雙力偶成分相對(duì)較高的解。因此，判斷本次地震的最佳震源深度在11～14 km；這與重新定位結(jié)果的震源深度14.7 km具有一定可比性。

圖3 約化方差VR隨反演深度的變化Fig.3 Variance of variance reduction versus inversion depth

3.4 參數(shù)穩(wěn)定性分析

影響地震矩張量反演結(jié)果穩(wěn)定性的因素較多，如區(qū)域地殼速度結(jié)構(gòu)、地震定位誤差(特別是震源深度誤差)、參與反演臺(tái)站的個(gè)數(shù)、臺(tái)站與震中的方位角誤差，等等[17]。為進(jìn)一步分析與驗(yàn)證反演參數(shù)的穩(wěn)定性，在反演過(guò)程中我們首先對(duì)參與反演的8個(gè)臺(tái)站，按每次增減一個(gè)臺(tái)站進(jìn)行反演；同時(shí)，基于Hypo2000重新定位結(jié)果，利用本區(qū)域最新的地殼速度模型進(jìn)行了多次反演。具體反演結(jié)果詳見(jiàn)圖4.圖4顯示，在5、6、7、8個(gè)臺(tái)站參與反演的情況下，最佳震源機(jī)制解節(jié)面I的走向最大偏差15°，傾角偏差8°，滑動(dòng)角偏差10°.在參與反演臺(tái)站數(shù)量減少的情況下，約化方差VR有所增高，但其余反演參數(shù)并未出現(xiàn)太大的變化。結(jié)合圖3中5～16 km深度范圍內(nèi)震源機(jī)制性質(zhì)基本穩(wěn)定，及11～14 km范圍內(nèi)VR值穩(wěn)定且取得最大值的結(jié)果，認(rèn)為本文基于8個(gè)地震臺(tái)站反演得到的參數(shù)基本穩(wěn)定。

圖4 不同數(shù)量臺(tái)站參與矩張量反演時(shí)最佳震源機(jī)制解及其他參數(shù)的取值Fig.4 Focal mechanisms and other parameters from moment tensor inversion based on different stations

3.5 可靠性與誤差分析

矩張量反演結(jié)果的合理性與可靠性通?？衫酶饔涗浥_(tái)站垂直分量上的P波初動(dòng)極性在震源球上的分布位置來(lái)分析判定。一般的做法是：首先利用區(qū)域地殼模型、震中距及震源深度計(jì)算記錄臺(tái)站的離源角；然后根據(jù)臺(tái)站方位與離源角將其投影到震源球上，位于壓縮區(qū)的P波初動(dòng)向上，反之向下；通過(guò)考察P波極性與反演得到的震源機(jī)制解的一致性來(lái)驗(yàn)證結(jié)果的合理性與可靠性[18]。

Strike=113°,Dip=53°,Rake=-74°,Source=Snoke圖5 原平ML4.7地震的矩張量反演結(jié)果(a)以及基于P波初動(dòng)的Snoke方法求解的結(jié)果(b)Fig.5 Focal mechanism of Yuanping ML4.7 earthquake by moment tensor inversion (a) and the focal mechanism result by Snoke method (b)

對(duì)于本次原平ML4.7地震，山西測(cè)震臺(tái)網(wǎng)中共有38個(gè)地震臺(tái)站垂直分量記錄到清晰的P波初動(dòng)，各個(gè)臺(tái)站的P波初動(dòng)極性及在震源球上的分布位置情況如圖5(a)所示。除了少數(shù)幾個(gè)臺(tái)站外，其余32個(gè)臺(tái)站的P波極性與基于8個(gè)臺(tái)站波形數(shù)據(jù)反演得到的最佳震源機(jī)制解基本一致。

此外，為進(jìn)一步分析結(jié)果的可靠性與可能誤差，基于38個(gè)地震臺(tái)站垂直分量記錄的清晰的P波初動(dòng)，利用Snoke方法求解了本次地震的震源機(jī)制斷面解(圖5(b))，并將矩張量反演得到的震源機(jī)制解與Snoke方法結(jié)果，及文獻(xiàn)[3]利用CAP方法得到的解進(jìn)行比較，具體震源機(jī)制解參數(shù)詳見(jiàn)表3.可以看到：3種方法得到的震源機(jī)制解存在一定的差異，但震源機(jī)制性質(zhì)一致，節(jié)面參數(shù)基本趨勢(shì)一致，總體上具有可比性；存在差異的主要原因在于，3種方法的理論原理不同，邏輯算法及利用的數(shù)據(jù)類型不同，等等[33]。

表3 原平ML4.7地震震源機(jī)制參數(shù)Table 3 Parameters of Yuanping ML4.7 earthquake focal mechanism solution

4 討論與結(jié)論

本文利用山西地震臺(tái)網(wǎng)寬頻帶數(shù)字波形資料，在Hypo2000重新定位的基礎(chǔ)上，采用TDMT_INV時(shí)間域矩張量反演方法，對(duì)2016年4月7日山西原平ML4.7地震進(jìn)行了矩張量反演計(jì)算與參數(shù)的穩(wěn)定性評(píng)價(jià)，并將獲得的最佳雙力偶解與利用Snoke方法及CAP方法得到的結(jié)果進(jìn)行對(duì)比，驗(yàn)證了結(jié)果的合理性與可靠性。初步得到以下認(rèn)識(shí)和結(jié)論：

1) 基于震中距120 km范圍內(nèi)的8個(gè)觀測(cè)臺(tái)站的波形資料反演得到的地震矩張量結(jié)果顯示：最佳雙力偶分量Percent DC=93，線性補(bǔ)償分量Percent CLVD=7，震源機(jī)制符合雙力偶點(diǎn)源的假設(shè)條件；節(jié)面Ⅰ參數(shù)為走向143°，傾向63°，滑動(dòng)角-66°；節(jié)面Ⅱ走向278°，傾向36°，滑動(dòng)角-128°；最佳震源深度11～14 km.主壓應(yīng)力P軸方位角94°，傾角64°；主張應(yīng)力T軸方位角216°，傾角15°；標(biāo)量地震矩M0=2.492 64×1015N·m，矩震級(jí)結(jié)果為MW=4.2，與MS震級(jí)基本一致。

2) 影響地震矩張量反演結(jié)果穩(wěn)定性的因素較多，如區(qū)域地殼速度結(jié)構(gòu)、地震定位誤差(特別是震源深度誤差)、參與反演臺(tái)站的個(gè)數(shù)、臺(tái)站與震中的方位角誤差，等等。本文基于研究區(qū)最新的地殼速度結(jié)構(gòu)，在Hypo2000重新定位的基礎(chǔ)上，按逐次增減一個(gè)臺(tái)站的方式進(jìn)行多次反演，分析與驗(yàn)證了結(jié)果的穩(wěn)定性，認(rèn)為本文基于8個(gè)地震臺(tái)站反演得到的參數(shù)基本穩(wěn)定。此外，各臺(tái)站垂向分量上清晰的P波極性在震源球上的分布與反演得到的最佳震源機(jī)制解基本一致。

3) 本文矩張量反演得到的最佳震源機(jī)制解與基于Snoke方法和CAP方法的結(jié)果存在一定的差異，但總體具有可比性。結(jié)果間存在誤差的原因在于，不同方法的理論基礎(chǔ)不同、利用的數(shù)據(jù)類型不同等?；赑波初動(dòng)的斷層面解反映的是初始破裂情況，而矩張量解能更好地詮釋整個(gè)地震破裂過(guò)程。

致謝感謝山西地震臺(tái)網(wǎng)中心為本研究提供地震波形數(shù)據(jù)。