施其斌， 韋生吉*

1 南洋理工大學(xué)亞洲環(huán)境學(xué)院， 新加坡 639798 2 南洋理工大學(xué)新加坡地球觀測(cè)中心， 新加坡 639798

0 引言

地震震源的運(yùn)動(dòng)學(xué)過(guò)程是地震學(xué)研究的一個(gè)基本內(nèi)容.隨著地震波觀測(cè)數(shù)據(jù)的質(zhì)量的提高——包括臺(tái)站密度的增加和觀測(cè)頻段范圍的拓展等——震源運(yùn)動(dòng)學(xué)的發(fā)展經(jīng)歷了從簡(jiǎn)單點(diǎn)源(例如地震的位置、大小、震級(jí)和矩張量解(Aki and Richards, 2002))到更接近真實(shí)地震的復(fù)雜的有限破裂模型(姚振興和紀(jì)晨, 1997; Ji et al., 2002; Xu and Chen, 2002; 許力生和陳運(yùn)泰, 2004; 張勇等, 2008; Zhang et al., 2009; Yue et al., 2012; 岳漢等, 2019)的進(jìn)步.點(diǎn)源對(duì)于中強(qiáng)度(MW6-7)地震的描述往往過(guò)于簡(jiǎn)單，然而有限破裂模型又常常用過(guò)多的參數(shù)來(lái)非唯一地表示地震的破裂過(guò)程.對(duì)后者而言，研究者一般很難估計(jì)斷層幾何等參數(shù)的誤差.地震的多點(diǎn)源解剛好填補(bǔ)了點(diǎn)源和有限破裂模型之間的空缺，彌補(bǔ)了上述兩者的不足之處.通常情況下，多點(diǎn)源解所適用的反演頻段也介于單點(diǎn)源和有限破裂模型的頻段之間.在此頻段，多點(diǎn)源解通常只需要用幾個(gè)(例如3～8個(gè))雙力偶位錯(cuò)點(diǎn)源來(lái)描述地震的破裂過(guò)程，要遠(yuǎn)比單矩張量解更能代表地震的復(fù)雜性.同時(shí)，多點(diǎn)源解的參數(shù)空間遠(yuǎn)小于有限破裂模型，因而不會(huì)過(guò)度擬合觀測(cè)數(shù)據(jù)，更易于對(duì)震源參數(shù)進(jìn)行誤差估計(jì).近年來(lái)，多點(diǎn)源反演已經(jīng)被多次應(yīng)用于大地震研究中(Chen et al., 2014; Duputel and Rivera, 2017; Yue and Lay, 2020; Jia et al., 2020)，此類(lèi)研究通常使用遠(yuǎn)震體波和長(zhǎng)周期的面波.對(duì)于那些密集強(qiáng)震記錄到的中強(qiáng)度地震，多點(diǎn)源反演也獲得了很好的效果 (Shi et al., 2018; Shi and Wei, 2020).然而，對(duì)中強(qiáng)度地震的多點(diǎn)源研究還沒(méi)有充分利用區(qū)域距離(<1000 km)的臺(tái)站記錄.實(shí)際上，除了少數(shù)地震頻發(fā)的發(fā)達(dá)地區(qū)(例如日本和美國(guó)西海岸)，絕大部分中強(qiáng)度地震頻發(fā)的地區(qū)并沒(méi)有密集的強(qiáng)震臺(tái)網(wǎng).相比之下，區(qū)域距離范圍內(nèi)臺(tái)站數(shù)量對(duì)全球絕大部分地區(qū)都更加可觀.因此，是否能利用此類(lèi)數(shù)據(jù)對(duì)中強(qiáng)度地震的破裂過(guò)程進(jìn)行多點(diǎn)源反演非常值得探究.為此，本文通過(guò)對(duì)2016MW7.1熊本地震的研究來(lái)測(cè)試?yán)脜^(qū)域?qū)掝l帶地震數(shù)據(jù)進(jìn)行多點(diǎn)源地震反演.我們將利用前/余震來(lái)對(duì)區(qū)域臺(tái)站的不同分量進(jìn)行路徑校正，以得到適用于主震反演的頻段.然后對(duì)比加入了更為密集的強(qiáng)震和高頻GPS數(shù)據(jù)得到的多點(diǎn)源解和現(xiàn)有的有限破裂模型，來(lái)進(jìn)一步驗(yàn)證多點(diǎn)源解.

我們改進(jìn)了Shi等(2018)所提出的基于馬爾可夫鏈蒙特卡洛(Markov Chain Monte Carlo，MCMC)的反演方法.該方法借鑒了前人所提出的“剪切黏貼”(Cut-And-Paste，CAP)理念 (Zhao and Helmberger, 1994; Zhu and Helmberger, 1996)，將每個(gè)臺(tái)站的三分量地震波數(shù)據(jù)截成Pnl波和面波部分，在使用理論地震波擬合數(shù)據(jù)時(shí)允許對(duì)不同分量使用不同的時(shí)移進(jìn)行擬合，然后計(jì)算時(shí)移后的波形擬合誤差.Shi等(2018)采用該方法，利用強(qiáng)震數(shù)據(jù)研究了2016熊本地震序列的MW6.2前震，結(jié)果揭示了發(fā)震區(qū)斷層幾何形態(tài)的復(fù)雜性.在本文中，我們將該多點(diǎn)源反演方法應(yīng)用到包含高頻GPS數(shù)據(jù)的密集近震觀測(cè)，并拓展到對(duì)區(qū)域地震數(shù)據(jù)的分析.在此基礎(chǔ)上，提出了使用正則化的Yoffe震源時(shí)間函數(shù)來(lái)代替之前的三角形函數(shù)，并且還提出了基于臺(tái)站和波形分量來(lái)選擇濾波參數(shù)，從而最大化地利用地震波形所包含的信息.為了對(duì)每個(gè)分量選擇合適的濾波頻段，我們利用比研究地震震級(jí)更小的地震來(lái)進(jìn)行路徑校正.這些地震的震級(jí)應(yīng)當(dāng)小到在主震反演的濾波頻段可以被合理地簡(jiǎn)化成一個(gè)點(diǎn)源，但同時(shí)也要大到可以被足夠多的臺(tái)站記錄到，不管是強(qiáng)震儀還是寬頻帶地震儀.在后文將給出更多相關(guān)的技術(shù)細(xì)節(jié).

1 2016年熊本地震序列

2016年的MW7.1熊本地震序列發(fā)生在日本九州島的Beppu-Shimabara地塹地區(qū)，起始于Futagawa斷層和Hinagu斷層相交處.該斷層系統(tǒng)的形成主要是由菲律賓板塊向西北方向略微斜向俯沖所導(dǎo)致，板塊俯沖的速度5～6 cm·a-1，斜沖投影在Futagawa斷層和Hinagu斷層上的滑動(dòng)速率為7～8 mm·a-1(Nishimura and Hashimoto, 2006)(圖1a).

序列始于4月14日發(fā)生的MW6.2前震，約2.5 h后，發(fā)生了一個(gè)MW6.0的余震，MW7.1主震發(fā)生在當(dāng)?shù)貢r(shí)間4月16日的凌晨，震中位于前震的東北方向(Kato et al., 2016).該地震序列的余震極其活躍，產(chǎn)生了至少30000個(gè)Mj(JMA(Japan Meteorological Agency)震級(jí)，即日本氣象局測(cè)定的地震震級(jí))大于2.5的余震 (Kato et al., 2016).地震序列還造成了復(fù)雜的地表破裂，以右旋走滑為主但有明顯的正斷層滑動(dòng) (Shirahama et al., 2016)，尤其在兩個(gè)斷層相交處和破裂的東北段(圖1b).

該地震序列有著高質(zhì)量的密集強(qiáng)震觀測(cè)、衛(wèi)星影像觀測(cè)以及實(shí)地考察記錄，為震源運(yùn)動(dòng)學(xué)研究提供了寶貴的數(shù)據(jù).前人利用地震學(xué)數(shù)據(jù)或者聯(lián)合測(cè)地學(xué)數(shù)據(jù)對(duì)MW6.2 前震(Asano and Iwata, 2016)和MW7.1主震(Asano and Iwata, 2016; Yagi et al., 2016; Hao et al., 2017; Kobayashi, 2017; Yue et al., 2017)進(jìn)行了有限破裂模型的反演.但是不同研究者所得到的主震破裂模型之間仍存在不小的差異.例如，Asano和Iwata(2016)的結(jié)果顯示主震的破裂主要集中在8～10 km 深度，位于震中東北約20 km 處；而Yue等(2017)和Kobayashi(2017)的結(jié)果則顯示主震主要的位錯(cuò)比較淺，已經(jīng)破裂到了地表.這些有限破裂模型使用2～4個(gè)斷層面來(lái)近似主震較為復(fù)雜的斷層幾何.相對(duì)的，前人對(duì)MW6.2前震進(jìn)行的有限破裂模型反演采用了相對(duì)簡(jiǎn)單的斷層幾何形態(tài)，僅為一個(gè)矩形的斷層面(Asano and Iwata, 2016).不同于有限破裂模型，Shi等(2018)提出了一種基于MCMC的地震多點(diǎn)源反演方法，用于在相對(duì)高頻反演MW6.2前震所產(chǎn)生強(qiáng)震波形.他們使用了一個(gè)MW5.4 地震的波形記錄來(lái)對(duì)所使用的波形分量進(jìn)行篩選并確定最佳濾波頻段，在此基礎(chǔ)上的多點(diǎn)源反演結(jié)果揭示了MW6.2前震較為復(fù)雜的斷層幾何，即發(fā)現(xiàn)該地震發(fā)生在沿相反方向傾斜的多個(gè)斷層上.

2 方法和數(shù)據(jù)

如前所述，Shi等(2018)所提出的多點(diǎn)源反演方法已經(jīng)在MW6.2前震中得到了成功的應(yīng)用.該方法盡管對(duì)強(qiáng)震數(shù)據(jù)反演非常有效，但還沒(méi)有在高頻GPS數(shù)據(jù)和區(qū)域?qū)掝l帶地震波的反演中嘗試過(guò).為此，本文將該方法拓展到高頻GPS和區(qū)域?qū)掝l帶地震記錄來(lái)對(duì)MW7.1 熊本主震的破裂過(guò)程進(jìn)行研究，并將結(jié)果和強(qiáng)震+高頻GPS記錄反演的多點(diǎn)源解以及已經(jīng)發(fā)表的有限破裂模型進(jìn)行對(duì)比，旨在了解該方法在中強(qiáng)度地震研究中的優(yōu)勢(shì)和局限性.我

們還對(duì)該方法做了進(jìn)一步的改進(jìn).Shi等(2018)的方法中對(duì)所有觀測(cè)臺(tái)站相同的波形分量采用相同的濾波頻段(例如用0.02～0.2 Hz對(duì)Pnl統(tǒng)一進(jìn)行濾波)，但是由于每個(gè)臺(tái)站所對(duì)應(yīng)的速度結(jié)構(gòu)不相同，其波形分量對(duì)應(yīng)的最佳擬合頻段也有所不同.這里我們采用了臺(tái)站相關(guān)的濾波參數(shù)，以最大限度地從波形數(shù)據(jù)中提取地震破裂過(guò)程的信息.另外，還改進(jìn)了反演所用的震源時(shí)間函數(shù)，不同于之前的等腰三角形函數(shù)，這里使用正則化的Yoffe函數(shù)(Tinti et al., 2005)來(lái)作為震源時(shí)間函數(shù)(圖2)，以測(cè)試所用的數(shù)據(jù)是否能進(jìn)一步分辨地震位錯(cuò)發(fā)生的時(shí)間過(guò)程.實(shí)際上Yoffe 函數(shù)更符合滑動(dòng)弱化震源動(dòng)力學(xué)破裂模型中斷層上應(yīng)力的演化過(guò)程(Yoffe, 1951; Tinti et al., 2005)，這對(duì)于理解震源破裂物理過(guò)程具有重要的意義.新的震源時(shí)間函數(shù)用α和β兩個(gè)參數(shù)來(lái)定義函數(shù)的形狀(圖2)，數(shù)學(xué)表達(dá)如下：

(1)

圖2 正則化Yoffe震源時(shí)間函數(shù)的數(shù)學(xué)定義該震源時(shí)間函數(shù)表示為黑色，由三個(gè)對(duì)稱(chēng)軸相同但周期不同的cosine/sine函數(shù)組成，分別顯示為紅色、藍(lán)色和綠色.Fig.2 The mathematical composition of the reguralized Yoffe-type source time functionThe Yoffe function is shown in black, composed of the three cosine/sine functions with the same symmetric axis but different periods and shown in red, blue and green, respectively.

其中f1(t),f2(t)和f3(t)分別用α和β來(lái)表示：

(2)

(3)

(4)

表達(dá)式(1)具有更多的靈活性，既可以表示一個(gè)對(duì)稱(chēng)的函數(shù)(α=0，β=1)，也可以表示非對(duì)稱(chēng)的“長(zhǎng)尾”的函數(shù)(例如α=0.5，β=0.2).

所使用的波形數(shù)據(jù)包括強(qiáng)震、高頻GPS和寬頻帶地震數(shù)據(jù)(圖3).三分量強(qiáng)震數(shù)據(jù)是從KiK-Net和K-net網(wǎng)站上(https:∥www.kyoshin.bosai.go.jp)直接下載，其原始波形為加速度記錄，將其積分至速度和位移波形用于反演.寬頻帶波形數(shù)據(jù)是從NIED(https:∥www.fnet.bosai.go.jp/top.php)上下載，去儀器響應(yīng)后轉(zhuǎn)化成位移和速度波形兩種記錄.我們還從日本的Nippon GPS Data Service Corporation購(gòu)買(mǎi)了1 Hz采樣的主震高頻GPS數(shù)據(jù).除了主震(MW7.1)，還從這些地震臺(tái)網(wǎng)上下載了兩個(gè)校正地震事件的波形數(shù)據(jù).這兩個(gè)校正事件分別發(fā)生在Futagawa斷層(MW5.4，32.776°N，130.815°E，2016-04-14-13∶07∶35.290)和Hinagu斷層上(MW6.0，32.7007°N，130.778°E，2016-04-14-15∶03∶46.450).這兩個(gè)事件均為MW6.2前震的余震，其破裂過(guò)程相對(duì)簡(jiǎn)單，在相對(duì)較高(例如0.3 Hz)的頻段可以看作點(diǎn)源.利用這兩個(gè)校正事件，我們確定了所使用的臺(tái)站和波形分量，以及臺(tái)站相關(guān)的濾波參數(shù)(例如以80%互相關(guān)系數(shù)來(lái)選擇濾波頻段和分量).最終選取了110個(gè)強(qiáng)震臺(tái)，包含276個(gè)Pnl波形分量和251個(gè)面波波形分量，36個(gè)寬頻帶臺(tái)站，包含90個(gè)Pnl波形分量和83個(gè)面波波形分量，以及26個(gè)高頻GPS臺(tái)站.圖3給出了這三個(gè)地震以及MW6.2前震在幾個(gè)代表性臺(tái)站上的0.2 Hz低通濾波濾波波形.可以看到兩個(gè)校正事件(MW5.4，MW6.0)的波形相對(duì)簡(jiǎn)單，而MW6.2前震則要復(fù)雜得多，說(shuō)明其破裂過(guò)程更加復(fù)雜.值得指出的是，在KMM005臺(tái)站上，MW6.2前震的波形甚至比MW7.1主震還要復(fù)雜.這一方面是因?yàn)橹髡鹗窍驏|北方向的單側(cè)破裂，而KMM005臺(tái)站(方位角38°)正好位于破裂的傳播方向上，因此多普勒效應(yīng)將波形壓縮顯得更為簡(jiǎn)單；另一方面是因?yàn)橹髡鸬恼鹪磿r(shí)間函數(shù)更長(zhǎng)、更平滑，導(dǎo)致了主震的波形的頻率成分更低一些，看起來(lái)沒(méi)那么復(fù)雜.而背離主震破裂方向的KMM008臺(tái)站其波形振幅要低很多，但是持續(xù)時(shí)間更長(zhǎng)，波形也更為復(fù)雜.這些都符合主震單側(cè)破裂的特征.同時(shí)寬頻帶臺(tái)站的波形數(shù)據(jù)則顯示MW7.1主震的波形要比其他地震復(fù)雜得多(圖3c，d)，因此將這些數(shù)據(jù)用于多點(diǎn)源反演可以對(duì)其破裂過(guò)程進(jìn)行很好的約束.相對(duì)于強(qiáng)震臺(tái)站，我們需要選擇更低頻段的濾波才能對(duì)區(qū)域臺(tái)站的校正事件波形實(shí)現(xiàn)較高互相關(guān)系數(shù)的擬合.區(qū)域臺(tái)站的路徑校正主要是通過(guò)MW6.0前震實(shí)現(xiàn)的.

圖3 近Futagawa斷層的熊本主震及前震的強(qiáng)震和寬頻數(shù)據(jù)對(duì)比(a)和(b)為強(qiáng)震儀KMM006和KMM005的觀測(cè)數(shù)據(jù)低通濾波后的波形； (c)和(d)為寬頻地震儀TKO和SIB記錄低通濾波后的波形.7.1級(jí)主震的波形由藍(lán)色特別顯示.由上至下分別為豎直、南北和東西分量的原始地震波形.所有波形均通過(guò)了角頻率為0.2 Hz的低通濾波.Fig.3 Strong-motion and broadband observations of the Kumamoto mainshock and foreshocks near the Futagawa fault(a) and (b) show the low-pass filtered waveforms of the strong-motion stations KMM006 and KMM005, respectively. (c) and (d) show the low-pass filtered waveforms of the broadband stations TKO and SIB, respectively. The mainshock waveforms are highlighted in blue. The vertical, north-south and east-west waveforms are arranged from the top to the bottom. All waveforms are low-pass filtered with the corner frequency of 0.2 Hz.

3 反演結(jié)果

我們使用經(jīng)過(guò)路徑校正得到的波形分量和濾波頻段，來(lái)對(duì)主震的多點(diǎn)源解進(jìn)行反演.從單個(gè)點(diǎn)源反演開(kāi)始逐漸增加點(diǎn)源的數(shù)量，直到反演得到的波形互相關(guān)系數(shù)的統(tǒng)計(jì)分布非常接近校正事件的波形擬合效果，這樣可以避免對(duì)數(shù)據(jù)的過(guò)度擬合.圖4給出了聯(lián)合使用強(qiáng)震數(shù)據(jù)和區(qū)域?qū)掝l帶數(shù)據(jù)，和只使用區(qū)域?qū)掝l帶數(shù)據(jù)得到的多點(diǎn)源反演結(jié)果，包括用速度波形和位移波形得到的結(jié)果.包含強(qiáng)震數(shù)據(jù)的反演結(jié)果顯示我們需要至少4個(gè)子事件來(lái)足夠準(zhǔn)確的表示主震的運(yùn)動(dòng)學(xué)破裂過(guò)程.我們也將三點(diǎn)源反演的結(jié)果和四點(diǎn)源解進(jìn)行比對(duì)(圖5)，發(fā)現(xiàn)四點(diǎn)源解包含更多細(xì)節(jié).如圖6所示，當(dāng)只用3個(gè)子事件時(shí)，得到的互相關(guān)系數(shù)的統(tǒng)計(jì)直方圖和用校正事件得到的直方圖還有明顯的差距，而當(dāng)用到4個(gè)子事件時(shí)統(tǒng)計(jì)的結(jié)果和校正事件更為接近.我們注意到3個(gè)和4個(gè)子事件的波形擬合結(jié)果最明顯的差別集中在離地震最近的幾個(gè)強(qiáng)震臺(tái)站上，把幾個(gè)代表性的臺(tái)站在圖7中標(biāo)注了出來(lái).對(duì)這幾個(gè)臺(tái)站，3個(gè)子事件明顯對(duì)波形和振幅的擬合都差于4個(gè)子事件的解.總體而言，這兩種數(shù)據(jù)組合的反演結(jié)果非常穩(wěn)定和可靠，這點(diǎn)從反演的誤差分布可以看出來(lái)(圖8).包含強(qiáng)震數(shù)據(jù)的反演所對(duì)應(yīng)的誤差分析和統(tǒng)計(jì)顯示，各個(gè)子事件的走向、傾角和滑動(dòng)角的誤差都在5°以?xún)?nèi)，震級(jí)的誤差小于0.1，深度和水平位置的誤差則控制在2 km以?xún)?nèi)(圖8).而區(qū)域?qū)掝l帶數(shù)據(jù)所對(duì)應(yīng)的誤差水平要相對(duì)更大一些.同時(shí)我們還看到用區(qū)域?qū)掝l帶速度波形得到的反演結(jié)果和強(qiáng)震數(shù)據(jù)得到的結(jié)果更一致，也更穩(wěn)定.這樣的差別是合理的，因?yàn)槭紫葏^(qū)域數(shù)據(jù)所對(duì)應(yīng)的濾波頻段要更低，而且臺(tái)站距離地震更遠(yuǎn)，對(duì)地震的時(shí)空分辨率不可避免的要更低；其次在相同的濾波頻段速度記錄所對(duì)應(yīng)的波形要比位移波形更復(fù)雜，因此對(duì)各子事件約束更為可靠.總之，區(qū)域?qū)掝l帶數(shù)據(jù)的結(jié)果和強(qiáng)震數(shù)據(jù)結(jié)果非常一致，所對(duì)應(yīng)的誤差范圍也足以分辨斷層幾何的變化和破裂過(guò)程的時(shí)空演化，驗(yàn)證了可以將該方法擴(kuò)展到區(qū)域臺(tái)站的反演中的有效性.

圖4 矩震級(jí)7.1級(jí)的主震破裂過(guò)程以及早期余震分布(a)強(qiáng)震和寬頻數(shù)據(jù)聯(lián)合反演所得的四點(diǎn)源解將主震分為四個(gè)由西南向東北傳播的子事件，其水平位置用較大圓圈表示，相應(yīng)的震源機(jī)制解繪制成紅色沙灘球形狀.余震則由較小圓圈表示，并且震源深度用不同顏色顯示.黑色和紅色線段分別表示地震前所測(cè)繪的活動(dòng)斷層和此次地震所產(chǎn)生的地表破裂 (Shirahama et al., 2016).近地震臺(tái)站用黑色三角形表示.(b)該四點(diǎn)源解(同(a))與測(cè)地學(xué)資料導(dǎo)出的有限斷層模型的對(duì)比.此次地震序列的主要事件由紅色五角星標(biāo)記.右下角地震矩張量投影依次為四點(diǎn)源解的合成解、Global Centroid-Moment-Tensor (CMT)反演解和美國(guó)地質(zhì)調(diào)查局的公布結(jié)果.(c)聯(lián)合高采樣GPS數(shù)據(jù)，強(qiáng)震地震數(shù)據(jù)和寬頻地震數(shù)據(jù)的反演結(jié)果.左上角為該四點(diǎn)源解的子事件震源時(shí)間函數(shù)(黑色)以及累加震源時(shí)間函數(shù)(灰色陰影).地震儀和GPS臺(tái)站分別由黑色三角形和白色方塊表示.(d)只利用寬頻地震儀反演的四點(diǎn)源解.其中黑色和灰色分別表示位移波形和速度波形的反演結(jié)果.左上角為位移波形反演得到的震源時(shí)間函數(shù)，而右下角為速度波形得到的震源時(shí)間函數(shù).Fig.4 Rupture process of the MW7.1 mainshock and early aftershocks(a) The four-point-source solution of the joint inversion of the strong-motion and broadband data divides the mainshock into four subevents propagating from south-west to north-east. Four large circles denote their locations and four red beachballs indicate their focal mechanisms. Aftershocks are indicated by small circles that are color-coded by event depth. Black and red lines represent the previously mapped active faults and the surface ruptures (Shirahama et al., 2016) of the MW7.1 earthquake, respectively. Near-fault strong-motion stations are denoted by black triangles. (b) The comparison between the four-point-source solution (same as (a)) and the finite-fault model derived from the geodetic data. The red stars show the epicenters of the major earthquakes that occurred in this sequence. The lower-right inset shows the total moment tensor of the four-point-source solution in comparison with the Global Centroid-Moment-Tensor (CMT) and United States Geological Survey (USGS) solutions. (c) The four-point-source solution obtained from the inversion incorporating the high-rate GPS, strong-motion and broadband seismic data. Top-left inset shows the source time functions of subevents (black) and the total sum (gray-shaded area). The seismic and GPS stations are denoted by black triangles and white squares, respectively. (d) Four-point-source solutions obtained from the inversions with only F-net broadband data. Black and gray color indicates inversion results obtained from the displacement and velocity waveforms, respectively. The top-left shows the source time functions corresponding to the solution obtained from the displacement waveform inversion. The lower-right is corresponding to the solution obtained from the velocity waveform inversion.

圖5 強(qiáng)震和寬頻數(shù)據(jù)聯(lián)合反演所得的四點(diǎn)源解和三點(diǎn)源解的對(duì)比(a)四點(diǎn)源解的震源參數(shù)和地圖展示.空心圓圈表示子事件的水平位置，與其相對(duì)應(yīng)的震源機(jī)制解相連.空心和實(shí)心的三角形分別表示朝向和遠(yuǎn)離破裂方向的地震臺(tái)站.底部小圖顯示四個(gè)子事件震源時(shí)間函數(shù)(黑色)以及累加震源時(shí)間函數(shù)(灰色區(qū)域).(b)與(a)相似，但為三點(diǎn)源解的震源參數(shù)和地圖表示.底部小圖中綠色虛線為(a)中的累加震源時(shí)間函數(shù).Fig.5 Comparison between the four-point-source and three-point-source solutions obtained from the joint inversions of the strong-motion and broadband data(a) Source parameters and the map view of the four-point-source solution. The subevent locations are denoted by open circles, connected with the beachballs that show the corresponding focal mechanisms. The open and solid triangles indicate the strong-motion stations that are toward and away from the rupture directivity. The bottom inset shows the source time functions of four subevents (black) and the sum of them (gray-shaded area). (b) Similar to (a), but for the source parameters and the map view of the three-point-source solution. The green dashed line in the bottom inset shows the cummulative source time function in (a).

圖6 路徑校正事件和主震三點(diǎn)源以及四點(diǎn)源反演的波形擬合效果比較此處只對(duì)比了利用強(qiáng)震和寬頻數(shù)據(jù)聯(lián)合反演的情況.四幅直方圖分別統(tǒng)計(jì)不同波形分量的互相關(guān)系數(shù).Fig.6 Comparison of the quality of the waveform fits between the path calibration event, three-point-source inversion and four-point-source inversionHere the strong-motion and broadband waveforms are jointly used. The four panels of histograms show the statistics of the cross-correlation coefficients of different waveform components.

圖7 部分波形擬合效果(a)和(b)分別為圖5中四點(diǎn)源解和三點(diǎn)源解的波形擬合.帶通濾波后的觀測(cè)數(shù)據(jù)和合成數(shù)據(jù)分別顯示為黑色和紅色，其擬合所需時(shí)移以及表示擬合效果的互相關(guān)系數(shù)分別標(biāo)記在波形下方.臺(tái)站名、臺(tái)站震中距、臺(tái)站方位角標(biāo)記在波形左側(cè)，波形分量名標(biāo)記在頂部.濾波頻段因臺(tái)站不同而有所區(qū)別，Pnl波的頻段總體在0.01～0.30 Hz，面波的頻段總體在0.01～0.22 Hz.Fig.7 Selected waveform fits(a) and (b) are waveform fits of the four-point-source and three-point-source shown in Fig.5. The bandpass filtered data and synthetic waveforms are shown in black and red, respectively. Numbers below each set of waveforms are the time shift used for alignment and cross-correlation coefficient of the two waveforms. Station names are shown on the left with the distance and azimuth indicated above and below respectively. The names of component are shown on the top. The filtering frequency is station-dependent, within 0.01～0.30 Hz for Pnl wave and 0.01～0.22 Hz for surface waves.

圖8 多組數(shù)據(jù)反演所得四點(diǎn)源解的不確定性分析(a)表示強(qiáng)震儀波形和寬頻地震儀波形數(shù)據(jù)的聯(lián)合反演不確定性.(b)表示GPS，強(qiáng)震儀波形和寬頻地震儀波形數(shù)據(jù)的聯(lián)合反演不確定性.(c)和(d)各代表只使用寬頻地震儀位移波形和只用寬頻地震以速度波形的反演的不確定性分析.以上四幅圖通過(guò)相似方式呈現(xiàn)MCMC采樣過(guò)程來(lái)反映其不確定性.四個(gè)子事件按照發(fā)震順序分別由紅、藍(lán)、紫、綠四色顯示.地圖上主要繪制了四個(gè)子事件水平位置采樣和反演得到的最佳震源機(jī)制.左上角繪制四個(gè)子事件斷層解抽樣.左側(cè)中間小圖為子事件震源時(shí)間函數(shù)的抽樣，其中粗黑線為反演最終震源時(shí)間函數(shù)，并且四點(diǎn)源累加震源時(shí)間函數(shù)由灰色填充區(qū)域顯示.右下角小圖為PP′豎直剖面上投影的子事件的深度抽樣，最佳深度顯示為黑色圓圈.Fig.8 Uncertainty analysis of the four-point-source inversions using multiple datasets(a) shows the uncertainty of the joint inversion of the strong-motion and broadband waveforms. (b) shows the uncertainty of the joint inversions incorporating the GPS, strong-motion and broadband data. (c) and (d) are for the inversions using only the broadband displacement and only the velocity waveforms, respectively. The uncertainty is represented by MCMC samples of the four subevents. The four subevents are shown in red, blue, purple and green in the order of occurrence. The sampled horizontal locations and the optimal focal mechanisms is shown in map view. The top-left inset shows the sampled fault-plane solutions of the subevents. The middle-left inset shows the sampled source time functions of the subevents with the final solution denoted by the thick black lines. The gray shaded area is the total source time function of the four-point-source solution. The bottom-right inset shows sampled depths projected to the PP′ cross-section with the optimal depths denoted by black circles.

此外，我們注意到當(dāng)把高頻GPS數(shù)據(jù)包括在反演中時(shí)，得到的結(jié)果和強(qiáng)震數(shù)據(jù)得到的結(jié)果幾乎完全相同(圖9)，這說(shuō)明數(shù)據(jù)之間的一致性非常好，也進(jìn)一步驗(yàn)證多點(diǎn)源反演的穩(wěn)定性和可靠性.

圖9 GPS、強(qiáng)震儀和寬頻地震數(shù)據(jù)的聯(lián)合反演所得到的五點(diǎn)源解(a)五點(diǎn)源解的地圖表示.最靠近主斷層的GPS臺(tái)站由白色方塊表示.右下角為GPS臺(tái)0701附近區(qū)域的放大視圖.左上角為五個(gè)子事件對(duì)應(yīng)的震源時(shí)間函數(shù)以及累積震源時(shí)間函數(shù).(b)部分GPS臺(tái)站位移波形擬合，臺(tái)站位置如(a)所示.帶通濾波后的觀測(cè)數(shù)據(jù)和合成數(shù)據(jù)分別顯示為黑色和紅色，其擬合所需時(shí)移以及表示擬合效果的互相關(guān)系數(shù)分別標(biāo)記在波形下方.臺(tái)站名、臺(tái)站震中距、臺(tái)站方位角標(biāo)記在波形左側(cè)，波形分量名標(biāo)記在頂部.濾波頻段因臺(tái)站不同而有所區(qū)別，總體在0.001～0.1 Hz.(c)0701觀測(cè)波形(黑色)與五點(diǎn)源解的合成波形(紅色)以及圖4a中四點(diǎn)源解的合成波形(淺藍(lán))的對(duì)比.Fig.9 Five-point-source solution obtained from the joint inversion of GPS, strong-motion and broadband data(a) The map view of the five-point-source solution. The nearest GPS stations are denoted by the white squares. The lower-right inset is a enlarged view near station 0701. The top-left inset shows the source time functions of the five subevents and the cumulative source time function. (b) Selected waveform fits of GPS stations shown in (a). The bandpass filtered data and synthetic waveforms are shown in black and red, respectively. The numbers below each set of waveforms are the time shift used for alignment and cross-correlation coefficient of the two waveforms. Station names are shown on the left with the distance and azimuth indicated above and below respectively. The names of component are shown on the top. The filtering frequency is station dependent and is within 0.001～0.1 Hz. (c) The comparison between the observation (black) of 0701, five-point-source synthetics (red), and four-point-source synthetics (cyan) that was calculated with the solution in Fig.4a.

近震(強(qiáng)震+高頻GPS)和區(qū)域數(shù)據(jù)聯(lián)合反演的結(jié)果顯示需要至少4個(gè)子事件來(lái)表示主震.第一個(gè)子事件(MW6.75)位于12 km 深度，非?？拷湔鹬?大約在Hinagu斷層和Futagawa斷層的交界處)，具有一個(gè)高傾角(78°)的純走滑震源機(jī)制解，其走向(208°)和Hinagu斷層的走向非常一致.第二個(gè)(MW6.67)和第三個(gè)(MW6.62)子事件的深度幾乎相同(7 km)，分別位于震中東北側(cè)12 km和20 km處，它們的震源機(jī)制解(圖5a)也非常相似，值得指出的是該機(jī)制解的斷層走向(230°±3°)和Futagawa斷層的走向(232°)更為接近.從位置和斷層幾何來(lái)看，這兩個(gè)子事件應(yīng)當(dāng)發(fā)生在Futagawa斷層上，顯示少量的正斷層分量.最后一個(gè)子事件(MW6.54)最淺(2 km)，其震源機(jī)制解顯示斜向的滑動(dòng)角，其走滑和正斷層分量幾乎相當(dāng).整個(gè)地震的震源持續(xù)時(shí)間大約20 s，最大的地震矩釋放發(fā)生在～8 s.我們得到的反演結(jié)果和Yue等(2017)以及Hao等(2017)的主震有限破裂模型具有很好的一致性，尤其是和Hao等(2017).我們的多點(diǎn)源解所揭示的斷層幾何形態(tài)，包括最北段的低傾角性質(zhì)也與Yue等(2017)通過(guò)余震分布約束得到的幾何狀態(tài)吻合.我們將多點(diǎn)源解相加得到的等效矩張量解和其他機(jī)構(gòu)給出的單個(gè)點(diǎn)源矩張量解也具有高度一致性(圖4b)，該一致性是對(duì)多點(diǎn)源解的有效驗(yàn)證(例如Wei et al.,2011; Avouac et al.,2014；Shi et al.,2018).

4 討論

反演的結(jié)果顯示，隨著地震破裂的發(fā)展，子事件(1～4)質(zhì)心深度不斷變淺并且越往后正斷滑移的分量越大，其水平位置處的余震的深度也隨之變淺.我們的子事件深度和Hao等(2017)和Yue等(2017)的有限破裂模型結(jié)果非常一致.Hao等(2017)將其解釋為斷層發(fā)震層深度(或者溫度)受附近的Aso火山影響所致，越靠近火山3～400 ℃等溫線的深度越淺，Lythgoe等(2021)對(duì)2018年Lombok 地震序列的解釋也類(lèi)似.第二和第三個(gè)子事件所在的位置(圖4a)正好對(duì)應(yīng)了余震的空區(qū)，余震則分布在更深和更淺處，這和其他地震(例如Wei et al., 2011)的觀測(cè)一致，可以很好地用應(yīng)力變化來(lái)解釋?zhuān)籝ue等(2017)也指出了這一點(diǎn).然而，第四個(gè)子事件附近幾乎沒(méi)有余震，而且在沿著斷層走向的方向約15 km處余震才重新出現(xiàn).Yue等(2017)提出了地形應(yīng)力加載對(duì)斷層上正應(yīng)力和剪應(yīng)力的影響，其模型計(jì)算發(fā)現(xiàn)，盡管設(shè)置了傾角為65°的正斷層，地形產(chǎn)生的正應(yīng)力卻顯著升高反而不利于正斷層滑動(dòng).除溫度和地形應(yīng)力的影響之外，不均勻分布的孔隙流體壓力的潛在影響也不能忽視.Aizawa等(2017)對(duì)于熊本地震區(qū)域電阻率研究結(jié)果顯示，F(xiàn)utagawa斷層淺部近Aso火山部分存在明顯高電導(dǎo)率，與高流體含量的特征相符.Matsumoto等(2018)利用Asano和Iwata(2016)的有限斷層模型，結(jié)合歷史地震活動(dòng)，計(jì)算了Futagawa斷層上的理論孔隙流體壓力，其結(jié)果表明Aso火山附近的淺部孔隙流體壓力超過(guò)了該斷層的理論強(qiáng)度.這對(duì)斷層黏滑特性的改變是比較可觀的，同時(shí)也促進(jìn)了主震破裂過(guò)程中該區(qū)域彈性勢(shì)能的充分釋放.值得注意的是，盡管本次地震序列中Aso火山口西北緣缺乏余震，但從2002年至2016年歷史地震活動(dòng)記錄中，Uchide等(2016)發(fā)現(xiàn)這一區(qū)域?qū)嶋H上發(fā)生了大量地震.因此，這一區(qū)域的地震活動(dòng)性有待未來(lái)進(jìn)一步研究.

我們通過(guò)多點(diǎn)源反演所得到的四個(gè)子事件顯示出Futagawa斷層相比于Hinagu斷層具有更多的正斷層成分，尤其是其北段靠近Aso火山口附近的最后一個(gè)子事件擁有幾乎同等大小的正斷層和走滑分量.Yue等(2017)則指出該部分主震滑移包含較大比例的正斷層成分.實(shí)際上，在相對(duì)遠(yuǎn)離Aso火山口的Futagawa斷層部分，第二個(gè)和第三個(gè)子事件由走滑主導(dǎo).如此差異可能是由于火山附近地形增加了正斷層性質(zhì)的相對(duì)剪切應(yīng)力，而未同時(shí)增加走滑方向的剪切應(yīng)力.精細(xì)的震后實(shí)地勘測(cè)(Shirahama et al., 2016)表明，F(xiàn)utagawa斷層地表破裂主要為走滑性質(zhì)，而其東南側(cè)大約2 km處的Idenokuchi斷層的地表破裂為正斷層性質(zhì).結(jié)合野外測(cè)量和衛(wèi)星遙感數(shù)據(jù)，Toda等(2016)提出了傾斜滑移分解模型，即Futagawa斷層主要釋放了長(zhǎng)期彈性積累的走滑分量，而Idenokuchi斷層主要釋放了正斷分量，從而形成小規(guī)模地塹.而靠近火山口的Futagawa斷層可能缺乏這種走滑和正斷層的分解機(jī)制，從而傾滑過(guò)程(正斷+走滑)主導(dǎo)了同震位錯(cuò).長(zhǎng)期的地殼尺度的地表形變觀測(cè)也可以為我們分析斷層性質(zhì)成因提供參考.通過(guò)對(duì)2000—2010年分布在九州島的GNSS數(shù)據(jù)的反演，Mochizuki和Mitsui(2016)發(fā)現(xiàn)Aso火山附近尤其是其西緣的GPS臺(tái)站(32.87°N，131.00°E)的下移分量會(huì)被均一結(jié)構(gòu)模型遠(yuǎn)遠(yuǎn)低估；而在Aso下方添加一塊與火山巖漿活動(dòng)相關(guān)的席狀形變?cè)磩t可以充分模擬出該臺(tái)站的下移.該結(jié)果表明，由于特殊的深部地質(zhì)構(gòu)造，Aso附近的Futagawa斷層可能在長(zhǎng)期地殼形變中積累比遠(yuǎn)離Aso區(qū)域更多的豎直方向的彈性應(yīng)變，從而提高其地震活動(dòng)釋放的正斷層分量.盡管從衛(wèi)星干涉影像中并未獲得最后一個(gè)子事件附近的正斷層所產(chǎn)生的明顯地表形變，但我們認(rèn)為最后一個(gè)子事件所揭示的淺部?jī)A滑的特性是可靠的，并且可能并未在地表產(chǎn)生容易被勘測(cè)到的正斷層成分.

我們注意到位于斷層上或者離斷層很近的臺(tái)站(GPS臺(tái)0701和強(qiáng)震儀KMMH16)的波形擬合不是很理想(圖9c，淺藍(lán)色波形).我們認(rèn)為這可能是多點(diǎn)源解對(duì)有限破裂過(guò)程近似所造成的影響.因?yàn)閿鄬由系呐_(tái)站更容易受到最近破裂的影響，而這些局部破裂振幅雖然不比其他破裂更大，但由于非?？拷_(tái)站，其近場(chǎng)項(xiàng)(隨距離的立方衰減)波場(chǎng)振幅在該臺(tái)站的記錄更大，因此沒(méi)有被反演.作為對(duì)比，我們參考了其他學(xué)者對(duì)于熊本地震震源過(guò)程的研究，發(fā)現(xiàn)大部分有限斷層模型都會(huì)明顯低估0701臺(tái)的總水平位移；而Zhang等(2018)的模型則是包含了靠近0701的近地表高滑移區(qū)，從而可以充分?jǐn)M合0701的靜態(tài)位移，盡管其他GPS臺(tái)站位移量擬合不是很好.為了更好的擬合在斷層上的波形數(shù)據(jù)，必須要考慮到有限破裂過(guò)程的影響，使用有限破裂模型反演對(duì)這些近臺(tái)數(shù)據(jù)的擬合將會(huì)有進(jìn)一步的提高.為了測(cè)試GPS臺(tái)站對(duì)于此次地震破裂究竟對(duì)破裂的矩心位置還是對(duì)最近的破裂更加敏感，我們?cè)黾恿巳M模擬測(cè)試(圖10).對(duì)于每組測(cè)試，先用六個(gè)點(diǎn)源構(gòu)建一個(gè)輸入震源模型.該輸入震源模型包含兩個(gè)滑移集中區(qū)，其中一個(gè)位于西北部，靠近熊本主震的震中，由一個(gè)MW6.7的點(diǎn)源表示；另一個(gè)位于東北部，靠近GPS觀測(cè)站0701，由線性排列的五個(gè)MW6.3的點(diǎn)源組成(如圖10中紅色圓圈所示).我們用構(gòu)建的輸入震源模型合成理論波形，并且加入相對(duì)強(qiáng)度為20%的高斯噪聲，用于反演測(cè)試.三組模擬測(cè)試分別針對(duì)了第二個(gè)滑移集中區(qū)的矩心位置相對(duì)臨近(圖10a)和遠(yuǎn)離GPS臺(tái)0701(圖10b和10c)的情況，其中圖10b和10c分別對(duì)應(yīng)該矩心位置位于0701的西側(cè)和東側(cè)的輸入震源模型.在接下來(lái)的反演測(cè)試中，我們只用兩個(gè)點(diǎn)源去恢復(fù)輸入震源模型，并且使用了和實(shí)際數(shù)據(jù)一樣的臺(tái)站選擇(包括近發(fā)震斷層的GPS臺(tái)站)和濾波頻段.三組反演結(jié)果(圖10)很好地恢復(fù)了輸入震源模型的矩心位置，(近似)震源時(shí)間函數(shù)以及震源機(jī)制解.因此，盡管兩個(gè)點(diǎn)源簡(jiǎn)化了震源模型，加入接近斷層的波形數(shù)據(jù)依然能有效恢復(fù)輸入模型.仍應(yīng)注意，在一些特殊情況下，即使使用了有限斷層模型，擬合非?？拷鼣鄬拥牡卣鹩^測(cè)也存在挑戰(zhàn).比如在對(duì)斷層結(jié)構(gòu)復(fù)雜性缺乏充分認(rèn)知的情況，或是臺(tái)站所在位置的地質(zhì)條件對(duì)波場(chǎng)傳播產(chǎn)生巨大影響的情況.例如，Lee等(2019)對(duì)2018年MW6.4中國(guó)臺(tái)灣花蓮地震的有限斷層反演中發(fā)現(xiàn)，盡管設(shè)置了三個(gè)斷層來(lái)模擬近震數(shù)據(jù)，但是該模型遠(yuǎn)遠(yuǎn)低估緊靠斷層的尤其是其中兩個(gè)斷層之間的臺(tái)站的速度振幅.他們將這解釋為該臺(tái)站所具有的較強(qiáng)的場(chǎng)地效應(yīng)以及厚沉積層的影響.對(duì)于MW7.1熊本地震，我們?cè)噲D添加一個(gè)額外的子事件來(lái)模擬0701的波形.如圖9所示，通過(guò)聯(lián)合反演包含0701在內(nèi)的GPS數(shù)據(jù)、強(qiáng)震儀波形和寬頻地震儀波形，獲得五點(diǎn)源解(圖9a).和四點(diǎn)源解類(lèi)似，五點(diǎn)源解依舊顯示了自西北向東南方向的連續(xù)破裂過(guò)程，并且主要破裂區(qū)域的斷面幾何、震級(jí)大小和震源深度等都非常相似.但主要區(qū)別在于，五點(diǎn)源解將四點(diǎn)源解的最后一個(gè)子事件分為兩個(gè)位于Aso火山口旁的更小的正斷層為主的子事件(MW6.2-6.3)，并且這兩個(gè)子事件代表不同的斷層走向(圖9a及其右下角小圖).考慮這兩個(gè)子事件的位置誤差在內(nèi)，它們反映出Futagawa斷層在Aso火山口附近形成了兩條分支斷層，其中一條保持主斷層的走向(大致232°，由第五個(gè)子事件表示)，另一條更加近東西走向(大致250°，由第四個(gè)子事件表示).此處的斷層分支很可能對(duì)地震破裂的停止起到了關(guān)鍵性作用.該斷層分支與Aso火山活動(dòng)以及火山特殊地質(zhì)結(jié)構(gòu)的關(guān)系值得進(jìn)一步探究，包括對(duì)長(zhǎng)期歷史地震數(shù)據(jù)、火山監(jiān)測(cè)數(shù)據(jù)、Aso區(qū)域地殼結(jié)構(gòu)以及巖石學(xué)特征等可以提供較為全面的研究角度.上述五點(diǎn)源解可以非常好地?cái)M合臺(tái)站0701的波形，并且合成圖和觀測(cè)波形的互相關(guān)系數(shù)達(dá)到100%(圖9c).值得注意的是，盡管先前所得到的四點(diǎn)源解無(wú)法理想地?cái)M合0701的觀測(cè)波形，但是可以與五點(diǎn)源解幾乎同等程度地?cái)M合其余GPS臺(tái)站(圖9b).總之，這印證了多點(diǎn)源反演方法的優(yōu)勢(shì)：在波形數(shù)據(jù)指導(dǎo)下，通過(guò)提高點(diǎn)源數(shù)量來(lái)測(cè)試觀測(cè)數(shù)據(jù)對(duì)地震破裂過(guò)程的敏感性.比如在不使用緊靠Futagawa斷層的0701GPS臺(tái)站，四點(diǎn)源解可以充分?jǐn)M合其余近震GPS、強(qiáng)震儀和區(qū)域?qū)掝l地震儀的波形數(shù)據(jù)；而更加反映局部斷層破裂過(guò)程的五點(diǎn)源解則需要利用非?？拷摼植繑鄬拥?701臺(tái)站的波形.這一經(jīng)驗(yàn)突出了緊靠潛在的發(fā)震斷層布置觀測(cè)臺(tái)站的重要性，這對(duì)普遍的地震破裂過(guò)程的研究非常關(guān)鍵(比如2018年MW6.4花蓮地震).

圖10 有限破裂效應(yīng)的模擬測(cè)試(a) 輸入震源模型由六個(gè)點(diǎn)源組成(紅色)，用來(lái)合成模擬波形，用于后續(xù)兩點(diǎn)源反演.反演結(jié)果(藍(lán)色)有效地恢復(fù)了輸入震源模型中的兩個(gè)滑移集中區(qū)的矩心位置、震源機(jī)制解以及震源時(shí)間函數(shù)(見(jiàn)右下小圖).輸入震源模型和反演結(jié)果所對(duì)應(yīng)的GPS臺(tái)站0701 的波形三分量分別由左上小圖中的紅線和藍(lán)線表示.濾波頻段為0.001～0.1 Hz. 圖(b)和(c)與(a)的繪制方法類(lèi)似，對(duì)應(yīng)另外兩組模擬測(cè)試.(b)中第二個(gè)滑移集中區(qū)相對(duì)于(a)朝西南方向平移了大約2 km.相反，(c)中第二個(gè)滑移集中區(qū)相對(duì)于(a)朝東北方向平移了大約2 km.Fig.10 Synthetic tests for the finite rupture effect(a) The input six-point-source model (red) are used to compute the synthetic waveforms for a two-point-source inversion. The inversion result (blue) well recovered the centroid locations, focal mechanisms and source-time functions (lower-right inset) of the two asperities of the input source model. The three-component waveforms computed for the input and inverted source models are shown in red and blue, respectively, for GPS station 0701 and filtered to 0.001～0.1 Hz. (b) and (c)Similar to (a), but for another two synthetic tests. The centroid locations of the second asperity in (b) is shifted to southwest by about 2 km relative to (a). The centroid locations of the second asperity in (c) is shifted to northeast by about 2 km relative to (a).

此外，0701臺(tái)的強(qiáng)振幅波形有可能包含了其他同震地質(zhì)事件的貢獻(xiàn)，比如該地震序列引發(fā)的最大的山體滑坡(Hung et al., 2018; Song et al., 2019)剛好發(fā)生在0701西北方約1 km處.Hung等(2018)認(rèn)為該滑坡附近的強(qiáng)震儀明顯記錄到了滑坡造成的位移波形，并由此確定滑坡對(duì)應(yīng)的波形起始時(shí)間為主震發(fā)震時(shí)刻之后18.45～21 s.這一結(jié)果支持該最大滑坡是由主震動(dòng)態(tài)觸發(fā).但是是否可以從0701臺(tái)站的波形數(shù)據(jù)中分離出明顯的對(duì)應(yīng)滑坡的成分，這需要對(duì)滑坡的運(yùn)動(dòng)學(xué)過(guò)程進(jìn)行進(jìn)一步模型設(shè)置，這超出了本文的探究范圍.我們用剪切位錯(cuò)(雙力偶)源對(duì)0701波形的擬合已接近完美，并且得到的斷層走向和地表破裂極其吻合，我們認(rèn)為滑坡對(duì)0701波形的貢獻(xiàn)應(yīng)當(dāng)比較有限.

熊本地震的寬頻帶臺(tái)站數(shù)據(jù)質(zhì)量相當(dāng)高，并且主要的傳播路徑是在大陸上，如果地質(zhì)條件更為復(fù)雜，例如地震波傳播要經(jīng)過(guò)大陸和海洋，尤其是經(jīng)過(guò)海溝，或者是發(fā)生在海溝附近的地震，其對(duì)應(yīng)的地震波記錄會(huì)更加復(fù)雜.在此情況下，一維格林函數(shù)就不再適用，如果仍舊使用一維格林函數(shù)則必須要使用更長(zhǎng)周期來(lái)對(duì)數(shù)據(jù)進(jìn)行濾波，那么對(duì)震源破裂過(guò)程的分辨率就會(huì)下降.而Qian等(2019)和Wu等(2018)所提出的方法和應(yīng)用表明，三維格林函數(shù)對(duì)源區(qū)結(jié)構(gòu)復(fù)雜的波場(chǎng)傳播在10～20 s的周期是不可忽略的，這些研究也表明在沒(méi)有很厚沉積層的情況下考慮地形和水層就可以對(duì)波形擬合大為改進(jìn).因此使用三維格林函數(shù)進(jìn)行多點(diǎn)源反演是將來(lái)一項(xiàng)重要的工作.我們注意到Duputel和Rivera(2017)在地震多點(diǎn)源反演中使用了全球SEM(Spectral Element Method，譜元法)程序(Komatitsch和Tromp, 2002)計(jì)算的三維格林函數(shù)，但是全球版SEM程序?qū)Ｋ褂昧藨?yīng)力近似邊界條件，并不能模擬水波及其多次波.可能也是由于這個(gè)原因，Duputel和Rivera(2017)盡管使用了三維格林函數(shù)，依然用了50 s以及更長(zhǎng)周期的地震波來(lái)對(duì)2016 年發(fā)生的MW7.8 Kaikouro地震進(jìn)行反演.

我們所使用的正則化的Yoffe函數(shù)在多點(diǎn)源地震反演更符合震源破裂動(dòng)力學(xué)(Tinti et al., 2005).在地震破裂過(guò)程的實(shí)際研究中，由于震源時(shí)間函數(shù)和其他震源參數(shù)之間存在取舍關(guān)系，并且高頻地震波的快速衰減不容忽視，因此破裂過(guò)程的震源時(shí)間函數(shù)往往極難分辨.顯然，在波形擬合可靠(經(jīng)過(guò)校正)的情況下，反演所用的頻段越高就越能分辨破裂過(guò)程的震源時(shí)間函數(shù).我們的結(jié)果顯示前三個(gè)子事件的震源時(shí)間函數(shù)更接近于高斯函數(shù)，在速度和位移波形反演中均是如此.最后一個(gè)子事件在速度波形反演結(jié)果中更接近于Yoffe函數(shù).對(duì)比于離該子事件最近的幾個(gè)臺(tái)站波形(圖7)，我們認(rèn)為這是一個(gè)可靠的結(jié)果.這說(shuō)明淺部的地震破裂很可能受到前面子事件造成的動(dòng)態(tài)應(yīng)力傳播的影響.在反演所用的頻段，我們不需要具有強(qiáng)不對(duì)稱(chēng)性的Yoffe函數(shù)來(lái)擬合數(shù)據(jù)，這有兩種解釋?zhuān)?1)熊本地震的主震的震源時(shí)間函數(shù)確實(shí)更接近于高斯函數(shù)；(2)需要用更高頻段的地震波來(lái)進(jìn)一步分析.我們提出的正則化的Yoffe函數(shù)震源時(shí)間函數(shù)可以在其他地震的破裂過(guò)程研究中使用，其靈活性將有助于我們進(jìn)一步了解破裂形成和發(fā)展的物理機(jī)理.

5 結(jié)論

我們將一個(gè)基于MCMC的地震多點(diǎn)源波形反演方法擴(kuò)展到了區(qū)域距離的地震波形數(shù)據(jù)，并對(duì)2016年MW7.1熊本地震進(jìn)行了應(yīng)用.結(jié)果顯示區(qū)域地震數(shù)據(jù)可以較好地分辨此類(lèi)中強(qiáng)度(M6-7)地震的時(shí)空破裂過(guò)程.利用離斷層很近的0701高頻GPS數(shù)據(jù)可以分辨出地震破裂最后階段更為精細(xì)的破裂過(guò)程.我們可以將該方法進(jìn)一步推廣到全球中強(qiáng)度地震的運(yùn)動(dòng)學(xué)過(guò)程研究中，來(lái)更準(zhǔn)確地確定基本的震源物理參數(shù)，以理解震源物理過(guò)程.

致謝謹(jǐn)此祝賀陳颙先生從事地球物理教學(xué)科研工作60周年.感謝兩位匿名審稿人寶貴的建設(shè)性意見(jiàn)，感謝日本Kik-net，K-net，F(xiàn)-net為本文研究所提供的波形數(shù)據(jù).