連續(xù)重力信號的非潮汐信息自動提取*

楊光亮 周莉娟申重陽韋進(jìn)

1)中國科學(xué)院研究生院地球科學(xué)學(xué)院地球動力學(xué)實驗室，北京1000392)中國地震局地震研究所，武漢4300713)地殼運(yùn)動實驗室，武漢430071

連續(xù)重力信號的非潮汐信息自動提取*

楊光亮1，2，3)周莉娟2，3)申重陽2，3)韋進(jìn)2，3)

1)中國科學(xué)院研究生院地球科學(xué)學(xué)院地球動力學(xué)實驗室，北京100039

2)中國地震局地震研究所，武漢430071
3)地殼運(yùn)動實驗室，武漢430071

采用基于EMD自動提取算法及地震觸發(fā)判斷算法實現(xiàn)對臺站連續(xù)重力數(shù)據(jù)的潮汐、尖峰、臺階、儀器漂移和趨勢性變化的消除，達(dá)到連續(xù)重力非潮汐信息的自動分離提取。研究實例表明，該算法可實現(xiàn)對重力臺站資料的自動處理，適于對大量臺站連續(xù)重力觀測資料進(jìn)行批量化處理，快捷地提取非潮汐變化信息。

連續(xù)重力;潮汐;非潮汐;自動提取;批量化處理

AbstractThe effect of tidal，peak，step，instrument drift and trend from the continuous tidal gravity data was removed by using the automatic extraction algorithm based on EMD and earthquake triggering algorithm，and achieved successfully the automatic extraction of tidal and non－tidal gravitational information.The examples show that the algorithm greatly reduces the manual work and it may be used for batch processing and quickly extracting the information of non－tidal changes of a large number of stations.

Key words:continuous gravity observation;tidal;non－tidal;automatic extraction;batch processing

1 引言

重力固體潮是由引潮力沿垂直方向的分量引起的地球重力場變化，其主要成分為長周期潮、月潮、半月潮、日潮、半日潮、1/3日潮等等。連續(xù)重力臺站數(shù)據(jù)處理技術(shù)主要集中在潮汐部分的研究。PRETERNA和TSOFT可去除連續(xù)重力數(shù)據(jù)中的尖峰(突跳)、臺階和間斷等數(shù)據(jù)異常。各種調(diào)和分析方法雖在算法過程上略有差異，但功能類似，均能實現(xiàn)長周期漂移與各種潮汐波信號的分離提取。

早期，連續(xù)重力觀測主要研究觀測數(shù)據(jù)的潮汐部分，如天體作用等，非潮汐部分通常作為干擾信息。隨著重力臺站連續(xù)重力數(shù)據(jù)由分鐘采樣到秒采樣過渡，重力固體潮汐不僅記錄了包括地球自身、日月及其他天體、地球外部大氣層等的引力作用，而且也記錄了地震發(fā)生前的巖石微破裂、錯動等地球內(nèi)部物理構(gòu)造變動造成的密度分布變化及其他小區(qū)域地質(zhì)作用過程。因此，重力固體潮觀測不僅是研究地球內(nèi)部物質(zhì)分布、大氣負(fù)荷、地球液核動力學(xué)、地殼運(yùn)動、地球自轉(zhuǎn)等地球動力學(xué)過程的利器，同時也可捕捉到反映地殼局部變化的高頻信息［1］。因而，連續(xù)重力觀測研究將在潮汐研究的基礎(chǔ)上，逐步向非潮汐研究延伸和擴(kuò)展。在充分利用現(xiàn)有觀測資料的前提下，進(jìn)一步發(fā)展新的連續(xù)重力非潮汐數(shù)據(jù)處理技術(shù)將是連續(xù)重力未來研究工作的重點。本文將對連續(xù)重力觀測的非潮汐信息快速、自動提取技術(shù)進(jìn)行探討。

2 連續(xù)重力信號特征

設(shè)t時刻的重力觀測值y(t)為:

其中s(t)為重力潮汐，d(t)為漂移，p(t)為非潮汐信息，ε(t)為觀測誤差。

連續(xù)重力的潮汐分析主要是提取各種潮波，這些潮波的周期從幾小時到幾千小時間變化［2－6］，在頻譜上表現(xiàn)為甚低頻特征［3］，連續(xù)重力原始記錄中的尖峰、間斷、臺階、儀器漂移及趨勢項也基本表現(xiàn)為低頻特征，在頻譜上與潮汐信息存在交叉，因而很難通過頻率濾波的方式去除，因此，連續(xù)重力數(shù)據(jù)中的尖峰、間斷、臺階通常是以人工操作的方式去除，趨勢項和儀器漂移等通過多項式擬合去除，潮波主要通過理論計算來提取。該過程通常比較繁瑣，也是連續(xù)重力數(shù)據(jù)處理過程中需要人工完成的部分。

隨著連續(xù)重力研究向局部、動態(tài)地震地質(zhì)變化研究的延伸，在連續(xù)重力處理中通常認(rèn)為是干擾的高頻信號的提取和分析也成為連續(xù)重力分析的一個重要方向。現(xiàn)有研究理論認(rèn)為，較大地震發(fā)生前都伴隨著巖石的微小破裂，因而連續(xù)重力信號可能包含有臨震或同震信息，可以通過對大量臺站的連續(xù)重力記錄進(jìn)行分析，提取出具有共同特征的信息，捕捉這種臨震或同震的信號。本文基于EMD的自動提取算法［7］對潮汐信號進(jìn)行分離提取，以期對大量臺站的連續(xù)重力進(jìn)行快速處理。

3 非潮汐信號自動提取

模態(tài)分解(EMD)技術(shù)是一種信號處理技術(shù)，它可以根據(jù)信號本身的尺度特征對連續(xù)時間序列進(jìn)行模態(tài)提取［8－12］。

地震信號與其相應(yīng)的干擾背景在頻譜與能量上比較接近，而從連續(xù)重力主要潮波的周期上看，其潮汐與非潮汐信號在頻譜和能量上相差較大，理論上更容易實現(xiàn)分離提取。一般，24小時或48小時的連續(xù)重力包含了主要能量的潮汐波，因而連續(xù)重力數(shù)據(jù)的分析長度通常是24或48小時。該時間長度的秒采樣數(shù)據(jù)量較大，小波分析由于涉及大量傅里葉變換技術(shù)及要求數(shù)據(jù)長度是2的n次方，而受到諸多限制。EMD算法不涉及譜域的計算，對數(shù)據(jù)長度也不作限制，在計算效率上優(yōu)勢明顯。采用HHT算法對北京地球物理觀象臺24小時連續(xù)重力數(shù)據(jù)進(jìn)行進(jìn)行了分離提取(圖1)，共提取出自動提取出12個模態(tài)，其中模態(tài)1是儀器漂移、長周期潮汐及趨勢項信息，模態(tài)2是主要的潮汐信息，其余各模態(tài)是高頻的非潮汐信息。理論計算的潮汐與分離出的主要潮汐(圖1(e))得到的潮汐因子與通過潮汐分析軟件(Tsoft)得出的潮汐因子結(jié)果相近。通過該方法分離出了潮汐與非潮汐信號(圖1(f))。從分離結(jié)果看，分離出的非潮汐信號基本不再含有潮汐信息及長周期變化，整個處理過程均可自動完成，無需人工干預(yù)。因此，該方法適合對大量臺站連續(xù)重力數(shù)據(jù)進(jìn)行批量、自動處理，從海量數(shù)據(jù)中初步提取非潮汐信息。

采用文獻(xiàn)［7］的算法，分別對成都臺連續(xù)重力進(jìn)行了自動潮汐與非潮汐分離計算，較好地分離出潮汐與非潮汐信息(圖2)。

4 非潮汐信號尖峰消除

STA/LTA觸發(fā)算法是檢測瞬時突變信號的經(jīng)典算法，對文獻(xiàn)［7］中設(shè)計的地震觸發(fā)算法作適當(dāng)調(diào)整，可自動判斷非潮汐時間序列中是否存在突變點，從而判斷尖峰位置。計算表明，上述比值設(shè)置為4時，可以基本識別出非潮汐信號的尖峰(圖3(a)、(b))。通過這種方式每次可去掉一部分尖峰，再對剩余時間序列重復(fù)上述步驟(1)和(2)，直到完全達(dá)到去除非潮汐時間序列中的尖峰為止。

圖3(c)和圖3(d)分別為迭代循環(huán)兩次后的結(jié)果，該非潮汐時間序列的峰值從600退化到100。通過人機(jī)交互運(yùn)算，尖峰可逐漸消除。

由此可見，使用該方法可使這一過程變得簡便。

5 結(jié)論

采用基于EMD模態(tài)分解的自動分離提取算法對連續(xù)重力信號進(jìn)行處理，能自動或交互完成連續(xù)記錄中的潮汐(固體潮、海潮)、尖峰、臺階、漂移及趨勢性變化等信息，成功提取到重力非潮汐變化信息。通過實例處理結(jié)果表明，該方法適合對大量連續(xù)重力數(shù)據(jù)進(jìn)行批量處理。

但本工作只是從海量連續(xù)重力數(shù)據(jù)中初步自動分離出了非潮汐信息，目前其作用主要體現(xiàn)在對海量連續(xù)重力數(shù)據(jù)的特征信息進(jìn)行初步篩選，更詳細(xì)、精確的分析還應(yīng)結(jié)合Tsoft、Baytap－G、gotic2等專業(yè)潮汐分析處理軟件。本文是實現(xiàn)自動非潮汐信號提取的嘗試，后續(xù)還需對非潮汐信息進(jìn)行更詳細(xì)的分解，分析其各頻帶所包含的物理意義，特別是對多臺站地震前后的信息進(jìn)行提取分析，以期能捕捉到一些有意義的地震前兆觀測結(jié)果。

圖1 連續(xù)重力數(shù)據(jù)的HHT提取Fig.1HHT extraction of continuous gravity data

圖2 連續(xù)重力潮汐與非潮汐信息自動提取Fig.2Automatic separation and extraction of tidal and nontidal information from continuous gravity data

圖3 連續(xù)重力非潮汐尖峰消除Fig.3Elimination of non－tidal peaks of gravity

1周江存，等.中國大陸精密重力潮汐改正模型.地球物理學(xué)報，2009，52(6):1 474－1 482.(Zhou Jiangcun，et al.Accurate correction models for tidal gravity in Chinese continent［J］.Chinese J.Geophys.，2009，52(6):1 474－1 482)

2周擎，等.固體潮的地震預(yù)測研究與地球動力學(xué)研究之分析比較［J］.地球物理學(xué)進(jìn)展，2005，20(1):118－122.(Zhou Qing，et al.Comparing the earthquake forecast of the earth tide with the geodynamics of the earth tide［J］.Progress in Geophysics.，2005，20(1):118－122)

3蔣駿，等.固體潮潮汐因子的特征頻譜及其應(yīng)用［J］.地殼形變與地震，1994，(2):88－95.(Jiang Jun，et al.Characteristics spectrum of the Earth tidal factor and its application［J］.Crustal deformation and Earthquake，1994，(2):88－95)

4Doodson A T.The Harmonic development of the tide－generating potential［J］.Proceedings of the Royal Society of London.Series A，1921，100(704):305－329.

5Casotto and Biscani F.A fully analytical approach to the harmonic development of the tide－generating potential accounting for precession，nutation，and perturbations due to figure and planetary terms［J］.AAS Division on Dynamical Astronomy，2004，36(2):67.

6Cartwright D E.Tides:a scientific history［M］.Cambridge University Press，2001.

7楊光亮，等.基于HHT的地震信號自動去噪算法［J］.大地測量與地球動力學(xué)，2010，(3):39－42.(Yang Guangliang，et al.Automatic de－nosing algorithm of earthquake signal based on HHT decomposition［J］.Journal of Geodesy and Geodynamics，2010，(3):39－42)

8Norden E，et al.The empirical mode decomposition and the Hilbert spectrum for nonlinear and nonstationary time series analysis［J］.Proceedings of the Royal Society of London，1998，454:903－995.

9Norden E，et al.A new view of nonlinear water waves——the Hilbert spectrum［J］.Annual Review of Fluid Mechanics，1999，31:417－457.

10Flandrin P，Rilling G and Gon?alves P.Empirical mode decomposition as a filterbank［J］.IEEE Signal Proc Lett.，2004，11:112－114.

11于德介，程圣軍，楊宇.希爾伯特黃變換在齒輪故障診斷中的應(yīng)用［J］.機(jī)械工程學(xué)報，2005，41(6):102－107.(Yu Dejie，Cheng Shengjun and Yang Yu.Hilbert－Huang Transform in gear fault diagnosis［J］.Mechanical Engineering，2005，41(6):102－107.)

AUTOMATIC SEPARATION AND EXTRACTION OF NON-TIDAL INFORMATION FROM CONTINUOUS GRAVITY OBSERVATION

Yang Guangliang1，2，3)，Zhou Lijuan2，3)，Shen Chongyang2，3)and Wei Jin2，3)
1)Geodynamics Laboratory，College of Earth Sciences，Graduate University of CAS，Beijing100039
2)Institute of Seismology，China Earthquake Administration，Wuhan430071
3)Crustal Movement Laboratory，Wuhan 430071

P207

A

1671-5942(2011)03-0075-04

2011－02－13

中國地震局地震研究所重點基金(IS200916004)

楊光亮，男，博士，主要從事重力動態(tài)變化與重力殼幔結(jié)構(gòu)反演研究工作.E－mail:vforyang@gmail.com