GM4磁通门磁力仪观测数据预处理分析

王秀敏 畅国平 王志敏 张建国

1）中国河北 054000 河北省地震局红山基准地震台

2）中国河北 075001 河北省地震局张家口地震监测中心站

3）中国河北 056001 河北省地震局邯郸地震监测中心站

0 引言

地磁台站的基本任务是取得连续、完整、准确、可靠的地磁观测资料，为地磁学及地震预报等相关学科研究和发展提供服务。然而，随着我国国民经济的快速发展，铁路、公路、轨道交通、国家电网等基础设施的大规模建设和投入使用，地磁台站观测环境受到不同程度的干扰（国家地震局科技监测司，1995；谢凡等，2011；刘敏等，2012），严重影响地磁观测数据质量，为干扰数据的预处理工作带来挑战。探索地磁干扰数据预处理方法，不仅有助于提高地磁观测数据质量，发挥地磁观测数据在防震减灾和其他科学研究领域的作用，而且可以促进地磁学科的应用发展，对地震地磁观测工作具有重要的实用价值和研究意义（王秀敏等，2016）。

文中选取北京、昌黎、红山地震台2017—2018 年GM4 磁通门磁力仪（下文简称GM4 磁力仪）观测资料，采用常规和小波变换方法进行数据预处理，剔除地铁轻轨、地电阻率、高压直流输电对地磁观测数据的影响，分析数据预处理效果，以便提高数据预处理效率，为地震研究提供连续、完整、可靠的基础数据。

1 地磁干扰及资料选取

1.1 干扰统计

表1 典型干扰对地磁观测的影响统计Table 1 Statistics of the effects of typical interference on geomagnetic observations

1.2 资料选取

北京、昌黎和红山地磁台站均位于首都圈地区，多年来受到台站周边地铁轻轨、地电阻率、高压直流输电等典型干扰因素的影响。其中：北京台自20 世纪90 年代后期起受到地铁运行干扰，干扰形态表现为毛刺状，地磁垂直分量Z受干扰较严重，幅度达十几nT，磁偏角D和水平分量H所受干扰幅度较小，在地铁停运时段，地磁观测数据恢复正常；昌黎台自2010 年10 月29 日起，地磁观测受到地电阻率观测的影响，干扰频次为1 次/h，干扰形态表现为台阶，地磁垂直分量Z所受干扰幅度最大，达3.5 nT，磁偏角D和水平分量H所受干扰幅度较小；红山台自2010 年10 月28 日起，地磁观测受到高压直流输电干扰，干扰形态表现为台阶，具有缓变、急始特点，地磁垂直分量Z所受干扰幅度最大，达几十nT，磁偏角D和水平分量H所受干扰幅度较小。

选取2017 年1 月1 日至2018 年11 月30 日红山、北京和昌黎地震台地磁观测数据，采用常规和小波变换分析方法进行数据预处理。

2 数据预处理方法

随着城市化进程的加速，地磁观测所受干扰因素逐渐增多，严重影响观测数据的完整性。目前，多采用秒转分方法对地磁观测数据进行预处理，对于受地铁轻轨、高压直流输电干扰数据的预处理效果明显，既提高了工作效率，又提高了数据质量。

2.1 常规处理方法

（1）分转秒数据。利用中国地震前兆数据处理系统（V2020.2 版本），对分数据进行干扰预处理，然后据此对秒数据进行自动化反处理，若有干扰数据（尖峰、台阶）遗漏，进行人工分转秒处理即可。

（2）秒转分数据。利用中国地震地磁前兆处理系统，对秒数据进行干扰预处理，使用高斯滤波算法对分数据进行反处理，且处理过程均不计入预处理日志，若预处理秒数据不存在，将跳过计算下一分量。

设计算第i分钟（imin）的分钟值，取imin00 s 及其前后各45 s 共91 s 的秒采样数据进行高斯滤波计算。公式如下

式中，Bi为第i分钟的分钟值；bi,j为第i分钟第j秒的秒采样数据，其中i、j的取值范围为00—59；Cn=C-n，为高斯系数（共91 个）。

计算00 h00 min 的分数据时，需要调用前一天的后45 个秒数据参与计算。当1 min的60 个秒数据中缺数≥10 时，对应的分数据为缺数。

“假作真时真亦假，无为有处有还无。黄梁村于万花谷固然是进入关键之一，万花谷在秦岭之中，秦岭在九州之中，九州在天地之中，又岂是偶然，说不定，它也会是天地之中的一子解双征。”

2.2 小波变换

小波变换也称多分辨分析，就是将信号或函数分解为不同频率的分量，依每个分的尺度（频率高低），按相应分辨率进行分析的方法。小波变换被广泛应用于各领域的研究工作，并得出不少有意义的成果，如：吴利辉等（2009）利用小波变换分析处理南京台地磁观测数据，发现地铁干扰被有效剔除；张明东等（2015）利用小波变换对天津地区地磁台站数据进行噪声频谱分析；张秀玲等（2018）对北京台地磁场环境干扰进小波变换分析，得出北京地磁场干扰特征。该方法的原理是，将基本小波（mother wavelet）函数位移τ，在不同尺度α下，与待分析信号χ（t）作内积，即

式中，α＞ 0，称为尺度因子，其作用是，对基本小波φ（t）函数作伸缩，τ反映了位移大小，其值可正可负，α和τ均为连续变量。式（2）又称为连续小波变换（continue wavelet transform，简记Ct-WT）。在不同尺度下，小波持续时间随α和τ值的加大而增宽，幅度则与呈反比，但波形保持不变。

文中利用Symlet 小波函数系列中的Sym4 小波，对北京地震台GM4 磁力仪原始观测数据进行分解，考虑到地铁轨道交通对地磁观测数据的干扰，频谱主要集中在0.04 Hz 及更高频率范围内，采用Sym4 小波对数据进行8 层分解，第一阶至第七阶细节系数置为零，进行信号重构得到滤波数据。

3 数据预处理分析

统计发现，北京地震台地磁观测受地铁干扰严重，昌黎地震台则受到地电阻率观测影响，红山地震台受到高压直流输电影响，选取3 个地震台站GM4 磁力仪观测数据，分别进行常规预处理，并选取北京地震台地磁观测数据进行小波变换，与不存在地铁干扰的红山地震台GM4 磁力仪观测数据进行频谱分析，以便为选择合适的数据预处理方法提供参考。

3.1 常规处理

3.1.1 地铁轻轨干扰分析。选取2017 年8 月8 日和8 月9 日北京台GM4（3）磁力仪观测数据，分别进行秒转分、分转秒处理，分析地铁轻轨对地磁观测数据的干扰特征。数据处理结果见图1、图2。由图1 中原始数据曲线可见，16 时20 分至20 时12 分数据变化正常，其他时段数据受到地铁轻轨干扰的影响，其中Z分量所受干扰较严重，幅度达11.3 nT。

图1 2017 年8 月8 日北京台GM4 磁力仪分转秒数据曲线Fig.1 The converted second-sampled data from minute-sampled Data of GM4 magnetometer observation of Beijing Station on Aug.8,2017

图2 2017 年8 月9 日北京台GM4 磁力仪秒转分数据对比曲线Fig.2 The converted minute-sampled data from second-sampled data of GM4 magnetometer observation of Beijing Station on Aug.9,2017

对2017 年8 月8 日地磁观测数据进行分转秒处理，结果见图1，可见秒数据出现缺记现象，且Z分量数据缺记明显，其中D分量秒数据完整率为34.03%，H分量为32.15%，Z分量为17.57%。对2017 年8 月9 日地磁观测数据进行秒转分处理，结果见图2，可见分数据D、H、Z分量数据完整率均为100%。

3.1.2 地电阻率干扰。选取2018 年3 月12 日和3 月13 日昌黎台GM4（2）磁力仪观测数据进行分析，发现地电阻率对地磁数据造成的干扰每小时一次，干扰持续时长180 s，其中Z分量所受干扰明显，干扰幅度达3.5 nT。对2018 年3 月12 日地磁观测数据进行分转秒处理，发现出现秒数据缺记现象。对秒数据完整率进行统计，其中Z分量为90.90%，H分量为92.22%，D分量为92.05%。对2018 年3 月13 日地磁观测数据进行秒转分处理，分数据3 个分量数据完整率均达100%。限于篇幅，文中仅列出干扰明显的Z分量预处理曲线，结果图3。

图3 2018 年3 月12 日和3 月13 日昌黎台GM4磁力仪数据预处理结果曲线Fig.3 The preprocessed curves of GM4 magnetometer observation at Changli Station on March 12 and 13,2018

3.1.3 高压直流输电干扰。据统计，目前我国有31 条高压直流输电线路对地磁台站造成影响（干扰类型有缓始型和急始型），给地磁观测数据的判别和预处理造成困难。

红山台GM4（1）磁力仪受到晋北至南京高压直流输电线路的干扰，选取2017 年5 月12 日该台地磁观测数据，发现垂直分量Z干扰幅度最大为12.9 nT。对数据进行分转秒处理，造成秒数据缺记。对秒数据完整率进行统计，可知Z分量为93.59%，H分量为97.50%，D分量为97.78%。对数据进行秒转分处理，可知分数据3 个分量完整率均达100%。限于篇幅，文中仅列出影响较严重的Z分量预处理曲线，结果见图4。

图4 2017 年5 月12 日红山台GM4磁力仪数据预处理结果曲线Fig.4 The preprocessed curves of GM4 magnetometer observation at Hongshan Station on May 12,2017

3.2 小波变换

北京台GM4（3）磁力仪观测受地铁轻轨干扰，选取该台2018 年12 月16 日地磁观测数据，进行小波变换和频谱分析。图5 给出北京台2018 年12 月16 日D、H、Z分量原始数据和小波变换处理数据对比曲线。选取不存在地铁轻轨干扰的红山台（LYH）Z分量原始数据，与北京台（BJI）Z分量原始数据及小波变换（CWT）数据，进行频谱分析，将频段分别设定为0—0.001 Hz、0.001—0.01 Hz、0.01—0.1 Hz 和0.1—0.5 Hz，频谱对比结果见图6。由图5 可见，数据经小波变换处理，在保持地磁日变形态基本不变的同时，降低了地铁干扰幅度。由图6 可知，在0—0.001 Hz 频带，北京台地磁观测数据滤波信号（小波变换处理）较好保留了原始信号的低频成分，未改变其低频成分频谱结构；在其他频带范围内，经小波变换后的滤波信号与红山台数据频谱一致性更高，表明小波变换可有效滤除地铁干扰。

图5 北京台GM4 仪原始数据和小波变换数据曲线（a）原始数据；（b）小波变换数据Fig.5 The curves of original and wavelet transform processed data of GM4 magnetometer observations at Beijing Station

图6 北京台（BJI）、红山台（LYH）原始波形与北京台（BJI）小波变换（CWT）Z 分量功率谱对比（a）0—0.001 Hz；（b）0.001—0.01 Hz；（c）0.01—0.1 Hz；（d）0.1—0.5 HzFig.6 The power spectrum of Z component original waveforms at BJI and LYH stations and the wavelet transform（CWT）waveforms at BJI station

4 讨论与结论

通过对GM4 磁通门磁力仪3 种典型干扰数据的预处理分析，可得出以下认识。

（1）采用常规方法，对受到地铁轻轨、地电阻率、高压直流输电干扰的GM4 磁力仪数据进行预处理时，分转秒均造成数据缺记，且完整率多低于90%，尤其是地铁轻轨干扰，数据完整率仅17.57%；而秒转分预处理数据无缺记，且完整率均达到100%。

（2）利用小波变换方法，对受到地铁轻轨干扰的GM4 磁力仪数据进行滤波处理，预处理数据具有正常日变形态，干扰幅度明显减小。

（3）对于受地铁轻轨干扰的数据，小波频谱分析显示，滤波信号未改变原始信号中低频部分的频谱结构，较好保留了原始信号的低频成分。因此，小波变换方法可有效抑制地铁轻轨对GM4 磁通门磁力仪观测的影响。

综上所述，受到地铁轻轨、地电阻率、高压直流输电的干扰，GM4 磁通门磁力仪观测数据可采用秒转分预处理，不但节省秒数据预处理时间，而且保证了数据的正确性和完整率，而小波变换方法对地铁轻轨干扰具有有效的抑制作用。该研究结果可为地磁观测数据预处理提供借鉴与参考。