曲浩鑫，王 怡，堵偉鵬

(內(nèi)蒙古自治區(qū)地震局海拉爾地震臺(tái)，內(nèi)蒙古 海拉爾 021000)

0 引言

小波分析是一種包含尺度伸縮和時(shí)間平移的雙參數(shù)函數(shù)分析方法，能將一個(gè)信號(hào)分解成對(duì)空間和尺度的獨(dú)立結(jié)果, 同時(shí)又保留原信號(hào)所包含的信息[1]?；谛〔ㄗ儞Q能分離出不同頻帶的信號(hào)，目前，小波分析法廣泛應(yīng)用于天然地震震相識(shí)別，重力、大面積形變測(cè)量和定點(diǎn)形變測(cè)量等地震異常分析中[2-7]。宋治平等研究表明，小波分析法對(duì)形變數(shù)字化資料中的干擾識(shí)別與消除、不同頻率潮汐波的分離、趨勢(shì)異常與短期異常的提取等方面都具有較強(qiáng)的功能[8]；張軍等應(yīng)用小波分析法對(duì)地下流體資料進(jìn)行研究，得出較好抑制流體數(shù)據(jù)中隨機(jī)噪聲的方法，并通過(guò)重構(gòu)小波系數(shù)達(dá)到消除降雨干擾的目的[9]；劉建明等研究得出小波分析可作為提取定點(diǎn)形變異常的一種有效方法，通過(guò)將定點(diǎn)形變觀測(cè)數(shù)據(jù)分解成不同頻段，對(duì)提取異常信號(hào)有較好的效果，提高了識(shí)別地震信息的能力[6]；倪友忠等通過(guò)小波變換法對(duì)面應(yīng)變?nèi)站颠M(jìn)行分析，發(fā)現(xiàn)在海域大震前，小波變換第三階細(xì)節(jié)部分在震前出現(xiàn)持續(xù)數(shù)月以上的高頻脈沖，且面應(yīng)變有鼓包現(xiàn)象[10]。

因此，在小波變換理論基礎(chǔ)上，運(yùn)用小波分析法對(duì)海拉爾地震臺(tái)(以下簡(jiǎn)稱(chēng)海拉爾臺(tái))前兆數(shù)據(jù)進(jìn)行處理，總結(jié)小波分析在前兆資料中的應(yīng)用經(jīng)驗(yàn)，以期為識(shí)別、提取及去除前兆異常信息提供新思路。

1 小波分析理論基礎(chǔ)

此時(shí)稱(chēng)Ψ(t)為一個(gè)基本母小波，將母函數(shù)Ψ(t)經(jīng)伸縮和平移后得到如下公式：

式中：a為伸縮因子；b為平移因子。稱(chēng)其為一個(gè)小波序列。在實(shí)際運(yùn)用中，連續(xù)小波必須加以離散化。任意函數(shù)f(t)的離散小波變換可表示為：

為區(qū)分不同分量信息，引入構(gòu)造正交小波的方法和小波變換的快速算法—Mallat。應(yīng)用Mallat算法不必知道具體的小波函數(shù)，同時(shí)，也使小波變換在計(jì)算上變得可行。Mallat一維分解算法可表示為：

式中：Cj+1f與Dj+1f是在分辨率為j時(shí)的近似部分和細(xì)節(jié)部分。對(duì)于Mallat一維重構(gòu)算法可表示為：

Cj-1[n]=∑nh(n-2k)Cj[k]+∑ng(n-2k)Dj[k],

式中：Cj[k]為低頻系數(shù)；Dj[k]為高頻系數(shù)；h，g分別為低通和高通濾波器。

2 小波基函數(shù)的理論選取

在Mallat算法中小波基函數(shù)的選擇不是唯一的，所有滿(mǎn)足小波條件的函數(shù)都可作為小波函數(shù)，不同的小波基函數(shù)使信號(hào)的分解和重構(gòu)有不同的結(jié)果。因此,選取合適的小波基函數(shù)就顯得較重要。目前，常用的小波基函數(shù)有：Haar小波、Daubechies(dbN)小波系、Biorthogonal(biorNr.Nd)小波系、Coiflet(coifN)小波系、SymletsA(symN)小波系、Morlet(morl)小波、Meyer小波等。常用小波基函數(shù)的主要特點(diǎn)如表1所示，選取的小波基函數(shù)既能突出異常信號(hào)又能最大限度地保留前兆資料的固體潮信息。

表1 常用小波基的主要特點(diǎn)Table 1 Main characteristics of commonly used wavelet bases

根據(jù)小波函數(shù)的定義和常用小波基的特點(diǎn)，結(jié)合前兆信號(hào)資料的特征及研究目的，綜合考慮小波基函數(shù)的對(duì)稱(chēng)性、正則性、緊支撐性及消失矩的階數(shù)[6,12-16]。理論上，Daubechies(dbN)小波系中的db4小波基函數(shù)能滿(mǎn)足前兆資料處理的要求。db4小波的尺度函數(shù)和小波函數(shù)如圖1所示，較好地兼顧了緊支撐性和正則性，同時(shí)，消失矩的階數(shù)也較合適，能在保留固體潮信息的基礎(chǔ)上較好地反應(yīng)異常信號(hào)。

為檢驗(yàn)db4小波在實(shí)際前兆數(shù)字化資料處理中的效果，選取2020年1月24日海拉爾臺(tái)垂直擺NS分量數(shù)據(jù)進(jìn)行小波變換，圖2為應(yīng)用db4小波基函數(shù)進(jìn)行離散小波變換的結(jié)果。可以看出，2020年1月24日15點(diǎn)09分塔吉克斯坦5.6級(jí)地震的信號(hào)在高頻分

圖1 db4小波尺度函數(shù)及小波函數(shù)Fig.1 db4 wavelet scaling function and wavelet function

解中變得突出，同時(shí)，應(yīng)用近似5階小波重構(gòu)消除了該地震波對(duì)垂直擺正常波形的影響，垂直擺的固體潮信息被較好地保留。由此可見(jiàn)，應(yīng)用db4小波基函數(shù)進(jìn)行離散小波變換對(duì)海拉爾臺(tái)前兆異常信號(hào)的識(shí)別、提取及消除是適用的。

圖2 垂直擺NS分量曲線Fig.2 NS component curve of vertical pendulum

3 小波分析剔除干擾

海拉爾臺(tái)目前共有5套數(shù)字化前兆觀測(cè)儀器，前兆資料信息量大且不同儀器會(huì)受到不同程度干擾，在原始曲線中不易識(shí)別和區(qū)分。小波變換理論具有較強(qiáng)的干擾識(shí)別和排查能力，因此，選取海拉爾臺(tái)受氣壓干擾和場(chǎng)地環(huán)境影響的典型事件，采用小波分析法進(jìn)行識(shí)別和消除干擾。

3.1 氣象干擾

第22頁(yè)圖3a為2020年8月2日至12日伸縮儀NS分量原始數(shù)據(jù)曲線，第22頁(yè)圖4a為2020年7月1日至3日水管儀NS分量原始數(shù)據(jù)曲線，從圖中難以識(shí)別出具體的干擾時(shí)段。應(yīng)用小波變換法對(duì)圖3a、圖4a原始數(shù)據(jù)曲線進(jìn)行細(xì)節(jié)8階分解，分解細(xì)節(jié)部分如圖3d、圖4d所示?？梢钥闯觯琩1-d2未見(jiàn)明顯異常，從d3開(kāi)始曲線出現(xiàn)明顯異常(圖中方框部分)，上下波動(dòng)幅度超過(guò)2倍中誤差。通過(guò)與同時(shí)段氣壓變化曲線進(jìn)行對(duì)比發(fā)現(xiàn)，d3-d6明顯的上下階躍部分與氣壓異常波動(dòng)時(shí)段一一對(duì)應(yīng)，因此，可以判定伸縮儀和水管儀的異常變化由氣壓突變引起。

在小波分解d1-d8中，d1-d4主要是高頻信號(hào)段部分，d5-d8為低頻信號(hào)段部分，氣壓影響在d3-d6部分清晰顯示，d8頻帶部分明顯過(guò)濾了異常變化，顯示的是固體潮信息。綜上所述，海拉爾臺(tái)形變儀受氣壓干擾在小波分解d3-d6部分較明顯。

圖3c為通過(guò)小波重構(gòu)剔除受氣壓影響后的曲線，與原始曲線圖3a相比，在保留原始信息的情況下，較好地去除了氣壓影響造成的固體潮畸變，證明小波分析法對(duì)排除氣壓干擾具有較好效果。

圖3 伸縮儀2020年8月2日至12日NS分量曲線Fig.3 NS component curve of extensometer from August 2 to 12, 2020

圖4 水管儀2020年7月1日至3日NS分量曲線Fig.4 NS component curve of water tube instrument from July 1 to 3, 2020

3.2 場(chǎng)地環(huán)境干擾

由第23頁(yè)圖5a看出,從2018年10月25日10時(shí)30分左右至27日8時(shí)30分左右，伸縮儀受干擾出現(xiàn)固體潮畸變，經(jīng)核查為體應(yīng)變鉆新井工程影響所致，屬于典型的場(chǎng)地環(huán)境干擾。通過(guò)小波細(xì)節(jié)8階分解顯示有3個(gè)異常時(shí)間段(見(jiàn)圖5c)，分別為25日10時(shí)至16時(shí)、26日 13時(shí)至22時(shí)、27日7時(shí)至8時(shí)，三個(gè)時(shí)間段與三次鉆井時(shí)間相吻合。從圖5a原始曲線中較難判斷出具體受干擾的時(shí)段，小波分解d3-d6部分壓制了其余信息，將場(chǎng)地環(huán)境干擾凸顯出來(lái)，其中細(xì)節(jié)5階最為明顯。將場(chǎng)地環(huán)境干擾信息提取之后對(duì)原始曲線進(jìn)行干擾剔除，剔除結(jié)果如圖5b所示。與原始曲線相比，曲線光滑度高，高頻干擾信號(hào)基本被抑制。

圖5 伸縮儀2018年10月25日至27日NS分量曲線Fig.5 NS component curve of extensometer from October 25 to 27, 2018

4 小波分析地震與爆破

4.1 小波分析近震與遠(yuǎn)震

海拉爾臺(tái)不同前兆儀器識(shí)別地震的能力不同，多種手段互補(bǔ)觀測(cè)可以清晰地記錄到遠(yuǎn)震、近震以及地方震的同震效應(yīng)[8，10]。對(duì)于前兆分析來(lái)講，這種同震效應(yīng)幅度大，夾雜在觀測(cè)曲線中，過(guò)大的振動(dòng)幅度會(huì)壓制前兆儀器正常的固體潮信息，也屬于一種干擾。因此，采用小波分析法對(duì)同震效應(yīng)造成的突跳進(jìn)行識(shí)別和處理是必要的。

圖6a為海拉爾臺(tái)2017年5月12日重力潮汐觀測(cè)值原始曲線，可以看出，呼倫貝爾市鄂倫春旗3.7級(jí)地震明顯，中美洲沿岸遠(yuǎn)海6.2級(jí)地震難以辨別。從圖6c中看出，d1-d8中有兩處明顯異常，是由19時(shí)11分的中美洲沿岸遠(yuǎn)海6.2級(jí)遠(yuǎn)震和22時(shí)19分的呼倫貝爾市鄂倫春旗3.7級(jí)近震造成的，在原始曲線上難以辨別的中美洲沿岸遠(yuǎn)海地震經(jīng)過(guò)小波變換后清晰可見(jiàn)。第24頁(yè)圖7a為海拉爾臺(tái)2019年9月30日垂直擺EW分量原始曲線，其中，智利中部沿岸近海地震為遠(yuǎn)震，呼倫貝爾市新巴爾虎右旗地震為近震，近震在原始曲線中不明顯，經(jīng)過(guò)小波變換后可輕易識(shí)別。

由圖6c、圖7c看出,近震和遠(yuǎn)震在小波變換細(xì)節(jié)分解中具有完全不同的特征，遠(yuǎn)震在d3-d5階中異常明顯；近震在d1-d3、d6-d8階中較明顯，d1-d3階異常幅度逐漸減弱，d6-d8階異常范圍的寬度逐漸增加。圖6b、圖7b為去除干擾后的重力潮汐觀測(cè)曲線，比原始數(shù)據(jù)曲線清晰，日變形態(tài)正常顯現(xiàn)。通過(guò)以上分析可以得出，小波分析對(duì)同震效應(yīng)的識(shí)別、提取及干擾消除均具有較好的效果。

4.2 小波分析爆破

第24頁(yè)圖8a為海拉爾臺(tái)垂直擺NS分量2019年6月10日原始數(shù)據(jù)曲線，可以看出，11時(shí)左右北南分量數(shù)據(jù)形態(tài)表現(xiàn)為上下階躍，難以判斷此階躍是由地震還是爆破引起。因此，對(duì)曲線進(jìn)行小波細(xì)節(jié)8階分解，結(jié)果如圖8b所示。對(duì)比圖8b、8c發(fā)現(xiàn)，爆破與地震兩種信號(hào)在小波細(xì)節(jié)分解上有明顯區(qū)別。爆破信號(hào)在小波分解d1-d8階中均清晰可見(jiàn)，異常的寬度和幅度變化不大；地震信號(hào)在d1-d8階中的某幾部分異常明顯，同時(shí)，地震類(lèi)型的不同，異常的寬度和幅度變化呈現(xiàn)不同的態(tài)勢(shì)。這是由于地震信號(hào)比爆破信號(hào)復(fù)雜，信號(hào)的頻帶范圍更寬所致。第25頁(yè)圖9為海拉爾臺(tái)2020年4月7日垂直擺EW分量記錄到的爆破，通過(guò)做小波分解發(fā)現(xiàn)與圖8具有相同的規(guī)律，即在小波分解d1-d8階中爆破信號(hào)的寬度和幅度變化不大。由此可見(jiàn)，利用小波細(xì)節(jié)分解辨別原始信號(hào)異常是由地震還是爆破引起具有較好的效果。

圖6 2017年5月12日重力儀曲線Fig.6 Gravimeter curve of May 12, 2017

圖7 2019年9月30日垂直擺曲線Fig.7 Vertical pendulum curve on September 30, 2019

5 結(jié)論與討論

通過(guò)以上分析，得出如下認(rèn)識(shí)：

(1) 小波變換可將不同頻帶的信號(hào)較好地分離，實(shí)現(xiàn)對(duì)各頻段小波系數(shù)的逐一分析，從而達(dá)到提取前兆異常、去除特定頻段干擾的目的。

(2) 通過(guò)對(duì)不同小波基函數(shù)特征的歸納分析,并對(duì)前兆資料進(jìn)行處理，驗(yàn)證了db4小波是海拉爾臺(tái)前兆數(shù)字化資料處理較適合的小波基函數(shù)，能滿(mǎn)足在保留固體潮信息的基礎(chǔ)上檢測(cè)前兆異常。

圖8 垂直擺2019年6月10日NS分量曲線Fig.8 NS component curve of vertical pendulum on June 10, 2019

圖9 2020年4月7日垂直擺EW分量曲線Fig.9 EW component curve of vertical pendulum on April 7, 2020

(3) 由海拉爾臺(tái)形變儀器小波分析結(jié)果可知，小波分析法對(duì)氣壓、場(chǎng)地環(huán)境干擾的識(shí)別和剔除具有較好效果。氣壓、場(chǎng)地環(huán)境干擾均在小波細(xì)節(jié)分解d3～d6部分異常顯著，其中d5階最明顯，d8階表征的基本是固體潮信息。應(yīng)用小波變換重構(gòu)信號(hào)基本去掉了氣壓、場(chǎng)地環(huán)境對(duì)原始信號(hào)的影響，去噪效果顯著，信號(hào)曲線光滑度高，基本抑制了噪聲信號(hào)。

(4) 小波變換法在海拉爾臺(tái)同震效應(yīng)的識(shí)別、提取和去除具有較好效果。通過(guò)對(duì)比遠(yuǎn)震和近震小波變換結(jié)果得出，遠(yuǎn)震在小波細(xì)節(jié)分解d3-d5部分異常顯著，近震在小波細(xì)節(jié)分解d1-d3、d6-d8部分明顯，且d1-d3部分異常幅度逐漸減弱，d6-d8部分異常范圍的寬度逐漸增加。

(5) 利用小波細(xì)節(jié)分解法可對(duì)原始異常相似的地震和爆破信號(hào)進(jìn)行辨別。地震信號(hào)只在小波細(xì)節(jié)分解d1-d8階中的某幾部分異常明顯，同時(shí)，地震類(lèi)型的不同，異常的寬度和幅度變化呈現(xiàn)不同的態(tài)勢(shì)；爆破信號(hào)在小波細(xì)節(jié)分解d1-d8階中均清晰可見(jiàn)，異常的寬度和幅度變化不大。