何先龍, 佘天莉, 高峰
中國地震局工程力學(xué)研究所地震工程與工程振動重點實驗室, 哈爾濱 150080
?
一種地震P波和S波初至?xí)r間自動拾取的新方法
何先龍, 佘天莉*, 高峰
中國地震局工程力學(xué)研究所地震工程與工程振動重點實驗室, 哈爾濱150080
摘要地震P波、S波初至?xí)r間的拾取是地震波分析的一項基礎(chǔ)性工作.本文提出了一種新的地震波初至?xí)r間自動拾取的方法:首先,把地震波的三分量時程曲線變換為一組空間向的能量變化率時程曲線;然后對能量變化率時程曲線進(jìn)行STA/LTA(Short Time Average/Long Time Average,短時間的均值/長時間的均值)處理,拾取地震P波和S波的大致初至?xí)r間;最后提出采用一種二次方自回歸模型對初至附近的能量變化率曲線進(jìn)行二次方自回歸處理,精確拾取出P波和S波的初至?xí)r間.本文采用了10組蘆山地震的記錄數(shù)據(jù)和150組汶川地震的記錄數(shù)據(jù)對此方法的可靠性進(jìn)行了檢驗.以人工拾取結(jié)果為參考,此方法具有很高的準(zhǔn)確率和穩(wěn)定性,同時,相比于常用的STA/LTA方法和AIC(Akaike Information Criterion,Akaike信息準(zhǔn)則)方法,此方法在計算時間效率方面稍微遜色,但是對S波初至?xí)r間的拾取精度和可靠性更高.此方法豐富了地震P波、S波初至?xí)r間的自動拾取方法.
關(guān)鍵詞地震波初至;能量變化率;STA/LTA方法和AIC方法;二次方自回歸模型
1引言
地震P波、S波初至?xí)r間是精確計算震源位置的必要參數(shù),因此地震P波、S波初至?xí)r間的拾取是地震數(shù)據(jù)處理的一項必要工作(周寶峰,2012).在地震發(fā)生后,大量地震記錄需要快速處理來實現(xiàn)P波、S波初至?xí)r間的精確和快速拾取,特別是在地震預(yù)警方面,地震P波、S波初至?xí)r間的拾取速度尤其重要.因此,基于傳統(tǒng)的人工拾取P波、S波初至?xí)r間的方法和手段已經(jīng)難以滿足社會發(fā)展需求,促使相關(guān)科研人員研究地震P波、S波初至?xí)r間的自動拾取技術(shù)(劉勁松等,2013;葉根喜等,2008;馬強(qiáng),2008).
地震P波、S波初至?xí)r間的自動拾取算法一直在不斷改進(jìn),迄今為止已有多種較成熟的方法,如相關(guān)法(Bai and Kennett,2000;Allen,1978)、能量比法(Allen,1982;Maeda,1985;王繼等,2006;Saragiotis et al.,2002)、最大振幅法、分形維法及神經(jīng)網(wǎng)絡(luò)法(劉伊克等,2001;Paige and Saunders,1982;Boschetti et al.,1996;Chang et al.,1999;Cao and Greenhalgh,1993)等.目前最常用的方法主要有以下幾種(馬強(qiáng),2008):
(1) 長短時間平均(STA/LTA)方法:Allen(1978,1982)提出的STA/LTA方法,是目前廣泛使用的一種地震P波、S波初至?xí)r間的自動拾取方法,其主要原理是根據(jù)地震波形特征函數(shù)的長短時均比值等特征拾取初至.
(2) AIC方法:該方法的基本原理是求取地震信號AIC函數(shù)局部最小值的位置(Maeda,1985;王繼等,2006).
(3) 基于高階統(tǒng)計量偏斜度和峰度的PAI-S/K方法:Saragiotis等(2002)提出了基于地震波形峰度和偏斜度的拾取初至方法,稱為PAI-S/K方法.
但是由于不同地震記錄的特點存在較大差異,目前存在的地震P波、S波初至?xí)r間的自動拾取方法沒有哪一種可以拾取出所有地震記錄的P波、S波初至?xí)r間,即單種方法不能保證所有地震記錄的P波、S波初至?xí)r間的100%的準(zhǔn)確拾取(葉根喜等,2008;馬強(qiáng),2008).因此,為了提高地震P波、S波初至?xí)r間的自動拾取結(jié)果的可靠性,馬強(qiáng)(2008)綜合多種自動拾取方法的特點,提出了地震P波、S波初至?xí)r間的多步自動拾取法.
相比于地震P波初至?xí)r間的拾取,地震S波初至?xí)r間的拾取更難,也是地震初至?xí)r間自動拾取技術(shù)的主要研究內(nèi)容和研究難點(劉勁松等,2013;葉根喜等,2008).P波初至波先于S波初至波到達(dá),但是由于S波初至波常在P波未完時就到達(dá),使S波初至波疊加在P波中而難以被區(qū)分.這也是基于現(xiàn)有的地震P波、S波初至?xí)r間的自動拾取方法難以精確拾取S波初至?xí)r間的主要原因(Bai and Kennett,2000;Allen,1978).因此,如何從P波、S波混疊部分識別S波初至波的到達(dá)時刻,是地震波初至?xí)r間的自動拾取技術(shù)的研究熱點(劉勁松等,2013;葉根喜等,2008;馬強(qiáng),2008).
本文提出一種新的自動拾取P波和S波初至?xí)r間的方法,首先把地震波的三分量時程曲線變換為一組空間向的能量變化率時程曲線,然后對能量變化率時程曲線進(jìn)行STA/LTA處理,自動拾取地震波的P波和S波大致初至?xí)r間,最后建立了一種二次方自回歸模型并基于此模型對地震P波、S波初至波附近的能量變化率曲線進(jìn)行二次方自回歸處理,獲得P波、S波初至波的二次方自回歸時程曲線,最后自動拾取P波和S波的初至波的二次方自回歸時程曲線的最大值所對應(yīng)的時間,所拾取的兩個時間值即為P波、S波的精確初至?xí)r間.
2基于信號能量變化率的二次方自回歸模型拾取初至原理介紹
對某個地震信號x(t)進(jìn)行等間隔Δt離散化后得序列{xi},假設(shè)序列長度為N,可定義此信號的震動瞬態(tài)幅值變化差為
(1)
則{pi}為一個等間隔時間為Δt的離散序列,它在物理意義上表示物體震動幅值的變化量值.
同時定義
(2)
假設(shè)由三分量傳感器所記錄到的某個地震信號的兩個水平向和豎向震動可分別表示為{xi},{yi},{zi},0≤i Qi=(xi+1-xi)2/Δt+(yi+1-yi)2/Δt +(zi+1-zi)2/Δt, 0≤i≤N-2 (3) 則序列{Qi}描述了地震信號在三維空間內(nèi)震動能量變化率的時程特征. 當(dāng)?shù)卣鸬竭_(dá)時,地震記錄波形在三維空間內(nèi)的振幅會突然增大;同樣在地震到達(dá)時刻,地震記錄波形在三維空間內(nèi)的震動能量變化率也會突然增大.因此可根據(jù)序列{Qi}的幅值時程曲線來分析出地震到達(dá)時刻.為了方便描述,本文定義序列{Qi}為地震記錄波形的震動能量變化率譜. (4)N4為STA/LTA的短尺度,N3為STA/LTA的長尺度.由公式(4)可見序列{ei}相對序列{Qi}的起始時間延遲了(N3-N4+1)Δt.序列{ei}的長度為N2=N-2-N3. 令 (5) en所在位置對應(yīng)地震記錄波上的時間為tn=(N3-N4+1+n)Δt.tn為地震波的大致初至?xí)r間,可以此值為中心確定一個包含了真實初至?xí)r間的狹窄時間段為[tn-sΔt,tn+sΔt],則此時間段對應(yīng)序列{Qi}的子序列為[QN3-N4+1+n-s,QN3-N4+1+n+s].如果定義此子序列為{qi},則可得 (6) 假設(shè)地震信號的準(zhǔn)確到達(dá)時刻為T1=T0+T2,T0為開始記錄波形的絕對時間,T2為地震信號準(zhǔn)確到達(dá)時間與開始記錄波形時刻的相對時間,則可得T2∈[tn-sΔt,tn+sΔt]. 為了進(jìn)一步精確自動拾取到達(dá)時刻,本文提出第二步拾取法:對序列{qi}基于二次方自回歸模型法,精確地自動拾取初至. 此模型的理想二次方程假設(shè)為: ci=w+gi+fi2, 0≤i<2s. (7) (8) 然后基于最小二乘法求解系數(shù)w,g,f.令 (9) (10) 由式(10)可得: (11) 由式(11)可解得系數(shù)w,g,f.然后基于式(7)求取序列{ci}最大值對應(yīng)的i,即為地震信號精確初至?xí)r間,得到地震信號到達(dá)的絕對精確時刻為: (12) 3實驗仿真分析 首先基于青川某斜坡地震監(jiān)測站的一組實測場地脈動和某個水平向余震信號進(jìn)行仿真實驗來驗證本文提出的地震初至自動拾取新方法的可靠性.余震信號為地震P波已經(jīng)到達(dá)后截取的某段地震波,采樣頻率為256 Hz.為了方便分析,本文人為地在場地脈動信號的3.4 s時刻疊加余震信號,得到信號S3.因此,由信號S3可以看出地震到達(dá)時刻為3.4 s,即P波到達(dá)時刻為3.4 s,對S2進(jìn)行人工拾取S波到達(dá)時刻為2.6006 s,則S3里的S波到達(dá)時刻實際值為6.006 s. 由圖2可得:對S3的絕對幅值按照長的平均取1 s、短的平均取0.05 s的參數(shù)進(jìn)行STA/LTA分析可得到P波的到達(dá)時間大概在3.4123 s附近,而S波的到達(dá)時間比較難確定,可初步確定在4.9436、5.3195、6.0456、6.2547這四個時刻的某個時刻的附近.對S3的震動能量變化率按照相同參數(shù)進(jìn)行STA/LTA分析,可得到P波的到達(dá)時間大概在3.4032 s附近,而S波的到達(dá)時間可確定在6.0456 s附近.因此,基于能量變化率的STA/LTA處理和基于絕對幅值的STA/LTA處理都能比較準(zhǔn)確地分析出P波的到達(dá)時間,但是基于能量變化率的STA/LTA處理更能準(zhǔn)確分析出S波的到達(dá)時刻. 圖1 模擬地震信號(a表示振動加速度)Fig.1 Simulated seismic signals 圖2 S3的能量變化率和長短時間平均處理Fig.2 Transforming S3 into its energy gradients and STA/LTA curves 為了更準(zhǔn)確地確定P波和S波的到達(dá)時間,對S3的能量變化率的時程波形分別以3.4032、6.0456 s所對應(yīng)的數(shù)據(jù)為中心,截取長度為200個的兩個序列,采用本文建立的二次方自回歸模型法來進(jìn)一步確定P波和S波的初至,即初步確定P波初至在(3.0116,3.7948)s內(nèi),S波初至在(5.6450,6.4372)s內(nèi).兩個序列的二次方自回歸模型曲線如圖3所示.讀取兩條回歸曲線最大幅值對應(yīng)的時間分別為0.3892 s、0.3720 s.因此,P波到達(dá)時刻為3.4008 s,S波到達(dá)時刻為6.0170 s. 圖3 S3的能量變化率的二次方的自回歸曲線Fig.3 Quadratic auto-regressive curves of energy gradients of S3 基于以上實驗仿真分析可初步得出以下結(jié)論: (1) 絕對幅值的STA/LTA處理結(jié)果的最大峰值對應(yīng)P波的初至,第2個最大峰值不一定對應(yīng)S波的初至.因此,可拾取絕對幅值的STA/LTA處理結(jié)果的最大峰值來確定P波的初至,但是難以依據(jù)此方法來拾取到S波的初至?xí)r間. (2) 能量變化率的STA/LTA處理結(jié)果的最大峰值對應(yīng)P波的初至,第2個最大峰值一定對應(yīng)S波的初至.因此,可拾取能量變化率的STA/LTA處理結(jié)果的最大峰值來確定P波的初至,可拾取能量變化率的STA/LTA處理結(jié)果的第2個最大峰值來確定S波的初至?xí)r間. (3) 基于能量變化率的STA/LTA處理結(jié)果和拾取的大致初至?xí)r間,分別對P波和S波初至附近的小段能量變化率曲線進(jìn)行二次方自回歸處理,分析出P波和S波的自回歸曲線,然后分別拾取自回歸曲線的最大幅值來確定P波和S波的精確初至?xí)r間. 此仿真實驗初步證明了本文提出的基于地震信號的能量變化率分別進(jìn)行STA/LTA處理和二次方自回歸模型分析,可準(zhǔn)確拾取出P波和S波初至?xí)r間.其拾取流程如圖4所示. 圖4 基于能量變化率的初至自動拾取流程Fig.4 Flow of picking up the arrival times of P and S waves automatically 為了進(jìn)一步論證噪聲強(qiáng)度對此方法的拾取結(jié)果精度的影響,本文對前面仿真分析的噪聲S1分別按照不同比例進(jìn)行放大,直到其最大值達(dá)到地震波S2最大值的50%為止,然后把地震波S2在3.4 s時刻加入到這些經(jīng)過放大后的噪聲信號S1里,最后編寫相關(guān)程序,采用本文提出的方法和傳統(tǒng)的AIC方法及STA/LTA方法,共計三種方法對疊加后的信號進(jìn)行地震波初至?xí)r間自動拾取.三種方法的分析結(jié)果如表1所示. 由表1可見:三種方法對P波和S波初至?xí)r間的拾取誤差隨著噪聲放大倍數(shù)的增加而逐步增大,如果以0.1 s為判斷拾取失效的閾值,則STA/LTA方法在放大倍數(shù)達(dá)到8倍時的拾取結(jié)果失效,AIC方法在放大倍數(shù)達(dá)到12倍時的拾取結(jié)果失效,新方法在放大倍數(shù)達(dá)到12倍時的P波初至?xí)r間拾取結(jié)果失效,但是S波初至?xí)r間拾取結(jié)果還未失效.在失效前,新方法和AIC方法的拾取誤差都明顯小于STA/LTA方法的拾取誤差;新方法對P波的拾取誤差與AIC方法的拾取誤差很接近,但是對S波初至?xí)r間的拾取誤差要小于AIC方法的拾取誤差. 表1 三種方法的自動拾取誤差統(tǒng)計 通過以上仿真實驗表明:本文提出的新方法和傳統(tǒng)的AIC方法,二者相比于STA/LTA方法,具有較強(qiáng)的抗噪聲的性能,自動拾取P波和S波的初至?xí)r間的誤差也更小些,但是新方法對S波初至?xí)r間的拾取精度要更優(yōu)于AIC方法. 如果以三種方法的計算時間長度來衡量計算效率,則STA/LTA方法的運(yùn)行效率要明顯優(yōu)于新方法和AIC方法,所需時間不到其他兩種方法的1/4;而新方法與AIC方法相比雖然稍微慢些,但是相差不到5%. 4實測地震信號的初至分析 為了進(jìn)一步檢驗基于能量變化率和二次方自回歸模型法自動拾取P波和S波初至的可靠性,本文采用此新方法、STA/LTA方法和AIC方法,分別對2013年4月20日發(fā)生的蘆山地震的10組余震記錄數(shù)據(jù)進(jìn)行自動拾取P波和S波的相對初至,同時與手工拾取結(jié)果進(jìn)行對比.圖5和圖6為第1組余震記錄數(shù)據(jù)的詳細(xì)分析結(jié)果.由圖5可得:手工拾取的P波相對初至為1.1835s、S波相對初至為2.5076 s;自動拾取能量變化率的STA/LTA分析結(jié)果的前兩個最大峰值得到P波相對初至為1.0663 s、S波相對初至為2.5702 s,且根據(jù)此拾取結(jié)果,截取(0.6751,1.4575)s和(2.1790,2.9614)s時間段內(nèi)的能量變化率曲線進(jìn)行二次方自回歸分析. 圖5 蘆山地震的余震記錄數(shù)據(jù)的初至分析Fig.5 Picking-up approximate arrival times of Lushan earthquake aftershocks 由圖6可得:P波的二次方自回歸曲線的最大值對應(yīng)時刻為0.4531 s、S波的二次方自回歸曲線的最大值對應(yīng)時刻為0.3571 s.因此,綜合兩步自動拾取結(jié)果可得P波的精確相對初至為: 圖6 能量變化率的二次方自回歸曲線Fig.6 Quadratic auto-regressive curves of energy gradient tP=0.6751+0.4531=1.1281 s S波的精確相對初至為: tS=2.1790+0.3571=2.5361 s 基于本文提出的方法對這10組余震記錄數(shù)據(jù)的相對初至的自動拾取結(jié)果如表2所示. 由表2可得: (1)如以手工拾取的結(jié)果為基準(zhǔn),則對10組地震記錄單獨(dú)采用STA/LTA法拾取的10個P波的初至?xí)r間都接近手工拾取結(jié)果;對S波的拾取結(jié)果里有9組的拾取結(jié)果接近手工拾取結(jié)果,但對第8組地震記錄S波的拾取結(jié)果與手工拾取結(jié)果相差了0.5424 s,存在較大偏差. 表2 對10組余震記錄數(shù)據(jù)初至的自動拾取結(jié)果的統(tǒng)計 (2) 如以手工拾取的結(jié)果為基準(zhǔn),分別采用AIC方法和新方法對這10組地震記錄自動拾取的初至?xí)r間拾取結(jié)果,無論是P波還是S波的初至?xí)r間的拾取結(jié)果,都接近手工拾取的結(jié)果,因此,這兩種方法的拾取結(jié)果精度都優(yōu)于STA/LTA方法. (3) 如以手工拾取的結(jié)果為基準(zhǔn),則新方法和AIC方法對P波的初至?xí)r間的拾取結(jié)果誤差基本一致,但是新方法對S波的初至?xí)r間的拾取結(jié)果誤差要小于AIC方法. 因此,通過以上比較表明了本文提出的新的地震波初至?xí)r間多步拾取方法具有較高精度和可靠性,特別是在S波的初至?xí)r間拾取上,其精度要優(yōu)于AIC方法. 為了進(jìn)一步論證本文提出的地震P波、S波自動多步拾取的新方法的可靠性,本文選擇了150條汶川地震記錄并編寫了相關(guān)程序分別采用新方法、STA/LTA方法和AIC方法,進(jìn)行了P波、S波初至?xí)r間的自動拾取,同時進(jìn)行了手工拾取.如以手工拾取結(jié)果為基準(zhǔn),則三種方法的自動拾取結(jié)果情況如表3所示. 表3 三種拾取方法對150組地震記錄的 由表3可見:STA/LTA方法拾取精度比較低,新方法和AIC方法的拾取精度都比較高.如果以0.1 s的誤差作為拾取結(jié)果失效的閾值,則AIC和新方法對P波的初至?xí)r間拾取準(zhǔn)確率都為99.333%,對S波的初至?xí)r間的拾取準(zhǔn)確率分別為95.333%和98%.STA/LTA方法的計算效率明顯優(yōu)于AIC方法和新方法,而新方法在此方面稍遜于AIC方法. 因此,通過汶川地震的150條記錄數(shù)據(jù)的P波、S波初至?xí)r間的拾取結(jié)果誤差的統(tǒng)計,進(jìn)一步檢驗了本文提出的多步自動拾取地震P波S波初至?xí)r間的新方法,雖然在計算效率方面稍遜于AIC方法,但是其對S波的初至?xí)r間的拾取精度和準(zhǔn)確率要稍微高于AIC方法,是一種較可靠的地震波初至?xí)r間自動拾取方法. 5結(jié)論 本文首先總結(jié)了常見的幾種地震P波和S波初至?xí)r間的自動拾取方法,指出了S波的初至?xí)r間自動拾取的難點,然后提出了一種新的地震P波和S波初至?xí)r間的自動拾取方法.此方法主要分三步來實現(xiàn). 第1步:把三個方向的地震記錄時程曲線合成為一組反映地震能量突變情況的能量變化率時程曲線. 第2步:地震P波和S波初至?xí)r間的初步拾取.對能量變化率時程曲線進(jìn)行STA/LTA處理,自動拾取能量變化率的STA/LTA處理結(jié)果的前兩個最大峰值,然后把這兩個最大峰值按照對應(yīng)時間的先后進(jìn)行排序.首先出現(xiàn)的峰值對應(yīng)的時間為P波初至,后面出現(xiàn)的為S波初至.依據(jù)拾取結(jié)果確定第3步分析數(shù)據(jù)范圍. 第3步:地震P波和S波初至?xí)r間的精確拾取.對由第2步拾取結(jié)果確定的兩小段能量變化率時程曲線進(jìn)行二次方自回歸分析,分別得到地震P波和S波到達(dá)時的二次方自回歸時程曲線,然后自動拾取兩條二次自回歸曲線最大值對應(yīng)的時間,分別作為地震P波和S波的精確初至?xí)r間. 基于實測地震數(shù)據(jù)的分析表明,能量變化率時程曲線的STA/LTA處理結(jié)果具有典型的特征:前兩個最大峰值分別對應(yīng)地震P波和S波的初至?xí)r間.地震P波和S波到達(dá)時的能量變化率時程曲線的二次方自回歸曲線也具有典型特征:曲線的最大值對應(yīng)P波或S波的精確初至?xí)r間. 基于此新方法、AIC方法和STA/LTA方法,分別對蘆山地震的10組地震記錄和汶川地震的150組地震記錄數(shù)據(jù)進(jìn)行了初至?xí)r間的自動拾取,且同時進(jìn)行了人工拾取獲得準(zhǔn)確的參考初至?xí)r間值.分析結(jié)果和誤差統(tǒng)計表明,新方法和AIC方法的拾取精度和可靠性明顯高于STA/LTA方法.新方法與AIC方法相比,雖然計算效率稍微差些,但是對S波的拾取精度和可靠性都稍微高些. 由于本文提出的方法的實現(xiàn)步驟是首先對地震波能量變化率曲線采用STA/LTA方法初步拾取出P波、S波的初至?xí)r間,然后再采用二次自回歸法進(jìn)一步精確拾取,因此理論上可以理解為是對現(xiàn)有STA/LTA法的補(bǔ)充和提升. 致謝感謝中國地震局工程力學(xué)研究所“國家強(qiáng)震臺網(wǎng)中心”對本文提供的數(shù)據(jù)支持. References Allen R V. 1978. Automatic earthquake recognition and timing from single traces. Bull. Seismol. Soc. Amer., 68(5): 1521-1532. Allen R V. 1982. Automatic phase pickers: their present use and future prospects. Bull. Seismol. Soc. Amer., 72(6B): S225-S242. Bai C, Kennett B L N. 2000. Automatic phase-detection and identification by full use of a single three-component broadband seismogram. Bull. Seismol. Soc. Amer., 90(1): 187-198.Boschetti F, Dentith M D, List R D. 1996. A fractal-based algorithm for detecting first arrivals on seismic traces. Geophysics, 61(4): 1095-1102. Cao S H, Greenhalgh S. 1993. Calculation of the seismic first-break time field and its ray path distribution using a minimum traveltime tree algorithm. Geophysical Journal International, 114(3): 593-600. Chang X, Liu Y K, Ashida Y. 1999. Hausdorff fractal algorithm for picking first break in seismic traces. Butsuri Tansa, 52(4): 316-322. Liu J S, Wang Y, Yao Z X. 2013. On micro-seismic first arrival identification: A case study. Chinese J. Geophys. (in Chinese), 56(5): 1660-1666, doi: 10.6038/cjg20130523. Liu Y K, Chang X, Wang H, et al. 2001. Estimation of near-surface velocity and seismic tomographic static corrections. Chinese J. Geophys. (in Chinese), 44(2): 272-278. Ma Q. 2008. Study and application on earthquake early warning [Ph. D. thesis] (in Chinese). Harbin: Institute of Engineering Mechanics, China Earthquake Administration. Maeda N. 1985. A method for reading and checking phase time in auto-processing system of seismic wave data. Zisin-JIsin, 38(3): 365-379. Paige C C, Saunders M A. 1982. LSQR: Sparse linear equations and least squares problems. ACM Transactions on Mathematical Software, 8(2): 195-209. Saragiotis C D, Hadjileontiadis L J, Panas S M. 2002. PAI-S/K: A robust automatic seismic P phase arrival identification scheme. IEEE Transactions on Geoscience and Remote Sensing, 40(6): 1395-1404. Wang J, Chen J H, Liu Q Y, et al. 2006. Automatic onset phase picking for portable seismic array observation. Acta Seismologica Sinica (in Chinese), 28(1): 42-51. Ye G X, Jiang F X, Yang S H. 2008. Possibility of automatically picking first arrival of microseismic wave by energy eigenvalue method. Chinese J. Geophys. (in Chinese), 51(5): 1574-1581. Zhou B F. 2012. Some key issues on the strong motion observation [Ph. D. thesis] (in Chinese). Harbin: Institute of Engineering Mechanics, China Earthquake Administration. 附中文參考文獻(xiàn) 劉勁松, 王赟, 姚振興. 2013. 微地震信號到時自動拾取方法. 地球物理學(xué)報, 56(5): 1660-1666, doi: 10.6038/cjg20130523. 劉伊克, 常旭, 王輝等. 2001. 三維復(fù)雜地形近地表速度估算及地震層析靜校正. 地球物理學(xué)報, 44(2): 272-278. 馬強(qiáng). 2008. 地震預(yù)警技術(shù)研究及應(yīng)用[博士論文]. 哈爾濱: 中國地震局工程力學(xué)研究所. 王繼, 陳九輝, 劉啟元等. 2006. 流動地震臺陣觀測初至震相的自動檢測. 地震學(xué)報, 28(1): 42-51. 葉根喜, 姜福興, 楊淑華. 2008. 時窗能量特征法拾取微地震波初始到時的可行性研究. 地球物理學(xué)報, 51(5): 1574-1581. 周寶峰. 2012. 強(qiáng)震觀測中的關(guān)鍵技術(shù)研究[博士論文]. 哈爾濱: 中國地震局工程力學(xué)研究所. (本文編輯何燕) 基金項目國家青年自然科學(xué)基金(51508536)資助. 作者簡介何先龍,男,1981年生,2012年畢業(yè)于中國地震局工程力學(xué)研究所,博士,主要從事工程地震與振動監(jiān)測技術(shù)的研究. E-mail:524245186@qq.com E-mail: agatha_iem@163.com *通訊作者佘天莉,女,1971年生,2012年畢業(yè)于中國地震局工程力學(xué)研究所,博士,主要從事工程地震與振動監(jiān)測技術(shù)的研究. doi:10.6038/cjg20160717 中圖分類號P315 收稿日期2015-10-24,2016-06-12收修定稿 A new method for picking up arrival times of seismic P and S waves automatically HE Xian-Long, SHE Tian-Li*, GAO Feng KeyLaboratoryofEarthquakeEngineeringandEngineeringVibration,InstituteofEngineeringMechanics,ChinaEarthquakeAdministration,Harbin150080,China AbstractTo pick up the arrival times of seismic P-wave and S-wave is the foundational work of seismic analysis. This paper presents a new method to do it. Firstly, the time-history curves of three-component seismic waves are transformed into the corresponding spatial energy gradient curves. Secondly, the STA/LTA method is applied to the energy gradient curves to pick up the approximate arrival times of P-wave and S-wave. Lastly, based on a quadratic auto-regressive model, the quadratic auto-regression is performed to the energy gradient curve around the arrival time to pick up the accurate arrival times of P-wave and S-wave. 10 groups of Lushan seismic waves and 150 groups of Wenchuan seismic waves were analyzed. By reference to the manual method, this new method has higher accuracy and stability. Compared with the STA/LTA method and AIC method, this new method is a little less in computing efficiency, but has higher accuracy and reliability in picking up the arrival time of S-wave. This new method has enriched the methods of picking up the arrival times of seismic waves automatically.KeywordsArrival time of seismic P and S wave; Energy gradient; STA/LTA method and AIC method; Quadratic auto-regression model 何先龍, 佘天莉, 高峰. 2016. 一種地震P波和S波初至?xí)r間自動拾取的新方法. 地球物理學(xué)報,59(7):2519-2527,doi:10.6038/cjg20160717. He X L, She T L, Gao F. 2016. A new method for picking up arrival times of seismic P and S waves automatically. Chinese J. Geophys. (in Chinese),59(7):2519-2527,doi:10.6038/cjg20160717.