连续重力信号的非潮汐信息自动提取*

杨光亮 周莉娟申重阳韦进

1)中国科学院研究生院地球科学学院地球动力学实验室，北京1000392)中国地震局地震研究所，武汉4300713)地壳运动实验室，武汉430071

连续重力信号的非潮汐信息自动提取*

杨光亮1，2，3)周莉娟2，3)申重阳2，3)韦进2，3)

1)中国科学院研究生院地球科学学院地球动力学实验室，北京100039

2)中国地震局地震研究所，武汉430071
3)地壳运动实验室，武汉430071

采用基于EMD自动提取算法及地震触发判断算法实现对台站连续重力数据的潮汐、尖峰、台阶、仪器漂移和趋势性变化的消除，达到连续重力非潮汐信息的自动分离提取。研究实例表明，该算法可实现对重力台站资料的自动处理，适于对大量台站连续重力观测资料进行批量化处理，快捷地提取非潮汐变化信息。

连续重力;潮汐;非潮汐;自动提取;批量化处理

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日潮等等。连续重力台站数据处理技术主要集中在潮汐部分的研究。PRETERNA和TSOFT可去除连续重力数据中的尖峰(突跳)、台阶和间断等数据异常。各种调和分析方法虽在算法过程上略有差异，但功能类似，均能实现长周期漂移与各种潮汐波信号的分离提取。

早期，连续重力观测主要研究观测数据的潮汐部分，如天体作用等，非潮汐部分通常作为干扰信息。随着重力台站连续重力数据由分钟采样到秒采样过渡，重力固体潮汐不仅记录了包括地球自身、日月及其他天体、地球外部大气层等的引力作用，而且也记录了地震发生前的岩石微破裂、错动等地球内部物理构造变动造成的密度分布变化及其他小区域地质作用过程。因此，重力固体潮观测不仅是研究地球内部物质分布、大气负荷、地球液核动力学、地壳运动、地球自转等地球动力学过程的利器，同时也可捕捉到反映地壳局部变化的高频信息［1］。因而，连续重力观测研究将在潮汐研究的基础上，逐步向非潮汐研究延伸和扩展。在充分利用现有观测资料的前提下，进一步发展新的连续重力非潮汐数据处理技术将是连续重力未来研究工作的重点。本文将对连续重力观测的非潮汐信息快速、自动提取技术进行探讨。

2 连续重力信号特征

设t时刻的重力观测值y(t)为:

其中s(t)为重力潮汐，d(t)为漂移，p(t)为非潮汐信息，ε(t)为观测误差。

连续重力的潮汐分析主要是提取各种潮波，这些潮波的周期从几小时到几千小时间变化［2－6］，在频谱上表现为甚低频特征［3］，连续重力原始记录中的尖峰、间断、台阶、仪器漂移及趋势项也基本表现为低频特征，在频谱上与潮汐信息存在交叉，因而很难通过频率滤波的方式去除，因此，连续重力数据中的尖峰、间断、台阶通常是以人工操作的方式去除，趋势项和仪器漂移等通过多项式拟合去除，潮波主要通过理论计算来提取。该过程通常比较繁琐，也是连续重力数据处理过程中需要人工完成的部分。

随着连续重力研究向局部、动态地震地质变化研究的延伸，在连续重力处理中通常认为是干扰的高频信号的提取和分析也成为连续重力分析的一个重要方向。现有研究理论认为，较大地震发生前都伴随着岩石的微小破裂，因而连续重力信号可能包含有临震或同震信息，可以通过对大量台站的连续重力记录进行分析，提取出具有共同特征的信息，捕捉这种临震或同震的信号。本文基于EMD的自动提取算法［7］对潮汐信号进行分离提取，以期对大量台站的连续重力进行快速处理。

3 非潮汐信号自动提取

模态分解(EMD)技术是一种信号处理技术，它可以根据信号本身的尺度特征对连续时间序列进行模态提取［8－12］。

地震信号与其相应的干扰背景在频谱与能量上比较接近，而从连续重力主要潮波的周期上看，其潮汐与非潮汐信号在频谱和能量上相差较大，理论上更容易实现分离提取。一般，24小时或48小时的连续重力包含了主要能量的潮汐波，因而连续重力数据的分析长度通常是24或48小时。该时间长度的秒采样数据量较大，小波分析由于涉及大量傅里叶变换技术及要求数据长度是2的n次方，而受到诸多限制。EMD算法不涉及谱域的计算，对数据长度也不作限制，在计算效率上优势明显。采用HHT算法对北京地球物理观象台24小时连续重力数据进行进行了分离提取(图1)，共提取出自动提取出12个模态，其中模态1是仪器漂移、长周期潮汐及趋势项信息，模态2是主要的潮汐信息，其余各模态是高频的非潮汐信息。理论计算的潮汐与分离出的主要潮汐(图1(e))得到的潮汐因子与通过潮汐分析软件(Tsoft)得出的潮汐因子结果相近。通过该方法分离出了潮汐与非潮汐信号(图1(f))。从分离结果看，分离出的非潮汐信号基本不再含有潮汐信息及长周期变化，整个处理过程均可自动完成，无需人工干预。因此，该方法适合对大量台站连续重力数据进行批量、自动处理，从海量数据中初步提取非潮汐信息。

采用文献［7］的算法，分别对成都台连续重力进行了自动潮汐与非潮汐分离计算，较好地分离出潮汐与非潮汐信息(图2)。

4 非潮汐信号尖峰消除

STA/LTA触发算法是检测瞬时突变信号的经典算法，对文献［7］中设计的地震触发算法作适当调整，可自动判断非潮汐时间序列中是否存在突变点，从而判断尖峰位置。计算表明，上述比值设置为4时，可以基本识别出非潮汐信号的尖峰(图3(a)、(b))。通过这种方式每次可去掉一部分尖峰，再对剩余时间序列重复上述步骤(1)和(2)，直到完全达到去除非潮汐时间序列中的尖峰为止。

图3(c)和图3(d)分别为迭代循环两次后的结果，该非潮汐时间序列的峰值从600退化到100。通过人机交互运算，尖峰可逐渐消除。

由此可见，使用该方法可使这一过程变得简便。

5 结论

采用基于EMD模态分解的自动分离提取算法对连续重力信号进行处理，能自动或交互完成连续记录中的潮汐(固体潮、海潮)、尖峰、台阶、漂移及趋势性变化等信息，成功提取到重力非潮汐变化信息。通过实例处理结果表明，该方法适合对大量连续重力数据进行批量处理。

但本工作只是从海量连续重力数据中初步自动分离出了非潮汐信息，目前其作用主要体现在对海量连续重力数据的特征信息进行初步筛选，更详细、精确的分析还应结合Tsoft、Baytap－G、gotic2等专业潮汐分析处理软件。本文是实现自动非潮汐信号提取的尝试，后续还需对非潮汐信息进行更详细的分解，分析其各频带所包含的物理意义，特别是对多台站地震前后的信息进行提取分析，以期能捕捉到一些有意义的地震前兆观测结果。

图1 连续重力数据的HHT提取Fig.1HHT extraction of continuous gravity data

图2 连续重力潮汐与非潮汐信息自动提取Fig.2Automatic separation and extraction of tidal and nontidal information from continuous gravity data

图3 连续重力非潮汐尖峰消除Fig.3Elimination of non－tidal peaks of gravity

1周江存，等.中国大陆精密重力潮汐改正模型.地球物理学报，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周擎，等.固体潮的地震预测研究与地球动力学研究之分析比较［J］.地球物理学进展，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蒋骏，等.固体潮潮汐因子的特征频谱及其应用［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］.大地测量与地球动力学，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于德介，程圣军，杨宇.希尔伯特黄变换在齿轮故障诊断中的应用［J］.机械工程学报，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)

杨光亮，男，博士，主要从事重力动态变化与重力壳幔结构反演研究工作.E－mail:vforyang@gmail.com