CN110389391B - 一种基于空间域的重磁位场解析延拓方法 - Google Patents

一种基于空间域的重磁位场解析延拓方法 Download PDF

Info

Publication number
CN110389391B
CN110389391B CN201910708065.4A CN201910708065A CN110389391B CN 110389391 B CN110389391 B CN 110389391B CN 201910708065 A CN201910708065 A CN 201910708065A CN 110389391 B CN110389391 B CN 110389391B
Authority
CN
China
Prior art keywords
continuation
plane
curved surface
extension
observation
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.)
Active
Application number
CN201910708065.4A
Other languages
English (en)
Other versions
CN110389391A (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.)
Second Institute of Oceanography MNR
Original Assignee
Second Institute of Oceanography MNR
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 Second Institute of Oceanography MNR filed Critical Second Institute of Oceanography MNR
Priority to CN201910708065.4A priority Critical patent/CN110389391B/zh
Publication of CN110389391A publication Critical patent/CN110389391A/zh
Application granted granted Critical
Publication of CN110389391B publication Critical patent/CN110389391B/zh
Active 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
    • G01V11/00—Prospecting or detecting by methods combining techniques covered by two or more of main groups G01V1/00 - G01V9/00
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/38—Processing data, e.g. for analysis, for interpretation, for correction
    • 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

  • Life Sciences & Earth Sciences (AREA)
  • Physics & Mathematics (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Engineering & Computer Science (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • Remote Sensing (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明公开了一种基于空间域的重磁位场解析延拓方法，重磁位场指的是重力场或者磁场；解析延拓指的是通过观测面上的重力场或者磁场计算其它曲面(包括平面)上的重力场或磁场的方法。解析延拓的求解有空间域和频率域两种思路，本发明是针对空间域的一种解析延拓新方法，可以直接对航磁测量的磁异常数据向下延拓至地面增强磁异常信号，也可以将地面观测数据向上延拓至更高的平面上以突出深部场特征。本发明可以将任意曲面观测的异常数据延拓至另一个任意曲面，突破了频率域方法中延拓面一定是平面的局限，解决了空间域解析延拓方法中的奇异点问题。

Description

一种基于空间域的重磁位场解析延拓方法
技术领域
本发明属于重力、磁异常数据处理技术领域，涉及一种基于空间域的重磁位场解析延拓方法。
背景技术
重磁勘探方法在地面、航空、海面及海底进行重力场或磁场观测，然后将原始数据经过一定校正得到重力异常或磁异常，反应的是观测面以下的密度和磁性不均匀性，从而应用于地质填图、矿藏勘探、管线探测等。重力勘探和磁法勘探这两种方法具有一定的相似性，故通常将重力异常和磁异常一起讨论，称之为重磁异常。重力场和磁场均是位场，故又称之为重磁位场。而重磁位场的解析延拓是对重磁异常数据处理的重要方法，通过解析延拓可以反映出重磁异常的多尺度信息。向上延拓可以突出深部场源信息，压制浅部场源的干扰；相反地，向下延拓则是增强浅部场源信息。
传统的解析延拓方法是基于快速傅里叶变换在频率域进行计算，这种方法的优势是计算速度快，但是其缺点是只要求重磁异常的观测面是平面。空间域方法可以弥补此缺点，但传统的空间域延拓方法只适用于延拓高度大于一个点距的情况，而且不适用于计算曲面之间的延拓。传统空间域延拓方法在计算较小延拓高度的情况时，会出现奇异值。但是通常的观测面均是曲面，比如航磁观测的飞行器大多是定高飞行，换言之，航磁观测曲面近似于地形平行。这种情况下的解析延拓都会遇到延拓高度小于一个点距和曲面之间延拓的问题。
发明内容
为了解决上述一般情况下传统解析延拓方法的不适用问题，就需要从解析延拓的理论公式中寻找解决方法。本发明通过对解析延拓的二重积分公式进行分段积分，求解得到了空间域解析延拓的一个新的核函数。通过新的核函数可以方便地实现平面到平面的向上延拓、平面到曲面的向上延拓、平面到平面的向下延拓、曲面到平面的向下延拓以及场源外的两个曲面之间的延拓。
本发明通过以下技术方案得以实现：一种基于空间域的重磁位场解析延拓方法，该方法包括以下步骤：
(a)观测面S(x,y)为平面，其高程为z0是常数，在S(x,y)上观测的重力或磁异常记为U(α,β,z0)，其中α,β分别为观测点在x方向和y方向上的坐标。
(b)将观测的重力或磁异常数据U(α,β,z0)向上延拓到延拓面S(x,y,z1)，z1(x,y)表示延拓面S(x,y,z1)的高程；延拓高度表示为Δz(x,y)＝z1(x,y)-z0，则U(α,β,z0)向上延拓结果U(x,y,z1)用积分表示为：
Figure BDA0002152804690000021
其中，
Figure BDA0002152804690000022
表示观测点(α,β)与延拓点(x,y)之间的距离，延拓点是观测点在延拓面上的垂直投影。
(c)对观测面S(x,y)上的位场数据进行网格化，得到网格化数据，总共有M行N列，点距分别为Δx和Δy，第i行第j列节点上的位场异常值表示为Ui,j、坐标为(αi,βi)。
(d)由于观测数据的有限性，结合步骤(c)，将公式(1)近似为观测数据范围内的分段积分的累加：
Figure BDA0002152804690000023
(e)由于观测数据的分辨率限制，假设观测数据在一个点距范围内是常数，则将公式(2)表示为：
Figure BDA0002152804690000024
其中Um,n表示延拓结果的第m行第n列的值；
(f)将公式(3)中的二重积分表示为公式(4)并求解得到一个反正切形式的精确解Pm,n；i,j：
Figure BDA0002152804690000025
Figure BDA0002152804690000031
(g)结合公式(3)和(4)将向上延拓表示为离散观测数据的加权求和，第m行第n列延拓结果表示为：
Figure BDA0002152804690000032
(h)将公式(5)用矩阵乘法表示为：
U＝KU0 (6)
其中U0表示位场观测数据，U为延拓结果，K是向上延拓的核矩阵，其元素由公式(4)计算得到，表示为：
Figure BDA0002152804690000033
对于平面到平面的向上延拓：如果延拓面是平面，Δz(x,y)＝h为延拓高度，是一个常数；则矩阵K所有元素均可以用第一行元素计算得到，只需要存储第一行元素即可；计算公式(6)所示的矩阵乘法即可得到平面到平面的向上延拓结果；
对于平面到曲面的向上延拓：如果延拓面是曲面，其高程为z1(x,y)，则延拓高度Δz(x,y)可以由延拓点高程与对应的观测点高程相减得到，然后带入公式(4)即可计算得到矩阵K，计算公式(6)所示的矩阵乘法即可得到平面到曲面的向上延拓结果。
进一步地，该方法还可以用于实现向下延拓，向下延拓是向上延拓的反问题，对于平面到平面的向下延拓及曲面到平面的向下延拓：
根据延拓高度和公式(4)和公式(7)计算核矩阵K，然后利用Landweber迭代法和光滑-拟合曲线求解式(6)得到向下延拓结果。
进一步地，所述利用Landweber迭代法和光滑-拟合曲线求解式(6)具体为：
Landweber迭代法的最佳迭代次数通过光滑-拟合曲线的第一个极小值确定；光滑-拟合曲线中的“光滑”指的是向下延拓结果的水平梯度的二范数，作为纵坐标；“拟合”指的是向下延拓结果与观测数据之差的二范数与观测数据二范数的比值，作为横坐标；对不同迭代次数的结果均可计算出“光滑”和“拟合”值，然后绘制曲线，找出第一个极小值点对应的迭代次数就是最优的迭代次数，最优迭代次数对应的向下延拓结果就是所求解的向下延拓结果。
进一步地，该方法还可以用于实现曲面到曲面之间的解析延拓：
曲面到曲面之间的解析延拓包括曲面到曲面的向下延拓及曲面到曲面的向上延拓；如果观测面是曲面，延拓面也是曲面，则首先需要求解观测曲面到最低的一个平面上的延拓场；如果是曲面到曲面的向上延拓，这个最低的平面取观测曲面的最低点所在平面；如果是曲面到曲面的向下延拓，这个最低的平面取为延拓曲面的最低点所在平面；将观测的位场向下延拓至最低平面得到等效场Ue；然后计算最低平面上的等效场Ue到延拓曲面的向上延拓，最终实现曲面到曲面的解析延拓，这个过程的理论依据是重磁位场的等效源原理。
本发明的有益效果是：本发明方法可以直接对航磁测量的磁异常数据向下延拓至地面增强磁异常信号，也可以将地面观测数据向上延拓至更高的平面上以突出深部场特征。本发明可以将任意曲面观测的异常数据延拓至另一个任意曲面，突破了频率域方法中延拓面一定是平面的局限，解决了空间域解析延拓方法中的奇异点问题。
附图说明
图1是本发明方法的流程图。
图2是本发明方法的网格节点及积分示意图。
图3是本发明方法模型实验所用的理论模拟数据；(a)为模拟的地形，黑色方框为立方体模型在平面上的投影；(b)为高程z＝0平面上的磁异常；(c)为图(a)所示的模拟地形上的磁异常；(d)为高程z＝8m平面上的磁异常；(e)为这三个磁异常数据和地形以及模型的三维视图。
图4是本发明方法对磁异常从平面到平面的向上延拓实例结果；(a)～(c)是用频率域方法计算的结果；(d)～(f)是本发明方法的计算结果；(a)和(d)是对图3(b)所示的磁异常向上延拓至z＝8m的平面的结果；(b)和(e)分别为向上延拓结果(a)和(d)与图3(d)所示的理论数据之间的误差分布；(c)和(f)分别是对误差分布(b)和(e)进行统计的柱状图，误差均方根(RMSE)分别为2.13和3.40nT。
图5是本发明方法对磁异常从平面向上延拓至曲面的实例结果；(a)是对图3(a)所示的磁异常向上延拓至图3(a)所示的模拟地形上的结果；(b)为向上延拓结果(a)与图3(c)所示的理论模拟数据之间的误差分布；(c)为对误差分布(b)进行统计的柱状图，误差均方根(RMSE)为1.76nT。
图6是本发明方法对磁异常从平面向下延拓至平面的结果；(a)是对图3(d)所示的磁异常向下延拓至z＝0m平面的结果；(b)为向下延拓结果(a)与图3(b)所示的理论模拟磁异常之间的误差分布；(c)为对误差分布(b)进行统计的柱状图，误差均方根(RMSE)为3.4nT。
图7是本发明方法对磁异常从曲面向下延拓至平面的结果；(a)是对图3(c)所示的曲面上的磁异常向下延拓至z＝0m平面的结果；(b)向下延拓结果与图3(b)所示的理论模拟数据之间的误差分布；(c)对误差分布进行统计的柱状图，均方根为0.69nT。
图8是本发明方法对含有噪声的磁异常从曲面向下延拓至平面的结果；(a)对图3(c)所示的曲面上的磁异常加入高斯白噪声(噪声均值为0nT，标准差为1nT)后的异常分布；(b)白噪声的分布；(c)对此白噪声进行统计的柱状图；(d)将(a)所示的含有噪声的曲面上的磁异常向下延拓至z＝0m平面的结果；(e)向下延拓结果与图3(b)所示的理论模型数据之间的误差分布；(f)对误差分布进行统计的柱状图，均方差为5.29nT。
图9是本发明在计算磁异常向下延拓(对应图8)的最优迭代次数选择曲线。
图10是本发明方法对实测航磁异常进行解析延拓的结果。(a)航磁异常；(b)航磁观测的飞行高度；(c)对观测的航次异常向下延拓至高程为970m的平面上的结果；(d)对(c)所示的延拓延拓场向上延拓至高程z＝2200m的结果，表示对观测磁异常的一个曲化平(从曲面到平面的延拓)的结果；(e)将(c)所示延拓场向上延拓至观测曲面之上200m曲面(与观测面平行)上的结果，表示从曲面到曲面的向上延拓。
具体实施方式
下面结合实施例对本发明的可行性做进一步的详细说明，以下实施例是对本发明的解释而本发明并不局限于以下实施例的结果。
本发明提供的一种基于空间域的重磁位场解析延拓方法，可实现平面到平面的向上延拓、平面到曲面的向上延拓、平面到平面的向下延拓、曲面到平面的向下延拓以及场源外的两个曲面之间的延拓；
(1)向上延拓，包括如下步骤：
(a)观测面S(x,y)为平面，其高程为z0是常数，在S(x,y)上观测的重力或磁异常(以下简称为位场)记为U(α,β,z0)，其中α,β分别为观测点在x方向和y方向上的坐标；
(b)将观测的重力或磁异常数据U(α,β,z0)向上延拓到延拓面S(x,y,z1)，z1(x,y)表示延拓面S(x,y,z1)的高程。延拓高度可以表示为Δz(x,y)＝z1(x,y)-z0，则U(α,β,z0)向上延拓结果U(x,y,z1)可以用积分表示为：
Figure BDA0002152804690000061
其中，
Figure BDA0002152804690000062
表示观测点(α,β)与延拓点(x,y)之间的距离，延拓点是观测点在延拓面上的垂直投影；
(c)公式(1)为连续积分，但是实际观测数据是离散点，所以需要对观测面S(x,y)上的位场数据进行网格化，得到网格化数据，总共有M行N列，点距分别为Δx和Δy，第i行第j列节点上的位场异常值表示为Ui,j、坐标为(αi,βi)；
(d)由于观测数据的有限性，结合步骤(c)，则可以将公式(1)近似为观测数据范围内的分段积分的累加：
Figure BDA0002152804690000063
(e)由于观测数据的分辨率限制，假设观测数据在一个点距范围内是常数，则公式(2)可以表示为：
Figure BDA0002152804690000071
其中Um,n表示延拓结果的第m行第n列的值；
(f)公式(3)中的二重积分是第i行第j列观测点和第m行第n列延拓点的坐标的函数，因此可以将此二重积分表示为公式(4)并求解得到一个反正切形式的精确解Pm,n；i,j：
Figure BDA0002152804690000072
(g)结合公式(3)和(4)可以将向上延拓表示为离散观测数据的加权求和，第m行第n列延拓结果可以表示为：
Figure BDA0002152804690000073
其中Pm,n；i,j可视为加权系数；
(h)将公式(5)可以进一步地用矩阵乘法表示为：
U＝KU0 (6)
其中U0表示位场观测数据，是一个M×N的列向量；U为延拓结果，也是一个M×N的列向量；K是向上延拓的核矩阵，是一个(M×N)×(M×N)的方阵，其元素由公式(4)计算得到，表示为：
Figure BDA0002152804690000074
(1.1)平面到平面的向上延拓
延拓面是平面的情况下，Δz(x,y)＝h为延拓高度，是一个常数。则矩阵K(公式7)所有元素都可以用第一行元素计算得到，所以只需要存储第一行元素即可，减少计算机内存消耗。计算公式(6)所示的矩阵乘法即可得到平面到平面的向上延拓结果；
(1.2)平面到曲面的向上延拓
如果延拓面是曲面，其高程为z1(x,y)，则延拓高度Δz(x,y)可以由延拓点高程与对应的观测点高程相减得到，然后带入公式(4)即可计算得到矩阵K，然后计算公式(6)所示的矩阵乘法即可得到平面到曲面的向上延拓结果；
(2)向下延拓，包含如下步骤：
(A)向下延拓是向上延拓的反问题，上述的向上延拓原理是向下延拓的基础。如果U是观测数据，则求解式(6)表示的线性方程组，可以得到向下延拓结果U0；
(B)采用Landweber迭代法求解线性方程组(6)，最佳迭代次数可以通过光滑-拟合曲线的第一个极小值确定。光滑-拟合曲线中的“光滑”指的是向下延拓结果的水平梯度的二范数，作为纵坐标；“拟合”指的是向下延拓结果与观测数据之差的二范数与观测数据二范数的比值，作为横坐标。对不同迭代次数的结果均可计算出“光滑”和“拟合”值，然后绘制曲线，找出第一个极小值点对应的迭代次数就是最优的迭代次数，最优迭代次数对应的向下延拓结果就是所求解的向下延拓结果。
平面到平面的向下延拓及曲面到平面的向下延拓：
根据延拓高度和公式(4)和公式(7)计算核矩阵K，然后利用Landweber迭代法和光滑-拟合曲线求解线性方程组(6)得到最优向下延拓结果；
(3)曲面到曲面之间的解析延拓
曲面到曲面之间的解析延拓包括曲面到曲面的向下延拓及曲面到曲面的向上延拓。如果观测面是曲面，延拓面也是曲面，则首先需要求解观测曲面到最低的一个平面上的延拓场。如果是曲面到曲面的向上延拓，这个最低的平面取观测曲面的最低点所在平面；如果是曲面到曲面的向下延拓，这个最低的平面取为延拓曲面的最低点所在平面。与“平面到平面的向下延拓及曲面到平面的向下延拓”相同的方法，将观测的位场向下延拓至最低平面得到等效场Ue；然后利用与(1.2)相同的方法计算最低平面上的等效场Ue到延拓曲面的向上延拓，最终实现曲面到曲面的解析延拓。
实施示例1：对模拟计算的磁异常进行解析延拓
利用立方体模型的磁异常正演公式分别计算得到一个立方体(图3e)在平面z＝0m，平面z＝8m以及一个曲面(图3a)上的磁异常，三个磁异常如图3(b)、(d)和(c)所示。利用上述的重磁位场解析延拓方法分别对模拟磁异常数据进行向上和向下延拓，然后与精确解进行对比，验证方法的可行性和准确性。
(1)磁异常从一个平面向上延拓至另一平面
对图3b所示的z＝0m平面上的磁异常数据向上延拓8m，结果如图4a所示，将向上延拓结果与z＝8m平面上的真实磁异常作差计算其误差，如图4b所示；对误差进行统计，如图4c所示。均方差为2.13nT。图4d-f为频率域计算结果，对比可发现空间域的计算结果更精确。
(2)磁异常从平面向上延拓至曲面
对图3b所示的z＝0m平面上的磁异常向上延拓至图3a所示的曲面上，结果如图5a所示，将向上延拓结果与曲面上的精确磁异常进行作差，计算误差分布如图5b所示。对误差统计得到均方差为1.76nT，如图5c所示。
(3)磁异常从平面向下延拓至平面
对图3d所示的z＝8m平面上的磁异常向下延拓8m，结果如图6a所示。将向下延拓结果与z＝0m平面上的精确磁异常进行作差，计算其误差分布，如图6b所示。对误差统计的到均方差为2.4nT，如图6c所示。
(4)磁异常从曲面向下延拓至平面
对图3c所示的曲面上的磁异常向下延拓至z＝0m平面上，延拓结果如图7a所示。将延拓结果与z＝0m平面精确磁异常作差，计算其误差分布，如图7b所示。对误差统计得到均方差为0.69nT，如图7c所示。
(5)对含有噪声的磁异常进行向下延拓
在图3c所示的磁异常中加入均值为0nT，标准差为1nT的高斯白噪声(如图8b所示，噪声的统计如图8c所示，加入白噪声的磁异常结果如图8a所示。对图8a所示的含噪声磁异常向下延拓至z＝0m平面，结果如图8d所示。对结果和z＝0m平面精确磁异常进行作差计算其误差分布，如图8e所示。对误差进行统计得到均方差为5.29nT，如图8f所示。
(6)最优迭代次数的选择
本发明在计算向下延拓时，需要进行迭代，最优迭代次数的选择依据提出的拟合-光滑曲线的第一个极小值点。以含噪声数据的向下延拓(图8)为例，其向下延拓最优迭代次数的选择所依据的拟合-光滑曲线如图9所示，最优迭代次数为277。
实施示例2：实测航磁数据的解析延拓
将本发明应用到航磁异常的解析延拓处理中，航磁异常如图10a所示，飞行曲面的高程如图10b所示。对于曲面上观测的磁异常，进一步处理之前，需要将其向下延拓至一个平面上。这个平面一般选择为飞行高度曲面的最低点所在平面，这里选择为z＝970m。将航磁异常延拓到z＝970m平面的结果如图10c所示，此向下延拓的高度并非一个常数，而是0～1430m变化的，传统的空间域方法无法实现。然后利用此中间结果进行其他平面之间的延拓处理，对图10c所示的磁异常向上延拓至z＝2200m平面得到航磁异常的一个曲化平的结果，如图10d所示。有些情况下，可能并不需要曲化平，而只需要将曲面上的磁异常延拓至另一个曲面上，比如延拓至于观测曲面平行的一个曲面。将图10c所示的结果向上延拓至与飞行曲面之上(平行)200m的结果如图10e所示。
以上两个实施示例，证明了新的空间域延拓方法在实际应用中的可行性和准确性。对重磁位场异常进行解析延拓或者曲化平，对异常进一步解释和处理具有重要实用价值。

Claims (4)

1.一种基于空间域的重磁位场解析延拓方法，其特征在于，该方法包括以下步骤：
(a)观测面S(x,y)为平面，其高程为z0是常数，在S(x,y)上观测的重力或磁异常记为U(α,β,z0)，其中α,β分别为观测点在x方向和y方向上的坐标；
(b)将观测的重力或磁异常数据U(α,β,z0)向上延拓到延拓面S(x,y,z1)，z1(x,y)表示延拓面S(x,y,z1)的高程；延拓高度表示为Δz(x,y)＝z1(x,y)-z0，则U(α,β,z0)向上延拓结果U(x,y,z1)用积分表示为：
Figure FDA0002616101430000011
其中，
Figure FDA0002616101430000012
表示观测点(α,β)与延拓点(x,y)之间的距离，延拓点是观测点在延拓面上的垂直投影；
(c)对观测面S(x,y)上的位场数据进行网格化，得到网格化数据，总共有M行N列，点距分别为Δx和Δy，第i行第j列节点上的位场异常值表示为Ui,j、坐标为(αi,βi)；
(d)由于观测数据的有限性，结合步骤(c)，将公式(1)近似为观测数据范围内的分段积分的累加：
Figure FDA0002616101430000013
(e)由于观测数据的分辨率限制，假设观测数据在一个点距范围内是常数，则将公式(2)表示为：
Figure FDA0002616101430000014
其中Um,n表示延拓结果的第m行第n列的值；
(f)将公式(3)中的二重积分表示为公式(4)并求解得到一个反正切形式的精确解Pm,n；i,j：
Figure FDA0002616101430000015
Figure FDA0002616101430000021
(g)结合公式(3)和(4)将向上延拓表示为离散观测数据的加权求和，第m行第n列延拓结果表示为：
Figure FDA0002616101430000022
(h)将公式(5)用矩阵乘法表示为：
U＝KU0 (6)
其中U0表示位场观测数据，U为延拓结果，K是向上延拓的核矩阵，其元素由公式(4)计算得到，表示为：
Figure FDA0002616101430000023
对于平面到平面的向上延拓：如果延拓面是平面，Δz(x,y)＝h为延拓高度，是一个常数；则矩阵K所有元素均用第一行元素计算得到，只需要存储第一行元素即可；计算公式(6)所示的矩阵乘法即可得到平面到平面的向上延拓结果；
对于平面到曲面的向上延拓：如果延拓面是曲面，其高程为z1(x,y)，则延拓高度Δz(x,y)由延拓点高程与对应的观测点高程相减得到，然后带入公式(4)即可计算得到矩阵K，计算公式(6)所示的矩阵乘法即可得到平面到曲面的向上延拓结果。
2.根据权利要求1所述的一种基于空间域的重磁位场解析延拓方法，其特征在于，该方法还用于实现向下延拓，向下延拓是向上延拓的反问题，对于平面到平面的向下延拓及曲面到平面的向下延拓：
根据延拓高度和公式(4)和公式(7)计算核矩阵K，然后利用Landweber迭代法和光滑-拟合曲线求解式(6)得到向下延拓结果。
3.根据权利要求2所述的一种基于空间域的重磁位场解析延拓方法，其特征在于，所述利用Landweber迭代法和光滑-拟合曲线求解式(6)具体为：
Landweber迭代法的最佳迭代次数通过光滑-拟合曲线的第一个极小值确定；光滑-拟合曲线中的“光滑”指的是向下延拓结果的水平梯度的二范数，作为纵坐标；“拟合”指的是向下延拓结果与观测数据之差的二范数与观测数据二范数的比值，作为横坐标；对不同迭代次数的结果均可计算出“光滑”和“拟合”值，然后绘制曲线，找出第一个极小值点对应的迭代次数就是最优的迭代次数，最优迭代次数对应的向下延拓结果就是所求解的向下延拓结果。
4.根据权利要求1所述的一种基于空间域的重磁位场解析延拓方法，其特征在于，该方法还用于实现曲面到曲面之间的解析延拓：
曲面到曲面之间的解析延拓包括曲面到曲面的向下延拓及曲面到曲面的向上延拓；如果观测面是曲面，延拓面也是曲面，则首先需要求解观测曲面到最低的一个平面上的延拓场；如果是曲面到曲面的向上延拓，这个最低的平面取观测曲面的最低点所在平面；如果是曲面到曲面的向下延拓，这个最低的平面取为延拓曲面的最低点所在平面；将观测的位场向下延拓至最低平面得到等效场Ue；然后计算最低平面上的等效场Ue到延拓曲面的向上延拓，最终实现曲面到曲面的解析延拓。
CN201910708065.4A 2019-08-01 2019-08-01 一种基于空间域的重磁位场解析延拓方法 Active CN110389391B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910708065.4A CN110389391B (zh) 2019-08-01 2019-08-01 一种基于空间域的重磁位场解析延拓方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910708065.4A CN110389391B (zh) 2019-08-01 2019-08-01 一种基于空间域的重磁位场解析延拓方法

Publications (2)

Publication Number Publication Date
CN110389391A CN110389391A (zh) 2019-10-29
CN110389391B true CN110389391B (zh) 2020-12-15

Family

ID=68288271

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910708065.4A Active CN110389391B (zh) 2019-08-01 2019-08-01 一种基于空间域的重磁位场解析延拓方法

Country Status (1)

Country Link
CN (1) CN110389391B (zh)

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111859251B (zh) * 2020-06-29 2021-06-15 中国地质大学（武汉） 一种基于pde的磁测数据等效源上延拓与下延拓方法
CN111856599B (zh) * 2020-06-29 2021-08-10 中国地质大学（武汉） 一种基于pde的磁测数据等效源化极与类型转换方法
CN111856598B (zh) * 2020-06-29 2021-06-15 中国地质大学（武汉） 一种磁测数据多层等效源上延拓与下延拓方法
CN111880236B (zh) * 2020-06-29 2022-02-18 中国地质大学（武汉） 一种构建多层等效源模型计算化极与数据类型转换的方法
CN111856615B (zh) * 2020-07-30 2023-01-24 中国人民解放军61540部队 重磁异常场重构方法及场源点物性参数解算方法、系统
CN112363236B (zh) * 2020-10-15 2022-04-01 中国地质大学（武汉） 一种基于pde的重力场数据等效源延拓与数据类型转换方法
CN112748471B (zh) * 2020-12-29 2022-06-14 中国地质大学（武汉） 一种非结构化等效源的重磁数据延拓与转换方法
CN113111500A (zh) * 2021-03-31 2021-07-13 中国地质大学（北京） 一种基于重磁物理场使用二维经验模态分解异常分析方法
CN114167511B (zh) * 2021-11-26 2023-08-11 兰州大学 一种基于连分式展开向下延拓的位场数据快速反演方法
CN115356784A (zh) * 2022-08-29 2022-11-18 西南交通大学 自适应阻尼系数的广义极小残差大深度位场向下延拓方法

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102419453A (zh) * 2011-07-15 2012-04-18 中国科学院地质与地球物理研究所 长导线源瞬变电磁地空探测方法
CN106772581A (zh) * 2016-12-13 2017-05-31 西京学院 一种基于重构技术的三维起伏地表物理模拟采集方法
CN107632964A (zh) * 2017-09-06 2018-01-26 哈尔滨工程大学 一种平面地磁异常场向下延拓递归余弦变换法
CN108594319A (zh) * 2018-05-11 2018-09-28 中国人民解放军61540部队 一种航空重力测量数据向下延拓方法及系统
CN109085652A (zh) * 2018-08-03 2018-12-25 吉林大学 基于改进迭代法的地空时间域电磁系统高精度下延拓方法
CN109471190A (zh) * 2018-11-12 2019-03-15 吉林大学 一种不同高度重力数据和井中重力数据联合反演方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8093893B2 (en) * 2004-03-18 2012-01-10 Baker Hughes Incorporated Rock and fluid properties prediction from downhole measurements using linear and nonlinear regression
US20100195434A1 (en) * 2009-01-30 2010-08-05 Conocophillips Company Heterodyned Seismic Source

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102419453A (zh) * 2011-07-15 2012-04-18 中国科学院地质与地球物理研究所 长导线源瞬变电磁地空探测方法
CN106772581A (zh) * 2016-12-13 2017-05-31 西京学院 一种基于重构技术的三维起伏地表物理模拟采集方法
CN107632964A (zh) * 2017-09-06 2018-01-26 哈尔滨工程大学 一种平面地磁异常场向下延拓递归余弦变换法
CN108594319A (zh) * 2018-05-11 2018-09-28 中国人民解放军61540部队 一种航空重力测量数据向下延拓方法及系统
CN109085652A (zh) * 2018-08-03 2018-12-25 吉林大学 基于改进迭代法的地空时间域电磁系统高精度下延拓方法
CN109471190A (zh) * 2018-11-12 2019-03-15 吉林大学 一种不同高度重力数据和井中重力数据联合反演方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
多层等效源曲面磁异常转换方法;李端 等;《地球物理学报》;20180731;第61卷(第7期);第3055-3073页 *
航磁异常深部弱信号提取技术研究;郭志馗 等;《地质与勘探》;20151130;第51卷(第6期);第1007-1014页 *

Also Published As

Publication number Publication date
CN110389391A (zh) 2019-10-29

Similar Documents

Publication Publication Date Title
CN110389391B (zh) 一种基于空间域的重磁位场解析延拓方法
JP6896180B2 (ja) 風流検知システム及び風流の速度場を求める方法
Wang et al. GPU-based, parallel-line, omni-directional integration of measured pressure gradient field to obtain the 3D pressure distribution
KR20180091372A (ko) 레이더의 표적 위치 추적 방법
Liu et al. A novel entropy-based method to quantify forest canopy structural complexity from multiplatform lidar point clouds
Nekkanti et al. Gappy spectral proper orthogonal decomposition
CN109458994A (zh) 一种空间非合作目标激光点云icp位姿匹配正确性判别方法及系统
WO2022232572A1 (en) Method and system for high resolution least-squares reverse time migration
CN110686610B (zh) 基于自适应网格的光学变形测量方法及电子设备
CN111665556B (zh) 地层声波传播速度模型构建方法
CN114896564A (zh) 采用自适应泰森多边形参数化的瞬变电磁二维贝叶斯反演方法
CN107590347A (zh) 一种基于设计模型匹配孤立点识别与删除方法及系统
CN114167511A (zh) 一种基于连分式展开向下延拓的位场数据快速反演方法
CN113419280A (zh) 基于改进椭圆拟合的叠前裂缝密度估算方法
CN110108271B (zh) 一种气动光学效应引起的星光偏折补偿方法
CN112346139A (zh) 一种重力数据多层等效源延拓与数据转换方法
CN111650570A (zh) 一种地基干涉雷达三维大气校正方法及系统
CN109975869B (zh) 一种沿地层走向光滑约束的反射波波形反演方法
CN114417565B (zh) 一种基于地面气象数据的对流层折射率梯度剖面建模方法
US20110196657A1 (en) Solving a Solute Lubrication Equation for 3D Droplet Evaporation on a Complicated OLED Bank Structure
Zhou et al. Application of Empirical Orthogonal Function Interpolation to Reconstruct Hourly Fine Particulate Matter Concentration Data in Tianjin, China
Han et al. Comparison of methods for curvature estimation from volume fractions
CN116125401A (zh) 基于电磁散射计算的动态海面回波仿真方法
CN109425858B (zh) 基于目标空间分布信息的GB-InSAR系统高程精度分析方法
Franklin et al. A high-order fast marching scheme for the linearized eikonal equation

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
CB03 Change of inventor or designer information
CB03 Change of inventor or designer information

Inventor after: Tao Chunhui

Inventor after: Guo Zhikui

Inventor before: Guo Zhikui

Inventor before: Tao Chunhui

GR01 Patent grant
GR01 Patent grant