徐文海，丁 勇

（四川省地震局，四川 成都 611730）

在數(shù)字化觀測系統(tǒng)中，系統(tǒng)的瞬態(tài)變化、儀器毛刺、數(shù)據(jù)記錄的階躍和尖峰以及由系統(tǒng)故障引起的信號失真和重大環(huán)境干擾源的出現(xiàn)等，都是影響數(shù)據(jù)記錄質(zhì)量的重要因素。除了真實的地震信號外，還包含了由于系統(tǒng)故障引起的各種失真信號以及臺站周圍的環(huán)境噪聲干擾等信息[1]。特別是地震烈度速報與預警工程的大背景下，如何應對地震臺站已存在及未來可能存在的各種類型的干擾，是我們需要長期面對的問題。這也就要求地震臺站或中心站的監(jiān)測人員需要依托自身熟悉自己臺站或轄區(qū)內(nèi)臺站觀測環(huán)境的優(yōu)勢，及時對可疑的干擾進行分析研究。本文認為在臺站干擾分析排查這項工作中，臺站監(jiān)測人員應該盡可能地將工作做在前面，才能有效避免潛在的負面影響。為了區(qū)分以及消除干擾，首先就需要對地震波進行分析。為此，本文采用了傅里葉變換、希爾伯特-黃變換進行頻譜分析，再結(jié)合功率譜對地震波干擾信息進行區(qū)分，根據(jù)干擾的頻率區(qū)間設計相應的濾波器來消除干擾。

1 成都地震監(jiān)測中心站職責

成都地震監(jiān)測中心站的主責主業(yè)分為3部分：①圍繞地震危險區(qū)地震活動態(tài)勢及未來地震形勢的嚴峻性，不斷提升區(qū)域地震監(jiān)測預報綜合能力。將日常監(jiān)測預報工作與重點監(jiān)視防御區(qū)的強化短臨跟蹤工作相結(jié)合，開展地震震情監(jiān)視跟蹤和年度地震趨勢會商工作。②制定針對性監(jiān)測業(yè)務拓展方案。以專業(yè)臺與市州局共促區(qū)域性地震監(jiān)測事業(yè)發(fā)展為基點，與市縣防震減災部門協(xié)同合作，共建共享，利用成都地震監(jiān)測中心站技術(shù)力量優(yōu)勢和地方良好的監(jiān)測基礎資源，優(yōu)化區(qū)域監(jiān)測臺網(wǎng)布局，提升區(qū)域地震監(jiān)測能力。③推動地震預警項目開展。組織人員開展預警項目建設實施指導服務和驗收工作，目前已完成轄區(qū)60余個站點的驗收工作。待地震速報與預警項目建設全面完成后，成都站將承擔所轄片區(qū)內(nèi)284個預警臺共計305臺預警儀器的運行維護工作。

2 測震觀測系統(tǒng)

成都站是龍門山地震帶和四川中北部地區(qū)唯一的多學科綜合型觀測臺站。近幾年，隨著“中國數(shù)字地震觀測網(wǎng)絡建設項目”“東半球空間環(huán)境地基綜合監(jiān)測子午鏈項目”“中國地震局背景場項目”“汶川地震災后恢復重建項目”“蘆山地震災后恢復重建項目”以及臺站優(yōu)化改造、標準化建設項目的相繼完成，成都站已建成一個基礎設施完善、觀測設備現(xiàn)代化、科研技術(shù)力量強勁的綜合臺站。成都站測震數(shù)據(jù)以及地磁數(shù)據(jù)均參與國際資料交換，服務于全球地震定位及地磁場空間變化研究。成都和松潘測震觀測系統(tǒng)運行穩(wěn)定，觀測資料記錄連續(xù)，震相分析準確，后續(xù)震相豐富，其觀測資料精度和分析水平在國內(nèi)始終處于前列，為中國地震臺網(wǎng)以及地震數(shù)據(jù)庫提供了大量、可靠的波形資料和震相分析資料，為地震監(jiān)測、科學研究及防震減災事業(yè)作出了應有的貢獻。成都臺作為離5·12汶川地震最近的國家級臺站，擁有一套周期360 s的超寬頻帶地震計，所記錄的寶貴地震資料被廣泛地應用于國內(nèi)外對地震的科學研究。

3 方法的MATLAB實現(xiàn)

第一步，讀取指定文件名的地震波形數(shù)據(jù)，確定參數(shù)，包括采樣頻率、地震計周期、數(shù)采轉(zhuǎn)換因子。

第二步，調(diào)用函數(shù)y=fft（x）實現(xiàn)快速傅里葉變換。

第三步，調(diào)用函數(shù)z=emd（x）對數(shù)據(jù)進行EMD分解，求得邊際譜和每個IMF分量及最后一個剩余分量residual的圖形。

第四步，用xt=x×C/2 000將count值轉(zhuǎn)化為速度值，其中C為數(shù)采轉(zhuǎn)化因子。然后調(diào)用mean()函數(shù)，用公式xt=xt-mean（xt）去均值。

第五步，調(diào)用pwelch()函數(shù)實現(xiàn)Welch算法，即：

式（1）中：pxx為功率譜估計值；f為得到的頻率點；x為進行功率譜估計的有限長序列，即上述已去均值的xt；window指所用窗函數(shù)；noverlap指重疊的點數(shù)；Nfft指采用的FFT算法長度；Fs為采樣頻率。

第六步，用函數(shù)sys=tf（[100]，[12×0.707×pi/T/2pi×pi/（T×T）/4]）求出頻率響應，再由函數(shù)[m p]=bode（sys，2×pi×f）求出幅值m，經(jīng)過矩陣轉(zhuǎn)置后，由公式pxx=pxx/abs（mag^2）扣除儀器響應，再用繪出速度功率譜。其中T為地震計周期。

第七步，用db4小波對原始信號進行5層分解并提取小波系數(shù)：[c，l]=wavedec（x，5，'db4'），再用appcoef()和detcoef()求得低頻與高頻部分，隨后進行小波消噪，之后再進行高通濾波，本文所選擇的是IIR高通濾波器。

4 應用實例

本文截取了2016-07-15T02：00記錄于成都地震監(jiān)測中心站JCZ-1T的垂直向波形數(shù)據(jù)，時間長度為300 s，包含了一個地震事件。

4.1 分析部分

由于原始波形受到嚴重的干擾，地震信息被干擾掩蓋難以分辨，如圖1所示?？焖俑道锶~變換（FFT）分析后的結(jié)果如圖2所示，由圖2可知，波形數(shù)據(jù)中存在低頻、高頻的干擾，2～4 Hz的頻率區(qū)間幅度較大。在Hilbert邊際譜上，同樣存在2～4 Hz頻率區(qū)間的較大幅度反應。這是由于2～4 Hz頻率范圍的干擾已成為成都地震臺（現(xiàn)成都地震監(jiān)測中心站）測震資料中的固有干擾[2]。

為了驗證分析的結(jié)果，設計了FIR帶阻濾波器對2～4 Hz的頻率區(qū)間進行濾波，毛刺干擾被壓制，如圖3所示，證實了2～4 Hz頻率區(qū)間的較大幅度反應是波形中存在毛刺的原因。

4.2 處理部分

db4小波分解的結(jié)果如圖4所示，由圖4可知，低頻信號振幅與原信號更接近，而高頻信號比較集中且振幅較小，與原波形相差較大，去噪便是對高頻信號進行處理。將去噪后進行IIR高通濾波，得到最終處理結(jié)果，如圖5所示。由處理結(jié)果可知，整個處理過程中原始波形中的干擾被有效地壓制。

圖1 原始波形

圖2 FFT幅頻圖與Hilbert邊際譜

圖3 消除2～4 Hz頻率區(qū)間后的波形

圖4 db4小波分解后各層低頻、高頻信號

圖5 小波去噪及IIR高通濾波后的波形

5 認識與討論

本文所提出的測震資料分析處理方法在地震臺站的實際應用中，取得了良好的效果，對干擾進行了有效的壓制。文中詳細地介紹了測震數(shù)據(jù)分析處理方法的設計原理、工作流程，核心是干擾類別的識別。測震工作人員不僅局限于地震震相的劃分，還需要對數(shù)據(jù)中存在的干擾類型進行判斷并區(qū)分是固定干擾還是隨機干擾。

本文所提出的分析處理方法的初衷是加強對干擾類別的判斷與識別，以傅里葉變換、希爾伯特-黃變換以及速度功率譜作為分析判別的理論依據(jù)，小波變換和濾波器作為分析結(jié)果的驗證手段。未來將對方法進行軟件化，并希望在更多地震臺站進行試用。