基于實時觀測數據快速判定地震破裂方向的方法研究

席楠 李小軍13) 楊天青

1)中國地震局地球物理研究所,北京 100081 2)中國地震臺網中心,北京 100045 3)北京工業大學建筑工程學院,北京 100124

0 引言

大震的破裂尺度一般為幾十千米到幾百千米,其震源不是一個點而是一個面。如四川汶川8.0級地震破裂長度約200km,其烈度等震線呈狹長橢圓狀,具體表現為沿地震破裂方向的烈度衰減緩慢,等震線的長軸尺度遠遠超過了點源烈度衰減關系所預測的尺度；垂直于地震破裂方向的烈度衰減較快,等震線的短軸尺度與烈度衰減關系所預測的尺度基本一致(帥向華等,2020)。因此,基于線源的烈度衰減模型能夠更準確地評估烈度分布,對地震破裂方向和破裂方式的快速準確判定十分重要。在目前常用的地震烈度評估模型中，地震破裂方向的判定方法相對簡單,即將震中位置與地震活動斷裂帶數據進行疊加分析,計算震中到周圍斷層的最小距離,并選擇距離最近的斷層走向作為地震破裂方向(地震影響場長軸方向),且破裂方式均為對稱雙側破裂。這種方法的優勢在于時效性高,可快速獲得評估結果。但該方法在實際工作中存在一定的局限性：一是地震活動斷裂帶數據不完備,盡管經過地質專家多年的調查,已經積累了大量的斷裂帶數據,但仍然不排除有未發現的隱伏斷裂,如2014年8月3日云南魯甸6.5級地震的發震斷層實際為NW向的包谷垴-小河斷裂,而當時調查記錄的距離震中最近的為NE向的昭通-魯甸斷裂；二是地震有時會發生在2條或多條斷裂交匯的地方,此時震中到幾條斷裂的距離十分接近,距離最近的斷裂不一定為發震斷層；三是大震的破裂方向往往是單側破裂或者不對稱的雙側破裂,如四川汶川8.0級地震以NE向為主、四川蘆山7.0級地震以SW向為主。因此,迫切需要一種快速且可靠性較高的地震破裂方向的判定方法,提高烈度速報的效率和精度,為地震應急救援及災害損失快速評估提供更加快速、準確的決策支撐和服務保障。

地震斷層破裂傳播的方向對地震波(體波和面波)的影響可以用有限移動源模式下的地震波譜理論來解釋。有限移動源模式是指地震斷層面上的各點不可能同時發生錯動,通常為先從某一點(或部分)開始破裂,然后以有限速度向斷層的其他部分傳播。因而在地震記錄中必定帶有某些震源的相關信息,如震源破裂速度等其他參數的信息。相關研究結果表明,震源的有限性和有限的破裂速度對地震波的影響主要體現在2個方面：一方面,地震波的振幅在方位上受到有效的調制,也就是所謂的“方位調制”作用,即地震波振幅在沿著破裂的傳播方向上加強,而在相反方向上減弱；另一方面,地震波卓越周期隨方位有類似于多普勒效應的變化(胡進軍等,2011)。因此,根據體波振幅和卓越周期在不同方位上的變化即可反推出地震斷層的破裂方向。

趨勢面分析是一種多元統計分析方法,能有效地分析某一種屬性數據在空間上的分布規律與變化趨勢。該方法最初應用于地理科學中,利用數學趨勢面模擬地理系統要素在空間上的分布規律及變化趨勢(王江萍等,2009),后來其被廣泛地應用于模擬環境、資源、疾病、人口、農業、經濟等要素在空間上的分布規律(Hota,2014)。

本文根據有限移動源理論,提出一種利用趨勢面分析和震中一定范圍內臺站初始P波地震動參數判斷地震斷層破裂方向的方法,并基于斷層距模型和最小二乘法統計回歸確定地震斷層破裂方式。通過實際震例,對比分析根據本文所提方法評估的地震斷層破裂方向與實際地震影響場長軸方向的吻合度,評估采用該方法進行地震烈度速報的可靠性,為利用地震監測臺站實時數據進行地震斷層破裂方向判定提供一種新思路,為地震烈度速報結果的動態修正提供參考。

1 地震斷層破裂方向研究現狀

確定地震斷層破裂的方向主要有現場地質探查、大地測量觀測以及基于地震學資料的反演分析3大類方法(何驍慧等,2015)。現場地質勘察方法為地質學專家在震區考察識別地震破裂帶、測定斷層位置及走向,是判定斷層面直接有效的方法,其得到的結果最精細,但是所需時間較長,不適合地震烈度快速評估需求。大地測量觀測方法采用合成孔徑雷達干涉測量技術(InSAR)測定地表形變,進而測定斷層破裂方向,該方法僅在地表裸露、植被不茂密的地區較為有效,且因衛星重復軌道數據往往無法快速獲取,不能用來快速(幾個小時或者當天)測定破裂斷層。地震學資料反演的方法主要是根據余震的空間分布、強震動臺觀測的地震動信息(峰值速度PGV、峰值加速度PGA等)或多個測震臺上記錄到的波形信息推斷地震斷層破裂的方向,從而確定發震斷層。其中余震分布法一般需要震后2h內的余震精定位信息,時效性較差,而我國強震動臺站數量較少的現狀會影響破裂方向評估的準確度,因此,基于測震地震波形判定地震斷層破裂方向亟待研究。

經驗格林函數法是研究地震破裂方向性的常用方法之一,該方法使用小震作為經驗格林函數校正地震波在三維地球介質中的傳播效應,從而對震源的特性進行更準確的估計(Mozaffari et al，1998)。Tan等(2010)在此基礎上,利用不同方位臺站記錄P波的形狀和振幅信息,發展了正演模擬技術,從而測定破裂方向。對于大震(M>7),簡單的破裂方向參數已不足以描述破裂的復雜性,一般使用有限斷層反演方法研究其破裂過程的細節(王衛民等,2008、2013)。另外,通過地震質心位置與破裂起始位置之間的差異,也可測定主震的破裂方向以及發震斷層的破裂方向性(秦劉冰等,2014)。上述研究中基于震源參數反演方法所需處理的資料在主震發生后很難快速得到。馮蔚等(2015)利用測震臺站波形數據,計算LQT坐標系下得到的S波振幅值與單力偶源S波輻射圖案匹配,判斷魯甸地震的發震斷層方向與調查結果吻合,但是對于大震,近場測震臺站S波記錄受限幅影響較大(馮蔚等,2015)。在地震破裂傳播方向上,地震波的周期變短；在其相反方向上,地震波的周期變長。對一些較大地震,在寬頻數字化記錄中能夠比較清楚地觀測到地震多普勒效應(Douglas et al，1981),也有學者研究了利用P波的多普勒效應來判斷地震破裂的方向(周云好等,2002)。

一次大震的斷層破裂過程非常復雜,會發生多次破裂,主震的首次破裂一般釋放的能量最大,對震源參數的確定起主導作用。隨著地震預警技術研究和發展,許多學者發現初始P波能夠較好地反映主震的震源特征。研究表明臺站記錄的P波前3s最大位移與相應臺站的PGV具有較好的相關性(Wu et al，2005)。Colombelli等(2015)在此基礎上利用日本 KiK-net 臺網記錄的本州島海域9.0級地震及其余震波形記錄的P波初始3s幅值來估算臺站的PGV,再由PGV與修正的麥卡利烈度IMM的關系來計算儀器烈度,最后通過內插來評估地震烈度分布。

2 數據準備

本文選取了2008—2014年8個破壞性地震作為實驗震例(表1),其中MS8.0地震1個,MS7.0地震1個,MS6.0～7.0地震6個,震源深度從最淺5km到最深的20km,分布在四川、云南、甘肅和新疆4個多震省份。

表1 實驗震例

實驗數據選取了以上震例中距震中300km范圍內的測震臺站和強震動臺站連續波形記錄,并對波形記錄進行預處理(圖1),經過基線校正和去線性化處理,對波形記錄的誤差進行修正,采用了0.075～10Hz的中通Butterworth濾波器進行濾波,然后利用STA/LTA方法進行P波震相的自動拾取,以1s為間隔截取2～10s的波形記錄,分別對記錄進行積分、水平向合成、三分量合成，即

(1)

其中,v(ti)為三分量合成速度,v(ti)E-W為EW向速度,v(ti)N-S為SN向速度,v(ti)U-D為垂直向速度。

圖1 數據預處理流程

之后取其極值得到不同時段垂直向、水平向和三分量合成的峰值速度(Peak Velocity,以下簡稱PV)、峰值位移(Peak Displacement,以下簡稱PD)以及對不同時段的記錄進行FFT變換得到卓越周期(Dominant Period,以下簡稱DP)。最后將計算數據寫入表格,形成數據文件,共得到240組臺站記錄。

3 方法

3.1 破裂方向分析方法

破裂方向分析采用趨勢面分析法。趨勢面分析法是研究一定區域內空間數據發展或演變趨勢的一種方法,該方法的實質是數據擬合(查文婷等,2015)。將局部區域內某一事件的觀測值分為趨勢值和殘差值,趨勢值反應了該事件空間要素的宏觀分布規律,可以排除一些偶然因素引起的變異影響。自變量為觀測點的經緯度坐標或平面坐標,因變量為觀測值,用最小二乘法獲得殘差平方和最小估計值,以反映觀測值變化的整體趨勢。本文采用的是多項式擬合,即利用線性模型對某一變量的觀察值進行多元回歸擬合,以產生該區域的趨勢面。趨勢面分析的多項式函數表達為

Z=A1+A2x+A3y

(2)

其中,x、y為觀測點的平面橫、縱坐標,Z為觀測值。

圖2 斷層破裂方式示意圖

3.2 破裂方式分析方法

對于本文研究的6級以上強烈地震,震源可以近似為一條地震沿斷層破裂的線段。破裂斷層一般分布于震中的一側或兩側,取決于地震破裂方式。如圖2 所示,設破裂斷層一端到震中的長度ao占斷層總長度ab的比例為per。在震中位置、破裂方向和破裂長度確定的情況下,變換per比例即可得到所有可能破裂方式的斷層分布。

觀測點的震中距應取其到破裂斷層線段最近的距離D(圖3)。當觀測點位于破裂斷層兩端之外時(如圖3中a點、d點),取其到斷層兩個端點的距離；當觀測點位于斷層兩個端點之間時(如圖3中b點、c點),取其到斷層的垂直距離。

圖3 斷層距示意圖

本文采用常用的地震動衰減關系模型為(Kanai,1961)

lgY=C1+C2M+C3lg(R+R0)

(3)

其中,Y為地震動參數,本文選取的是三分量合成的峰值速度(PGV)；M為震級；R為震中距,本文采用斷層距表示；R0控制近場Y為一個有限值,由于近場地震動的飽和特性,可以假定震中一定范圍內的地震動飽和,所以本文設定R0為10km,即距震中10km內地震動飽和。

將觀測點的斷層距與PGV代入式(3),利用最小二乘法進行回歸分析,每種斷層破裂方式可以擬合出一個結果,這些結果對應的標準差越接近0,說明擬合的誤差越小,也就是說最小標準差對應的破裂方式與觀測值的匹配度最高。本文設定per比例為0到100%之間,變換步長為5%,通過迭代回歸分析21種破裂方式的衰減模型,以擬合標準差最小的per比例作為斷層破裂方式評估結果。

4 結果及分析

4.1 破裂方向分析結果

本文利用ArcGIS軟件對一定分布范圍內臺站的初始P波參數進行插值分析以獲得地震斷層破裂的趨勢方向,評估流程如圖4 所示。其中,在選取插值范圍時,應保證震中的各個方向上均選取一定數量臺站,且臺站分布較均勻。

圖4 基于初始P波的地震斷層破裂方向評估流程

破裂方向分析以四川蘆山7.0級、云南景谷6.6級和四川汶川8.0級地震為例,采用了距震中200km范圍內的測震臺站波形記錄,選取了P波初始3s的最大位移和初始5s的卓越周期,再利用反向距離加權和趨勢面插值方法分析地震影響場長軸方向,研究結果表明：蘆山地震P波幅值和卓越周期的反向距離插值IDW結果呈現出一定的破裂方向性,P波初始3s峰值沿發震斷層NE方向較低,SW方向較高,5s內的卓越周期沿發震斷層NE方向較高,SW方向較低,表現出蘆山地震不對稱破裂方式的特點,2種參數的趨勢面插值結果得到的趨勢方向均與實際調查烈度影響場的長軸方向基本一致(圖5、圖6)；景谷地震P波幅值的插值結果與該地震較對稱的破裂方式特點一致,即P波初始3s峰值沿發震斷層的2個方向上沒有明顯差別,但是5s內的卓越周期的插值結果呈現出延發震斷層SE向的破裂趨勢,2種參數的趨勢面插值結果得到的趨勢方向均與實際調查烈度影響場的長軸方向基本一致(圖7、圖8)；汶川地震P波幅值的IDW插值結果具有比較明顯的破裂方向性,即P波初始3s幅值沿發震斷層NE方向較高,SW方向較低,而且有顯著的不對稱破裂特點,趨勢面插值得到的趨勢方向與實際調查烈度影響場的長軸方向一致(圖9),汶川地震5s內的卓越周期不具有顯著方向性。

圖5 蘆山地震P波幅值反向距離插值IDW和趨勢面插值分析紅色影響場為實際地震烈度分布；底部為插值結果；黑色箭頭為由趨勢面插值得到的斷層破裂方向

圖6 蘆山地震P波卓越周期反向距離插值IDW和趨勢面插值分析

圖7 景谷地震P波幅值反向距離插值IDW和趨勢面插值分析

圖8 景谷地震P波卓越周期反向距離插值IDW和趨勢面插值分析

圖9 汶川地震P波幅值反向距離插值IDW和趨勢面插值分析

在實驗過程中發現,初始時刻和P波參數的選擇對分析結果有較大影響。因此,本文利用實驗震例數據,分析不同初始時段PD、PV和DP進行趨勢面分析的可靠性。具體方法如下,分別利用初始2～10s的PD、PV和DP進行趨勢面分析,將分析得到的趨勢方向與影響場方向進行比較,完全一致計3分,基本一致計2分,不一致計0分,以打分之和代表某一時段對應參數預測的準確率。當趨勢方向與影響場方向一致時可判斷為正確,正確震例占所有震例的比例即為正確率。

分析結果表明，PV的準確率總體最高,其次是PD,DP最低,即參數的可靠性 PV>PD>DP(圖10)。PD在4s和5s時正確率最高(83%)；PV在4s、9s和10s時正確率最高,均達到了100%；DP在7s時正確率最高,也達到了100%。所有參數的總體正確率達到70%以上；DP在6s前的正確率較低,原因可能是近場臺站密度較低,6s內的卓越周期在近場差異性不夠顯著。

圖10 不同初始時段PD、PV和DP趨勢分析的準確性

圖11 不同斷層破裂方式下衰減模型擬合標準差

4.2 破裂方式分析結果

破裂方式分析以汶川8.0級地震為例,利用震中300km范圍內臺站初始P波3s位移峰值的趨勢面分析得出斷層破裂走向為NE∠50°,破裂長度計算采用Wells等(1994)關于破裂長度L與震級M之間的統計關系

lgL=0.69M-3.22

(4)

汶川地震斷層破裂長度由式(4)計算約為200km。在確定方向和破裂長度情況下,變換per比例可得到21個可能的斷層分布。根據線源斷層距衰減模型(式(3)),對觀測臺站lgPGV擬合結果的標準差如圖11 所示,最小標準差對應的per為0,即評估破裂方式為NE∠50°的單側破裂。

將本文利用臺站觀測數據和模型評估計算得到的汶川地震破裂斷層分布繪制到地圖上,與汶川地震現場調查烈度分布進行對比(圖12)。評估破裂斷層從震中汶川縣沿NE向經北川縣延伸至青川縣,評估的斷層破裂方向和破裂方式與現場實際調查結果基本一致。

圖12 汶川地震評估破裂斷層分布與烈度分布對比

5 結論

在地震烈度速報中,對地震斷層破裂方向效應的考慮十分重要,如何快速評定地震斷層破裂方向和破裂方式成為相關研究的關注點。本文提出了一種利用趨勢面分析法和臺站觀測數據快速評估地震破裂方向和破裂方式的方法,并基于近年來我國西部發生的一些強震實例分析了該方法的可靠性和實用性。主要研究結果如下：

(1)提出的基于趨勢面分析和觀測臺站初始P波的地震斷層破裂方向快速判斷方法具有較高的可行性,對于不同區域、震源深度及破裂類型的中強地震均具有較好的適應性。但方法的準確率受地震動參數(速度峰值、位移峰值和卓越周期)的選取、P波初始時段選取及臺站選擇等影響,對于震中附近臺站非常稀疏的地震無法適用。

(2)在初始P波信息處理中,PD的最優初始時段選擇為3～5s和9～10s,PV的最優初始時段選擇為4～5s和9～10s,DP的最優初始時段選擇為6s以后；該方法使用的地震動參數可靠性為PV>PD>DP,但總體正確率較高,達到了70%以上。

(3)提出利用臺站觀測數據和最小二乘法擬合不同破裂方式的PGV衰減關系,以擬合結果最優的破裂方式來確定地震破裂斷層分布的方法具有可行性,在震后能夠利用實時獲取的臺站數據快速評估破裂斷層分布,能夠較大地提高地震破裂斷層評估的效率。

(4)實驗震例中7級以上地震較少,因此該方法對于大震破裂方向判定的可靠性還有待檢驗。另外,本文選取的測震臺站分布密度較低,需要進一步利用密度較高的地震監測臺網來進行記錄,以檢驗該方法的適用性。

(5)本文快速判定地震破裂方向的整個處理過程包括波形數據預處理、P波震相拾取、一定時窗內P波峰值速度或峰值位移參數提取,然后通過趨勢面分析法評估地震破裂方向,再利用最小二乘法分析地震破裂方式,最終得到破裂斷層分布。整個過程可以通過編程自動計算完成,能夠構成一個系統，處理計算時間不超過1min,震后10min內可以產出判定結果。