CN111399074B - 一种重力和重力梯度模量联合三维反演方法 - Google Patents

一种重力和重力梯度模量联合三维反演方法 Download PDF

Info

Publication number
CN111399074B
CN111399074B CN202010347476.8A CN202010347476A CN111399074B CN 111399074 B CN111399074 B CN 111399074B CN 202010347476 A CN202010347476 A CN 202010347476A CN 111399074 B CN111399074 B CN 111399074B
Authority
CN
China
Prior art keywords
gravity
model
inversion
data
derivative
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN202010347476.8A
Other languages
English (en)
Other versions
CN111399074A (zh
Inventor
熊盛青
高秀鹤
于长春
张迪硕
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Aero Geophysical Survey and Remote Sensing Center for Natural Resources
Original Assignee
China Aero Geophysical Survey and Remote Sensing Center for Natural Resources
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China Aero Geophysical Survey and Remote Sensing Center for Natural Resources filed Critical China Aero Geophysical Survey and Remote Sensing Center for Natural Resources
Priority to CN202010347476.8A priority Critical patent/CN111399074B/zh
Publication of CN111399074A publication Critical patent/CN111399074A/zh
Application granted granted Critical
Publication of CN111399074B publication Critical patent/CN111399074B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V7/00—Measuring gravitational fields or waves; Gravimetric prospecting or detecting
    • G01V7/02—Details
    • G01V7/06—Analysis or interpretation of gravimetric records

Landscapes

  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明公开了一种重力和重力梯度模量联合三维反演方法，首先，通过位场转换得到重力梯度模量，将重力梯度模量和重力数据一起加入到反演中，与单独重力反演相比，能够提高反演的分辨率，与单独重力梯度模量反演相比，无需二次判断密度差的正负。同时，本发明利用重力梯度模量构建水平加权函数，对模型进行约束，进一步提高了水平分辨率。在求解最优解过程中，使用重力数据和先验地质信息，建立优化的初始搜索模型，加快了收敛速度，并提高了反演的准确性。

Description

一种重力和重力梯度模量联合三维反演方法
技术领域
本发明属于重力勘探技术领域，涉及一种重力和重力梯度模量联合三维反演方法。
背景技术
重力数据三维物性反演技术能够提供有关地质体的形状、规模以及物性大小、分布特征等信息，可为地质能源勘查提供明确、可靠的参考，因此，重力数据三维反演技术一直是研究的热点问题。与波场在一个观测点可采集多个频率的数据不同，重力场属于位场，位场在一个观测点仅能采集到一个数据，这一特征决定了观测重力数据中含有的信息量有限。因此，仅使用重力场进行三维反演建模，将导致多解性严重和分辨率低的问题。为了降低多解性、提高分辨率，最常用的做法是：根据先验地质信息引入模型约束条件，例如地质体倾向约束(Guillen A&Menichetti V,1984.Gravity and magnetic inversion withminimization of a specific functional.Geophysics,49(8):1354-1360)，这种情况下，反演结果准确性依赖先验地质信息的可靠性和丰富性。但是，在勘探程度较低、先验地质信息较少的地区，很难提取合适的先验地质信息。如何在不依赖先验地质信息的情况下，仅依靠观测重力数据，提高反演准确性和分辨率，一直是重力数据反演的技术难题。重力梯度模量数据与重力数据相比，具有更高的水平分辨率，但是，重力梯度模量反演得到的有意义的物理量是密度差的幅值，后续还需结合重力资料或先验地质资料进一步判断密度差的正负，这在一定程度上增加了反演建模的复杂性。
发明内容
本发明的目的是针对上述现有技术的不足，提出一种重力与重力梯度模量联合三维反演方法。
发明思想：充分开发观测重力数据中的有效信息，通过位场转换得到重力梯度模量，与重力数据相比，重力梯度模量具有水平分辨率高的特点。本发明充分利用重力梯度模量水平分辨率高的特点，一是将重力梯度模量视为观测数据，和原始观测重力数据联合反演，这在保证拟合原始重力数据的前提下，提高了反演分辨率，同时，反演得到的物理量是密度差，无需再判断正负；二是利用重力梯度模量构建水平加权函数，对模型加权，进一步提高了反演的水平分辨率；三是利用重力数据与先验地质信息，共同优化初始搜索模型，提高反演计算的收敛速度和准确性。
为实现上述目的，本发明是通过以下技术方案实现的：
一种重力和重力梯度模量联合三维反演方法，包括以下步骤：
1、对观测重力数据g进行位场转换，得到重力梯度模量M，其元素M与对应位置的重力g的关系为：
Figure BDA0002470671450000021
其中，gxz为重力在x方向导数，gyz为重力在y方向导数，gzz为重力在z方向导数。
2、利用重力梯度模量构建水平加权函数Wh，其对角线元素表达式为：
Figure BDA0002470671450000022
其中，τ为水平加权强度因子，τ越大横向加权强度越大。通过上式可以看出，同一水平位置的块体的wh数值相同，与深度z无关。
3、重力g和重力梯度模量M联合反演的目标函数为：
Figure BDA0002470671450000031
其中，
Figure BDA0002470671450000032
Figure BDA0002470671450000033
g为重力数据，gxz为重力在x方向导数，gyz为重力在y方向导数，gzz重力在z方向导数，A、Axz、Ayz、Azz是对应g、gxz、gyz、gzz的正演算子，ρ为待恢复的密度模型，Wv为模型垂向加权矩阵，Wh为模型水平加权矩阵，正则化因子
Figure BDA0002470671450000034
平衡因子
Figure BDA0002470671450000035
符号
Figure BDA0002470671450000036
表示哈达玛乘积，即两个向量对应元素相乘。
4、利用重力数据g和先验地质信息构建优化的初始搜索模型ρ0，构建方法为：将重力数据g垂直映射到地下多层模型空间，得到无量纲模型，再根据先验地质资料，将无量纲模型线性映射到密度空间，得到初始模型。再根据地质体深度的先验地质信息，对不同深度层的块体再加权，得到优化的初始模型ρ0。
5、在加权密度域求解目标函数的最优解，目标函数φ对加权密度ρW的导数：
Figure BDA0002470671450000037
其中，ρW＝WvWhρ，AW＝AWh -1Wv -1，AxzW＝AxzWh -1Wv -1，AyzW＝AyzWh -1Wv -1，AzzW＝AzzWh -1Wv -1，diag()表示将向量变换为对角矩阵。
6、利用最优化算法以及优化的初始模型即可求解目标函数最优解，即得到反演结果。
本发明的有益效果：
常规三维反演仅将原始的观测重力数据加入到三维反演中，在没有可靠先验地质信息约束的情况下，反演结果的分辨率低，本发明充分挖掘观测数据中的信息，通过位场转换得到重力梯度模量，并充分利用其分辨率高的特点，提升反演效果。
1、联合反演重力梯度模量和原始重力数据，与仅反演重力数据相比，本发明方法的反演结果的分辨率高；与仅反演重力梯度模量相比，本发明的反演结果是密度差，无需二次判断密度差的正负。
2、利用重力梯度模量构建水平加权函数，对模型加权，进一步提升了反演结果的水平分辨率。
3、利用重力数据与先验地质信息，设计优化的初始搜索模型，提高了反演计算的收敛速度快和准确性。
附图说明
图1是本发明的流程图。
图2是本发明的理论模型及其正演数据，(a)理论密度模型，(b)理论密度模型正演计算的地面重力异常。
图3是本发明的重力梯度模量，(a)由理论模型正演的重力梯度模量，(b)由重力数据位场转换计算的重力梯度模量，(c)转换重力梯度模量和正演的重力梯度模量的残差。
图4本发明由重力梯度模型构建的水平加权函数的切片，(a)水平切片z＝0.4km，(b)竖直切片y＝12km，(c)竖直切片x＝12km。
图5本发明优化的初始搜索模型的建立流程。
图6本发明重力和重力梯度模量联合三维反演方法结果；(a)水平切片z＝0.4km；(b)竖直切片y＝12km；(c)竖直切片x＝12km。
具体实施方式
一种重力和重力梯度模量联合三维反演方法的流程如图1所示。为结合实施例对本发明进行说明，建立理论模型，如图2(a)所示，双立方体模型边长为400m，其顶部埋深都为200m，密度差分别为1g/cm3、-1g/cm3。由理论模型产生的地面重力数据如图2(b)所示，共有31×24＝744个采样点，采样间距为100m。针对这个模型实例，对本发明进行详细说明。
1、对观测重力数据g进行位场转换，得到重力梯度模量M，如图3(b)所示，其元素M与对应位置的重力g的关系为：
Figure BDA0002470671450000051
模型正演计算得到的理论重力梯度模量如图3(a)所示，位场转换重力梯度模量图3(b)和正演计算的理论重力梯度模量图3(a)的残差如图3(c)所示，可以看出，位场转换得到的重力梯度模量准确性较高，可与观测重力数据一起用于反演。
2、利用重力梯度模量M构建水平加权函数Wh，其对角线元素表达式为：
Figure BDA0002470671450000052
其中，水平加权强度因子
Figure BDA0002470671450000053
由重力梯度模型构建的水平加权函数构建的水平加权函数切片如图4所示，由图可以看出，同一水平位置的块体的wh数值相同，与深度z无关。
3、重力g(如图2(b)所示)和重力梯度模量M(如图3(b)所示)联合反演的目标函数为：
Figure BDA0002470671450000061
其中，
Figure BDA0002470671450000062
Figure BDA0002470671450000063
g为重力数据，gxz为重力在x方向导数，gyz为重力在y方向导数，gzz重力在z方向导数，A、Axz、Ayz、Azz是对应g、gxz、gyz、gzz的正演算子，ρ为待恢复的密度模型，Wv为模型垂向加权矩阵，Wh为模型水平加权矩阵，正则化因子
Figure BDA0002470671450000064
平衡因子
Figure BDA0002470671450000065
符号
Figure BDA0002470671450000066
表示哈达玛乘积，即两个向量对应元素相乘。
4、利用数据g和先验地质信息构建优化的初始搜索模型ρ0，流程如图5所示，构建方法为：将重力g垂直映射到地下多层模型空间，得到无量纲模型，再根据先验地质资料，将无量纲模型线性映射到密度空间，得到密度模型。再根据地质体深度的先验地质信息，对不同深度层的块体再加权，得到优化的初始模型ρ0。
5、在加权密度域求解目标函数的最优解，目标函数φ对加权密度ρW的导数：
Figure BDA0002470671450000067
其中，ρW＝WvWhρ，AW＝AWh -1Wv -1，AxzW＝AxzWh -1Wv -1，AyzW＝AyzWh -1Wv -1，AzzW＝AzzWh -1Wv -1，diag()表示将向量变换为对角矩阵。
6、利用最优化算法以及优化的初始模型即可求解目标函数最优解，即得到反演结果，如图6所示，异常体位置识别准确，边界较为清晰，最大密度差为0.8668g/cm3，最小密度差为-0.8668g/cm3。

Claims (1)

1.一种重力和重力梯度模量联合三维反演方法，其特征在于：包括以下步骤：
1)、对观测重力数据g进行位场转换，得到重力梯度模量M，其元素M与对应位置的重力g的关系为：
Figure FDA0002470671440000011
其中，gxz为重力在x方向导数，gyz为重力在y方向导数，gzz为重力在z方向导数；
2)、利用重力梯度模量构建水平加权函数Wh，其对角线元素表达式为：
Figure FDA0002470671440000012
其中，τ为水平加权强度因子，τ越大横向加权强度越大；通过上式可以看出，同一水平位置的块体的wh数值相同，与深度z无关；
3)、重力g和重力梯度模量M联合反演的目标函数为：
Figure FDA0002470671440000013
其中，
Figure FDA0002470671440000014
Figure FDA0002470671440000015
g为重力数据，gxz为重力在x方向导数，gyz为重力在y方向导数，gzz重力在z方向导数，A、Axz、Ayz、Azz是对应g、gxz、gyz、gzz的正演算子，ρ为待恢复的密度模型，Wv为模型垂向加权矩阵，Wh为模型水平加权矩阵，正则化因子
Figure FDA0002470671440000016
平衡因子
Figure FDA0002470671440000017
符号
Figure FDA0002470671440000018
表示哈达玛乘积，即两个向量对应元素相乘；
4)、利用重力数据g和先验地质信息构建优化的初始搜索模型ρ0，构建方法为：将重力数据g垂直映射到地下多层模型空间，得到无量纲模型，再根据先验地质资料，将无量纲模型线性映射到密度空间，得到初始模型；再根据地质体深度的先验地质信息，对不同深度层的块体再加权，得到优化的初始模型ρ0；
5)、在加权密度域求解目标函数的最优解，目标函数φ对加权密度ρW的导数：
Figure FDA0002470671440000021
其中，ρW＝WvWhρ，AW＝AWh -1Wv -1，AxzW＝AxzWh -1Wv -1，AyzW＝AyzWh -1Wv -1，AzzW＝AzzWh .1Wv -1，diag()表示将向量变换为对角矩阵；
6)、利用最优化算法以及优化的初始模型即可求解目标函数最优解，即得到反演结果。
CN202010347476.8A 2020-04-28 2020-04-28 一种重力和重力梯度模量联合三维反演方法 Expired - Fee Related CN111399074B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010347476.8A CN111399074B (zh) 2020-04-28 2020-04-28 一种重力和重力梯度模量联合三维反演方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010347476.8A CN111399074B (zh) 2020-04-28 2020-04-28 一种重力和重力梯度模量联合三维反演方法

Publications (2)

Publication Number Publication Date
CN111399074A CN111399074A (zh) 2020-07-10
CN111399074B true CN111399074B (zh) 2022-08-12

Family

ID=71435371

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010347476.8A Expired - Fee Related CN111399074B (zh) 2020-04-28 2020-04-28 一种重力和重力梯度模量联合三维反演方法

Country Status (1)

Country Link
CN (1) CN111399074B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112147709B (zh) * 2020-08-03 2022-07-29 中国海洋石油集团有限公司 一种基于部分光滑约束的重力梯度数据三维反演方法
CN113514900B (zh) * 2021-07-12 2022-05-17 吉林大学 基于密度约束的球坐标系重力和重力梯度联合反演方法
CN113552637B (zh) * 2021-07-30 2023-11-14 中国自然资源航空物探遥感中心 一种航空-地面-井中磁异常数据协同三维反演方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2003023447A2 (en) * 2001-09-07 2003-03-20 Conocophillips Company A nonlinear constrained inversion method to determine base of salt interface from gravity and gravity tensor data
CN109471190A (zh) * 2018-11-12 2019-03-15 吉林大学 一种不同高度重力数据和井中重力数据联合反演方法
CN110133716A (zh) * 2019-06-03 2019-08-16 吉林大学 基于组合模型加权函数的磁异常数据三维反演方法
CN110398782A (zh) * 2019-07-17 2019-11-01 广州海洋地质调查局 一种重力数据和重力梯度数据联合正则化反演方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2003023447A2 (en) * 2001-09-07 2003-03-20 Conocophillips Company A nonlinear constrained inversion method to determine base of salt interface from gravity and gravity tensor data
CN109471190A (zh) * 2018-11-12 2019-03-15 吉林大学 一种不同高度重力数据和井中重力数据联合反演方法
CN110133716A (zh) * 2019-06-03 2019-08-16 吉林大学 基于组合模型加权函数的磁异常数据三维反演方法
CN110398782A (zh) * 2019-07-17 2019-11-01 广州海洋地质调查局 一种重力数据和重力梯度数据联合正则化反演方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于阈值约束的协克里金法联合反演重力与重力梯度数据;高秀鹤;《地球物理学报》;20190331;第62卷(第3期);全文 *

Also Published As

Publication number Publication date
CN111399074A (zh) 2020-07-10

Similar Documents

Publication Publication Date Title
CN111399074B (zh) 一种重力和重力梯度模量联合三维反演方法
CN105158808B (zh) 一种浅海瞬变电磁海空探测及其解释方法
Mallick et al. A finite element method for computation of the regional gravity anomaly
CN110133716B (zh) 基于组合模型加权函数的磁异常数据三维反演方法
CN109657314B (zh) 基于资料循环同化和模式循环积分的台风涡旋初始化方法
CN102901985B (zh) 一种适用于起伏地表的深度域层速度修正方法
Jia et al. Deep learning for 3-D magnetic inversion
CN102841375A (zh) 一种复杂条件下基于角度域共成像点道集的层析速度反演方法
CN110261902B (zh) 一种基于多谱能量合成的地下浅层震源定位方法
CN106886047A (zh) 一种接收函数和重力联合反演地壳厚度和波速比的方法
CN110261903B (zh) 一种基于逆时能量聚焦的地下震源被动定位方法
CN110007365A (zh) 一种基于信号数据稀疏空间快速计算的联合反演方法
CN106846246B (zh) 一种基于对象的遥感影像超分辨率制图方法
CN108845353B (zh) 一种重震联合反演的方法及装置
Li et al. Self-supervised deep learning for 3D gravity inversion
CN115586577A (zh) 一种定源瞬变电磁非中心点观测数据全时转换方法
CN108008456B (zh) 一种圈定热液型铀矿深部三维重点铀成矿有利靶区的方法
CN109100798B (zh) 实现折射多次波层析反演的方法、装置及处理终端
CN104155690B (zh) 基于椭球展开的三维地震数据叠加速度求取方法
Yang et al. Crustal structure of the Dabie orogenic belt (eastern China) inferred from gravity and magnetic data
CN109752762A (zh) 单发多收观测装置瞬变电场数据动校正方法及装置
CN110490354B (zh) 一种基于roms模拟结果计算潮汐锋位置的优化方法
CN117950068A (zh) 一种稀疏域可控源电磁反演方法
CN113777654B (zh) 一种基于伴随状态法初至波走时层析的海水速度建模方法
CN117933049A (zh) 基于深度学习的航空瞬变电磁三维反演方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20220812

CF01 Termination of patent right due to non-payment of annual fee