地球物理学报  2012, Vol. 55 Issue (12): 4147-4159   PDF    
庐枞矿集区大地电磁测深强噪声的影响规律
汤井田1 , 徐志敏2 , 肖晓1 , 李晋1     
1. 中南大学地球科学与信息物理学院, 长沙 410083;
2. 西北综合勘察设计研究院, 西安 710003
摘要: 天然大地电磁场信号微弱,极化方向随机,极易受电磁噪声污染.张量阻抗分析、远参考技术、Robust估计 等对随机噪声和不相关噪声有比较好的压制效果,但对强的相关噪声目前还没有有效的压制方法.庐枞矿集区人 烟稠密、工业发达,是我国著名的铁、铜、硫等矿产基地之一,区内强烈的工业、通讯、矿山、民用等电磁干扰严重污 染了大地电磁测深数据.本文首先根据实测的电磁场时间域波形和卡尼亚电阻率测深曲线形态,挑选出基本未受 噪声污染的测点(Y1650).然后利用数学形态学从受严重污染的电磁场时间序列中提取出类方波、三角波、阶跃、脉 冲和充放电5种典型噪声的波形,并以不同的方式将这些噪声波形与Y1650的电磁场波形叠加,对比分析加噪后 Y1650点电阻率和相位测深曲线的变化,进而归纳出典型强噪声对庐枞大地电磁测深资料的影响规律,为进一步 压制强噪声和资料处理提供依据.
关键词: 大地电磁测深      矿集区      强噪声      波形      相干度      信噪比     
Eefect rules of strong noise on magnetotelluric (MT) sounding in the Luzong ore cluster area
TANG Jing-Tian1, XU Zhi-Min2, XIAO Xiao1, LI Jin1     
1. Institute of Geosciences and Info-Physics,Central South University,Changsha 410083,China;
2. Northwest Research IssLiLuLe of Engineering Investigations and Design,Xian 710003,China
Abstract: Natural earth electromagnetic field signal is weak,with random polarization directions, highly susceptible to electromagnetic noise. Tensor impedance analysis, the remote reference technology,and Robust estimation can suppress the effect of random noise and uncorrelated noise, but currently no effective methods can suppress strong noise. The Luzong ore concentration area is densely populated, with highly developed industry, also one of China7 s famous bases of iron, copper, sulfur and other minerals. In this region, strong noise for industry, communication, mining, civll and other electromagnetic interference is a serious pollution to magnetotelluric sounding data. Firstly, according to the measured electromagnetic tield time domain waveform and Cagniard resistivity sounding curve shape,we pick out the basic non-noise pollution measurement points (Y1650). Then we use mathematical morphology to extract the class of square wave,triangle wave,step,impulse, and charging and discharging of five kinds of typical noise waveforms. And we stack these noise waveforms in a different way and Y1650 of the electromagnetic field waveform from the electromagnetic field time series of highly contaminated, and make contrast analysis added noise Y1650 point resistivity and phase sounding curve changes. Finally we summarize the effect rules of the typical strong noises on magnetotelluric sounding data in Luzong, and provide the basis for the further suppression of strong noise and data processing..
Key words: Magnetotelluric sounding      Ore concentration area      Strong noise      Waveform      Coherence      Signal-to-noise ration (SNR)     
1 引 言

大地电磁测深法(MT)假设来自遥远的天然电磁场以均匀平面电磁波垂直入射到地球表面.天然电磁场频率范围很宽,信号微弱,极化方向随机,因此大地电磁信号极易受到各种噪声干扰[1-2].影响大地电磁测深资料的噪声可以分为工频干扰、地质噪声和其它外界和观测系统不稳定引起的随机干扰,它们分别来自不同的源[3-5].其中工频干扰基本上产生于测点周围的人工电磁系统与环境特征,为近场干扰;地质噪声主要包括静位移畸变和地形影响,通常是全频域的;而大功率供电系统的不稳定则主要影响低频.

抑制噪声干扰,最初甘布尔等提出了一种互功率谱法[6-7],由于自功率谱中包含有噪声的成分,因此低信噪比将会不可避免的降低,采取在计算阻抗时完全使用互功率谱而不用自功率谱的方法就会大大降低噪声的影响,具有一定的抑制噪声干扰的作用;而在对大地电磁资料误差分布规律的分析与研究中许多学者提出,MT 资料的误差分布并不完全遵循高斯误差分布[8-11],所以建立在高斯误差分布之上的最小二乘阻抗张量估算法,将会导致处理结果的分散或偏离.为此,EgbertG D、BookerJR、林长佑等将Robust统计法引入到大地电磁测深法的阻抗张量估算中,提出了大地电磁Robust处理方法;在压制人文干扰对MT 数据资料的影响以及抑制局部和区域性电磁噪声上,远参考大地电磁测深法(Remote Reference MT,RRMT),以远参考点与测点之间噪声的不相关性特征为依据,通过阻抗估计进而压制了不相关噪声的干扰[12-15].近地表电性不均匀体和地形起伏的影响使大地电磁测深数据的卡尼亚电阻率-相位曲线发生畸变,称为静态效应或静位移,国内外大地电磁学者对静态效应的校正进行了广泛的研究[16-18],其成果可归纳为直接校正法、联合反演法、电磁阵列剖面(EMAP)法和阻抗张量分解技术等.综合以上调研,目前存在的问题是对强的相关噪声还没有有效的压制方法.

庐江—枞阳矿集区位于安徽省境内,经过三县两市辖区,涉及约40个乡镇.测区内人口稠密、水系发育、交通网密布、通信电力网发达,另有较多的矿山正在开采,错综复杂的干扰源,为MT 数据的采集和处理带来了许多困难.按深部探测实验研究专项课题(SinoProbe-03)的要求,在区内部署了5 条综合地球物理测线(图 1).测线总长约325km,设计MT 测点655个,频率范围320Hz~2000s.由于高山地形、水域、城镇、矿山及电力干扰等影响,实际完成MT 测点523个.数据采集共投入6套加拿大凤凰公司的V5-2000 宽频带大地电磁系统,配备MTC-50磁传感器.数据采集由于高频采样率较高,如果全时间段采集,数据量将会很大,因此高频采集采取抽样采集的办法,采集的起止时间段与低频起止时间段相同,采用1-8-5模式,即每5 min采集一次高频或中频的样(高频和中频交替采集),其中有1s的高频数据(采样率2560个样/s),连续8s的中频数据(采样率320 个样/s).低频数据(采样率24个样/s)为全时间段采集.滤波频率设为50 Hz.通过测量AC 和DC 电位差,观察饱和数据的比例,设置合理的增益.通过试验,每个MT 测点数据采集时间不低于20个小时.

图 1 庐枞矿集区MT测线布置 Fig. 1 MT survey line layout in the ore cluster area of Luzong
2 参考点(Y1650)时间序列及测深曲线特征

图 2为参考点Y1650时间域波形片段,分析该点电磁场时间域数据,可以发现很多脉冲状波形,且幅值大(比通常的大地电磁信号大1 个量级以上),能量较强,但整体时间域数据无其它非平稳信号,符合平稳随机的天然大地电磁场特征.

图 2 参考点Y1650时间序列曲线 (a)2560 Hz采样率;(b)320 Hz采样率;(c)24 Hz采样率. Fig. 2 Curve of time series at reference point Y1650 (a)2560 Hz sampling rate; (b)320 Hz sampling rate; (c)24 Hz sampling rate.

图 3是选自云南某地数据点Y1650 的MT 卡尼亚电阻率及相位测深曲线,及相干度和信噪比曲线.在等间隔双对数坐标中,Y1650 点(图 3a)的卡尼亚电阻率曲线全频段基本是连续光滑变化的,误差棒小,符合理论预期.320~0.1 Hz频段两种模式轻度分离,应是受到“静态效应"影响;相干度仅低频段最后四个频点较差,其余频段均在0.8以上,320~0.2Hz频段达到了0.9以上;信噪比仅在0.1 Hz、0.004Hz、0.0007 Hz、0.0005 Hz这四个频点处低于0.8,其它频点均在0.9左右,从卡尼亚电阻率及相位测深曲线、相干度和信噪比的分析表明,该点数据质量较好,符合实验要求.

图 3 Y1650卡尼亚电阻率及相位测深曲线(a)、相干度曲线(b)以及信噪比曲线(c) Fig. 3 (a) Cagniard resistivity and phase sounding curves, (b) degree of coherence, and (c) SNR at Y1650
3 强噪声波形特征及其分布规律

对庐枞矿集区523 个MT 测点的时间域波形进行了详细的分析统计,得到了矿集区典型强噪声的波形特征及其分布规律.

图 4所示,方波噪声在大地电磁测深数据的时间序列中表现为非正弦曲线的波形,呈类方波形态,其幅值很大,通常是正常信号的几个数量级,多出现于24Hz采样率电道数据中,通过对方波噪声的频谱分析,我们发现其影响频带范围为10 Hz以后的中低频段.

图 4 方波噪声(a)时间序列曲线及其(b)Ex频谱,(c)Ey频谱 Fig. 4 (a) Time series graph and (b)Ex spectrogram, (c) Eyspectrogram for square-wave noise

图 5所示,三角波噪声在大地电磁测深数据的时间序列中表现为非正弦曲线的锯齿波形,呈类三角形态,其幅值很大,多出现于24Hz采样率磁道数据中.通过对三角波噪声的频谱分析,我们发现其影响频带范围在10Hz以后的中低频段,其中0.1~0.01Hz的低频段最为严重.

图 5 三角波噪声(a)时间序列曲线及其(b)Hi频谱,(c)Hy频谱 Fig. 5 (a) Time series graph and (b) Hi spectrogram, (c) Hy spectrogram for triangular wave noise

图 6所示,阶跃噪声在大地电磁测深数据的时间序列中表现为大地电磁信号的突然抬升(或下降)然后向下(或向上)逐渐趋于正常大地电磁信号幅值的波形形态,其幅值可以是正常信号的若干倍甚至几个数量级,阶跃噪声存在于大地电磁测深数据的24Hz、2560 Hz采样率的电道ExEy或者磁道HxHy中,相应的磁道或者电道表现为脉冲噪声.通过对阶跃噪声的频谱分析,我们发现阶跃噪声提取于高频2560 Hz采样率中时,其频谱能量在10~1Hz频段达到最大值,当阶跃噪声提取于低频24Hz采样率中时,其频谱能量在0.1 Hz以后的低频段达到最大值,因此其影响频带范围主要为10~1Hz、0.1Hz以后的低频段.

图 6图 5,但为阶跃噪声 Fig. 6 Same as Fig. 5,but for step wave noise

图 7所示,脉冲噪声在大地电磁测深数据的时间序列中表现为尖峰形态,其幅值可以是正常信号的若干倍甚至几个数量级,脉冲噪声存在于大地电磁数据的所有采样率中,可见其影响范围之广.通过对脉冲噪声的频谱分析,我们发现脉冲噪声提取于高频2560Hz采样率中时,其频谱能量在100~1Hz频段均匀分布,当脉冲噪声提取于低频320 Hz采样率中时,其频谱能量在10~0.1Hz频段均匀分布,当脉冲噪声提取于低频24Hz采样率中时,其频谱能量在1~0.001 Hz频段均匀分布,而脉冲噪声一般在大地电磁原始数据中全频段均有出现,因此其影响频带范围为大地电磁数据的全频段.

图 7 脉冲噪声(a)时间序列曲线及其(b) Ex频谱,(c) 频谱,(d)Hx频谱,(e) Hy频谱 Fig. 7 (a) Time series graph and (b)Ex spectrogram, (c) Hxspectrogram, (d) Hyspectrogram,(e) Hy spectrogram for impulse wave noise

图 8所示,充放电噪声在大地电磁测深数据的时间序列中表现为充电、放电形态,其幅值也较大,可以是正常信号的若干倍甚至几个数量级,充放电噪声通常存在于大地电磁数据320 Hz采样率电道和磁道数据中.通过对充放电噪声的频谱分析,我们发现充放电噪声仅出现于320 Hz采样率的原始数据中,其频谱能量在10~0.1 Hz频段达到最大值,因此其影响频带范围为10~0.1Hz.

图 8图 7,但为充放电噪声 Fig. 8 Same as Fig. 7,but for charge and discharge wave noise

对庐枞强噪声的分析发现,类方波噪声是影响强度最大的一类噪声,一般出现在24 Hz采样率的电场信号中,在磁场信号中同时伴随有类三角波噪声.类阶跃波噪声也多出现在电道或者磁道的24Hz采样率信号中,2560Hz采样率信号中出现的多是尖峰干扰.类充放电噪声只存在于320Hz采样率信号中,一般同时出现在电道和磁道中.频谱分析表明,这些明显由人工源电磁场引起的强噪声一般在10~0.1Hz频率范围内,对卡尼亚电阻率和相位的影响严重.

4 强噪声对MT测深曲线的影响规律

为进一步分析这类强噪声对测深曲线的影响规律,我们选择基本没受干扰影响的Y1650点作为参考点,根据前文对庐枞强噪声的归纳结果,应用广义形态滤波对实测时间序列进行噪声分离处理,结果对于大尺度强干扰噪声具有较好的效果[19].以采样率、幅值和频谱的差异对提取出的噪声进行统计分析,并将其中典型的强干扰时间域波形叠加到Y1650的电磁场时间波形中,然后以合成后的时间域数据进行张量阻抗分析,计算卡尼亚视电阻率、相位、相干度以及信噪比等参数,对比合成前后参考点视电阻率及相位等参数的变化,总结出强噪声的影响规律.

以24 Hz采样率为例,图 9a是参考点Y1650在14:00—22:00 这段时间内任意10 min 的电场Ex时间域波形,图 9b是采用广义形态滤波从电场Ex数据中提取的典型的类方波噪声波形,时长也是10min,将其与图 9a波形数据逐点相加,得到加噪后波形如图 9c所示.对其它电磁场分量及不同采样率数据采用同样的处理方法,可得到完整的加噪的5个电磁场分量的时间域数据,用V5-2000 自带的软件即可进行阻抗分析,计算各种参数.

图 9 噪声添加过程示意图 (a)参考点Y1650原始时间序列曲线;(b)方波噪声电道Ex时间序列曲线;(c)参考点Y1650电道Ex添加方波噪声后时间序列曲线. Fig. 9 Schematic diagram of noise adding process (a) Original time series graph at reference point Y1650; (b) Time series graph for square-wave noise power road Ex ;(c) Time series graph for power road Ex add the square-wave noise at reference point Y1650.
4.1 强噪声对测深曲线的影响规律

图 10是对参考点Y1650时间域数据添加不同强噪声后计算的视电阻率和相位曲线.与图 10a没有加噪的原始曲线对比,可得出如下规律:

图 10 参考点Y1650添加噪声前后卡尼亚电阻率-相位测深曲线对比 (a)加噪前;(b)加入方波噪声后;(c)加入三角波噪声后;(d)加入阶跃噪声后; (e)加入脉冲噪声后; (f) 加入充放电噪声后. Fig. 10 Comparison of Cagniard resistivity and phase sounding curves at reference point Y1650 before and after adding noise (a) Before adding noise; (b) After adding square-wave noise; (c) After adding triangular wave noise; (d) After adding step wave noise; (e) After adding impulse wave noise; (f) After adding charge and discharge wave noise.

(1) 添加类方波、三角波和充放电噪声后,在等间隔双对数坐标下,视电阻率曲线在10~0.01 Hz频段呈45°(或大于45°)上升,相位接近或等于0.

(2) 类方波、三角波和阶跃噪声对10Hz以下的中低频段影响强烈,但基本不影响10Hz以上频段.

(3) 添加类阶跃噪声,1~0.01 Hz频段的视电阻率曲线呈45°(或大于45°)上升,相位趋于0;相比于原始曲线,视电阻率曲线在10~1Hz频段下降.

(4) 脉冲噪声严重影响全频段阻抗分析的稳定性,使视电阻率-相位曲线在全频段跳变不连续,且在1~0.01Hz频段曲线抬升.

4.2 强噪声对相干度与对信噪比的影响规律

相干度是用来度量两个场之间的相干程度,其定义为:

式中f是频率,T是时间序列的长度,PxxPyy是各自的自功率谱,Pxy是它们的互功率谱,相干度可以为0~1之间的实数值.当Coh 为0 时两个场为不相关的序列,当为1时两个场为线性相关序列.当大地电磁测深数据采集中不包含噪声时,大地电磁场的两对正交分量Ex-HyEy-Hx是线性相关的,即相干度为1,反之当含有的噪声越多,相干度就越差,数值越接近0.

信噪比是研究电、磁场分量资料受到干扰程度的一种方法.因为庐枞大地电磁测深数据采用的是以本地磁场为参考信号计算信噪比,因此磁场的信噪比恒为1,电场的信噪比为单道电场与两道磁场相关系数的平方值(图 11 中字母FB、SJB、JY、MC和CFD 分别为添加方波噪声、三角波噪声、阶跃噪声、脉冲噪声和充放电噪声后数据的简称).

图 11 参考点Y1650添加噪声前后相干度(a — e)对比、信噪比(f一j)对比 (a,f)加入方波噪声前后;(b,g)加入三角波噪声前后;(c,h)加入阶跃噪声前后;(d,i)加入脉冲噪声前后;(e,j)加入充放电噪声前后. Fig. 11 Comparisons of (a一e) coherence degree and (f一j) SNR at reference point Y1650 before and after adding noise (a,f ) Before and after adding square-wave noise; (b,g ) Before and after adding triangle-wave noise; (c,h ) Before and after adding step-wave noise; (d,i) Before and after adding impulse noise; (e,j) Before and after adding charge-discharge noise.

我们通过研究不同种类型噪声对大地电磁数据相干度的影响规律发现,方波、三角波噪声对大地电磁10Hz以下低频数据相干度影响巨大,且可以分为两个区间,一个是0.1~10 Hz频段,该频段内相干度与信噪比均受到了强烈干扰;一个是0.1Hz以下的低频段,该频段由于大地电磁信号微弱,极易受到干扰,而我们所添加的方波噪声强,因此造成该频段相干度与信噪比普遍较低;阶跃噪声对大地电磁10~0.01Hz频段数据影响巨大,0.01 Hz以下的低频段由于原始信号信噪比低导致我们添加阶跃噪声后相干度与信噪比均偏低;脉冲噪声对大地电磁全频段数据相干度均有影响,其中0.1Hz附近干扰最为严重;充放电噪声对大地电磁10~1 Hz和0.1Hz以下低频数据相干度影响巨大.本文通过对相干度和信噪比的研究进一步证明了各种强噪声的影响频带范围.

5 结 论

庐枞矿集区内强烈的工业、通讯、矿山、民用等电磁干扰严重污染了大地电磁测深数据,类方波、三角波、充放电波形、尖峰脉冲等强噪声是主要的干扰信号,通常方波、三角波仅存在于24Hz采样率信号中,充放电波形仅存在于320Hz采样率信号中.

通过将强噪声添加到未受干扰的大地电磁信号中进行仿真,结果表明,类方波、三角波、阶跃波和充放电噪声对大地电磁10Hz以下的中低频段数据影响剧烈,视电阻率曲线呈近似45°抬升,相位趋于零,与可控源音频大地电磁法(CSAMT)近区效应一致.尖峰脉冲噪声使视电阻率曲线在1~0.01 Hz频段抬升,且全频段均有不同程度的飞点,曲线形态不明确.这些结论与矿集区实测MT 测深曲线的结果是一致的,具有一定的普遍性,也为进一步压制这类强干扰提供了依据.

对此类强干扰噪声可采用广义形态滤波的方法进行压制,汤井田等一文[19]中有较详细论述,为节省篇幅,本文不再讨论.

参考文献
[1] 考夫曼, 凯勒. 大地电磁测深法. 北京: 地震出版社, 1987 . Kaufman A A, Keller G V. Magnetotelluric Sounding Method (in Chinese). Beijing: Seismological Press, 1987 .
[2] 陈乐寿, 王光锷. 大地电磁测深法. 北京: 地质出版社, 1990 . Chen L S, Wang G E. Magnetotelluric Sounding Method (in Chinese). Beijing: Geological Publishing House, 1990 .
[3] 徐志敏, 汤井田, 强建科. 矿集区大地电磁强干扰类型分析. 物探与化探 , 2012, 36(2): 214–219. Xu Z M, Tang J T, Qiang J K. An analysis of the magnetotelluric strong interference types in ore concentration areas. Geophysical and Geochemical Exploration (in Chinese) , 2012, 36(2): 214-219.
[4] 孙洁, 普光文, 白登海, 等. 大地电磁测深资料的噪声干扰. 物探与化探 , 2000, 24(2): 119–127. Sun J, Jin G W, Bai D H, et al. The noise interference of magnetotelluric sounding data. Geophysical and Geochemical Exploration (in Chinese) , 2000, 24(2): 119-127.
[5] 胡家华, 陈清礼, 严良俊, 等. MT资料的噪声源分析及减小观测噪声的措施. 江汉石油学院学报 , 1999, 21(4): 69–71. Hu J H, Chen Q L, Yan L J, et al. Analyzing noise sources of MT data and minimizing measurement noise. Journal of Jianghan Petroleum Institute (in Chinese) , 1999, 21(4): 69-71.
[6] Gamble T D, Goubau W M, Clarke J. Magnetotellurics with a remote magnetic reference. Geophysics , 1979, 44(1): 53-68. DOI:10.1190/1.1440923
[7] 邓前辉, 白改先. 互功率法在大地电磁阻抗张量估算中的应用. 石油地球物理勘探 , 1982(4): 57–64. Deng Q H, Bai G X. Crosspower method applied in the estimation of magnetotelluric impendance tensor. Oil Geophysical Prospecting (in Chinese) , 1982(4): 57-64.
[8] 杨生. 大地电磁测深法环境噪声抑制研究及其应用. 长沙: 中南大学 , 2004. Yang S. The study of restraining environmental noise and its application in magnetotelluric sounding (in Chinese). Changsha: Central South University (in Chinese) , 2004.
[9] Egbert G D, Booker J R. Robust estimation of geomagnetic transfer functions. Geophys. J. Roy. Astr. Soc. , 1986, 87(1): 173-194. DOI:10.1111/gji.1986.87.issue-1
[10] Sutarno D, Vozoff K. Robust M-estimation of magnetotelluric impedance tensors. Expl. Geophys. , 1989, 20(3): 383-398. DOI:10.1071/EG989383
[11] 林长佑, 武玉霞, 刘晓玲. 大地电磁响应函数的除偏估算和误差的研究. 西北地震学报 , 1988, 10(3): 25–38. Lin C Y, Wu Y X, Liu X L. On the unbiased estimation and error research of magnetotelluric response functions. Northwestern Seismological Journal (in Chinese) , 1988, 10(3): 25-38.
[12] 严良俊, 胡文宝, 陈清礼, 等. 远参考MT方法及其在南方强干扰地区的应用. 江汉石油学院学报 , 1998, 20(4): 34–38. Yan L J, Hu W B, Chen Q L, et al. Application of remote reference MT to noisy area. Journal of Jianghan Petroleum Institute (in Chinese) , 1998, 20(4): 34-38.
[13] 陈清礼, 胡文宝, 苏朱刘, 等. 长距离远参考大地电磁测深实验研究. 石油地球物理勘探 , 2002, 37(2): 145–148. Chen Q L, Hu W B, Su Z L, et al. Study for long-distant and far-referential MT. Oil Geophysical Prospecting (in Chinese) , 2002, 37(2): 145-148.
[14] Ritter O, Junge A, Dawes G. New equipment and processing for magnetotelluric remote reference observations. Geophys. J. Int , 1998, 132(2): 535-548.
[15] 杨生, 鲍光淑, 张全胜. 远参考大地电磁测深法应用研究. 物探与化探 , 2002, 26(1): 27–31. Yang S, Bao G S, Zhang Q S. A study on the application of remote reference magnetotelluric sounding technique. Geophysical and Geochemical Exploration (in Chinese) , 2002, 26(1): 27-31.
[16] 王家映. 电磁阵列剖面法的基本原理. 地球科学 , 1990, 15(Suppl): 1–11. Wang J Y. The basic principle of the electromagnetic array profiling. Earth Science (in Chinese) , 1990, 15(Suppl): 1-11.
[17] Bahr K. Interpretation of the magnetotelluric impedance tensor: regional induction and local telluric distortion. J. Geophys. , 1989, 62: 119-129.
[18] 杨生, 鲍光淑, 李爱勇. MT法中静态效应及阻抗张量静态校正法. 中南工业大学学报(自然科学版) , 2002, 33(1): 8–13. Yang S, Bao G S, Li A Y. The static migration to MT data and the impedance tensor static correction method. Journal of Central South University of Technology (Natural Science) (in Chinese) , 2002, 33(1): 8-13.
[19] 汤井田, 李晋, 肖晓, 等. 数学形态滤波与大地电磁噪声压制. 地球物理学报 , 2012, 55(5): 1784–1793. Tang J T, Li J, Xiao X, et al. Mathematical morphology filtering and noise suppression of magnetotelluric sounding data. Chinese J. Geophys. (in Chinese) , 2012, 55(5): 1784-1793.