地震初至震相自動識別方法研究

李恩來，王承偉，安祥宇

（遼寧省地震局，遼寧 沈陽110034）

0 引言

在傳統(tǒng)地震學(xué)中，地震震相的拾取一般采用人工分析的方式，雖然精度較高，但效率低。隨著地震臺網(wǎng)的不斷發(fā)展和臺站密度的不斷增加，記錄的數(shù)據(jù)量也不斷增加，特別是實時地震學(xué)的發(fā)展，簡單的人工拾取震相已不能滿足要求，對震相的自動精確拾取成為地震學(xué)家日益關(guān)注的問題，同時也發(fā)展了多種地震震相識別方法。

目前，絕大多數(shù)震相自動識別方法是提取信號和噪聲的不同特征來作為提取震相到時的判據(jù)，比如用波形的幅值變化、頻率變化、地震波形的相關(guān)性和運動學(xué)特征等來判別地震震相類型。歸納起來，常用的判別方法有：能量變化方法、偏振分析方法、自回歸方法（ARAIC）以及神經(jīng)網(wǎng)絡(luò)方法等。能量方法中最常用的為短長時間平均（STA/LTA）方法［1，3-6］，短長時間平均方法反應(yīng)了幅值的瞬時變化，短時平均代表能量的瞬時變化，長時平均代表能量的穩(wěn)態(tài)變化。震相到來時，質(zhì)點運動的偏振方向會發(fā)生改變，地震波的偏振分析由于提取出了地震波的偏振特征經(jīng)常被用在震相撿拾和判別震相類型［2，7-9］。自回歸方法是假定可以把地震波劃分為局部平穩(wěn)段，并且在觸發(fā)點前后是不同的平穩(wěn)過程［10］。根據(jù)這個假定AR-AIC通常被用作P 波和S波的自動撿拾中［11-12］。不同于AR-AIC算法，Maeda（1985）提出了一種不用求AR系數(shù)的時域AIC方法［14］，Zhang 等（2003）應(yīng)用基于小波變換的AIC方法，直接對小波系數(shù)應(yīng)用AIC方法［13］。Baiand Kennett（2000）指出目前存在的震相自動識別方法中沒有哪一種可以撿拾出所有的震相到時，各種方法都有其局限性，特別是信噪比的影響［15］。比如，在噪聲較高或者事件淹沒在前面事件的尾波中，最常用的能量方法就不能有效識別震相，或者拾取的震相誤差較大。本文使用最常用的短長時窗能量比（STA/LTA）加自回歸方法（AIC）的方式，采用兩步法進(jìn)行初至震相識別，并分析特征函數(shù)對能量變化的敏感度，得到了比較好的結(jié)果。

1 震相識別的基本原理

1.1 長短時窗能量比（STA/LTA）方法

我們給出一定長度的滑動長時窗，再在這個長窗中取一個固定時長的短時窗，兩窗口的起點或者終點重合，然后求取短時間窗內(nèi)信號的平均值（STA）和長時間窗內(nèi)信號的平均值的比值，用這個比值來反映信號幅值或者能量的變化。其中，STA主要反映的是信號瞬時變化的平均值，LTA主要反映的是背景噪聲變化的平均值。在地震信號到達(dá)時，STA的變化要比LTA的變化大的多，STA/LTA值會有明顯的增加。我們設(shè)定一個閾值，當(dāng)STA/LTA值大于該閾值時，即認(rèn)為有震相到達(dá)， 從而達(dá)到自動識別震相的目的。

其中，i為采樣時刻；ns為短時窗長度；nl是長時窗長度；λ為設(shè)定的觸發(fā)閾值；CF（j）為在j時刻的關(guān)于信號的特征函數(shù)值，表示信號的振幅、能量或其變化。

1.2 自回歸方法（AIC）

七十年代，日本學(xué)者Akaike提出一個基本信息量的定階準(zhǔn)則，即AIC準(zhǔn)則。地震震相到時識別的自回歸（AR）技術(shù)是假設(shè)震相到時前后的地震記錄是兩個不同的穩(wěn)態(tài)過程。不同于AR-AIC方法，Maeda（Maeda N.1985）建議可由地震波形數(shù)據(jù)在時域直接計算AIC函數(shù)，而不需要求出AR 系數(shù)，對地震記錄x（i）（i=1，2，…， L）來說，AIC檢測器定義為：

AIC（k）=k·log｛var（x［1，k］）｝+（L-k-1）·log｛var（x［k+1，L］）｝ （2）

其中，k的范圍為地震圖某窗口內(nèi)所有的采樣點；var 表示方差。震相到時對應(yīng)于AIC函數(shù)的最小值。

本文中，我們不在地震圖上直接應(yīng)用AIC準(zhǔn)則，而是在用STA/LTA方法粗略撿拾到P 波位置后，在固定窗內(nèi)用地震記錄的垂直向特征函數(shù)來進(jìn)行精確識別震相，即在粗略撿拾點向前和向后截取一定長度的時間窗，在窗內(nèi)應(yīng)用AIC準(zhǔn)則，窗內(nèi)AIC的最小值，即認(rèn)為是精確到時點。因為在對P 波進(jìn)行STA/LTA 粗略撿拾時，其得到的觸發(fā)點一般滯后于真實觸發(fā)點，所以在本文中我們對P 波撿拾后，在拾取到的P 波位置處向前推100個數(shù)據(jù)點，向后推10個數(shù)據(jù)點，即對采樣頻率為100Hz的地震記錄，前推時間為1s，后推時間為0.1s。

2 特征函數(shù)分析

在STA/LTA方法中，特征函數(shù)的選擇對震相識別至關(guān)重要，選擇特征函數(shù)應(yīng)遵循如下原則：能靈敏地反映地震信號到達(dá)時其振幅和頻率的變化特征，并能增強(qiáng)這些變化。目前常用的特征函數(shù)有如下五個：

本文分別測試了上述五個特征函數(shù)對頻率和振幅變化的響應(yīng)，并作了分析。

2.1 特征函數(shù)對振幅變化的響應(yīng)

為了測試上述特征函數(shù)對振幅變化的響應(yīng)，生成一個余弦信號，頻率為40Hz，在信號的先后兩段，振幅相差2 倍。

圖1 特征函數(shù)對振幅變化的響應(yīng)（a 為原始波形，b-f分別為各特征函數(shù)對信號的響應(yīng)）Fig.1Responseof eigenfunctionstoamplitudechanges

由測試結(jié)果可見，五個特征函數(shù)對振幅變化都有響應(yīng)，但響應(yīng)的結(jié)果明顯不同，CF1、CF3、CF5對振幅的響應(yīng)幅度較小，且響應(yīng)處兩邊毛刺較大。CF2對振幅的響應(yīng)較大，但響應(yīng)處兩邊同樣存在較大的毛刺。CF4對振幅的響應(yīng)最為明顯，響應(yīng)處振幅較大，光滑，兩邊無毛刺。

2.2 特征函數(shù)對頻率變化的響應(yīng)

為了測試特征函數(shù)（公式（3））對頻率變化的響應(yīng)，生成一個余弦信號，在信號的先后兩段，頻率相差2 倍。

圖2 特征函數(shù)對振幅變化的影響（a 為原始波形，b-f分別為各特征函數(shù)對改信號的響應(yīng)）Fig.2 Responseof eigenfunctions tofrequencychanges

從測試結(jié)果可見，5個特征函數(shù)隨頻率變化都有不同程度的響應(yīng)，其中，CF1、CF2對頻率變化的響應(yīng)不明顯。CF3和CF5對頻率的變化有響應(yīng)，但響應(yīng)處兩側(cè)的毛刺較多。CF4對頻率的響應(yīng)最強(qiáng)，響應(yīng)光滑且兩邊無毛刺現(xiàn)象。

總之，經(jīng)過上述測試，發(fā)現(xiàn)CF4（第四個特征函數(shù)）對振幅和頻率的變化均表現(xiàn)出非常好的敏感性。因此，在對P 波初至震相進(jìn)行拾取時，選擇CF4做為特征函數(shù)為最佳。

3 資料的選取

本研究共使用遼寧地區(qū)2010—2015年記錄到的地震521個，因震級較小或者震中距較大都會造成記錄的信號較弱，信噪比較低，影響震相的拾取精度。故選擇臺站記錄的垂直分向（UD向），臺站的震中距＜150 km，地震震級ML>2.0 級。所有事件的震中分布如下圖所示。

臺網(wǎng)保存的事件數(shù)據(jù)格式有兩種，分別為SEED格式和EVT 格式。本文使用SEED格式，使用rdseed軟件解壓*.seed 波形事件文件，去除儀器響應(yīng)，去均值、去傾斜校正，并進(jìn)行帶通濾波，濾波頻率0.1～20Hz。

4 結(jié)果分析

本次計算共使用地震記錄4662 條，其中檢測到的地震記錄4597個，未檢測到地震記錄65個，檢測到的初至震相比例達(dá)到96.5%。通過與觀測報告中人工拾取的初至震相到時對比，對拾取的震相到時進(jìn)行分析。由圖5可見，在檢測到的所有記錄中，大部分記錄拾取的震相誤差都0.5秒以內(nèi)，達(dá)到3193條，誤差1秒到2 秒的記錄一秒的記錄363個，其余的為誤差1秒以上記錄。從檢測的結(jié)果總可以看到，該方法可以很好的對地震記錄的初至震相進(jìn)行拾取。

在所有使用的全部記錄中，有65個未被檢測到，還有個別拾取到時誤差較大的記錄，分析原因，主要為信噪比較低，干擾較大，或者有連續(xù)地震事件發(fā)生，且兩地震的時間間隔很小，第二個事件的震相到時拾取誤差較大，不能應(yīng)用本文方法，需進(jìn)一步分析研究。除此以外，通過分析發(fā)現(xiàn)，還有一些觀測報告中震相拾取時間的明顯錯誤。

5 結(jié)論和討論

本文通過利用STA/LTA加AIC方法這種兩步拾取地震波初至震相的方式，對遼寧臺網(wǎng)近年來記錄的地震初至震相進(jìn)行了分析，認(rèn)為本方法可以有效的自動識別較近震中距（150km）、信噪比較好的記錄的初至震相。在分析的所有記錄中，誤差小于0.5秒的記錄占到70%以上。同時對該方法中的特征函數(shù)（公式3）進(jìn)行了分析，認(rèn)為CF4無論對頻率變化，還是對振幅變化，均表現(xiàn)出清晰的響應(yīng)，所以在選擇使用該方法時，應(yīng)優(yōu)先選擇CF4作為特征函數(shù)。

該方法原理簡單，操作便捷，能較好的識別初至震相，但也有其使用條件。在震中距較遠(yuǎn)、信噪比較小的情況下，可能會導(dǎo)致識別震相到時的誤差較大，或者不能有效的識別震相；在有連續(xù)地震事件發(fā)生，特別是震群發(fā)生的情況下，由于事件波形間隔時間較短，或者后續(xù)事件淹沒在前一事件當(dāng)中時，該方法不適用。