基于模拟退火法的磁共振测深多源谐波噪声压制方法

doi:10.11720/wtyht.2022.1158

基于模拟退火法的磁共振测深多源谐波噪声压制方法

陈亮^,¹, 付立恒¹, 蔡冻¹, 李凡², 李振宇², 鲁恺²

1.中国电建集团江西省电力设计院有限公司,江西南昌 330096

2.中国地质大学(武汉) 地球物理与空间信息学院,湖北武汉 430074

Suppression method of multi-source harmonic noise in magnetic resonance sounding based on simulated annealing method

CHEN Liang^,¹, FU Li-Heng¹, CAI Dong¹, LI Fan², LI Zhen-Yu², LU Kai²

1. POWERCHINA Jiangxi Electric Power Engineering Co., Ltd., Nanchang 330096, China

2. Institute of Geophysics and Geomatics, China University of Geosciences (Wuhan), Wuhan 430074, China

责任编辑: 沈效群

收稿日期: 2021-03-23 修回日期: 2021-09-1

基金资助:

国家重点基础研究发展计划(973计划)项目(2011CB710606)

Received: 2021-03-23 Revised: 2021-09-1

作者简介 About authors

陈亮(1994-),男,助理工程师,硕士,主要从事电磁法勘探数据处理及应用研究工作。Email: chenliang1994@cug.edu.cn

摘要

磁共振测深找水方法在高电磁噪声环境下应用时,往往会因为电磁噪声的干扰导致实测数据信噪比降低,反演结果难以准确地确定实际地下水体的分布状况,降低方法的应用效果。针对磁共振测深找水方法在实际野外数据采集工作中常面临的多源谐波噪声干扰问题,在模型去噪的基础上推导出网格搜索同步删除法,进一步提出了效率更高的模拟退火同步删除法。仿真结果表明,2种方法均能有效压制多源谐波噪声,而模拟退火同步删除法能够更快速地搜索谐波基频值,相较网格搜索同步删除法在双基频谐波情况下工作效率提高2.35倍,大大降低了去噪流程的时间成本,且对多源谐波噪声也能取得较好的去噪效果。将该去噪算法应用于实例中,反演结果和钻孔资料的对比表明,模拟退火同步删除法可以有效压制磁共振测深方法实测数据中的多源谐波噪声,显著提高磁共振测深方法的应用效果。

关键词： 模拟退火法; 压制噪声; 磁共振测深; 多源谐波

Abstract

When the magnetic resonance sounding (MRS) method is applied in an environment with high electromagnetic noise, the signal-to-noise ratio of the measured data is often reduced due to the interference of electromagnetic noise. As a result, it is difficult to accurately determine the aquifer distribution using the inversion results, thus reducing the application effects of the method. In this paper, aiming at the common problem of multi-source harmonic noise interference in the field data acquisition using the MRS method, this paper derives the grid search simultaneous removal method based on the model denoising and further proposes the more efficient simulated annealing simultaneous removal method. The simulation results show that both methods can effectively suppress multi-source harmonic noise. Compared with the grid search simultaneous removal method, the efficiency of the simulated annealing simultaneous removal method is improved by 2.35 times in the case of double fundamental frequency harmonics, which greatly reduces the time cost of the denoising process. Meanwhile, the simulated annealing simultaneous removal method allows for great denoising effects of multi-source harmonic noise. Finally, the proposed denoising algorithm was applied to a field example. The comparison of the inversion results and borehole data shows that the simulated annealing simultaneous removal method can effectively suppress the multi-source harmonic noise in the measured data obtained using the MRS method and can significantly improve the application effects of the method.

Keywords： simulated annealing; magnetic resonance sounding; multi-source harmonic; denoising

PDF (4080KB) 元数据多维度评价相关文章导出 EndNote| Ris| Bibtex 收藏本文

本文引用格式

陈亮, 付立恒, 蔡冻, 李凡, 李振宇, 鲁恺. 基于模拟退火法的磁共振测深多源谐波噪声压制方法[J]. 物探与化探, 2022, 46(1): 141-149 doi:10.11720/wtyht.2022.1158

CHEN Liang, FU Li-Heng, CAI Dong, LI Fan, LI Zhen-Yu, LU Kai. Suppression method of multi-source harmonic noise in magnetic resonance sounding based on simulated annealing method[J]. Geophysical and Geochemical Exploration, 2022, 46(1): 141-149 doi:10.11720/wtyht.2022.1158

0 引言

磁共振测深(magnetic resonance sounding,MRS)方法是一种基于水中氢核的核磁共振特性实现地下水探查的地球物理勘探方法,地下水中氢核在地磁场的作用下形成宏观磁矩,当被拉莫尔频率的交变电流脉冲激发时,该宏观磁矩将受到激发而偏离地磁场方向,撤去交变电流脉冲,磁矩发生弛豫,产生核磁共振(nuclear magnetic resonance,NMR)信号^[1]。由于具有直接找水、勘探尺度大、反演信息丰富等优势,近年来,该方法已广泛应用于水资源勘察、水文地质调查、文物保护^[2]、堤坝检测、滑坡监测^[3]等领域。MRS方法接收的NMR信号一般非常微弱(nV级别),因而极易受到周围环境电磁噪声的干扰。在城镇人居区域、工业厂区以及矿区、隧道等地区谐波电磁噪声干扰很复杂,成为制约方法应用效果的主要因素^[4]。MRS方法所面临的环境电磁噪声一般可分类为尖峰噪声、谐波噪声和随机噪声。常规的MRS方法的信号处理流程为:① 去除尖峰噪声;② 去除谐波噪声;③ 叠加;④ 包络线求取和参数提取^[5]。在提取到NMR信号的3个关键参数(初始振幅E₀、初始相位φ₀和平均弛豫时间 $T_{2}^{*}$ )后,通过反演即可得到地下各含水层的深度、厚度、单位体积含水量、孔隙特征等水文地质信息^[6]。目前的研究结果表明,尖峰噪声主要是由太阳磁暴、雷暴或物体的突然放电引起,针对这一类噪声的去除方法主要有中位数绝对偏差法、尖峰建模方法^[7]、能量运算方法^[8]等。对于随机噪声,一般采用多次重复观测然后叠加的方式,Lin T T等^[9]提出的时频峰值滤波和Li F等^[10]提出的最大似然估计都取得了一定的效果。谐波噪声的来源主要为发电设备、输配电设备以及用电设备的非线性负载等,其在MRS探测中影响范围最广,影响程度最大^[11]。Legchenko等^[12]基于谐波的频率周期性提出陷波滤波方法,然而该方法在拉莫尔频率与谐波频率相近时会导致NMR信号失真。Dalgaard等^[13]提出的自适应算法与Larsen等^[14]采用的多道维纳算法都是利用谐波在主通道和参考道中的相关性进行去噪,但是,在实际工作中相关性与多道线圈的摆放位置以及噪声源的个数等因素有很大关系。Li F等^[15]进一步提出了基于分帧运算的参考道噪声压制方法,提高了多道维纳滤波在变化谐波噪声中的压制效果。Larsen等^[14]提出了模型去噪方法,该方法构建谐波模型,利用最小二乘法匹配出最接近实测含噪信号的模型参数,进而从含噪信号中减去构建的谐波模型。该方法在谐波基频单一且稳定的条件下可以取得较好的去噪效果,然而,在多源谐波情况下去噪效果会明显受到制约。Kremer等^[16]证明针对多源谐波噪声连续使用模型去噪的方法,计算时间成本较大且只有在谐波基频频率差很小或很大时才有一定的效果。为更好地压制实际MRS方法采集的信号中常出现的多源谐波噪声,本文提出了一种基于模拟退火法的压制方法,并验证了其有效性和优越性。

1 MRS方法的基本原理及多源谐波电磁噪声特征

1.1 MRS方法的基本原理

MRS方法以NMR理论为基础,即地下水中的氢核磁矩受到一定频率的交变电流脉冲激发而发生偏转,磁矩弛豫过程中产生核磁共振信号,位于地面的接收线圈接收到一个自由感应衰减信号,即NMR信号,其过程如图1所示。图中:0~t₁为仪器激发前的纯噪声采集时间;t₁~t₂为仪器发射拉莫尔频率的激发脉冲持续时间,脉冲包络线为矩形;t₂~t₃为发射转接收时间,称为死区时间;t₃~t₄为接收NMR信号时间,NMR信号的包络线按指数规律衰减。NMR信号E_NMR为:

(1)

E_{NMR} = E_{0} \exp (- t / T_{2}^{*}) \cos (2 π f_{L} t + φ_{0}),

式中:E₀为NMR信号的初始振幅,其大小与地下含水量有关; $T_{2}^{*}$ 为自旋—自旋弛豫时间(通常称为平均弛豫时间),与含水层的孔隙特征密切相关;f_L为拉莫尔频率,其与地磁场强度有关;φ₀为NMR信号的初始相位,是二次场相对激发电流的相位移,与地下岩层的导电性有关。

图1

新窗口打开| 下载原图ZIP| 生成PPT

图1 MRS方法激发/接收NMR信号时序图

Fig.1 Transmitting and receiving NMR signal in time domain of MRS

1.2 干扰噪声的种类及其特征

1.2.1 噪声种类

MRS方法接收到的NMR信号通常非常微弱,一般在10~1 000 nV之间,由于仪器接收灵敏度高易受到外部环境电磁噪声的干扰,故接收到的信号信噪比通常较低,提取NMR信号的难度较大^[17]。按照噪声成分及来源,电磁噪声分为尖峰噪声、谐波噪声和随机噪声等,可以将接收到的NMR信号和电磁噪声表示为

(2)

E_{obs} = E_{N MR} + E_{spike} + E_{harmonic} + E_{random},

式中:E_obs为仪器接收到的信号数据,E_NMR为磁共振信号,E_spike为尖峰噪声,E_harmonic为谐波噪声,E_random为随机噪声。

1.2.2 多源谐波电磁噪声特征

在经典谐波噪声模型条件下,2013年Larsen等^[14]提出假设谐波噪声的基频在信号采集持续期间内是单一且稳定的,由此写出如下谐波模型式:

(3)

E_{harmonic} (f_{0}, n) = \sum_{k} A_{k} \cos (2 π k \frac{f_{0}}{f_{s}} n + ϕ_{k}),

式中:E_harmonic(f₀,n)表示采样时间n处的谐波模型的振幅,A_k、ϕ_k分别为第k次谐波分量的振幅和相位,f_s为采样频率,f₀为谐波噪声的基频。

随着MRS勘查工作逐步向城市靠近,面对的电磁噪声情况也愈发复杂,在测量区域可能存在2个或2个以上的谐波干扰源。本文分析存在2个谐波源情形下的谐波噪声特征,多源谐波与之相通^[16]。

假设有2个基频不同的谐波干扰源,且在测量时间内不变,可以写出2个谐波源的谐波模型式:

(4)

\begin{array}{l} E_{harmonic} (f_{1}, f_{2}, n) = \sum_{k} A_{k} \cos (2 π k \frac{f_{1}}{f_{s}} n + ϕ_{k}) + \\ \sum_{k} B_{k} \cos (2 π k \frac{f_{2}}{f_{s}} n + φ_{k}) 。 \end{array}

式中:f₁和f₂分别表示要搜索的2个不同的基频,A_k和B_k分别表示每个基频的第k次谐波的振幅,ϕ_k和φ_k分别表示每个基频的第k次谐波的相位。

含有双基频的谐波噪声频谱相对于单一且稳定基频的谐波噪声频谱会出现2种特殊的峰值频谱特征:“弯曲”和“双峰”(图2),“双峰”峰值频谱特征出现的条件较“弯曲”峰值频谱要更为严苛,要求频率间隔大且幅值相近,因此“双峰”峰值频谱的出现概率要低于“弯曲”峰值频谱。

图2

新窗口打开| 下载原图ZIP| 生成PPT

图2 双基频的谐波噪声频谱特征示意

Fig.2 Spectrum characteristics of harmonic noise with double fundamental frequencies

2 同步删除法和模拟退火法原理

2.1 同步删除法原理

由双基频谐波噪声的特征可知^[17],主要问题是谐波的频率彼此接近但不完全相同,使得其频率内容重叠但不完全一致。为解决双基频谐波噪声的相互干扰问题,可以拓展模型去噪的搜索方式,在一个二维频率空间内同步搜索2个基频,求取观测信号与双基频谐波模型的最小值:

(5)

‖ E_{obs} - E_{harmonic} (f_{1}, f_{2}, n) ‖ \to \min 。

由于式(4)中的f₁、f₂、A_k、B_k、ϕ_k和φ_k的确定是一个非线性优化问题,可以通过标准化改写式(4)中的余弦项:

(6)

\{\begin{array}{l} A_{k} \cos (2 π k \frac{f_{1}}{f_{s}} n + ϕ_{k}) = α_{k} \cos (2 π k \frac{f_{1}}{f_{s}} n) + \\ β_{k} \sin (2 π k \frac{f_{1}}{f_{s}} n), \\ B_{k} \cos (2 π k \frac{f_{2}}{f_{s}} n + φ_{k}) = γ_{k} \cos (2 π k \frac{f_{2}}{f_{s}} n) + \\ δ_{k} \sin (2 π k \frac{f_{2}}{f_{s}} n) 。 \end{array}

式中变量α_k、β_k、γ_k和δ_k可通过以下关系与对应谐波源的振幅和相位相关联:

(7)

\{\begin{array}{l} A_{k} = \sqrt[]{α_{k}^{2} + β_{k}^{2}}, \\ B_{k} = \sqrt[]{γ_{k}^{2} + δ_{k}^{2}} 。 \end{array}

(8)

\{\begin{array}{l} \tan ϕ_{k} = - β_{k} / α_{k}, \\ \tan ϕ_{k} = - δ_{k} / γ_{k} 。 \end{array}

将式(6)中的余弦项展开变换,假设f₁和f₂已知,观测信号与谐波模型的差值最小值问题就变成线性,可以用矩阵表示法写成:

(9)

\begin{array}{l} [\begin{array}{l} \cos (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{0}) & \sin (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{0}) & \cos (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{0}) & \sin (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{0}) & \dots & \cos (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{0}) \\ ︙ & ︙ & ︙ & ⋱ & ︙ & ︙ \\ \cos (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{f}) & \sin (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{f}) & \cos (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{f}) & \sin (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{f}) & \dots & \cos (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{f}) \end{array} \\ \begin{array}{l} \sin (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{0}) & \cos (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{0}) & \sin (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{0}) \\ ︙ & ︙ & ︙ \\ \sin (k_{0} 2 π \frac{f_{1}}{f_{s}} n_{f}) & \cos (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{f}) & \sin (k_{0} 2 π \frac{f_{2}}{f_{s}} n_{f}) \end{array}] \cdot [\begin{array}{l} α_{k_{0}} \\ β_{k_{0}} \\ γ_{k_{0}} \\ δ_{k_{0}} \\ ︙ \\ α_{k_{N}} \\ β_{k_{N}} \\ γ_{k_{N}} \\ δ_{k_{N}} \end{array}] = [\begin{array}{l} E_{obs} (n_{0}) \\ E_{obs} (n_{1}) \\ ︙ \\ E_{obs} (n_{f - 1}) \\ E_{obs} (n_{f}) \end{array}], \end{array}

式中:[k₀,k₁,…,k_N]表示构建谐波模型的阶数,[n₀,n₁,…,n_f]为采样时间样本。同时,可将式(9)整理成矩阵形式:

(10)

Ax = b,

式中:x=[ $α_{k_{0}}$ , $β_{k_{0}}$ , $γ_{k_{0}}$ , $δ_{k_{0}}$ ,…, $α_{k_{N}}$ , $β_{k_{N}}$ , $γ_{k_{N}}$ , $δ_{k_{N}}$ ]^T,b=[E_obs(n₀),E_obs(n₁),…,E_obs(n_f_-1),E_obs(n_f)]^T。利用最小二乘法求解该标准线性方程可得:

(11)

x = (A^{T} {A)}^{- 1} A^{T} b 。

每次代入1对假设基频值f₁和f₂,可以得到1组谐波参数。利用式(5)寻找1组能使得该差值最小的假设的基频值作为构建谐波模型的基频值,再利用测量数据E_obs减去构建的E_harmonic,即可达到同步压制2个干扰源谐波噪声的效果。

2.2 模拟退火法原理

模拟退火法的发展受到金属物理退火的启发^[18]:当金属从红热状态缓慢冷却时,金属内部形成有序的最小能量晶体结构。最初,当金属处于红热状态时,原子的运动完全由随机的热涨落控制,但是,随着温度的缓慢降低,原子间的作用力变得越来越重要,最后,原子变成了一个代表最小能量结构的晶格^[19]。在模拟退火算法中,参数T代表模拟温度,误差E代表模拟能量,T值很大时算法的行为类似于Monte Carlo搜索,T值很小时算法的行为更有方向性。在物理退火中,一开始是一个大值的T,然后随着越来越多的试验被检验,T值慢慢降低;最初,大量的模型空间是随机采样的,但是随着搜索的进行,搜索变得越来越有方向性^[20]。

模拟退火法从一个初始模型m⁽^p⁾和其对应的误差E(m⁽^p⁾)开始。然后生成一个在m⁽^p⁾附近的测试模型m⁽^*⁾,计算其对应的误差E(m⁽^*⁾),m⁽^*⁾可以通过向m⁽^p⁾增加一个从高斯分布中提取的增量Δm。当E(m⁽^*⁾)≤E(m⁽^p⁾)时,m⁽^*⁾会替代m⁽^p⁾成为新的比较值。但是,有时E(m⁽^*⁾)>E(m⁽^p⁾),m⁽^*⁾也可以接受。为了决定后面的情况,这个测试参数t可通过

(12)

t = \frac{e^{- E (m^{*}) / T}}{e^{- E (m^{p}) / T}} = e^{- \frac{E (m^{*}) - E (m^{p})}{T}}

计算。在区间[0,1]生成一个均匀分布的随机数r,如果t>r,则接受m⁽^*⁾。当T值很大时,测试参数t趋近于1,故不管误差值是多少,m⁽^*⁾几乎总是被接受的,这对应着“热运动”情况下,模型参数在空间上以无向的方式探索;当T值很小时,测试参数t趋近于0,m⁽^*⁾几乎不被接受,这与定向搜索情况相对应,只有使误差E减小的测试模型才被接受。

3 仿真数据验证

3.1 网格搜索同步删除法

前文论述过,若要式(5)成立需要多次求解线性问题,同时优化二维网格搜索,每个维度对应不同的干扰源谐波基频。为了加快搜索速度,可以将搜索分2步进行:首先,覆盖相对宽且粗糙的网格,每次计算模型减去后的剩余信号;然后,从在第一步中获得的最小值开始探索更细的网格,并且通过第二步获得的最小剩余信号产生最佳频率估计值。

仿真条件:一号谐波基频f₁=49.985 Hz,二号谐波基频f₂=50.032 Hz,谐波从2 000 Hz开始,共10次谐波,每次谐波振幅A_k随机在0~100 nV,平均幅值50 nV,每次谐波相位ϕ_k相互独立且服从 U(-π,π)分布;初始振幅E₀=200 nV;拉莫尔频率f_L=2 260 Hz;初始相位φ₀=π/6 rad;平均衰减时间 $T_{2}^{*}$ =300 ms;采样频率f_s=19 200 Hz;采样时间t=1 s;随机噪声为5 nV的高斯白噪。

图3a是在一个相对粗糙的网格内搜索2个基频值,即图中2个相互对称的深蓝色色块;图3b是选取图3a中的一个色块范围进一步精细化搜索2个基频,可以明显看出找到的基频值与仿真设定的一致。搜索基频频率的精度要达到0.001 Hz需要进行600多次的搜索,处理时间为60.36±1.27 s。图4显示,同步删除法能够压制双基频谐波噪声,经过处理的信号与真实NMR信号相似程度较高,呈指数衰减趋势。图中纵坐标变量A表示振幅。

图3

新窗口打开| 下载原图ZIP| 生成PPT

图3 频率空间谐波基频网格搜索

Fig.3 Harmonic fundamental frequency grid search in frequency space

图4

新窗口打开| 下载原图ZIP| 生成PPT

图4 网格搜索同步删除法处理双基频谐波时域

Fig.4 Time domain diagram of double fundamental frequency harmonics processed by grid search simultaneous removal method

3.2 模拟退火同步删除法

为了加快基频的搜索和避免陷入局部极值点,利用模拟退火法来搜索干扰源谐波的基频值。仿真条件同上,迭代次数200次,搜索路径如图5所示。图5中黑色圆圈代表初始模型(50.000;50.000),白色圆圈代表每次迭代被接受的模型,绿色圆圈代表迭代结束最终模型(50.031 7;49.984 6),红色线段代表基频搜索路径。可以看出,每次被接受的模型都在逐步逼近真实模型,虽然模拟退火法最终能找到最小值,但在有限的迭代次数时搜索路径会存在一定的偏差。为了更好地评价去噪效果,特意引入均方根误差(root mean square error,RMSE)这个指标参数^[21]:

(13)

RMSE = \sqrt[]{\frac{1}{N} \overset{N}{\sum_{i = 1}} [E_{denoised} (i) - E_{NMR} {(i)]}^{2}},

式中:N为仪器采集到的信号长度,E_denoised(i)是去噪后的第i个数据,E_NMR(i)是真实NMR信号的第i个数据。

图5

新窗口打开| 下载原图ZIP| 生成PPT

图5 模拟退火同步删除法的基频搜索路径

Fig.5 Fundamental frequency search path of simulated annealing simulated annealing

为了选取更适合的迭代次数,特意设置不同的迭代次数情况下搜索的效果,结果见表1。由表可见:随着迭代次数的增加,RMSE和标准差均逐步变小,拟合参数E₀和 $T_{2}^{*}$ 更接近真实值,处理时间逐步增加;当迭代次数超过100次时,迭代次数的增加对压制效果提升增长不明显。为了兼顾处理时间和压制效果,建议设置迭代次数为100次。在同一双基频谐波仿真情况下,相较网格搜索,同步删除法的处理时间缩短了一半以上,计算效率提高2.35倍。

表1 不同迭代次数的模拟退火同步删除法去噪效果

Table 1 Denoising effect of simulated annealing simultaneous removal method with different number of iterations

迭代次数	RMSE/nV	E₀/nV	$T_{2}^{*}$ /ms	t/s
10	67.16±22.96	177.12±4.01	1038.44±522.84	3.19±0.29
30	37.89±23.20	184.84±8.65	489.32±175.95	8.89±0.68
50	26.72±21.20	186.83±8.86	465.99±78.27	13.76±0.74
75	22.42±5.56	188.93±3.68	366.05±28.97	20.05±0.86
100	15.69±3.77	193.06±2.04	334.47±13.16	25.68±0.88
150	14.60±3.71	193.64±2.25	330.44±11.83	38.45±1.78
200	13.14±3.54	194.33±1.71	325.34±11.36	50.40±1.89
500	11.46±2.61	195.10±1.75	320.59±7.54	123.76±4.18

新窗口打开| 下载CSV

复杂的电磁噪声可能存在3个或更多基频的多源谐波噪声。为了验证模拟退火同步删除法能够适用于多源谐波电磁噪声的情况,选取了3个干扰源进行仿真,谐波基频分别为49.950、50.010、50.050,迭代次数300次,初始基频模型(50.000;50.000;50.000),其他仿真条件同上,图6为仿真结果。从图6a可以看出模拟退火同步删除法在3个基频谐波的情况下能够较好地压制多源谐波噪声,去噪后的蓝色曲线与真实NMR信号的红色曲线基本吻合。图6b中黑色圆圈显示了迭代初始模型的频率空间位置,蓝色圆圈是300次迭代结束最终模型的频率空间位置(50.009 2;49.951 2;50.049 9),红色圆圈是3个真实谐波基频不同排列顺序的频率空间位置,灰色曲线代表基频搜索路径;经过一定迭代次数的搜索,最终结果会归于6个红色圆圈中任意一个的附近。

图6

新窗口打开| 下载原图ZIP| 生成PPT

图6 模拟退火同步删除法压制多源谐波噪声

Fig.6 Simulated annealing simultaneous removal method suppressing Multi- source harmonic noise

4 实测数据检验

为了验证本文提出的噪声压制方法在实测数据中的应用效果,选择天门市多宝镇汉江边的一处场地进行实验。该场地地下含水层已由场地内的钻孔所揭露,场地附近有多个人居村庄,符合多源谐波的噪声特点,实测数据可以很好地应用于去噪处理的实验工作中。数据采集使用法国IRIS公司的NUMIS^poly仪器。本次测量采集技术参数:方形线圈边长100 m;拉莫尔频率f_L=2 135.2 Hz;激发脉冲矩范围为0.12~10.866 A·s,分为16个;采样率19 200 Hz;采样时间1 s;死区时间0.04 s。

通过分析野外实际工作环境和噪声,发现2个较强的谐波干扰源。根据数据处理流程对噪声数据进行了尖峰噪声的处理后,采用100次迭代的模拟退火同步删除法来压制双谐波噪声(图7),再将去噪后的信号进行叠加、希尔伯特变换和非线性拟合,结果见图8。

图7

新窗口打开| 下载原图ZIP| 生成PPT

图7 模拟退火同步删除法在实测数据处理中的效果

Fig.7 Effect of simulated annealing simultaneous removal method in measured data processing

图8

新窗口打开| 下载原图ZIP| 生成PPT

图8 去噪后各脉冲矩的包络线及其拟合曲线

Fig.8 Envelope and fitting curve of each pulse moment after denoising

图7是利用模拟退火同步删除法来压制实测双谐波噪声。图7a中经过模拟退火同步删除法压制后的去噪曲线(红色曲线)有明显的指数衰减趋势,仅剩余随机噪声的干扰;图7b中相对于实测含噪信号(灰色曲线),去噪曲线在50 Hz整数附近的峰值被较好压制,拉莫尔频率的信号保存较为完整。

图8中灰色曲线是去噪后的信号经过叠加处理和希尔伯特变换得到的包络曲线,红色曲线是灰色曲线经过非线性拟合得到的拟合曲线。包络曲线和拟合曲线重合性较好,表明模拟退火同步删除法能够较好地压制多谐波噪声,也能看出残余谐波噪声对包络曲线的影响在可接受的范围,且对拟合曲线的参数提取影响较小。

通过对拟合曲线进行零时外延,可以获取NMR信号的初始振幅E₀和平均弛豫时间 $T_{2}^{*}$ ,再进行反演与已知钻孔进行对比得到图9。图9a是不同脉冲矩对应的初始振幅的关系曲线,即MRS方法一个测点的测深曲线,该曲线平滑无异常突起,说明谐波噪声被较好压制;图9b是不同脉冲矩对应的平均弛豫时间曲线,曲线平滑无异常突起。较为准确的初始振幅曲线和平均弛豫时间曲线,能够为后期反演地下真实的含水量和渗透率提供有力保证。利用根据测量参数构建的核函数和初始振幅曲线、平均弛豫时间曲线一起反演得到图9c,与该测点钻孔资料(图9d)有较好的对应。该测点表层为冲击砂层,有一定的植被,地下0~7.4 m含水量很低,对应位于潜水面上的砂质黏性土;7.4~15.9 m的含水层含水量约为17%,对应潜水面下的弱渗透性细砂层;15.9~48.3 m为一个较大的含水层,含水量最大为24%,平均18%,与地下渗透性较好的中砂层和中粗砂层相对应;48 m以下含水量很低,与隔水性的粉质黏土相对应。对比反演结果和钻孔资料可以证明,模拟退火同步删除法能够较好压制多源谐波噪声,改善地下含水层反演信息的准确度和稳定性。

图9

新窗口打开| 下载原图ZIP| 生成PPT

图9 天门河滩MRS方法某测深点实测数据及其反演结果与钻探结果对比

Fig.9 Measured data, inversion results and drilling results of a sounding point by MRS method in River Beach, Tianmen

5 结论

研究与实践表明,对于单一且规律性的谐波电磁噪声,利用模型去噪的方法能够有效地压制,而变化、多源的谐波噪声的压制是磁共振测深方法信号处理中的难点和研究热点,压制效果的好坏直接影响后续反演解释工作的准确性。根据多源谐波电磁噪声的特点出发,验证了网格搜索同步删除法的可行性,为了节约计算时间成本,提出了模拟退火同步删除法这一新的磁共振测深方法信号处理手段,通过仿真数据验证和实测数据检验,本文论述的信号处理手段取得了较好的噪声压制效果。

1)针对磁共振测深信号中模型去噪方法无法有效压制多源谐波噪声的问题,利用网格搜索同步删除法和模拟退火同步删除法均能有效找到较为准确的基频值,达到压制多源谐波的目的。

2)通过模拟仿真实验,证明网格搜索同步删除法和模拟退火同步删除法均能压制多谐波噪声。在双基频谐波情况下比较不同迭代次数下的压制效果,模拟退火同步删除法较网格搜索同步删除法的计算效率提高了2.35倍,大大降低了时间成本,且在多源谐波噪声压制方面也能有较好的效果,为开展更加复杂的谐波情况的压制提供了可能。

3)通过野外MRS方法的实测资料证明了模拟退火同步删除法压制多源谐波噪声的有效性,反演得到的含水量与钻探结果基本相符,极大地拓展了磁共振测深方法在人居区等电磁情况复杂地区作业的有效性,也为多通道磁共振多维探测提供新的数据处理思路。

仿真和实测结果表明,模拟退火同步删除法能够有效压制磁共振测深方法中的多源谐波噪声,为后期反演资料的解释提供了技术保证和选择方法。本文主要讨论单通道情况下多源谐波噪声的压制,未来将继续研究二维\三维多通道变化谐波噪声的压制方法。

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

Behroozmand A

, Keating

, Auken

A review of the principles and applications of the nmr technique for near-surface characterization

[J]. Surveys in Geophysics, 2014,36(1):27-85.

DOI:10.1007/s10712-014-9304-0 URL [本文引用: 1]

[2]

, Li Z

, Niu R

, et al.

Using surface nuclear magnetic resonance and spontaneous potential to investigate the source of water seepage in the Jindeng Temple grottoes, China

[J]. Journal of Cultural Heritage, 2020,45(1):42-51.