車正聲

（浙江省有色金屬地質(zhì)勘查局，浙江 紹興 312000）

1 研究動(dòng)機(jī)與目的

地質(zhì)變化的發(fā)展過(guò)程是一個(gè)復(fù)雜的行為，諸多地質(zhì)結(jié)構(gòu)源參數(shù)可以幫助了解地質(zhì)結(jié)構(gòu)的破裂行為，而其中地質(zhì)結(jié)構(gòu)發(fā)展時(shí)間函數(shù)是了解結(jié)構(gòu)變化機(jī)制的重要參數(shù)之一。地質(zhì)結(jié)構(gòu)變化時(shí)間函數(shù)除反應(yīng)了地震從發(fā)生到結(jié)束的時(shí)間，也能描繪整個(gè)破裂過(guò)程的復(fù)雜變化。

傳統(tǒng)上，利用水平(water-level)頻率域解回旋求得結(jié)構(gòu)變化時(shí)間函數(shù)，但經(jīng)常會(huì)因頻譜相除而造成不穩(wěn)定的現(xiàn)象，近年來(lái)，非負(fù)解(non-negative solution)時(shí)間域解回旋(time-domain deconvolution)被利用去調(diào)查地質(zhì)結(jié)構(gòu)變化及災(zāi)害時(shí)間函數(shù)(例如Ammon et al.,2006)[1]，能夠解出大地震的破裂特征。

本研究利用Hwang(2016)[2]所發(fā)展的非負(fù)解阻尼平滑化時(shí)間域解回旋技術(shù)調(diào)查是否能穩(wěn)定求得地質(zhì)結(jié)構(gòu)變化的歷時(shí)，希冀未來(lái)重新建立地質(zhì)結(jié)構(gòu)變化歷時(shí)與結(jié)構(gòu)矩陣的關(guān)系 (例如 Furumoto and Nakanishi,1983)[3]。

本研究先以四個(gè)地質(zhì)結(jié)構(gòu)變化為例(表1)，將其結(jié)果與先前的研究比較，探討非負(fù)解平滑化時(shí)間域解回旋方法的可用性。

表1 研究所分析地質(zhì)結(jié)構(gòu)變化歷時(shí)

2 遠(yuǎn)場(chǎng)P波資料

本研究使用的遠(yuǎn)場(chǎng)P波數(shù)據(jù)來(lái)自于美國(guó)IRIS地質(zhì)變化數(shù)據(jù)中心。

首先，挑選地質(zhì)變化中心距離30至90的測(cè)站數(shù)據(jù)，避免PcP波及PP波的干擾；接著，依據(jù)各個(gè)地震的余震分布，判斷此地震的破裂方向；最后，挑選垂直破裂方向且具清晰的初達(dá)P波進(jìn)行分析。而挑選垂直破裂方向的數(shù)據(jù)，主要在避免破裂方向性的影響。最后每筆數(shù)據(jù)在經(jīng)儀器修正后，再經(jīng)0.01至0.5 Hz濾波，而截切P波的時(shí)間長(zhǎng)度先依各地質(zhì)結(jié)構(gòu)變化的矩陣配合Furumoto and Nakanishi(1983)的經(jīng)驗(yàn)公式所得震源歷時(shí)的兩倍計(jì)算，最后再繪圖判讀合適的截切時(shí)間長(zhǎng)度。

3 方法與原理

本研究所探討的地質(zhì)結(jié)構(gòu)變化災(zāi)害皆為單向地質(zhì)破裂。對(duì)于單向地質(zhì)破裂的情況而言，在某一方向觀測(cè)到的結(jié)構(gòu)變化歷時(shí)可表示為：

TSPT為地質(zhì)結(jié)構(gòu)破裂歷時(shí)，L是破裂長(zhǎng)度，Vr為破裂速度，C為在地質(zhì)結(jié)構(gòu)災(zāi)害源區(qū)的表面波相速，θ為站方位與災(zāi)害源破裂方向的夾角。

為了移除破裂方向?qū)Y(jié)構(gòu)變化歷時(shí)的影響，選取與破裂方向近乎垂直的測(cè)站進(jìn)行研究。

首先對(duì)這些測(cè)站觀測(cè)到的P波進(jìn)行濾波得到測(cè)站的觀測(cè)波，然后依照Ammon et.al.(2006)的方式產(chǎn)生合成P波，合成P波up（t）可由下式求得:

其中M0為結(jié)構(gòu)矩陣αh、ρh、βh分別為P波的速度、密度以及S波的速度，g（△）為幾何擴(kuò)散因子（geometrical spreading factor），a為 地 球 半 徑，Rp、RpP、RsP分 別P波、pP波及sP波的地質(zhì)結(jié)構(gòu)變化范圍模式（radiation pattern），VP、VpP、VsP分 別 為P波、pP波 及sP波 與 地表反射的反射系數(shù)（coefficients of reflection），ih、jh分別為P波、S波，出射角（take-off angle）CP（ih）為接收站的地表反應(yīng)，f（t）為災(zāi)害源時(shí)間函數(shù)，tP、tpP、tsP分別為P波、pP波、sP波的走時(shí)（travel time），Q（t）為地質(zhì)結(jié)構(gòu)災(zāi)害衰減項(xiàng)，I（t）為儀器反應(yīng)相。其中關(guān)于與擴(kuò)散幾何因子、反射系數(shù)、輻射模式、接收站的地表反應(yīng)等等的計(jì)算，可參考Kanamori and Stewart(1976)[4]、Aki and Richard（1980）[5]等等，災(zāi)害波動(dòng)衰減項(xiàng)目可參考Yoshid（1988）[6]。速度模型采用IASP91速度模型（Kennett and Engdahl,1991）。最后進(jìn)行利用觀測(cè)波的數(shù)值與合成波的數(shù)值進(jìn)行解回旋以得到地質(zhì)結(jié)構(gòu)變化歷時(shí)[7，8]（如圖1所示）。

圖1 由上往下依次是測(cè)站BKE的觀測(cè)波圖、合成波圖以及觀測(cè)波與合成波解回旋后的圖

若g（t）是格林函數(shù)，是s（t）是災(zāi)害源時(shí)間函數(shù)，b（t）是觀測(cè)災(zāi)害波動(dòng)范圍，則g＊s=b，＊表示回旋，其積分型態(tài)為：

上述公式得出，L為基本形態(tài)。

在離散（disrtet）且有限的長(zhǎng)度下，其離散形式為：

△t是取樣率，m是災(zāi)害源時(shí)間函數(shù)s（t）的總點(diǎn)數(shù)。（2）式可改寫(xiě)成矩陣形式如下：

令G、S、b分別等于上式三項(xiàng)中的一項(xiàng)，表示為：

從（4）式中求出S，一般是利用最小平方法逆推（如Menke，1989），其觀念是使觀測(cè)值與預(yù)測(cè)值得誤差平方總和（，eT是e的轉(zhuǎn)置矩陣）最小，即使

圖2 觀測(cè)站解回旋后的波形圖

圖3 地質(zhì)結(jié)構(gòu)變化災(zāi)害分布四種分布位置

最后得到的解（Sest）可表示為：

GT為G的轉(zhuǎn)置矩陣。然而（5）式所得結(jié)果并不穩(wěn)定所以在這一過(guò)程中加入平滑因子（W）與阻尼值（λ2）最終的到：

Sest即為索要求的震源時(shí)間函數(shù)。

4 結(jié)果與討論

表2 震源歷時(shí)與WCMT、GCMT的結(jié)果比較[9]

5 結(jié)論

本研究利用非負(fù)解阻尼平滑化時(shí)間域解回旋法調(diào)查地質(zhì)結(jié)構(gòu)災(zāi)害變化時(shí)間函數(shù)，并在主要變化位置產(chǎn)生合成波用以代替地質(zhì)變化進(jìn)行時(shí)間域解回旋，解決了因在地質(zhì)結(jié)構(gòu)變化附近沒(méi)有機(jī)構(gòu)機(jī)制相似的災(zāi)害而無(wú)法求出地質(zhì)結(jié)構(gòu)變化源時(shí)間函數(shù)的問(wèn)題。通過(guò)加入非負(fù)解與阻尼值來(lái)使所得的結(jié)果變得穩(wěn)定，嘗試多個(gè)阻尼值最終確定當(dāng)阻尼值為100時(shí)的波形圖是較為理想的，說(shuō)明此法是相對(duì)較穩(wěn)定的。