張 琪, 王有學, 熊 彬, 蔣嬋君, 王心宇, 曾 成
(桂林理工大學 地球科學學院,桂林 541006)
在天然地震數(shù)據(jù)的處理中,拾取初至信號是地震波走時成像、研究地球內(nèi)部結(jié)構(gòu)的基礎。實際工作中,根據(jù)地震波的振幅和波形特征,利用人工拾取初至信號。然而,天然地震波含有各種干擾,人工提取存在的主觀偏見給提取初至信號帶來了極大的不唯一性。因此,以計算機為基礎的自動拾取方法毋容置疑地可以給我們提供更為準確的地震波旅行時間信息。
對于自動拾取地震波初至到時,地震科學家們提出了許多方法,諸如相關(guān)法、神經(jīng)網(wǎng)絡法、分形法等[1]。Gelchinsky等[2]提出了一種基于地震信號的特定相關(guān)屬性來確定初至的方法;王繼等[3]用多道互相關(guān)和最小二乘方法確定遠震震相初至時間;潘樹林等[4]利用優(yōu)化相關(guān)法拾取地震事件的初至信號,這些方法均基于地震道的整體特征,對噪聲有較好的抑制作用,但其拾取結(jié)果仍然受到續(xù)至波的影響較大。Murat等[5]提出了在低信噪比的信號中反向追蹤神經(jīng)網(wǎng)絡技術(shù)能在噪聲干擾背景下識別初至,但神經(jīng)網(wǎng)絡拾取初至學習速度慢,實現(xiàn)復雜,需要綜合考慮網(wǎng)格結(jié)構(gòu)、訓練方法等問題[6];美籍數(shù)學家曼德布羅特[7]提出分形思想并于1975年創(chuàng)建分形理論,分形盒維數(shù)(Fractal Box)作為分形維數(shù)的一種定義,它的大小及變化可以反映信號的不規(guī)則度和復雜度。這種方法容易理解,易于數(shù)學計算且直觀形象,所以分形理論在解決自動提取地震波初至到時的問題上發(fā)揮了長足的作用[8-10]。
對地震信號的抽樣序列來說,將其作無量綱圖形化處理可得到一條自仿射曲線,具有不規(guī)則性和近似的自相似性,可以用分數(shù)維來描述其不規(guī)則程度[11]。研究發(fā)現(xiàn),在地震波初至到達之前,地震監(jiān)測系統(tǒng)監(jiān)測到的地震信號時間序列為隨機噪聲的響應,地震波初至到達之后,地震信號時間序列為隨機噪聲與地震波動的復合響應[12]。由于地震信號的分段特征,在地震波初至到達前后地震波時間序列的分形特征會發(fā)生突變,根據(jù)其特征的變化可以定量地確定初至波的到達時間[13]。
筆者基于分形盒維算法(FBM- Fractal Box Method)提出了一種天然地震波初至拾取的方法,該方法能夠較好地解決自動拾取地震波初至的問題。
維的定義基于δ規(guī)模下的測量法的思想。設X是Rn中任一非空有界子集,對于一段離散的數(shù)字信號x(i) (i=1,2,…,N),x(n)∈X,將一系列尺度為δ(i) (i=1,2,…,M)的方網(wǎng)格對X進行覆蓋(圖1),得到有效覆蓋網(wǎng)格數(shù)Nδ(i)(X) (i=1,2,…,M),則X的上、下盒維定義[14]分別為
(1)
(2)
如果X的上、下盒維相等,則X盒維數(shù)為
(3)
離散的數(shù)字信號無法求得式(3)中的極限值,實際計算中,通過最小二乘法得到lnNδ(i)(X)~ lnδ(i)的擬合直線,其斜率就是的分形盒維數(shù)p。
圖1 數(shù)字信號分形維數(shù)的計算原理示意圖Fig.1 Calculation principle of fractal dimension on digital signal
對于地震波記錄來說,它是一個離散序列x(i) (i=1,2,…,N),我們在一個相對固定的時窗Wk(k=1,2,…,N)(長度為l)內(nèi),將時窗M等分,一般M取2的冪次方。以信號零值為零線,以相同的等分間隔畫平行于零線的平行線,使得最外側(cè)的平行線包括了整個時窗內(nèi)的地震信號,此時,在時窗內(nèi)形成一個邊長為l/(M-1)的正方形陣列。然后對時窗內(nèi)的地震信號以Δt=l/(M-1)進行重采樣,并對時窗內(nèi)重新采樣的地震信號計算與一系列尺度為δ(i)=l/(M-1)*(i-1)對應的Nδ(i),最后求得該時窗的分形盒維數(shù)pk。按照地震記錄的采樣間隔移動時窗,重復上述步驟,得到每一個時窗所對應的分形盒維數(shù),直到覆蓋整個信號序列。
通過對每個時間窗口下進行分形盒維數(shù)的計算,可得到地震波信號的分形盒維數(shù)軌跡。但是,時窗長度過小,計算量增大;時窗長度過大,降低了信號識別能力。針對時窗長度選取及它對分形維的影響,又分別選取時窗長度l=0.3 s、0.6 s、1 s、2 s這四個取值進行測試,其結(jié)果如圖2所示。
由圖2可以看出,當時窗長度取0.6 s時,分形盒維數(shù)曲線在有效信號附近下降趨勢明顯(圖2(c)),計算速度合適。
根據(jù)分形盒維方法,我們對廣西地區(qū)寬頻帶地震流動臺網(wǎng)記錄的地震事件 (表1、圖3(b)、圖4(b)、圖5(b)、圖6(b))進行盒維分析。地震波采樣間隔10 ms,時窗長度取0.6 s。為了求取分形盒維數(shù)的變化率,對分形盒維數(shù)曲線進行了平滑[15](圖3(c)、圖4(c)、圖5(c)、圖6(c)),其計算結(jié)果如圖3~圖6(b)~圖6(d)所示。
由圖3可以看出,在10 s以前分形盒維數(shù)p在1.5左右,而在10 s以后,分形盒維數(shù)p急劇降至1.2以下,由分形盒維數(shù)p的變化率(圖3(d))自動拾取初至時間為t=9.51 s;在圖4中,在3 s以前分形盒維數(shù)p在1.5左右,而在3 s以后,分形盒維數(shù)p明顯地降至1.4左右,由分形盒維數(shù)p的變化率(圖4(d))自動拾取初至時間為t=3.08 s。在圖5和圖6中,分形盒維數(shù)p具有明顯變化,由分形盒維數(shù)p的變化率(圖5(d)與圖6(d))自動拾取初至時間分別為t=7.70 s、t=18.91 s。
地震事件時間/UTC經(jīng)度/°E緯度/°N震源深度/km震級震中距/°馬里亞納群島地震2016/7/29 21:18:26145.5418.54207.627.733.47廣西梧州地震2016/7/31 9:18:13111.4724.1322.554.92.05俾斯麥群島地震2016/8/31 3:11:35152.79-3.69476.006.849.77日本本州地震2016/11/21 20:59:50141.3537.3914.696.929.57
圖3 馬里亞納群島地震記錄及其分形計算結(jié)果Fig.3 Seismic data from earthquake in Mariana islands and computed results by FBM(a)地震信號;(b)地震信號的分形盒維數(shù);(c)平滑后的分形盒維數(shù);(d)分形盒維數(shù)的變化率
圖4 廣西梧州地震記錄及其分形計算結(jié)果Fig.4 Record data from earthquake in Wuzhou and computed results by FBM(a)地震信號;(b)地震信號的分形盒維數(shù);(c)平滑后的分形盒維數(shù);(d)分形盒維數(shù)的變化率
圖5 俾斯麥群島地震記錄及其分形計算結(jié)果Fig.5 Seismic data from earthquake in Bismarck islands and computed results by FBM(a)地震信號;(b)地震信號的分形盒維數(shù);(c)平滑后的分形盒維數(shù);(d)分形盒維數(shù)的變化率
圖6 日本本州地震記錄及其分形計算結(jié)果Fig.6 Record data from earthquake in Japan's Honshu and computed results by FBM(a)地震信號;(b)地震信號的分形盒維數(shù);(c)平滑后的分形盒維數(shù);(d)分形盒維數(shù)的變化率
根據(jù)分形理論,對基于分形理論的地震波初至拾取的方法進行了研究,主要結(jié)論如下:
1)分形維理論可以根據(jù)地震波初至到達前后地震波時間序列的分形特征突變,定量地確定初至波的到達時間,能夠較好地解決自動拾取地震波初至的問題。
2)基于分形理論的地震波初至時間的拾取方法易于數(shù)學計算,其拾取精度可以達到0.01 s;拾取有唯一點,解決了手動拾取帶來的不唯一性。對于采樣間隔10.0 ms,記錄長度在30 s之內(nèi)的地震信號,耗時為3.42 s,以后可以通過改進算法更好的縮短拾取耗時。
3)試驗結(jié)果表明,當窗口大小為0.6 s時,計算速度合適,可以保證地震波波至時間的拾取精度。