地球物理资料异常处理和解释方法

如题所述

所谓常规数据处理是指在重、磁数据处理中经常使用到的位场转换或滤波处理,如上延、求导数、化极等。下面对本次研究中使用的常规数据处理方法做一简述。

空间域异常处理、转换的基本公式均可以写成如下的褶积形式

东北地球物理场与地壳演化

式中:fa、fb分别表示原始异常处理和转换前、后的异常;为权函数。不同的处理、转换间只是它们的权函数不同而已。

由于上述两式与电子电工学中描述滤波器滤波特性的卷积的形式完全相同,因此异常的处理和转换又称为异常的滤波。由卷积定理有

东北地球物理场与地壳演化

式中: 分别为原始异常和处理、转换后异常的波谱; 为权函数φ的波谱(也称为处理与转换因子、波数响应、滤波算子)。

将(4-5)、(4-6)、(4-7)式比较可知,空间域的处理和转换是褶积运算,而波数域是乘积运算。而且,波数谱的连乘可以完成连续的多种变换。因此数域的转换方法要简单得多。随着电子计算机的广泛应用,特别是快速傅里叶变换算法的问世,使区域重磁资料数据处理中的波数域方法成为主要的方法。

在本次研究工作中,数据处理的计算工作是在波数域中进行的。十分明显,为实现波数域的异常处理、转换,必须已知或设计出处理、转换因子 。

根据计算,重、磁异常波谱公式是一些独立因子的乘积,其通式为

东北地球物理场与地壳演化

式中:A为参数因子,只与地质体的剩余密度或剩余质量有关;B为参数因子,只与地质体的磁化强度或磁矩有关;H(ϖ,h)为深度因子,仅与地质体的深度有关;S(ϖ,a,b)为水平尺寸因子,仅与地质体的水平尺寸有关;L(ls,ms,ns,mt,nt)为方向因子,只与磁化强度和地磁场的方向有关;D(ϖ,ξ,η)为位移因子,它是由于坐标原点任选而增加的因子。

(一)解析延拓

根据异常波谱特征的计算,可知无限延深直角棱柱体异常波谱深度因子H=e(-hϖ)。若原始异常体的深度为h1,解析延拓后异常体的深度为h2,则有

东北地球物理场与地壳演化

式中:Δh=h2-h1

于是 因子延拓因子应为

东北地球物理场与地壳演化

式中:ϖx、ϖy分别为与x轴、y轴的对应的圆波数, 为径向圆波数。

上延时Δh>0,下延时Δh<0。上延可压制高波数成分(即突出低波数成分),属低通滤波;下延可放大高波数成分(即突出高波数成分),属高通滤波,但对低波数成分无压制作用。

(二)导数计算

重、磁异常的导数计算广泛应用于异常的处理和解释。原因在于:①异常的导数在不同形状的地质体上有不同的特征,有助于对异常的解释和分类;②异常导数可以突出反映浅部地质因素,而压制区域深部地质因素的影响,在一定程度上可以划分不同深度和大小异常源产生的叠加异常;③在利用欧拉反褶积对重、磁异常进行构造反演计算时,要用到异常水平一阶导数和垂向一阶导数。

1.垂向导数

由于重、磁异常函数f(x,y,z)的n阶垂向导数可以用下列公式表示

nf(x,y,z)/∂n=-∂nf(x,y,h)/∂hn,因此由深度因子H=e(-hϖ)可以求出异常的n阶垂向导数因子。

东北地球物理场与地壳演化

显然一、二阶垂向导数因子分别为:

东北地球物理场与地壳演化

重、磁异常垂直导数可放大高波数成分(即突出高波数成分),但对低波数成分有压制作用。

2.水平导数

由傅里叶变换(FT)的微分性质可知,沿x和y方向的异常水平导数因子分别为

东北地球物理场与地壳演化

式中: 当n取1,即为一阶水平导数的转换因子。如果异常f(x,y,h)对任意水平方向l的导数为: (其中,α为l与x轴的夹角)。依FT的微分性质可得到异常的方向导数因子。

由上式可知,异常的水平导数可突出某一方向上异常特征(或构造线),如α=45°时,能突出135°方向的构造线。

(三)化地磁极

球体总强度磁异常T的谱Δ 为:

东北地球物理场与地壳演化

式中:qs=j(lscosθ+mssinθ)+ns,qt=j(ltcosθ+mtsinθ)+nt;而ls,ms,ns为磁化强度M的方向余弦;lt,mt,nt为地磁场T0的方向余弦;θ为径向圆波数的极角。 为引力位的谱;G为万有引力常数;ρ为密度;ϖ为圆频率;μ0为真空导有磁率。

令qt=nt=1时,即垂直磁异常的谱为:

东北地球物理场与地壳演化

再令qs=ns=1时,即化极后垂直磁异常的谱为:

东北地球物理场与地壳演化

比较以上各式,便可得到化极转换因子为:

东北地球物理场与地壳演化

从(4-18)式看出,化极因子与ϖ无关,因此磁异常化极无滤波作用。另外,从(4-18)式可知,磁异常化极需已知磁化强度方向。然而剩磁与感磁方向不一致时,磁化强度的方向是难以确定的,尤其在大面积磁测资料处理时,区内磁体很多,更无法了解它们的磁化强度方向,因此往往假设磁化强度的方向与地磁场方向一致。另外,实践中还认为,在研究区范围内的地磁场方向是相同的。这种假设在测区不大的情况下对结果的影响较小。图4-2显示出化磁时磁化倾角对结果的影响。根据这一计算结果,若测区的南北纬度差在10°之内,用统一的地磁场方向余弦来作化磁极运算,其结果受影响不大。因此,对一般成矿预测为目的的区域性资料研究问题并不突出。但当作深部地质研究或大面积区域地质研究时,就要注意这种情况,有些学者已开始了测点地磁场方向余弦各异时的化磁研究工作。

图4-2 不同磁化方向化磁极后的曲线对比

(四)谱分析方法

谱分析方法作为重、磁异常数据处理、转换的重要方法,有着广泛的应用。利用径向平均对数能谱分析可以估算重、磁场源的平均深度,为进一步的处理和解释提供基础信息。

下面对径向平均对数能谱分析和平均深度估算的原理进行简介:

经计算可知,球体重力异常的波谱为:

东北地球物理场与地壳演化

则球体重力异常功率谱为:

东北地球物理场与地壳演化

对数径向功率频谱为:

东北地球物理场与地壳演化

式中A=2πGρv=2πGm。上式表明,球体重力异常对数径向功率谱与径向圆波数中呈线性关系,见图4-3a,故可利用lnE(ϖ)的拟合直线斜率求解出球体中心深度:

东北地球物理场与地壳演化

式中:f为波数,而ϖ=2πf。

另外,由球体磁异常功率谱也可计算深度。把引力位谱( =2πGme-hϖ/ϖ)代入(4-15)式,经整理可求得球体垂直磁化(qs=1),垂直磁异常(qt=1)的功率谱为:

东北地球物理场与地壳演化

式中B=2πμ0VM/4π=2πμ0m/4π;m为磁矩;V为球体体积;M为磁化强度。

东北地球物理场与地壳演化

上式说明,球体垂直磁化垂直磁异常的对数径向功率谱与圆波数呈非线性关系(图4-3b)。但是,在高波数段近于线性关系,可用(4-21)式计算深度。

图4-3 对数径向功率谱

(五)重、磁对应分析

基于泊松定理发展起来的重磁异常对应分析方法,是重磁数据综合解释的重要方法,能对重磁异常的相关性进行定量研究,有效地将重磁信息进行综合,对重磁资料定性地赋予地质意义,并突出地质目标的反映,为重磁资料的地质解释提供有用的信息,特别是在强磁性火山岩解释中具有重要的作用。

重磁异常对应分析方法的基本理论如下:

由同一场源引起的重力异常和磁异常间的关系可以简单地用泊松方程描述。当垂直磁化时,泊松方程可表示为:

东北地球物理场与地壳演化

式中:Δz为垂直磁化的垂直磁异常;M为场源磁化强度;G为万有引力常数;Δρ为场源剩余密度;Δg为重力异常; 为重力异常的垂向一阶导数。

上式表明垂直磁化的垂直磁异常与重力场的垂向一阶导数满足线性关系,而且拟合直线的截距为零。

由于原始资料不可避免地存在某些干扰因素,通常进行重磁异常的线性回归分析时,选用如下稍加推广的泊松方程:

东北地球物理场与地壳演化

式中:b为斜率;A为截距。

将Δz与 作线性回归分析则可得到斜率b与截距A的估计值。两个离散序列的相关导数可以由下式求得:

东北地球物理场与地壳演化

式中:Cxy(k) 为两个离散序列x(t)={x1,x2,…,xn}和y(t)={y1,y2,…,yn}的相关函数,k为延迟时间。当x(t)=y(t)时,称为自相关函数Cxx(k)或Cyy(k)。

计算处理时,给定适当大小的分析窗口,将窗口内各点垂直磁化磁异常和重力异常的垂向一阶导数进行最小二乘线性回归,求得中心点的相关系数R、斜率b和截距A。

相关系数R反映了在给定窗口内重磁异常的线性相关程度,即宏观地反映了重磁异常的“同源性程度”。相关系数绝对值接近于1的窗口区间重磁异常的“同源性”较好,它们或者同源、或者都离场源较远、或者同处异常的拐点等。其中R接近+1时,重磁异常正相关;R接近-1时,重磁异常负相关。当R绝对值较小时,重磁异常相关性差,重磁异常可能不同源,或存在邻近异常干扰,或是存在方向不同于地磁场的强剩磁磁性体等。

斜率b反映了所有场源泊松比的加权平均值,称为广义泊松比。只有在重磁异常同源的前提下,回归所得的斜率b才有意义。仅由b不能直接确定M和Δσ,但若在解释中结合其他地质、地球物理信息,就能从中获得关于物性分布的有用信息,从而为进一步的定量解释提供依据。

截距A反映了实测资料中的长波长成分,它主要反映重磁异常数据的背景变化。在重磁异常完全同源的理想情况下,A=0。

由于重磁异常对应分析是对场源之间的相关系数进行定性和半定量研究的方法,它能分离和鉴别不同类型的异常,从而勾画出与异常场源相对一致的地质单元和构造分区,不相关区说明重磁异常不同源或存在邻近异常干扰。

(六)欧拉反褶积与构造反演

欧拉反褶积方法使用欧拉(Euler)齐次关系,对经方向谱分析过的数据快速估计重、磁场源的位置和深度,是一种既能够利用重磁网格数据,又对剖面数据有效地确定地质体位置(边界)和深度的定量反演方法(Reid等,1990)。这种方法并不需要已知地质信息(密度、磁化率等)的控制。使用该方法可以将位场及其梯度以及场源位置之间的关系用欧拉齐次方程表示,而场源的不同形状即地质构造的差异则表现为方程的齐次程度,就是所谓的地质构造指数,地质构造指数或齐次程度实质上表现了场随离开场源距离的衰减率。模型研究和应用实例表明,这个方法对确定断层、磁性接触带、岩脉、喷出岩体等构造位置或勾绘它们的轮廓有较高的精度。

位场的欧拉方程是由Thompson推导的。首先建立一个直角坐标系,取观测平面为z=0,z轴向下为正,x轴指北,y轴指东。考虑在此坐标系中的任一函数f(x,y,z),如果

东北地球物理场与地壳演化

则称函数是n阶齐次的。此外可证明,如果f(x,y,z)是n阶齐次的,则满足下列方程

东北地球物理场与地壳演化

此偏微分方程称为欧拉齐次方程,或称欧拉方程。

对于位于(x0,y0,z0)的点磁源,在观测平面上任一点(x,y,z)处的总磁场强度具有如下形式:

东北地球物理场与地壳演化

式中 N=1,2,3,……。G不依赖于(x,y,z)。对于(3-23)这样的函数,其欧拉方程可写成

东北地球物理场与地壳演化

方程(4-29)是n=-N阶齐次的。三个坐标方向的梯度值可以利用空间域或波数域的一般位场变换计算出来。如果梯度值通过观测获得,直接用于方程(4-29)则更可取。

方程(4-29)虽然是根据磁源异常推导的,但对重力异常也同样适用。该方程用于平面网格重、磁异常数据的反演计算。如果假定方程中横向梯度∂ΔT/∂y为零,则可得到适用于剖面数据计算的方程。这对于众多走向方向不变的二度情况很显然就是这样。

齐次度N被定义为“构造指数”,它是重、磁异常场源深度变化“陡缓”的量度。特定的地质构造具有特定的衰减率(即:构造指数)。例如:倾斜断层的磁场、水平薄岩脉的磁场按线性的规律变化,构造指数就为1。表4-1列出了构造指数对应的地质构造。

利用不同坐标点(x,y,z)上的场值ΔT及其三个方向上的梯度值以及方程(4-29)组成的线性方程组,最后可以解出未知变量x0、y0、z0,进而确定构造形迹及位置。

表4-1 欧拉构造指数表

但是,直接用方程(4-29)及其变换的二度形式解决构造问题,会使解的精度极不可靠和不稳定。主要原因有如下几个方面:

(1)很难知道磁场ΔT的绝对水平,区域场或邻近磁异常的影响几乎总是存在的。

(2)根据线性方程组与系数的关系,较低的构造指数才会有较好的深度估计值。但大多数磁异常是偶极性的,有较高的构造指数。同时又有许多线性构造的指数接近于零而使反褶积发散。

(3)实测异常是多种构造指数特征的复杂叠加,很难用一些简单模型来模拟,亦很难将具有线性特征的构造识别与分离出来。

为克服上述三个方面的问题,釆用下面的一些办法:

从观测数据中消除偏差是通过网格数据进行窗口计算解决的。对网格数据假定异常在方程(4-29)求值的窗口范围内有一常量偏差,观测值为

东北地球物理场与地壳演化

这里B是常数。从方程(4-30)中解出ΔT,代入(4-29),整理得

东北地球物理场与地壳演化

如果构造指数小于0.5,即构造指数接近零时,这样可能造成对深度值的过低估计。为此需要提供一个补偿值A,使得欧拉齐次方程在构造指数较低时写为

东北地球物理场与地壳演化

式中:A是与场幅值有关的一个参数。不同的构造形体有不同的A,A可以通过将已知的某一构造的解析式代入欧拉方程(4-29)而求出。

图4-4 重力异常欧拉反褶积计算结果示意图

图4-4中的曲线为重力异常等值线,圆圈为反演解的构造位置。圆圈直径的大小代表了不同的构造深度。

(七)重、磁人机交互剖面正反演

该项技术的优点是便于将重、磁异常的处理、转换方法得出的结果和其他地质、地球物理方法获取的先验信息输入到模型里,形成初始模型。并且根据计算结果和实际重、磁异常的差异,随时方便地修改模型,直观地监督和指导正反演过程。重、磁人机交互剖面正反演流程见图4-5。

图4-5 人机交互正反演流程图

1.重力人机交互正反演技术

重力人机交互正反演技术(Gamble,1979)主要是依据A截面为多边形的二度体重力异常计算方法来实现的。通过对初始模型计算出的重力效应与测线上的布格重力异常进行对比,不断修正模型,直至达到计算出的重力效应与测线上的布格重力异常之差满足预定精度。重力人机交互正反演流程见图4-5。

图4-6 二度体A截面

A截面为多边形的二度体重力异常计算方法:

假设二度体的剩余密度为σ,以计算点作为坐标原点,x轴与二度体走向垂直,z轴铅垂向下(图4-6)。若n边形第k个顶点的坐标为(ξk,ζk),其中k=1,2,...,n。则(ξk,ζk)与(ξk+1,ζk+1)两个顶点连线上ξ与ζ有如下关系:

东北地球物理场与地壳演化

引用解Δg正问题的基本公式,首先对ξ求积分,得

东北地球物理场与地壳演化

式中:s为多边形的A截面积;l为A截面的周长。

将式(4-33)代入(4-34)得

东北地球物理场与地壳演化

对上式积分可得到如下结果

东北地球物理场与地壳演化

或写成下面的形式

东北地球物理场与地壳演化

(4-36)式或(4-37)可以编成计算机程序,用以计算A截面为任意多边形的二度体的重力异常,进而可以进行重力人机交互正反演。在具体编程计算时应注意以下几个问题:

(1)因多边形的边数为n,故ξn+11,ζn+11;

(2)(4-36)和(4-37)两式是假设计算点位于原点时导出的,因此,当任意计算点P(x,y)的重力异常时,式中的ξk和ζk应以ξk-x和ζk-z来代替;

(3)在(4-36)和(4-37)式中,反正切函数的取值范围应在-π到π之间,即当ξk+1k时,反正切函数在0到π之间取值;反之,则在-π到0之间取值。

2.磁法人机交互正反演技术

磁法人机交互正反演技术主要是依据A截面为多边形的二度体磁力异常计算方法来实现的。其基本思想同重力人机交互正反演技术相一致。

由于V2=Δg,所以根据公式(4-37)可以求出引力位的二阶导数

东北地球物理场与地壳演化

将(4-38)和(4-39)式代入下面二度体的磁异常公式,就可以利用该式进行磁力人机交互正反演计算。

东北地球物理场与地壳演化

式中:μ0为真空的磁导率;MS为有效磁化强度;is为有效磁化倾角;I0为地磁场倾角;A'为x轴与磁北的夹角。在具体编程序上机计算时应注意的问题等方面与重力人机交互方法相同(图4-5)。

温馨提示:答案为网友推荐,仅供参考
相似回答