内容摘要：用在勘探重力学领域和位场研究中广泛采用的反演和解释方法，即位场数据的三维欧拉反褶积方法，来研究重力场变化数据。该方法适合在先验信息缺乏的情况下开展位场反演和解释。三维欧拉反褶积方法不仅合适于重力异常反演，同样也适合磁异常的反演问题。
1、欧拉反褶积
欧拉反褶积方法适合在已知先验信息较少的情况下，开展自动或半自动的场源位置确定，圈定构造范围和解释异常源的特征参数等。该方法最早由Reid提出（Reid et al.,1990），主要是基于欧拉（Euler）齐次方程（公式1）建立的，通过求解欧拉方程可以确定场源位置的三维空间位置，可以用于对不同类型地质体给出参数估计。
欧拉反褶积方法适合在已知先验信息较少的情况下，开展自动或半自动的场源位置确定，圈定构造范围和解释异常源的特征参数等。该方法最早由Reid提出（Reid et al.,1990），主要是基于欧拉（Euler）齐次方程（公式1）建立的，通过求解欧拉方程可以确定场源位置的三维空间位置，可以用于对不同类型地质体给出参数估计。
公式1 欧拉公式
其中，(x,y,z) 为观测点的位置， (x0,y0,z0))为场源位置，dT/dx 、dT/dy 、dT/dz 为位场异常 在x、y、z三个方向的导数，B 为背景场，N 为构造指数，表示位场异常强度随着深度变化的衰减率，具体与场源的几何结构相关。对于构造指数 的选择问题，前人已经开展了大量研究，一般对于重力异常问题，对于不同的构造特征，N的范围可以参考图1的不同几何形体进行选择。
图1 不同几何体的构造指数N
在应用三维欧拉反褶积方法进行反演时，需要首先对离散数据进行网格化，然后计算位场异常在x、y、z三个方向的梯度，并选择适当的构造指数 。通过对观测数据网格指定一个大小的数据窗口，采用窗口滑动的方法进行反演。
图2 反褶积算法中的滑动窗口示意图
窗口每个观测点数据带入公式1可以组成线性方程组，选择合适的窗口大小，要求解1式所示方程组窗口大小至少应该大于等于3×3。每次移动一次窗口，都可以求解得到一个场源位置参数，最后可以根据这些解的汇聚位置，来确定最佳反演结果，或者调整反演参数。
在选择合适的窗口大小时，应该使窗口大小覆盖异常范围，在一个合适窗口内可以覆盖一个场源的异常范围。在实际的应用中，常常需要通过计算网格尺寸和估计异常规模形态，来进行反复调整，最终得到一个最佳的反演结果。
2、程序实现
Geoist的Euler模块中的EulerDeconvEW和EulerDeconvMW函数可以实现欧拉反演问题。这里面要注意的是欧拉反演不但要输入位场异常还需要输入位场异常在x，y，z三个方向的导数。还记的pftrans模块中的derivx/derivy/derivz函数吧，组合起来用即可。
(1)扩展窗模式
EulerDeconvEW函数采用扩展窗模式反演，这种模式可以指定窗口中心位置，然后以该中心扩展窗口大小，最后选取误差最小的结果为估计值。这个方法的好处是，一次反演只有一个结果输出。
from geoist.pfm import pftrans, euler
sol1 = euler.EulerDeconvEW(x, y, z, data, xderiv, yderiv, zderiv,
structural_index=3, center=(-2000, -2000),
sizes=np.linspace(300, 7000, 20))
sol1.fit()
print("The anomaly location:", sol1.estimate_)
图1 扩展窗Euler反演结果
图1中两个估算场源位置信息如下：
Lower right anomaly location: [-1072.04297339 -830.29615323 1428.85976886]
Upper left anomaly location: [1018.12838822 1576.71821498 1039.07466633]
(2)滑动窗模式
EulerDeconvMW函数采用滑动窗口反演，这种模式的窗口会根据设置进行滑动，windows参数可以设置ny和nx方向的窗口数量，每个方向滑动过程的窗口大小由size参数控制。最后，输出的结果数量不定，一般聚集程度最高的位置可以认为是场源点。
from geoist.pfm import pftrans, euler
sol2 = euler.EulerDeconvMW(x, y, z, data, xderiv, yderiv, zderiv,
structural_index=3, windows=(10, 10),
size=(1000, 1000))
sol2.fit()
print(sol2.estimate_)
图2 滑动窗Euler反演结果
图1和图2是两个球体磁异常的Euler反演结果，理论异常体中心点位置分别为：（x=-1000, y=-1000, z=1500）和（x=1000, y=1500, z=1000）。
相同的异常数据采用两个方法得到的计算结果不完全相同，具体怎么选择需要在实际中具体体会。
一句话总结：由于引起位场的场源本身十分复杂，很难估算场源的合理密度参数范围及深度分布，这就决定采用常规的三维密度结构反演技术，进行密度结构变化的计算存在很大技术难度。而欧拉方法的主要优点就是不需要知道太多先验信息，这就使得该方法在解释位场异常场源特征上具备独到的优点。
参考文献：
[1] Reid A B. 1995. Euler deconvolution: past, present and future, a review[J] . 65th Ann. Internat. Mtng., Soc. Expl. Geophys., Expanded Abstracts, 1995,272~273.
[2] Reid A B, Allsop J M, Granser, H, et al. 1990. Magnetic interpretation in three dimensions using Euler deconvolution[J]. Geophysics,55, 80~91.
网友评论