E-mail Alert Rss
 

物探与化探, 2019, 43(3): 589-597 doi: 10.11720/wtyht.2019.0043

方法研究·仪器研制

磁异常和梯度的频率域三维成像方法

崔亚彤, 郭良辉

中国地质大学(北京) 地球物理与信息技术学院,北京 100083

3-D imaging of magnetic anomalies and gradients in the frequency domain

CUI Ya-Tong, GUO Liang-Hui

School of Geophysics and Information Technology, China University of Geosciences (Beijing), Beijing 100083, China

责任编辑: 王萌

收稿日期: 2019-01-8   修回日期: 2019-03-12   网络出版日期: 2019-06-20

基金资助: 国家自然科学基金面上项目.  41774098
中央高校基本科研业务费专项资金.  2652018266

Received: 2019-01-8   Revised: 2019-03-12   Online: 2019-06-20

作者简介 About authors

崔亚彤(1993-),女,博士在读,主要研究方向为重磁数据处理和反演算法。Email:3010170005@cugb.edu.cn 。

摘要

三维反演在磁数据定量解释中具有重要作用。常用的空间域三维反演方法通常需要大量的正演和反演计算,因此对大规模数据的反演效率较低。三维成像是另一种定性和定量解释磁数据的重要方法。文中给出了一种磁异常与梯度三维成像的频率域迭代方法,该方法可以提高成像效率,适用于大规模数据的快速成像。笔者推导了磁总场异常和异常梯度频率域正演公式和成像公式,并将一种深度尺度因子引入成像公式中,提高了深度精度;笔者采用了迭代优化算法,减小了拟合误差,进一步提高了成像精度。通过理论模型数据试验和中国新疆某金属矿床实测数据,验证了本文方法的有效性、可行性。

关键词： 磁异常 ; 梯度 ; 三维成像 ; 正演 ; 频率域 ; 深度尺度因子

Abstract

3-D inversion plays an important role in the quantitative interpretation of magnetic data. However, the commonly used space-domain 3-D inversion algorithms usually require a large number of forward modeling and inversion calculations. Hence, the inversion based on a large-scale data is usually inefficient. 3-D imaging is another significant algorithm for the qualitative and quantitative interpretation of magnetic data. This paper implements a frequency-domain iterative approach for 3-D imaging of magnetic anomalies and gradients, which can improve imaging efficiency and is suitable for rapid imaging of large-scale data. The frequency-domain forward formulae and imaging formulae of magnetic total field anomaly and magnetic gradients are derived in this paper. A depth scaling factor is added to the imaging formulae to significantly improve the depth resolution. In order to reduce the fitting error and improve the imaging accuracy, this paper adopts an iterative optimization algorithm. The effectiveness and feasibility of the presented approach were verified by the synthetic data and real data from a metallic deposit area in Xinjiang.

Keywords： magnetic anomalies ; gradients ; 3-D imaging ; forward modeling ; frequency domain ; depth scaling factor

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

本文引用格式

崔亚彤, 郭良辉. 磁异常和梯度的频率域三维成像方法. 物探与化探[J], 2019, 43(3): 589-597 doi:10.11720/wtyht.2019.0043

CUI Ya-Tong, GUO Liang-Hui. 3-D imaging of magnetic anomalies and gradients in the frequency domain. Geophysical and Geochemical Exploration[J], 2019, 43(3): 589-597 doi:10.11720/wtyht.2019.0043

0 引言

磁力勘探方法是地球物理方法的主要分支之一,长期在资源勘查、工程勘察与区域构造研究中发挥着重要作用。传统的磁力勘探是测量总磁场或3个分量。近年来,磁梯度因其较高的分辨率,在磁力勘探领域拥有了举足轻重的地位。三维反演是对磁数据进行定量解释的重要工具,它可以反演出地下半空间磁化强度或磁化率的三维分布情况,从而圈定深部地质体的赋存状态[1,2,3,4,5]。该方法假设观测面为平面,地下模型即可分为多个水平层,每个水平层由一组大小相同、物理性质各异的规整排列的直立长方体组成。在空间域中,每个直立长方体在观测面上任意一点产生的磁场,可以应用简单的数学解析式进行正演计算得到[6,7],而观测面上每一个测点的理论磁场值则是由正演计算得到的地下所有长方体在该点引起的磁场的累加求和。而反演问题的求解通常采用基于Tikhonov正则化算法[8]的线性或非线性反演理论来实现,以达到观测磁场值与理论磁场值之间偏差的最小化,进而反演估计地下异常体的物理性质及其赋存状态。

然而,三维磁反演常常是一个病态的、不稳定的问题,而且反演结果的深度精度低,还会产生趋肤效应,不符合实际地质关系。因此,很多学者针对上述问题提出了相关约束算法和理论,如深度加权约束算法[1]、物性约束算法[9,10,11]、聚焦反演理论[12,13]、先验地质信息约束算法[14,15,16]、结构倾向约束算法[17]等等。传统三维反演方法另一个很重要的问题是其计算量大、计算机内存占用量大。因此,在处理大规模数据体时,灵敏度矩阵庞大,需要占用大量的计算机内存,使得计算效率极低。为此,很多学者给出了相关快速算法,如数据压缩算法[9, 13, 18]、随机子空间快速算法[19]、随机奇异值分解算法(RSVD)[20]等等。虽然基于快速算法的反演效率显著提高,但这些算法在迭代过程中仍需要进行大量的正演和反演计算,甚至对于大规模数据的反演也有难以克服的困难。

频率域方法因其计算简单、快速,在地球物理数据处理和反演中得到了广泛的应用。Parker[21]提出了一种获取由非均匀层引起的重磁异常的频率域方法。Oldenburg[22]在Parker[21]提出的理论基础上,给出了基于频率域的重力异常迭代反演算法。Cribb[23]提出了磁偶极子频率域的广义线性成像方法,并且其过程与频率域向上延拓算法相类似。Kobrunov和Varfolomeev[24]提出了一种基于重力异常的频率域成像反褶积方法。Pedersen[25]对磁场数据的功率谱进行分析,进而估计地下半空间磁化强度分布情况。根据Kobrunov和Varfolomeev[24]的理论,Priezzhev[26,27]简化了重力异常频率域的反褶积滤波算子,同时该算子可以降低吉布斯效应,使模型得到光滑解,根据泊松方程,其反褶积滤波算子可以应用于磁异常成像算法中。Kobrunov[28]进一步给出了显式频率域成像公式及其较为复杂的迭代算法。上述频率域成像方法可以快速获得地下三维半空间的等效物性分布结果,因此,它们属于三维成像方法,而不是传统意义上的三维反演方法。然而,针对磁异常及其梯度数据成像方法,上述方法很少涉及,或是简而言之,也很少提及相关约束方法及迭代优化算法。

笔者基于Priezzhev的频率域成像原理[26,27],给出了一种磁异常与梯度三维成像的频率域迭代方法:首先,推导了磁总场异常和梯度的正演和成像的频率域公式,同时,在成像过程中引入了一种深度尺度因子来提高深度精度;为了进一步提高成像精度,采用了迭代优化算法;最后,结合理论模型数据和中国新疆某金属矿区实测数据,对本文方法进行了验证。

1 磁异常和梯度的频率域三维成像方法

1.1 磁异常与梯度的频率域三维正演公式

假设观测面为平面,则三维网格模型可以从上到下划分为多个水平层。假设其中一个水平层的上界面深度为z1,下界面深度为z2,z2>z1。因此,在频率域中,该水平层引起的磁位U(p0)为[6]:

F2⁃D[U(p0)]=μ02|k|2ψme|k|z0(e-|k|z1-e-|k|z2)·F2⁃D[M(p)],

其中:F2-D[·]代表二维FFT变换;p0=(x0,y0,z0)为观测点位置,p=(x,y,z)代表场源内任意点;μ0=4π×10-7 Henry/m为真空磁导率;k=(kx,ky)是x,y 方向的波数;|k|= kx2+ky2为径向波数;z0是观测面高度,且有z2>z1,z1>z0;ψm是磁化强度方向单位矢量 M˙=(Mx,My,Mz)与波数k的函数,ψm=(iMxkx+iMyky+Mz|k|);M(p)为磁化强度模量。为了方便表达,因此将式(1)简化成[6]

F2-D[U(p0)]=F2-D[ξ]·F2-D[M(p)],其中:F2-D[ξ]=μ02|k|2ψme|k|z0(e-|k|z1-e-|k|z2)。

根据磁位与磁总场异常的函数关系,可以得到水平层引起的磁总场异常的频率域公式[6]:

F2-D[ΔT(p0)]=HΔT·F2-D[M(P)],

其中:HΔT为磁总场异常频率域正演滤波算子;ψf是地磁场方向单位矢量 F˙=(Fx,Fy,Fz)与波数k的函数,ψf=(iFxkx+iFyky+Fz|k|),即:

HΔT=ψf·F2-D[ξ]。

令α=(x,y,z), 根据傅里叶求导原理和位场理论,可以分别求得x、y和z方向磁异常梯度的频率域公式:

F2-D[ΔTα(p0)]=HΔTα·F2-D[M(p)],
HΔTα=φα·ψfF2-D[ξ],φα=φx=ikxφy=ikyφz=|k|

其中, HΔTα为磁异常梯度频率域正演滤波算子。因此,单一水平层的各参量异常频率域正演公式可以很容易地通过磁化强度频谱与相应参量的正演滤波算子相乘得到,如式(3)、(5)。然后,对上述异常频谱进行二维快速傅里叶逆变换(IFFT),即可在空间域中获得每个水平层的二维异常。对于由多个水平层组成的三维模型,通过对每个水平层产生的异常进行简单的累加求和,最终可计算出该三维模型在平面观测面上任意一点所产生的异常。在实际应用中,可以在模型x、y方向的边缘外补充多个点再进行正演计算,计算结束后再消去相应的补充点,这样可以压制由FFT产成的混叠效应、边缘效应问题对正演计算的影响。

1.2 磁异常与梯度的频率域三维成像公式

假设观测面为平面,则磁场U(p0)与磁源内一点的函数关系为[6]

$U(p_{0})=\frac{1}{4 \pi}\int_{R^{3}} M(p)\hat{M}·▽\lgroup \frac{1}{|r|} \rgroup dp$,

其中:R3将积分函数定义在一个有限的三维空间; ▽为场源点p的梯度表示符号,r为观测点p0到场源点p的矢量,r=p0-p。式(7)可以由三重积分化为二重积分:一重是沿x,y平面积分,另一重是沿z方向积分。因此,式(7)可进一步推导为

$U(p_{0})=\frac{1}{4 \pi}\int_{S} dS \int_{z=z_{0}}^{∞} M(p)\hat{M}·▽\lgroup \frac{1}{|p_{0}-p|} \rgroup dz$,

其中:S将积分函数定义在有限的二维空间。对式(8)两端同时进行二维FFT变换可得

F2-D[U(p0)]=ψm2|k|∫z=z0∞F2-D[M(p)]e-|k|zdz,

为了方便后续演算,将式(9)表示为

Ũ(kx,ky,z0)=ψm2|k|∫z=z0∞M̃(kx,ky,z)e-|k|zdz,

式(10)为三维磁化强度分布引起的磁势场频率域表达公式。根据磁势场与磁总场异常的频率域函数关系,即Δ T̃(kx,ky,z0)=μ0ψfŨ(kx,ky,z0),则由一个三维异常体引起的磁总场异常的频率域表达公式为

ΔT̃(kx,ky,z0)=μ0ψmψf2|k|∫z=z0∞M̃(kx,ky,z)e-|k|zdz。

基于Kobrunov和Varfolomeev[24]提出的理论,可以由式(11)推导出各水平层磁化强度的频率域成像公式

M̃j(kx,ky,zj)=2|k|μ0ψmψf·hj(kx,ky,zj)·ΔT̃(kx,ky,z0),
hjkx,ky,zj=Kkx,ky,zj∫z=z0∞Kkx,ky,ze-kzdz,

其中:hj(kx,ky,zj)为单一水平层的深度反褶积滤波算子,j=1,2,3,…;N为沿z方向各个水平层的顺序编号, 为水平层的总个数;zj为每个水平层j的掩埋深度;K(kx,ky,zj)是与成像深度相关的函数。Priezzhev[26,27]将式(13)简化为

hj(kx,ky,zj)=(n+1)n+1n!(zj)n|k|n+1e-n|k|zj,

其中:n为垂向导数幂次,且为正整数,n∈[1,10]。使用上述深度反褶积滤波算子可以有效地削弱吉普斯现象,使成像结果更加平滑。然而,仅仅使用上述深度反褶积滤波算子得到的成像结果往往容易出现明显的趋肤效应,使深度精度降低。为了压制趋肤效应,提高深度精度,本文在深度反褶积滤波算子中引入一种深度尺度因子λ,即,

hmod_j(kx,ky,zj)=(n+1)n+1n!zjλn|k|n+1e-n|k|zjλ,

其中:λ为深度尺度因子,λ≥1。式(15)为改进的深度反褶积滤波算子。基于先验信息(如已知地质和钻井资料,或其他地球物理资料),λ最优值可以通过试错法来获得,例如,将不同的λ值带入式(16)中,可以得到一系列不同的磁化强度三维成像模型;再根据先验信息或其他深度估计方法,得到异常体的大致埋深情况,然后将其与三维成像结果对比,即可选择最佳深度尺度因子λ。

M̃jkx,ky,zj=HΔT˙jkx,ky,zj·ΔT˙kx,ky,z0。

式(16)即为各水平层磁化强度改进的频率域磁总场异常成像公式。其中 HΔT˙_j(kx,ky,zj)= 2|k|μ0ψmψfhmod_j(kx,ky,zj)为磁总场异常频率域成像滤波算子。

基于傅立叶求导定理,可以得到磁总场异常3个方向梯度的频率域成像公式

M̃j(kx,ky,zj)=HΔT˙α_j(kx,ky,zj)·ΔT̃α(kx,ky,z0)。

其中: HΔT̃α_j(kx,ky,zj)= 2|k|μ0φαψmψfhmod_j(kx,ky,zj)为磁梯度频率域成像滤波算子,α=(x,y,z),φα与式(6)中表达方式一致: 综上所述,将各参量异常二维频谱与相应成像滤波算子相乘,即可得到各参量异常的频率域成像公式。因此,应用式(16)、(17),可以获得每一水平层的磁化强度的二维频谱分布情况,对上述磁化强度频谱进行IFFT变换,即可在空间域中获得每个水平层磁化强度的二维分布情况,然后,沿着z方向将获得的水平层磁化强度顺序排列,即可得到最终的三维磁化强度成像模型。然而,在实际应用中,观测面往往是起伏的,因此在进行成像流程之前,需要将观测面进行曲化平处理[29,30,31]。

1.3 迭代算法

为了进一步提高成像精度,可以将迭代的思想引入到频率域成像过程中。频率域迭代算法流程见图1。

图1

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

图1   频率域迭代算法流程

Fig.1   The flow chart for frequency-domain iterative approach


下述将以磁总场异常频率域三维成像为例,其他参数只需替换相应的频率域成像滤波算子即可。

步骤 1:对观测磁总场异常进行二维FFT变换,得到其二维异常频谱。

步骤 2:根据先验信息,对深度尺度因子进行筛选,选取最佳λ值。应用式(16)可计算成像出每个水平层的磁化强度二维频谱。对其进行IFFT变换,即可得到每个水平层的磁化强度分布情况。然后沿z方向顺序排列各个水平层,即可获得三维成像模型。此外,可以对三维成像模型进行物性范围硬约束。根据先验信息,如果研究区内的物性值在[m1,m2]范围之间,对于一个三维网格模型数据,遍历每个网格单元,当物性值小于m1时,修改其物性值,使其强制等于m1。同理,当物性值大于m2时,使其强制等于m2。进行上述步骤后,可获得基于物性约束的三维成像模型M1。

步骤 3:应用式(3)对上述三维成像模型M1进行正演计算,得到其理论磁总场异常F1。

步骤 4:计算理论异常F1与观测异常F0的偏差,即δF1=F0-F1。

步骤 5:如果偏差δF1的均方差(RMS)没有达到规定阈值,则继续对δF1进行三维成像计算,得到三维偏差模型δM1。随后更新三维成像模型,即M2=M1+δM1,并对其进行物性范围硬约束,进入下一次迭代计算。

步骤 6:重复步骤3~5,直到RMS达到规定阈值。结束迭代计算后,即可获得最终三维成像模型。

2 理论模型数据试验

假设一个理论模型由2个大小不同、掩埋深度不同的直立长方体组成(表1,图2a)。该模型是一个规则的三维网格模型,网格大小为101×101×31,x,y,z方向的间距均为50 m ;假设磁倾角为30°,偏角为45°,且磁化方向与地磁场方向一致;观测面高度为0 m。根据频率域三维正演公式,上述三维模型引起的磁总场异常和其垂直梯度可分别由式(3)和(5)计算得到,并加入了2%的高斯白噪声。如图2可见,小而浅的长方体产生的磁场波长较短,频率高,而大而深的长方体产生的磁场波长较长,频率低。利用这两种含噪异常数据验证前述频率域迭代方法,根据已知信息,利用试错法可确定成像相关参数为:垂向导数幂次n=10,深度尺度因子λ=2.5。

表1   直立长方体的几何和物性参数

Table 1  The parameters about geometry and magnetization of two cuboids

编号东向长度/m北向长度/m厚度/m东向中心位置/m北向中心位置/m纵向中心位置/m磁化强度/(A·m-1)
A300300200-4000-3001
B7007007006000-11001.5

新窗口打开| 下载CSV


图2

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

图2   理论模型及其含噪数据

a—理论模型三维示意;b—含噪磁总场异常;c—含噪垂直梯度异常

Fig.2   Synthetic model and noisy data

a—the 3-D diagram of synthetic model;b—the noisy total magnetic anomaly;c—the noisy vertical gradient


首先,应用含噪磁总场异常对本文方法进行验证。图3给出了迭代1次和迭代10次的磁总场异常成像结果。迭代1次后得到两个不同深度和尺寸的异常体,在空间位置上都与模型的真实长方体近似,但尺寸略大,磁化强度相对较小,最大值约为0.195 A/m(图3a)。随着迭代次数的增加,成像结果逐渐改善,异常偏差逐渐下降。经过10次迭代后,成像得到的两个异常体在形状和尺寸上都更接近真实长方体模型,磁化强度提高到0.451 A/m(图3b)。图4a为10次迭代后的理论磁异常与原始含噪异常的偏差,图4b为磁总场异常成像收敛曲线,随着迭代次数的增加,异常均方根值(RMS)快速下降,最后逐渐趋于一个最小值。显然,成像模型的理论异常与原始含噪异常相近,偏差较小(-2.64~0.81 nT),RMS为0.73 nT。

图3

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

图3   基于本文频率域方法的含噪磁总场异常三维成像结果

a—迭代1次成像结果;b—迭代10次成像结果

Fig.3   The 3-D imaging results of the noisy total magnetic anomaly by using the presented frequency-domain approach

a—after one iteration;b—after ten iterations


图4

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

图4   成像结果误差评价

a—理论磁总场异常与实际磁总场异常的偏差;b—收敛曲线

Fig.4   The evaluation of the 3-D imaging result

a—the deviation between the theoretical noisy total magnetic anomaly and the observed noisy ones;b—the convergence curve


其次,应用含噪磁异常垂直梯度对本文频率域迭代方法进行验证。经过10次迭代后,磁垂直梯度三维成像的结果及其偏差图分别如图5a和5b所示。垂直梯度成像结果(图5a)与磁总场异常成像结果(图3b)十分相似。其成像结果显示浅层和深部两个异常体都非常接近实际模型的两个长方体的掩埋位置,最大磁化强度为0.456 A/m,RMS为0.0013 nT/m。

图5

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

图5   基于本文频率域方法的含噪磁异常垂直梯度成像结果及误差

a—迭代10次三维成像结果;b—理论垂直梯度异常与实际异常的偏差

Fig.5   The imaging result and deviation of the noisy vertical gradient by using the presented frequency-domain approach

a—the 3-D imaging result after ten iterations;b—the deviation between the theoretical noisy vertical gradient and the observed noisy ones


最后,将空间域反演方法(UBC-地球物理数据处理软件)与由MATLAB(R2017b)编译的频率域迭代方法进行运行时间对比。电脑的中央处理器(CPU)为Intel Xeon E5-2620,随机存取存储器(RAM)为16 GB。图6为基于3种不同网格大小的三维网格模型的成像时间对比,网格大小分别为51×51×21,101×101×21,201×201×21。很明显,本文方法的运行时间远低于空间域方法(UBC),这意味着本文方法的计算效率明显高于空间域方法。

图6

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

图6   本文频率域方法与空间域反演方法(UBC)的运行时间对比

Fig.6   The calculation time comparison between the frequency-domain approach and space-domain approach (UBC) based on three size of the 3-D regular grid model


3 实际数据试验

实测数据来自中国新疆某金属矿区[32]。研究区内结构自上而下分别为无磁性的第四纪松散堆积物、弱磁性或中等磁性(约(0~100)×10-3 SI)的泥盆纪凝灰岩地层和弱磁性(约(0~16)×10-3 SI)的泥盆纪石英云闪长岩地层。该研究区内的金属矿主要赋存于凝灰岩地层内部,包含磁铁矿、黄铁矿和黄铜矿。这些金属矿表现为强磁性,且磁化率约为(30~1 000)×10-3 SI。同时,在该区域内存在两口不同深度的钻井zk1和zk3(图7)。zk1钻孔实测资料显示,浅层为强磁性金属矿,深部为弱磁性蚀变凝灰岩地层。zk3钻孔实测资料显示,浅层为弱磁性蚀变凝灰岩地层,深部为弱磁性角页岩凝灰岩地层。

图7

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

图7   新疆某金属矿区实测数据

a—磁总场异常;b—异常垂直梯度;黑色虚线A-B为图8所示成像结果剖面位置,zk1和zk3为两口深度不同的钻井

Fig.7   The real data from a metallic deposit area in Xinjiang

a—the total magnetic anomalies;b—magnetic vertical gradient;black dashed line A-B locates the profile shown in Figure 8, and zk1 and zk3 locate two boreholes


研究区磁总场异常(图7a)和垂直梯度异常(图7b)数据均为91×98的规则网格数据,E向、N向的网格间距均为10 m。磁倾角和磁偏角分别为64.3°和-0.4°。整个区域呈现低磁异常,而在钻孔zk1周围区域出现一个小而闭合的高磁异常、高磁梯度异常,推断该异常很有可能是由金属矿引起。

利用本文给出的频率域迭代方法对磁总场异常和垂直梯度异常分别进行了三维成像计算,得到了地下半空间磁化率的三维分布情况。成像深度范围为-400~0 m,深度步长为10 m。成像相关参数为:垂向导数幂次n=3,深度尺度因子λ=1。图8显示了沿剖面A-B深度切片的磁总场异常和垂直梯度的三维成像结果。两种结果在zk1浅部周围均呈现强磁性分布,成像异常形状大致为透镜状。磁总场异常三维成像结果(图8a)最大磁化率约为0.84 SI,磁异常垂直梯度成像结果(图8b) 最大磁化率约为0.612 SI,这种两种成像结果所显示的强磁性异常体与zk1钻孔揭示的小尺度强磁性金属矿相吻合,并且这种强磁异常体掩埋位置与实际钻井数据显示的金属矿顶底界面轮廓线也具有很好的一致性。同时,在zk1和zk3钻井深部均呈现弱中磁性分布,最大磁化率约为0.16 SI左右,而这种弱中磁性分布与两口钻井揭示的深部泥盆纪大型凝灰岩地层相对应。

图8

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

图8   基于本文频率域方法的三维成像结果

a—实际磁总场异常三维成像结果;b—磁垂直梯度三维成像结果;红色实线为两个钻孔位置,白线为沿zk1钻孔的金属矿轮廓线

Fig.8   The 3-D imaging results by using the presented frequency-domain approach

a—the 3-D imaging result of the real total magnetic anomalies;b—the 3-D imaging result of the real magnetic vertical gradient. Red solid lines show the two boreholes;white line outlines the metallic ores along the zk1 borehole


4 结论

文中给出了一种磁异常与梯度三维成像的频率域迭代方法,为了改善成像深度精度,在成像方法中加入了深度尺度因子。同时,为了进一步提高成像精度,采用了迭代优化算法。此外,该方法要求磁数据的观测面必须为平面,若观测面为起伏面,则在进行成像计算之前要对其进行曲化平处理。

本文频率域迭代方法,简单、快速、易于实现,适用于三维大尺度数据的正演和成像,并能够快速模拟出地下半空间物性大致分布情况。通过对理论模型数据和新疆某金属矿区实测数据试验,验证了该方法的可行性、有效性。

致谢:

衷心感谢王万银教授及物探与化探编辑部提出的宝贵意见。

The authors have declared that no competing interests exist.
作者已声明无竞争性利益关系。

参考文献

Li Y, Oldenburg D W .

3-D inversion of magnetic data

[J]. Geophysics, 1996,61(2):394-408.

DOI:10.1190/1.1443968      URL     [本文引用: 2]

Pilkington M .

3D magnetic data-space inversion with sparseness constraints

[J]. Geophysics, 2009,74(1):L7-L15.

DOI:10.1190/1.3026538      URL     [本文引用: 1]

杨文采 . 地球物理反演的理论与方法[M]. 北京: 地质出版社, 1997.

[本文引用: 1]

Yang W C. Theory and Method of Geophysical Inversion[M]. Beijing: Geological Publishing House, 1997.

[本文引用: 1]

管志宁, 侯俊胜, 黄临平 , 等.

重磁异常反演的拟BP神经网络方法及其应用

[J]. 地球物理学报, 1998,41(2):242-251.

DOI:      Magsci     [本文引用: 1]

把神经网络与重磁异常反演理论相结合，提出了用于重磁反演的一种拟BP神经网络方法．基于3层神经网络结构，把隐含层神经元设定为三维空间物性（磁化强度或密度）单元．对实测与理论重磁异常经S型函数变换，采用自动修改物性单元物性值的拟BP算法，反演三维空间的物性分布．利用该网络对理论模型数据和内蒙古某花岗岩体上的航磁资料进行了反演计算，取得了满意的反演效果．

Guan Z N, Hou J S, Huang L P , et al.

Inversion of gravity and magnetic anomalies using pseudo-BP neural network method and its application

[J]. Chinese Journal of Geophysics, 1998,41(2):242-251.

Magsci     [本文引用: 1]

王万银, 冯旭亮, 高玲举 , 等.

重磁方法在吐尔库班套铜镍矿区勘查中的应用

[J]. 物探与化探, 2014,38(3):423-429.

[本文引用: 1]

Wang W Y, Feng X L, Gao L J , et al.

The application of gravity and magnetic techniques to the prospecting for the Tuerkubantao copper-nickel ore district

[J]. Geophysical and Geochemical Exploration, 2014,38(3):423-429.

[本文引用: 1]

Blakely R J. Potential Theory in Gravity and Magnetic Applications[M]. Cambridge: Cambridge University Press, 1995.

[本文引用: 5]

李焓, 邱之云, 王万银 .

复杂形体重、磁异常正演问题综述

[J]. 物探与化探, 2008,32(1):36-43.

Magsci     [本文引用: 1]

<p>总结了复杂形体重、磁异常各种正演方法,并对各种方法所具备的优点和存在的缺点进行了分析说明。据此讨论了复杂形体重、磁异常正演方法将来的发展方向。</p>

Li H, Qiu Z Y, Wang W Y .

A review of the forward calculation of gravity and magnetic anomalies caused by irregular models

[J]. Geophysical and Geochemical Exploration, 2008,32(1):36-43.

Magsci     [本文引用: 1]

Tikhonov A N, Arsenin V Y .

Solutions of Ⅲ-posed problems

[M]. Washington D.C: V. H. Winston & Sons, 1977.

[本文引用: 1]

Li Y, Oldenburg D W .

Fast inversion of large-scale magnetic data using wavelet transforms and a logarithmic barrier method

[J]. Geophysics Journal International, 2003,152(2):251-265.

DOI:10.1046/j.1365-246X.2003.01766.x      URL     [本文引用: 2]

Commer M .

Three-dimensional gravity modelling and focusing inversion using rectangular meshes

[J]. Geophysical Prospecting, 2011,59(5):966-979.

[本文引用: 1]

姚长利, 郝天珧, 管志宁 .

重磁反演约束条件及三维物性反演技术策略

[J]. 物探与化探, 2002,26(4):253-257.

Magsci     [本文引用: 1]

<p>重磁资料反演与其他地球物理反演一样也存在严重的多解性,要想得到好的结果,必须附加约束条件,而且尽可能是各种约束的组合。三维反演中多解性更加严重,同时与约束的结合又更加艰难。非线性的广义随机算法使反演求解过程稳定,约束条件容易结合,但计算速度和维数困难同样制约其发挥作用,采取针对性措施后,使三维反演进入实用化阶段。</p>

Yao C L, Hao T Y, Guan Z N .

Restrictions in gravity and magnetic inversions and technical strategy of 3D properties inversion

[J]. Geophysical and Geochemical Exploration, 2002,26(4):253-257.

Magsci     [本文引用: 1]

Portniaguine O, Zhdanov M S .

Focusing geophysical inversion images

[J]. Geophysics, 1999,64(3):874-887.

DOI:10.1190/1.1444596      URL     [本文引用: 1]

Portniaguine O, Zhdanov M S .

3-D magnetic inversion with data compression and image focusing

[J]. Geophysics, 2002,67(5):1532-1541.

DOI:10.1190/1.1512749      URL     [本文引用: 2]

罗凡, 严加永, 付光明 .

基于已知信息约束的重磁三维反演在深部磁铁矿勘查中的应用——以安徽泥河铁矿为例

[J]. 物探与化探, 2018,42(1):50-60.

[本文引用: 1]

Luo F, Yan J Y, Fu G M .

The application of gravity and magnetic three-dimensional inversion based on known information constraint in deep magnetite exploration: A case study of the Nihe iron deposit in Anhui Province

[J]. Geophysical and Geochemical Exploration, 2018,42(1):50-60.

[本文引用: 1]

Lv Q, Qi G, Yan J .

3D geologic model of Shizishan ore field constrained by gravity and magnetic interactive modeling: A case history

[J]. Geophysics, 2012,78(1):B25-B35.

[本文引用: 1]

Zhang Y, Yan J, Li F , et al.

A new bound constraints method for 3-D potential field data inversion using Lagrangian multipliers

[J]. Geophysical Journal International, 2015,201(1):267-275.

DOI:10.1093/gji/ggv016      URL     [本文引用: 1]

Li Y, Oldenburg D W .

Joint inversion of surface and three-component borehole magnetic data

[J]. Geophysics, 2000,65(2):540-552.

DOI:10.1190/1.1444749      URL     [本文引用: 1]

姚长利, 郝天珧, 管志宁 , 等.

重磁遗传算法三维反演中高速计算及有效存储方法技术

[J]. 地球物理学报, 2003,46(2):252-258.

Magsci     [本文引用: 1]

将地下场源区域规则划分成很多小长方体单元，并且通过反演确定这些单元的物性变 化，勾画出场源的分布图像，这种方式逐步成为重磁反演，特别是三维反演的重要方向；遗 传算法等非线性技术进行该类反演将逐步成为发展趋势. 本文指出，在应用遗传算法进行该 类反演过程中，隐含着数据量较大时超常规的计算量，它已成为制约该类反演充分发挥作用 的瓶颈问题；同时，本文提出了针对性的分离并存储几何格架的计算策略、以及独特的几何 格架等效压缩存储技术，可以从根本上提高非线性反演计算速度，为该类反演的有效应用奠 定了坚实的基础.

Yao C L, Hao T Y, Guan Z N , et al.

High-speed computation and efficient storage in 3-D gravity and magnetic inversion based on genetic algorithms

[J]. Chinese Journal of Geophysics, 2003,46(2):252-258.

Magsci     [本文引用: 1]

姚长利, 郑元满, 张聿文 .

重磁异常三维物性反演随机子域法方法技术

[J]. 地球物理学报, 2007,50(5):1576-1583.

DOI:      Magsci     [本文引用: 1]

<FONT face=Verdana>本研究针对三维物性反演中存在的困难和问题，提出三维物性反演的随机子域方法技术，首先是将正反演中保持不变的几何格架分离计算并存储，避免重复计算，从而提高正反演计算速度；其次是利用对称性等实现等效计算，明显降低格架计算和存储要求；再通过随机子域方式，降低反演的维数问题；另外，通过概率方式控制子域生成的分布，实现约束新机制. 模型和实例计算表明了方法技术的效果，为大面积重磁数据的三维反演提供了有效的途径. </FONT>

Yao C L, Zheng Y M, Zhang Y W .

3-D gravity and magnetic inversion for physical properties using stochastic subspaces

[J]. Chinese Journal of Geophysics, 2007,50(5):1576-1583.

Magsci     [本文引用: 1]

Vatankhah S, Renaut R A, Ardestani V E .

A fast algorithm for regularized focused 3D inversion of gravity data using randomized singular-value decomposition

[J]. Geophysics, 2018,83(4):G25-G34.

DOI:10.1190/geo2017-0386.1      URL     [本文引用: 1]

Parker R L .

The rapid calculation of potential anomalies

[J]. Geophysical Journal International, 1973,31(4):447-455.

DOI:10.1111/j.1365-246X.1973.tb06513.x      URL     [本文引用: 2]

Oldenburg D W .

The inversion and interpretation of gravity anomalies

[J]. Geophysics, 1974,39(4):526-536.

DOI:10.1190/1.1440444      URL     [本文引用: 1]

Cribb J .

Application of the generalized linear inverse to the inversion of static potential data

[J]. Geophysics, 1976,41(6):1365-1369.

DOI:10.1190/1.1440686      URL     [本文引用: 1]

Kobrunov A I, Varfolomeev V A .

On one method of ε-equivalent redistribution and its practical application in the interpretation of gravity fields

[J]. Earth Physics USSR Academy of Science, 1981,10:25-44.

[本文引用: 3]

Pedersen L B .

Relations between potential fields and some equivalent sources

[J]. Geophysics, 1991,56(7):961-971.

DOI:10.1190/1.1443129      URL     [本文引用: 1]

Priezzhev I I .

Integrated interpretation technique of geophysical data for geological modeling

[D]. Moscow:State university of Sergo ordjonikidze, 2010.

[本文引用: 3]

Priezzhev I I, Scollard A, Lu Z, Schlumberger.

Regional production prediction technology based on gravity and magnetic data from the Eagle Ford formation, Texas, USA

[C]//SEG Technical Program Expanded Abstracts, 2014, 1354-1358.

[本文引用: 3]

Kobrunov A I .

The method of functional representations in the solution of inverse problems of gravimetry

[J]. Izvestiya Physics of the Solid Earth, 2015,51(4):459-468.

DOI:10.1134/S1069351315030076      URL     [本文引用: 1]

王万银, 刘金兰, 邱之云 , 等.

频率域偶层位曲面位场处理和转换方法研究

[J]. 地球物理学报, 2009,52(10):2652-2665.

DOI:10.3969/j.issn.0001-5733.2009.10.026      Magsci     [本文引用: 1]

<FONT face=Verdana>在空间域偶层位法的基础上，研究了完整的频率域偶层位曲面位场处理和转换方法.该法可应用于平面或曲面、规则网或非规则网的位场数据处理和转换.通过对偶层面<EM>z</EM>坐标和计算面<EM>z</EM>坐标平移不同的量来加速正演快速收敛和保证反演稳定、快速收敛；提出了适合于不规则网曲面处理和转换的核心算法——单点快速Fourier变换；提出了频率域不规则网曲面处理和转换方法技术.通过以上技术措施解决了大数据量特别是曲面不规则网的位场处理和转换问题，模型试算以及实际资料处理验证了该方法的应用效果.</FONT>

Wang W Y, Liu J L, Qiu Z Y , et al.

The research of the frequency domain dipole layer method for the processing and transformation of potential field on curved surface

[J]. Chinese Journal of Geophysics, 2009,52(10):2652-2665.

Magsci     [本文引用: 1]

程振炎 .

重磁场的有限元法曲化平

[J]. 物探与化探, 1981,5(3):153-158.

Magsci     [本文引用: 1]

<p>有限元法是解椭圆形体偏微分方程的有效数值方法。这种方法通用性较强,适合解比较复杂的边值问题。我们将它用于重磁场的曲化平,在观测面为任意曲面的情况下,依据实际观测场值,按三角形单元剖分求解,得到了较好的效果。特别是在场的边界处理上,体现了有限元法的灵活方便之处。整个计算过程可以在电子计算机上实现。</p>

Cheng Z Y .

Curved leveling of gravity and magnetic field by finite element method

[J]. Geophysical and Geochemical Exploration, 1981,5(3):153-158.

Magsci     [本文引用: 1]

刘天佑, 刘大为, 詹应林 , 等.

磁测资料处理新方法及在危机矿山挖潜中的应用

[J]. 物探与化探, 2006,30(5):377-381,396.

Magsci     [本文引用: 1]

<p>研究了危机矿山深部找矿的高精度磁测资料处理解释新方法技术:基于三度体数值积分法的三维可视化精细反演, 小波分析方法提取深部弱磁信号,起伏地形下高精度磁测资料快速化平处理的等效偶层法,以及旧资料 换算为 并与新资料对比的综合处理解释等方法。将这些新方法运用于大冶铁矿危机矿山的深部找矿,获得明显的地质效果。2005年根据解释结果在15线布置的钻孔打到了26.65 m厚的磁铁矿体,说明在危机矿山进一步挖潜中,新方法技术的运用是十分必要的。</p>

Liu T Y, Liu D W, Shan Y L , et al.

A new method of magnetic survey data processing and its application in tapping potential of crisis mines

[J]. Geophysical and Geochemical Exploration, 2006,30(5):377-381,396.

Magsci     [本文引用: 1]

Guo L, Yan J .

3-D wavelet-based fusion approach for comprehensively analyzing multiple physical-property voxel models inverted from magnetic data

[J]. Journal of Applied Geophysics, 2017,139:47-53.

DOI:10.1016/j.jappgeo.2017.02.006      URL     [本文引用: 1]

/

〈 〉

京ICP备05055290号-3
版权所有 © 2021《物探与化探》编辑部
通讯地址：北京市学院路29号航遥中心 邮编：100083
电话：010-62060192；62060193 E-mail：whtbjb@sina.com