内容摘要： 地质建模是在已知的地质和地球物理解释基础上，合理运用地质学方法，综合分析通过定量化的方法实现数字化的地质模型构建。通常可分为缺定性建模和随机性建模两大类。地球物理反演虽然具有多解性，这无法避免，但通过与地质建模工具相结合，在已有认识基础上开展三维反演，来提高反演结果的可解释性，这条路是值得深入研究的。
1、认识GemPy
GemPy是一套开源的三维地质建模工具包，使用python语言编写。Github的项目地址是https://github.com/cgre-aachen/gempy。
最新版本2.0.0，官方的网站地址是https://www.gempy.org/，该项目迭代很快，目前还在快速发展之中，有时候网上的一些例子可能接口变化，会出现运行错误。
该工具包具备创建复杂地质模型的能力，可以实现层序特征和结构特征的地质模型构建任务，例如：背斜、向斜类的褶皱构造、考虑相互接触关系的复杂断层和各种不整和的复杂模型。下图是该软件的典型建模结果。
2、用法简介
首先，肯定就是安装，在python虚拟环境下运行
$ pip install gempy
在使用Gempy过程中，大家只要掌握：如何计算模型、修改模型和在二维视图及三维视图中查看模型，这些基本过程就可以。对于主要的函数如何使用，建议通过笔记来帮助记忆。
地质建模核心技术还是三维空间插值算法，Gempy的核心支撑就是kriging和co-kriging等算法，可以在地表、剖面和定点的一些有限地质资料情况下，通过三维空间差值算法，实现最优地质概念三维模型的建立，模型可以你和剖面解释结果，层序地层学结果等。
地表的地质图/露头信息，概要地层划分，地质剖面解释，钻孔资料，断层系统的划分是实现这些建模的基础。要先整理好这些数据才能顺利完成建模过程。
第一个程序
步骤1：读入数据
import gempy as gp
#读入数据
geo_data = gp.create_data(extent = [0,1000,0,1000,0,1000],
resolution = [50, 50, 50],
path_o = 'model1_orientations.csv',
path_i = 'model1_surface_points.csv')
gp.map_series_to_surfaces(geo_data, {"Strat_Series": ('rock2','rock1'),"Basement_Series":('basement')})
读入数据后可以检查一下正确性，函数为：plot.plot_data
gp.plot.plot_data(geo_data, direction='y')
步骤2：计算模型
interp_data = gp.set_interpolation_data(geo_data, compile_theano=True,
theano_optimizer='fast_compile')
#计算模型，这块有点慢
sol = gp.compute_model(geo_data)
步骤3：查看结果
gp.plot.plot_section(geo_data, cell_number=25,
direction='x', show_data=True)
查看剖面模型结果的函数为：plot.plot_section，效果如下：
还可以修改参数后，再重新画图
gp.plot.plot_section(geo_data, cell_number=25,
direction='y', show_data=True)
大家如果手头没有数据，可以在我们团队的数据仓库下载
https://gitee.com/cea2020/geodataset
在geodataset的gempy目录下，可以找到更多测试用的CSV格式数据
下面对地球物理专业的小盆友补以下几个地质概念：
不整合与侵入岩。不整合是一种受到侵蚀的地质事件，由于没有受到沉积作用，但是还需要正确地知道其地质形成年代。上面的Series就需要被不整合分开。侵入岩是一种典型的跨地层切割方式。他们具有侵入和侵蚀加之后沉积之间的拓扑相似性。
理解各种地层之间的关系是使用该软件进行地质建模和解释的基础。这个综合例子就展示了三套地层体系。
有时候，在先后沉积的地层之间缺失了某一时期的地层，造成上、下地层时代的不连续。上、下地层之间的这种接触关系称为不整合接触（nonconformable contact）。不整合接触的上、下地层之间隔着一个大陆剥蚀面，这个面就叫不整合面。
3、注意事项
（1）第三方依赖
PyMC3：这个包提供了概率编程框架,使用Theano通过变分推理进行梯度计算,并使用了C实现加速运算，Gempy赖它实现高级随机性建模。PyMC3包很复杂，安装过程中，可能会出现本地C语言环境方面的错误提示。
Theano：这个包可以让用户去定义,优化,计算数学表达式,特别是numpy等多维数组。依赖它可以实现模型计算加速。也是与Caffe，Tensorflow齐名的包，但是现在遗憾的是该包已经不在提供升级，作者放弃了维护，如果是新项目建议选择其它替代方案吧
（2）一些概念
在使用软件过程中，注意Formations、Pile这两个最基本概念的定义。其中Formation可以理解为地层的概念，一共模型中反映出那些地层都需要事先在这里定义出来，并给出合适的颜色区分；
Pile的概念可以理解为这些地层的组合关系，Pile把这些地层的关系定义出来，如那些地层是什么年代的，并且给出地层之间的关系，共分两种Onlap，Erode；然后还要给出建模的参考面：Top或Bottom。
这两步做完后，就是建立Section的过程，因为该软件的建模过程主要参考剖面地质解释和钻孔资料等数据，最终目的就是能让三维地质模型结果与各种观测资料相拟合。
一句话结语：GemPy是基于Python语言的开源软件，与商业软件3D GeoModeller的功能相似，融入了基本的构造地质学概念，可以实现多个地质体之间的拓扑关系建模，适合与反演方法研究相互配合，提高反演效果。
参考文献
[1] Varga M., Schaaf A., Wellmann F., 2019. GemPy 1.0: open-source stochastic geological modeling and inversion, Geoscientific Model Development, 12, 1-32
[2] Calcagno, P., Chiles J.P., Courrioux, G., Guillen, A., 2008, Geological modelling from field data and geological knowledge Part I. Modelling method coupling 3D potential-field interpolation and geological rules, Physics of the Earth and Planetary Interiors, 171, 147-157.
