大地电磁法人文噪声干扰特点及处理方法综述

葛双超， 李 斌

(1.中北大学 仪器与电子学院，太原 300051;2.太原卫星发射中心技术部，太原 300051)

0 引言

我国矿产资源种类繁多，为我国经济发展提供了坚实的能源保障，但这些矿产资源的勘探程序较低，其中大部分资源的存储情况未得到准确探测，同时浅部资源开发殆尽，矿产资源日益短缺，为国家的安全发展带来隐患，因此亟需加强深部复杂矿藏勘探能力，提高矿产资源勘查的深度和效率。

在各种地球物理勘探方法中，电磁勘探是目前应用最多的矿产资源勘探手段之一。大地电磁法(MT)通过在地表测量正交电磁场分量的变化情况，来研究地下介质的电性结构，是被动源电磁勘探的一个重要分支。MT法无需大功率供电设备且探测深度大，在矿产资源勘查、地震前兆监测、环境保护及地下水资源勘查等方面得到广泛应用。但大地电磁测深要求以平面电磁波作为有效场源，场源噪声、地质噪声以及人文噪声等所有不满足平面电磁波的电磁信号，都会对MT法观测数据造成影响[1]。尤其是随着社会经济不断发展，高压输电网、通讯电缆、地下金属管网、广播电台、信号塔及各种交通工具等带来的电磁干扰，严重污染了大地电磁实测资料，极大地影响了地下物性结构和电性分析的可解释性及数据本身的可靠性，导致大地电磁阻抗估算偏差过大，测得的视电阻率曲线失真。而实际工作中又无法避开各类干扰源，故开发有效的大地电磁数据信噪分离方法，提高大地电磁测深数据质量，始终是大地电磁数据采集与处理的核心问题[2]。

笔者首先分析了MT观测资料中的常见人文噪声及特点，然后对国内、外近几年发展起来的重要MT法数据处理方法进行了综述。主要对Hilbert-Huang 变换、小波分析、统计分析和稀疏分解法等方法的原理、特点及其在MT数据处理中的发展和应用进行了讨论，最后对MT观测数据信噪分离方法做出总结与展望。

1 观测资料常见噪声及特征

MT法的观测对象是瞬息万变的天然大地电磁信号，频带范围约为104Hz～10-4Hz。MT观测数据由有效信号和干扰噪声两部分构成，其中前者为天然大地电磁信号，由太阳风与地球磁层、电离层之间复杂的相互作用，以及雷电等地球外部不同频率段的电磁场相叠加混合形成，属于非平稳信号，具有很强的随机性；而干扰噪声主要来自工频干扰、游散电流、电设备的开关、车辆噪声、风扰、水流等，通常具有较强的规律性。MT观测系统和各类强干扰源示意图见图1。

图1 MT观测系统和干扰噪声源

1.1 工频噪声

工频干扰来自观测点附近的高压输电线，主要体现在电通道，两个正交电通道的工频分量具有较好的相关性，磁道中偶尔也会出现工频成分。工频干扰时域波形(图2)表现为等振幅类正弦波干扰，频域表现为工频干扰基频(我国通常为50 Hz)及奇次谐波干扰。这类干扰通常强度很大，往往会淹没有用信号，给后续数据处理带来很大困难。

图2 实测MT观测数据中的工频干扰

1.2 游散电流

电动设备突然开、关或负荷发生突变时，接地电流导入大地会引起游散电流噪声干扰，是MT法中的一类常见干扰，其波形图如图3所示。该类噪声通常出现在各种采样率的电道信号和磁道信号中，在时间序列上通常呈正弦阻尼振荡，见图3(a)红色部分局部放大图，其幅值是正常有用信号的几个数量级。

图3 实测MT观测数据中的游散电流干扰

1.3 电子设备开关

电子设备开关瞬间会对MT观测资料带来强烈干扰，这类干扰类似方波噪声，是矿集区内影响强度最大的噪声之一，其时频波形如图4所示。该类噪声通常出现在中低频段的电场通道中，且两正交电通道数据时域波形相关性较好，其幅值一般较大，会淹没正常大地电磁有用信号，导致阻抗估算出现严重偏差。

图4 电子设备开关造成的方波干扰

1.4 电动机噪声

电动机调速和阀门控制产生的干扰在观测数据中表现为不规则的三角波形，一般出现在磁通道中，与磁场信号具有较好的相关性(图5)。电场信号很少受到三角波噪声的影响。

图5 电动机马达干扰

1.5 车辆干扰

大型矿山机械进行矿产开采作业时，会产生大范围高强度电磁干扰，噪声波形类充放电三角波。该类噪声强度较大，观测数据时域波形有明显跳变，各通道具有较好的相关性(图6)。由图6可以看出，观测数据中除车辆干扰外还混有工频干扰。

图6 大型矿山机械车辆干扰

1.6 其他干扰

除上述几种常见干扰外，MT观测资料中还包含周围人员或动物活动、风扰等一些随机噪声，在观测数据时域波形图上体现为短时不规则噪声(图7)。

图7 短时随机干扰

2 大地电磁观测数据信噪分离方法

进行大地电磁测深时，通常采集正交电磁场的时域信息，然后进行傅里叶变换并进行阻抗估计，由于傅里叶变换适用于确定性的平稳信号，所以采用这种方法进行MT信号处理时，通常假设MT信号具有高斯分布、地质模型是最小相位系统、电磁场是平面电磁波等。王书明等[3]通过对国内多个不同地区的MT观测资料进行统计分析，证实MT信号一般具有非高斯、非线性和非最小相位性。因此，直接采用傅里叶变换来分析和处理大地电磁信号具有一定的局限性。为提高阻抗估计准确度，地球物理学家不断改进和优化MT资料处理和参数计算方法，比如时域处理中的Hilbert-Huang变换、小波分析、统计分析、经验模态分解等。

2.1 Hilbert-Huang变换

Hilbert-Huang变换是20世纪90年代中期发展起来的一种处理非线性、非平稳信号的有效方法[4]。利用Hilbert-Huang变换进行MT观测资料处理时，首先对MT信号进行经验模态分解(EMD)，获得各阶固有模态函数；然后根据噪声的不同来源和特征，对各分量进行提取和重构，从而提高大地电磁信号信噪比。

蔡建华等[5]研究了Hilbert-Huang变换在大地电磁信号噪声识别、功率谱估计和阻抗估算等方面的应用，开发了相应的数据处理算法对大地电磁信号中几种典型干扰噪声进行了处理。结果表明，所述方法对基线漂移、高频噪声和工频干扰有较好的抑制效果。

图8 小波变换去噪流程

2.2 小波变换

小波变换是一种分析瞬时时变信号的有力工具，通过伸缩和平移等运算功能对函数或信号进行多尺度细化分析，能够成功处理傅立叶变换不能解决的问题[6]。因此小波变换成为大地电磁信号信噪分离的有效方法之一。

小波变换与傅立叶变换的原理类似，但小波变换利用有限长会衰减的小波基代替无限长的三角函数基。连续小波变换公式为式(1)。

(1)

其中：x(t)为原始信号；ψ(t)为小波基函数；α为尺度因子;τ为平移因子;Q控制小波基函数的伸缩;τ控制小波基函数的平移。不同的小波基函数在处理信号时各有不同，为保证最优信噪分离效果，需要事先优选合适的小波基函数。

小波变换结果Wt(α,τ)为分解系数，是信号在所选择的小波函数空间的投影，与α成反比。

由于规则信号的能量集中在少数小波系数上，而白噪声在小波变换域上仍然分散在大量小波系数之上。利用非线性阈值对小波变换系数进行处理，将较小的小波系数置零并保留较大的小波系数(硬阈值法)，或者对较大的小波系数向零收缩(软阈值法)。然后通过逆变换对信号进行重构，进而实现信噪分离。

小波变换信噪分离流程见图8。

小波变换是一个多分辨率分析过程，以三层小波分解为例，其分析树结构如图9所示，分解结果如下：

图9 多分辨分析树结构图

x=sss+dss+ds+d

(2)

随着分解层数的增加细节信号(即高频信号、随机信号)会逐渐突出。高层分解的小波系数对应低频信号分量，通常情况下低频部分主要由规则信号或干扰构成。因此分解层次越高，去掉的低频成分越多，去噪效果越明显，但失真度也增大。所以在实际应用中，应综合考虑信号的实际变化趋势及其理论变化规律来合理确定小波分解层数，在保证有效分离不同频率的信号成分的同时，又要避免过度分解导致信号内在变化规律和趋势降低。Daniel Trad等[7]采用离散小波变换和Robust阻抗估计相结合的方法处理MT观测数据，当干扰成分在时域和频域均表现为脉冲噪声时，该方法能够获得较好的观测结果。

凌振宝等[8]针对矿集区人文噪声中的方波、基线漂移的长周期干扰，以及周边大型机械突然开关产生的高频干扰对MT观测资料的影响，基于多分辨率分析算法，先采用db3小波基对观测信号进行10层分解来去除频率低于0.2 Hz长周期噪声；然后利用经验贝叶斯估计优化小波阈值函数去除频率高于0.375 Hz的高频干扰。实测结果证实所提方法可以有效抑制近源效应。Cai Jianhua等[9]为克服小波阈值选择不合理对MT观测数据处理的误差问题，提出了一种基于中点优化算法和多层分解Stein无偏似然估计的自适应小波阈值去噪技术，并用此方法对MT观测资料中的强脉冲干扰进行了抑制。曹小玲等[10]结合小波分析和盲源分离的相关理论，提出了一种改进的独立分量分析去噪方法。模拟仿真实验证明该方法相比传统小波阈值去噪方法性能更优。

2.3 统计分析法

局部均值分解(LMD)是一种新的自适应非平稳信号处理方法。该方法自适应地将一个复杂非平稳信号，分解为若干个瞬时频率具有物理意义的乘积函数(Product function,PF)，其中每一个PF分量由一个包络信号和一个纯调频信号相乘获得[11]。前者对应该PF分量的瞬时幅值，后者可用于计算PF分量的瞬时频率。该算法由两层循环构成，内循环用于求解某PF分量的包络和调频信号，外循环通过迭代处理，将所有PF分量分离出来，最终可获得原始信号的时频分布(图10)。

图10 LMD算法流程图

李晋等[12]根据大地电磁信号的特征，结合局域均值分解法(LMD)和小波分析法，提出了一种基于局域均值分解和小波阈值的大地电磁噪声压制方法。首先利用局域均值分解将大地电磁信号自适应地分解为一系列PF分量之和；然后保留PF1分量，仅对其余各阶PF分量选择合适的软阈值进行小波阈值去噪；最后叠加重构获得大地电磁有用信号。计算机模拟仿真结果表明,所述方法可用于压制矿集区典型大尺度方波和充放电三角波干扰。Zhang Gang等[13]为补偿相关检测技术在处理MT观测数据时磁场数据相关性问题，利用磁通道和电通道的相干系数分析法，对龙门山地区的远参考MT数据进行预处理，然后再用Robust法进行阻抗张力估计。实验结果证明,在信噪比较低情况下相干分析可以有效改进阻抗估计效果。

2.4 形态滤波

傅里叶变换和小波变换都属于线性滤波，形态滤波属于非线性滤波器，可以有效提取信号的边缘轮廓以及形状特征。

数学形态学由腐蚀和膨胀两个基本运算单元构成。对原始离散信号序列x(n)，其中n=1、2、…、N；定义结构元素g(m)，其中m=1、2、…、M，且M≤N。x(n)关于g(m)的腐蚀和膨胀操作为：

腐蚀(xΘg)(n)=min[x(n+m)-g(m)],m∈0,1,…,M-1

膨胀(x⊕g)(n)=max[x(n+m)+g(m)],m∈0,1,…,M-1

(3)

其中，min和max分别为求极小和极大值函数。

腐蚀和膨胀运算等价于离散函数在滑动滤波窗(g(m))内的最小值和最大值滤波。根据腐蚀和膨胀运算，x(n)关于g(m)的开运算和闭运算如下：

开运算(x°g)(n)=(xΘg)⊕g

闭运算(x·g)(n)=(x⊕g)Θg

(4)

形态学开闭运算是一对对偶变换，其中形态开运算可以抑制信号中的正脉冲噪声，而形态闭运算可以抑制信号中的负脉冲噪声[14]。

形态学滤波的效果由变换形式和结构元素共同决定。由于形态开、闭运算存在统计偏倚现象，导致单独使用时滤波效果不够好，因此通常采用开-闭或者闭-开组合滤波方法。结构元素的选取包括结构元素的形状、宽度和高度等[15]。为达到理想的信噪分离效果，进行形态学滤波时,需要根据信号的几何特征预先定义合适的结构元素对信号进行匹配。

李晋等[16]基于数学形态学的思想，研究了传统形态滤波和广义形态滤波的大地电磁强干扰分离方法，并在形态滤波的基础上，研究了Top-hat变换、中值滤波和信号子空间增强的大地电磁二次信噪分离方法。仿真和实验结果表明,基于数学形态学的大地电磁强干扰分离方法,对大地电磁强干扰中的大尺度干扰和基线漂移有较好的抑制作用，能够有效改善MT数据质量。

进行MT阻抗估计时，当噪声满足独立Gauss分布时最小二乘估计为最优无偏估计[17]。为了得到更可靠的阻抗估计，人们又相继提出了更优的Robust阻抗估计法和远参考法。远参考大地电磁法可消除非相关噪声，其基本思想就是在同步的情况下，利用与测点相隔一定距离的远参考点的磁信号作为测点处的磁信号估算张量阻抗。当远参考点和测点间距达到一定程度，且两点位的噪声特性非相关的情况下，远参考处理能够有效提高张量阻抗计算精度。这两种方法在目前的大地电磁阻抗估计法中相对比较成熟，得到了广泛地应用[18]。以上几种阻抗估计方法的使用前提均是假设其残差服从高斯分布，然而在实际电磁法勘探中，这一假设并不一定成立。

2.5 其他数学分析方法

Shalivahan等[19]采用相干加权估计、抗差M估计方法和视电阻率加权估计处理10-3Hz～103Hz的宽频MT观测数据，结果表明三种方法联合处理能够有效抑制死区频带内的异常干扰。Antoine Saucier等[20]针对谐波干扰，在Nyman-Gaiser估计的基础上改进谐波干扰模型,并利用频域迭代法估计谐波基频和各次谐波频率，然后利用最小二乘法估计各谐波分量的振幅和相位。相比传统最小二乘法，该方法的数据处理效率更高。但是该方法要求谐波基频预估值与真实值相差不大且基频波动范围较小。Kappler等[21]为消除MT观测资料中的脉冲和阶跃噪声，利用站间数据不一致性识别干扰成分,并利用维纳滤波进行误差补偿。该方法在多台站联合观测中表现出良好的应用效果。Ge Shuangchao等[22]利用全通滤波器减去梳状滤波器方式，设计了一种全通带谐波陷波器，可以一次性滤除观测数据频带内的工频干扰基频及所有奇次谐波干扰。该滤波器适合作MT观测系统的前端滤波。李晋等[23-24]为了分离出微弱的大地电磁有用信号，提出基于递归分析和聚类的大地电磁信噪辨识及分离方法。利用该方法对大地电磁信号和强于扰的本质特征进行了分析，并研究了基于信噪辨识的矿集区大地电磁噪声压制方法。仿真和实测实验表明该方法,对于低频段的缓变化信息识别具有较好的效果。燕欢等[25]针对矿集区典型的强干扰类型，引入子空间增强方法对大地电磁观测数据中的大尺度强干扰进行去噪处理。针对矿集区MT观测常见的大尺度方波和充放电三角波干扰，通过构造协方差矩阵和特征值分解等方法将噪声子空间和信号子空间分离，并重构原始信号来实现信噪分离。

2.6 稀疏分解法

Mallat等[26-27]在小波分析的基础上,提出信号在一定条件下可以用超完备字典进行稀疏表示，并在此思想上引入匹配追踪算法。Candès等[28]提出压缩感知技术后，信号的稀疏表示作为压缩感知的重要组成部分,得到了系统性的研究和发展。利用冗余字典对信号进行稀疏分解，将信号的特征信息在变换域直观呈现出来，然后采用合适的手段进行信噪分离，使得稀疏分解方法成为一种新的信号去噪技术[29]，在波达方向估计[30]、音频信号中风声[31]和随机噪声[32]等干扰处理、图像去噪[33]、数据采集系统的重尾脉冲干扰抑制[34]等多领域得到了广泛应用。

压缩感知和信号的稀疏分解在地球物理领域的应用,主要集中在地震勘探和合成孔径雷达成像等方面。Liu Wei等[35]综合利用压缩感知和小波分析处理地震数据，通过梯度投影稀疏重构算法,将噪声分解问题转换为L1范数正则匹配问题，进而对随机噪声进行了有效抑制;Liu Sicong等[36]针对车载通信中的窄带和脉冲噪声，搭建了基于压缩感知的时频测量正交频分复用框架，利用先验稀疏度自适应匹配追踪对上述两种非高斯噪声进行重构，提高了车载通信的效率;祁锐[37]针对块稀疏信号重构问题，提出基于BAOMP 算法的分布式压缩感知重构算法，该算法的原子选择准则类似于正交匹配追踪算法，不同之处在于每次选择固定个数的块指标，从而减少迭代次数，提高了算法的重构性能;Mohammad Amir Nazari Siahsar等[38]基于矩阵低秩和稀疏分解的方法处理地震数据中的非相干随机噪声，通过对原始观测数据进行稀疏时频表示将地震有效信号分解为低秩和稀疏分量，并利用匹配追踪算法进行降噪处理，仿真和实测结果验证了算法的有效性。

在电磁勘探领域，已经有学者尝试采用稀疏表示或匹配追踪技术对MT观测资料进行处理，汤井田等[39]构造了谐波、方波和尖冲干扰原子库，利用子空间追踪和正交匹配追踪相结合的方法处理MT观测数据中的强干扰成分。由于原子库冗余度较高，该算法的缺点是耗时较长，算法效率不够高。为了改进匹配追踪在大地电磁强干扰分离中的计算效率，李晋等[40]引入匹配追踪和遗传算法处理MT观测资料中的大尺度强干扰时间序列。实验结果证实,该方法能从过完备原子库中快速、自适应地搜寻最优匹配原子重构大地电磁有用信号。但是在重构过程中会造成有用信号的损失，另外该算法对于稀疏度的自动获取的研究还不够完善。李广等[41]针对音频大地电磁观测资料中的近场源干扰问题，首先利用形态滤波处理观测数据中的大尺度强干扰，然后利用压缩感知重构算法分离尖冲噪声。该方法在强干扰区应用效果良好，但是当干扰信号强度较弱时，算法性能会随干扰强度衰减而降低。另外为提高算法效率和信噪分离效果,亟需构建更加合理的冗余字典库和开发更智能的匹配算法。汤井田等[42]通过分别设计与方波信号、脉冲信号以及谐波信号匹配,而对有用的大地电磁信号不敏感的稀疏表示字典，利用IOMP算法进行重构，有效地抑制了非平稳信号中的方波噪声和脉冲噪声。汤井田等[42]在分析噪声的规律与特征基础上，提出一种新的AMT数据处理方法。利用移不变稀疏编码字典学习方法从观测数据中自主学习人文噪声的特征结构，构建冗余字典，然后利用学习到的冗余字典对观测数据进行信噪分离，最后用去噪后的数据替换原始数据计算视电阻率和相位曲线。野外试验结果表明，所述方法能够快速、准确地分离出AMT信号中的人文干扰，保留有用信号，显著改善AMT数据质量。叶华等[43]提出了基于稀疏表示与粒子群优化算法的非平稳信号去噪方法，以冲击原子作为稀疏表示基，构建了仅对人文噪声敏感的冗余字典；采用粒子群优化算法选择最优原子，提高了匹配追踪算法的效率。

3 结语

3.1 问题总结

Hilbert-Huang变换关键步骤是EMD分解，该过程需进行循环筛选，算法效率低，且易出现模态混叠现象。由于MT法观测数据量较大，故该方法难以满足实时性需求。所以为提高工作效率，需要对该方法进行优化处理。

小波变换是最经典的非平稳信号处理方法，在诸多领域取得较好的效果，但小波滤波过分依赖小波基函数的选取，如果小波基函数选择不合适或阈值设置不合理，可能会损失部分有用信号，且小波基函数选定后，在分解与重构过程中将无法改变，缺乏自适应性。

经验模态分解类的信号处理方法，对于谐波类的干扰具有良好的抑制效果，能够较好地处理工频电力谐波等干扰，但对于方波、充放电噪声以及脉冲噪声不甚理想，且存在模态混叠和端点效应问题。形态滤波法能够压制大尺度的强人文噪声，但是对于脉冲噪声效果不佳，经过形态滤波处理后，常常残留大量的脉冲噪声，此外，形态滤波容易损伤低频有用信号。然而，矿集区的噪声干扰情况复杂，当采集的大地电磁原始数据中的噪声干扰形态特征不规则时，信号与噪声出现混淆，信噪辨识的精度将有所下降。

信号稀疏表示常用算法属于贪婪算法，由于原子库冗余度较高，算法耗时较长，算法效率不够高、容易陷入局部最优原子。另外，由于实测数据中的噪声形态复杂多样，预先设定的冗余字典无法与各类噪声完全匹配，从而不可避免的割裂噪声的完整结构，并损失部分有用信号。

3.2 展望

近年来，仿生智能优化算法在图像降噪处理、机械振动噪声识别等多个研究领域得到广泛关注。该类算法是研究者通过分析生物集群活动机理、觅食现象等而提出的用以解决复杂问题的新方法，代表性算法包括蚁群算法、萤火虫算法、蝙蝠算法等。基于智能算法的信号特征提取、信噪分离等，成为各类信号处理领域的重要发展趋势。通过系统分析和建模，构建各类不同噪声模型的超完备冗余字典库，基于智能算法实现干扰噪声的快速自适应提取，是提高MT观测资料信噪比的一种新思路。在进行此类问题研究时，需要构造合理的寻优模型、研究合理的参数调整判据，保证在干扰信号精确重构的同时最大限度保留有效MT信号，最终实现强干扰下微弱MT信号的有效提取。