CN110501751A - 一种基于多分量梯度数据联合和深度加权的相关成像方法 - Google Patents

Abstract

本发明涉及地球物理数据解译技术领域，提供一种基于多分量梯度数据联合和深度加权的相关成像方法。首先，使用重力仪采集待成像地质区域的重力梯度数据；然后，对重力梯度数据进行预处理；接着，对预处理后的重力梯度数据的部分分量进行联合，选取成像效果最好的多分量梯度数据组合计算相关系数矩阵；再根据待成像地质区域的先验地质资料，得到待成像地质区域的顶面埋深和底面埋深，对待成像地质区域的地下空间进行分块处理，引入深度加权函数，计算加权后的相关系数矩阵；最后，利用加权后的相关系数矩阵对待成像地质区域进行三维成像。本发明能够提高成像结果在深度方向上的分辨率，更加精准地确定地质体的上下边界，且计算效率高。

Description

一种基于多分量梯度数据联合和深度加权的相关成像方法

技术领域

本发明涉及地球物理数据解译技术领域，特别是涉及一种基于多分量梯度数据联合和深度加权的相关成像方法。

背景技术

随着我国资源开采向地球深部进军，传统和现有的物探方法技术还存在着许多缺陷，不能满足要求，特别是重、磁方法要想在新一轮的资源勘探中发挥应有的作用，其方法技术的突破势在必行。地球物理勘探中，数据往往是海量的，如何能够更快、更精确解译出勘探结果是地球物理未来发展的趋势。

重、磁三维反演根据反演方式的不同大体上可分为两类，一类是基于反演理论、在最小二乘的意义下使目标函数达到极小的线性与非线性反演，在引入足够约束条件和相对简单的情况下，能够给出接近实际情况的物性分布和几何形态。由于重、磁三维反演在数学上是典型的病态问题，具有固有的多解性，要使病态问题转化成非病态的稳定问题，需要补充、增加信息，即引入先验约束条件,但约束条件的提取及与反演过程的结合，特别是一些地质信息如何转化成具有解析形式的约束条件并同反演过程融合，目前还存在着许多问题。另外，三维反演的另一个问题是计算维数巨大。

重、磁三维反演的另一类方法是快速自动解释，如欧拉反褶积、解析信号法和位场相关成像等。解析信号和欧拉褶积方法可对大面积的平面网格数据进行自动反演解释，具有较强的适应性和灵活性，近年来成为位场反演方法研究的一个热点。

综合以上两类反演方法，要想快速处理海量地球物理勘探数据，就需要用到上述方法中的位场相关成像方法。

位场相关成像方法的基本原理是，将下半空剖分成均匀网格，计算每一网格结点(扫描点)单位物性差所产生的异常(扫描函数)与实测场在一定窗口范围内的归一化互相关(场源发生的概率)，并逐点移动扫描点和窗口，计算所有网格结点上场源发生的概率，然后根据场源出现的概率情况勾画出地质体的分布情况。

传统的位场相关成像方法在横向方向上可以很好地确定地质体的空间位置和形态，但在纵向上的分辨率较低，存在趋肤效应以及异常体底部相连等情况。现有技术中针对该问题进行了改进，如将一种基于先验深度信息的加权函数应用于重力聚焦反演中，随后进一步提出空间梯度加权函数，能够提升反演效果，改善地质体在深度方向上分辨率较低的问题；将相关成像分别应用于磁法和重力勘探数据，提出基于异常分离的三维相关成像；磁总场异常垂直梯度的相关成像，在模型数据与实际磁测数据中验证了该方法的稳定性。然而，上述现有技术中，重力勘探中相关成像的纵向分辨率仍然较低。

发明内容

针对现有技术存在的问题，本发明提供一种基于多分量梯度数据联合和深度加权的相关成像方法，能够提高成像结果在深度方向上的分辨率，更加精准地确定地质体的上下边界，且计算效率高。

本发明的技术方案为：

一种基于多分量梯度数据联合和深度加权的相关成像方法，其特征在于，包括下述步骤：

步骤1：进行数据采集：使用重力仪采集待成像地质区域的重力梯度数据；所述重力梯度数据有6个分量{V_xx,V_xy,V_xz,V_yy,V_yz,V_zz}；

步骤2：进行数据预处理：对重力梯度数据进行网格化处理、插值处理、地形校正及高通滤波处理；

步骤3：对预处理后的重力梯度数据的部分分量进行联合，得到多种多分量梯度数据组合，选取成像效果最好的多分量梯度数据组合；

步骤4：根据成像效果最好的多分量梯度数据组合计算相关系数矩阵C；

步骤5：根据待成像地质区域的先验地质资料，得到待成像地质区域的顶面埋深和底面埋深，对待成像地质区域的地下空间进行分块处理，引入深度加权函数，计算加权后的相关系数矩阵C_later；

步骤6：利用加权后的相关系数矩阵C_later对待成像地质区域进行三维成像。

所述步骤2中，采用空间波长为2000-5000m的高通滤波器对重力梯度数据进行高通滤波处理以移除区域背景场。

所述步骤3包括下述步骤：

步骤3.1：对预处理后的重力梯度数据的部分分量进行联合，得到多种多分量梯度数据组合；其中，每种多分量梯度数据组合中第j个质量元的相关系数为该多分量梯度数据组合中第j个质量元的参与联合的各个分量的相关系数的乘积C_joint,j＝C_1,j×…×C_n,j，n∈{1,2,...,6}；其中，C_n,j为第j个质量元的第n个重力梯度分量的相关系数；

步骤3.2：利用每种多分量梯度数据组合对待成像地质区域进行三维成像，得到每种多分量梯度数据组合下的成像结果图；

步骤3.3：选取分辨率最高的成像结果图对应的多分量梯度数据组合作为成像效果最好的多分量梯度数据组合。

所述步骤4中，根据成像效果最好的多分量梯度数据组合计算相关系数矩阵C，包括：

根据成像效果最好的多分量梯度数据组合计算第j个质量元的相关系数为

其中，αβ∈A，A为成像效果最好的多分量梯度数据组合中参与联合的分量集合，V_αβ,j(x_i,y_i,z_i)为地下第j个质量元在第i个测点上引起的重力梯度的αβ分量，N为测点总数，B_αβ,j(x_i,y_i,z_i)为地下第j个质量元在第i个测点上引起的重力梯度的αβ分量的几何函数；

所有质量元的相关系数构成相关系数矩阵C。

所述步骤5包括下述步骤：

步骤5.1：根据待成像地质区域的先验地质资料，得到待成像地质区域的顶面埋深和底面埋深，对待成像地质区域的地下空间进行分块处理；

步骤5.2：引入深度加权函数为

其中，z为地下网格的深度，dz为待成像地质区域地下空间范围内的最大深度值；α为经验值，α＝0.001；r为缩放因子；z_c1为待成像地质区域的顶面埋深，z_c2为待成像地质区域的底面埋深；

步骤5.3：将相关系数矩阵C与深度加权函数W(z)结合，计算加权后第j个质量元的相关系数为C_later,j＝C_j×W(z_j)，得到加权后的相关系数矩阵C_later；其中，C_j为相关系数矩阵C中第j个质量元的相关系数，W(z_j)为第j个质量元的深度加权函数值，z_j为第j个质量元所在地下网格的深度。

所述步骤6中，利用加权后的相关系数矩阵C_later通过matlab中的slice函数对待成像地质区域进行三维成像。

本发明的有益效果为：

(1)本发明通过多分量梯度数据联合，得到多种多分量梯度数据组合，根据成像效果最好的多分量梯度数据组合计算相关系数矩阵，能够提高成像的分辨率。

(2)本发明引入基于先验信息的深度加权函数，通过对待成像地质区域的地下空间进行分块处理，能够进一步提高成像结果在深度方向上的分辨率，更加精准地确定地质体的上下边界。

(3)本发明避免了常规反演方法需要进行的多次迭代运算，当处理的数据量很大时，计算效率仍然很高。

附图说明

图1为本发明的基于多分量梯度数据联合和深度加权的相关成像方法的流程图。

图2为具体实施方式中未引入深度加权函数得到的垂直剖面x＝442.56km处的V_xz|V_yz|V_zz组合成像结果图。

图3为具体实施方式中未引入深度加权函数得到的垂直剖面y＝3334.44km处的V_xz|V_yz|V_zz组合成像结果图。

图4为具体实施方式中使用本发明的基于多分量梯度数据联合和深度加权的相关成像方法得到的垂直剖面x＝442.56km处的V_xz|V_yz|V_zz组合成像结果图。

图5为具体实施方式中使用本发明的基于多分量梯度数据联合和深度加权的相关成像方法得到的垂直剖面y＝3334.44km处的V_xz|V_yz|V_zz组合成像结果图。

具体实施方式

下面将结合附图和具体实施方式，对本发明作进一步描述。

如图1所示，本发明的基于多分量梯度数据联合和深度加权的相关成像方法，包括下述步骤：

步骤1：进行数据采集：使用重力仪采集待成像地质区域的重力梯度数据；所述重力梯度数据有6个分量{V_xx,V_xy,V_xz,V_yy,V_yz,V_zz}。

本实施例中，待成像地质区域为美国文顿盐丘，位于德克萨斯州与路易斯安那州的西南部相交地带，主要由一个大的盐岩和上方的盖岩组成，盖岩主要是由硬石膏和石灰岩组成的。

本实施例中，待成像地质区域的重力梯度数据为美国BellGeospace公司提供的数据，包含一个10000*6的矩阵，矩阵的每一行代表对应质量元的重力梯度数据的6个分量。

步骤2：进行数据预处理：对重力梯度数据进行网格化处理、插值处理、地形校正及高通滤波处理。

本实施例中，对文顿盐丘进行网格化处理，将文顿盐丘划分为100*100*20的网格区域；采用空间波长为2000-5000m的高通滤波器对重力梯度数据进行高通滤波处理以移除区域背景场。

步骤3：对预处理后的重力梯度数据的部分分量进行联合，得到多种多分量梯度数据组合，选取成像效果最好的多分量梯度数据组合，包括下述步骤：

步骤3.1：对预处理后的重力梯度数据的部分分量进行联合，得到多种多分量梯度数据组合；其中，每种多分量梯度数据组合中第j个质量元的相关系数为该多分量梯度数据组合中第j个质量元的参与联合的各个分量的相关系数的乘积C_joint,j＝C_1,j×...×C_n,j，n∈{1,2,...,6}；其中，C_n,j为第j个质量元的第n个重力梯度分量的相关系数；

步骤3.2：利用每种多分量梯度数据组合对待成像地质区域进行三维成像，得到每种多分量梯度数据组合下的成像结果图；

步骤3.3：选取分辨率最高的成像结果图对应的多分量梯度数据组合作为成像效果最好的多分量梯度数据组合。

其中，一般来讲，3-4个的分量组合效果最好，因为每个分量都包含不同的地质信息，分量太少，则信息量不全面，分量太多，则会密度堆积，识别出的部分模型范围减小，分辨率反而有所下降。本实施例中，文顿盐丘中，多分量梯度数据组合V_xz|V_yz|V_zz的成像结果最好。

步骤4：根据成像效果最好的多分量梯度数据组合计算相关系数矩阵C，包括：

根据成像效果最好的多分量梯度数据组合计算第j个质量元的相关系数为

其中，αβ∈A，A为成像效果最好的多分量梯度数据组合中参与联合的分量集合，V_αβ,j(x_i,y_i,z_i)为地下第j个质量元在第i个测点上引起的重力梯度的αβ分量，N为测点总数，B_αβ,j(x_i,y_i,z_i)为地下第j个质量元在第i个测点上引起的重力梯度的αβ分量的几何函数；

所有质量元的相关系数构成相关系数矩阵C。

其中，相关系数的取值范围是[-1,1]，相关系数可用来反映地下质量分布的情况。若相关系数为正，则表明该点质量盈余，反之则代表亏损。相关系数越接近1，表明该点质量剩余的可能性越大，反之相关系数越接近-1，表示该点质量亏损的可能性越大。本实施例中，相关系数矩阵C为100*100*20的三维矩阵。

步骤5：根据待成像地质区域的先验地质资料，得到待成像地质区域的顶面埋深和底面埋深，对待成像地质区域的地下空间进行分块处理，引入深度加权函数，计算加权后的相关系数矩阵C_later，包括下述步骤：

步骤5.1：根据待成像地质区域的先验地质资料，得到待成像地质区域的顶面埋深和底面埋深，对待成像地质区域的地下空间进行分块处理；

步骤5.2：引入深度加权函数为

其中，z为地下网格的深度，dz为待成像地质区域地下空间范围内的最大深度值；α为经验值，α决定近地表处权值的大小，为克服趋附效应，取α＝0.001；r为缩放因子；z_c1为待成像地质区域的顶面埋深，z_c2为待成像地质区域的底面埋深；

步骤5.3：将相关系数矩阵C与深度加权函数W(z)结合，计算加权后第j个质量元的相关系数为C_later,j＝C_j×W(z_j)，得到加权后的相关系数矩阵C_later；其中，C_j为相关系数矩阵C中第j个质量元的相关系数，W(z_j)为第j个质量元的深度加权函数值，z_j为第j个质量元所在地下网格的深度。

本实施例中，根据文顿盐丘的先验地质资料，包括地震与测井数据，得到地质体的顶面埋深与底面埋深的大概范围。人们已经用多种地球物理数据对文顿盐丘地区进行了研究，包括VSP方法、反射地震数据等多种手段。归纳这些实验数据，取顶面埋深z_c1＝160m、底面埋深z_c2＝450m，应用深度加权函数对相关系数矩阵C进行加权处理，得到加权后的相关系数矩阵C_later，新矩阵C_later的维数也是100*100*20。

步骤6：利用加权后的相关系数矩阵C_later对待成像地质区域进行三维成像。

本实施例中，对文顿盐丘进行了两组成像对比实验。第一组中，未应用深度加权函数对相关系数矩阵C进行加权处理，直接利用相关系数矩阵C对文顿盐丘进行三维成像，得到垂直剖面x＝442.56km处、y＝3334.44km处的V_xz|V_yz|V_zz组合成像结果图分别如图2、图3所示；第二组中，则应用了深度加权函数对相关系数矩阵C进行加权处理，利用加权后的相关系数矩阵C_later对文顿盐丘进行三维成像，得到垂直剖面x＝442.56km处、y＝3334.44km处的V_xz|V_yz|V_zz组合成像结果图分别如图4、图5所示。其中，两组均是通过matlab中的slice函数对文顿盐丘进行三维成像，得到的成像结果图能够显示地质体的空间位置与形态，从而可进一步分析地下存在的断裂以及构造边界。

由两组成像结果图的对比可以看出，本发明的引入了深度加权函数的相关成像方法得到的成像结果在深度方向上的分辨率较未引入深度加权函数的相关成像方法有明显提高，能够更加精准地确定地质体的上下边界。

显然，上述实施例仅仅是本发明的一部分实施例，而不是全部的实施例。上述实施例仅用于解释本发明，并不构成对本发明保护范围的限定。基于上述实施例，本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，也即凡在本申请的精神和原理之内所作的所有修改、等同替换和改进等，均落在本发明要求的保护范围内。

Claims (6)

1.一种基于多分量梯度数据联合和深度加权的相关成像方法，其特征在于，包括下述步骤：

步骤1：进行数据采集：使用重力仪采集待成像地质区域的重力梯度数据；所述重力梯度数据有6个分量{V_xx,V_xy,V_xz,V_yy,V_yz,V_zz}；

步骤2：进行数据预处理：对重力梯度数据进行网格化处理、插值处理、地形校正及高通滤波处理；

步骤3：对预处理后的重力梯度数据的部分分量进行联合，得到多种多分量梯度数据组合，选取成像效果最好的多分量梯度数据组合；

步骤4：根据成像效果最好的多分量梯度数据组合计算相关系数矩阵C；

步骤5：根据待成像地质区域的先验地质资料，得到待成像地质区域的顶面埋深和底面埋深，对待成像地质区域的地下空间进行分块处理，引入深度加权函数，计算加权后的相关系数矩阵C_later；

步骤6：利用加权后的相关系数矩阵C_later对待成像地质区域进行三维成像。

2.根据权利要求1所述的基于多分量梯度数据联合和深度加权的相关成像方法，其特征在于，所述步骤2中，采用空间波长为2000-5000m的高通滤波器对重力梯度数据进行高通滤波处理以移除区域背景场。

3.根据权利要求1所述的基于多分量梯度数据联合和深度加权的相关成像方法，其特征在于，所述步骤3包括下述步骤：

步骤3.1：对预处理后的重力梯度数据的部分分量进行联合，得到多种多分量梯度数据组合；其中，每种多分量梯度数据组合中第j个质量元的相关系数为该多分量梯度数据组合中第j个质量元的参与联合的各个分量的相关系数的乘积C_joint,j＝C_1,j×…×C_n,j，n∈{1,2,…,6}；其中，C_n,j为第j个质量元的第n个重力梯度分量的相关系数；

步骤3.2：利用每种多分量梯度数据组合对待成像地质区域进行三维成像，得到每种多分量梯度数据组合下的成像结果图；

步骤3.3：选取分辨率最高的成像结果图对应的多分量梯度数据组合作为成像效果最好的多分量梯度数据组合。

4.根据权利要求1所述的基于多分量梯度数据联合和深度加权的相关成像方法，其特征在于，所述步骤4中，根据成像效果最好的多分量梯度数据组合计算相关系数矩阵C，包括：

根据成像效果最好的多分量梯度数据组合计算第j个质量元的相关系数为

其中，αβ∈A，A为成像效果最好的多分量梯度数据组合中参与联合的分量集合，V_αβ,j(x_i,y_i,z_i)为地下第j个质量元在第i个测点上引起的重力梯度的αβ分量，N为测点总数，B_αβ,j(x_i,y_i,z_i)为地下第j个质量元在第i个测点上引起的重力梯度的αβ分量的几何函数；

所有质量元的相关系数构成相关系数矩阵C。

5.根据权利要求1所述的基于多分量梯度数据联合和深度加权的相关成像方法，其特征在于，所述步骤5包括下述步骤：

步骤5.1：根据待成像地质区域的先验地质资料，得到待成像地质区域的顶面埋深和底面埋深，对待成像地质区域的地下空间进行分块处理；

步骤5.2：引入深度加权函数为

其中，z为地下网格的深度，dz为待成像地质区域地下空间范围内的最大深度值；α为经验值，α＝0.001；r为缩放因子；z_c1为待成像地质区域的顶面埋深，z_c2为待成像地质区域的底面埋深；

步骤5.3：将相关系数矩阵C与深度加权函数W(z)结合，计算加权后第j个质量元的相关系数为C_later,j＝C_j×W(z_j)，得到加权后的相关系数矩阵C_later；其中，C_j为相关系数矩阵C中第j个质量元的相关系数，W(z_j)为第j个质量元的深度加权函数值，z_j为第j个质量元所在地下网格的深度。

6.根据权利要求1所述的基于多分量梯度数据联合和深度加权的相关成像方法，其特征在于，所述步骤6中，利用加权后的相关系数矩阵C_later通过matlab中的slice函数对待成像地质区域进行三维成像。