基于GOCAD软件的多源地质勘探数据接口开发
2023-01-03汪文刚李凯
汪文刚，李凯
(中铁第四勘察设计院集团有限公司 地质路基设计研究院，湖北 武汉 430074)
0 引言
自20世纪80年代以来，三维可视化技术被引入地学领域，三维地质建模技术取得了飞速发展，各种建模软件也不断涌现。以GOCAD为代表的先进地质建模软件大大提高了地质建模的效率和精度，许多学者以GOCAD为平台，围绕三维地质建模开展了大量的研究与应用工作[1-4]。最近几年，基于GOCAD平台的矿产资源开发、地质结构研究、水利工程建设、城市地下空间建设等三维地质建模应用越来越多，发展前景十分广阔[5-7]。
三维地质建模是建立在各种地质勘探等数据之上，其中较为基础的有钻探数据，例如井位、井深、岩性分层情况等，地球物理勘探数据(例如：高密度视电阻率数据、地震数据等)可以作为其补充，丰富地质模型的属性信息。GOCAD本身提供了通用的数据导入功能，在导入钻探数据方面，需要将钻孔位置、井深、岩性分层情况等数据按照指定的格式分类整理成多个文本文件，分别依次导入。最基础的为位置文件，包括钻孔编号、井口坐标(北坐标、东坐标、高程)、井深等5列数据；其次是岩性分层文件，包括钻孔编号、地层深度、地层名称3列数据；对于斜井还需要测斜文件，包括钻孔编号、方位角、倾角和深度4列数据[8-11]。在导入地球物理勘探数据时，GOCAD无法直接导入各种原始测量数据文件，只能将其以散点数据的形式整理为文本文件，包括测点坐标(北坐标、东坐标、高程)、测量值4列数据；在此之前，还需要将测点位置从测线坐标转换为大地坐标。此外，对于普遍采用的地球物理勘探成果数据文件—Surfer网格数据文件(*.grd)，GOCAD并未提供直接的导入功能。鉴于以上地质勘探数据导入的繁琐、复杂，无法满足实际生产的需求，有必要探索自动化的导入手段，使其既能提高数据导入的效率，又能避免人为编辑带来的误差。GOCAD提供了非常开放的开发套件，在此基础上，针对不同的数据来源开发接口，可以实现多源数据的快速导入，以解决上述问题。
1 GOCAD二次开发介绍
GOCAD是由法国Nancy大学开发的一款地质对象计算机辅助设计软件，具有强大的三维建模、可视化、地质解释和分析功能，既可以进行表面建模，也可以进行实体建模，既可以设计空间几何对象，也可以表现空间属性分布[12]。GOCAD是目前国际上公认的主流建模软件，在地质工程、地球物理勘探、矿业开发、水利工程中均有广泛的应用。
GOCAD之所以功能强大并不断发展，这与它插件式的开发模式是分不开的，在核心模块的基础上，来自全世界的上百名开发人员开发了500多款功能插件，为不同的行业应用提供了优秀的解决方案。商业化的GOCAD软件提供了非常开放的开发套件，通过开发套件，开发者可以进行以下操作：①以读写方式访问所有GOCAD对象(C++类，图1)；②访问所有 GOCAD算法(C++ 类或 API)；③创建新算法和新 CLI 命令；④在 CLI 命令中使用第三方库；⑤创建、编辑、绘制、操作新的GOCAD对象；⑥引入第三方Qt小部件；⑦创建新工作流。GOCAD的二次开发套件提供插件运行时开发环境，采用基于命令框架 (CLI) 的快速原型设计模式(图2)，能够自动生成 UI 和可编写的脚本命令，大大提高了编程效率。
图1 GOCAD类的继承关系图(部分)
图2 GOCAD开发快速原型设计模式
GOCAD提供的开发套件可以一键配置开发环境，自动生成解决方案工程，开发者只需在.cli文件中编写描述代码。其基本结构如下：
Command
function
do(
};
dialog(
};
}
开发者可以在do代码块中编写逻辑执行代码，其语法规则与C++十分类似；在dialog代码块自定义UI，通过简单的关键字设计对话框控件的排列组合。GOCAD命令框架自定义的数据类型可以自动生成对应的输入控件，省去了控件设计的麻烦。此外，开发者还可以自定义输入参数的外观显示、引入第三方QT小部件等。整体来说，GOCAD提供的二次开发框架既有严格的开发规则限制，又给予了相对充分的扩展空间，在模式化开发的同时，也可以方便开发者按需求定制插件。
2 钻探数据
钻探数据是三维地质建模的基础，GOCAD提供了通用的钻探数据导入功能，主要包括井位数据、地层数据以及测井曲线数据，用户需要按照指定格式将钻孔数据整理为文本文件，不同的信息分别依次导入。当钻孔很多时，这个工作是非常繁琐的。目前大多数生产单位均有自己的钻探数据管理系统，通过搭建平台，将钻探数据管理系统里的数据直接导入到GOCAD中，将会大大提高三维地质建模工作的效率[13]。下面以中铁第四勘察设计院集团有限公司的地质信息系统(以下简称“地质系统”)为例，介绍如何实现这一功能。
该套地质系统是专为钻探数据管理所开发的一套专业软件，其功能强大，已投入使用十多年，外业生产的钻探数据均通过此系统录入。地质系统后台数据存储采用的是微软的ACCESS数据库，该数据库包括了勘察点表、地层表、静力触探表等，囊括了三维地质建模所需要的基本钻探信息。
开发地质系统与GOCAD之间的数据接口可分为3步：①读取ACCESS数据库钻探数据；②根据数据生成GOCAD中的Well(继承自GObj类)对象；③将Well对象批量导入到GOCAD中。C++对ACCESS数据库进行操作，常用的方法有DAO(data access objects)和ADO(activeX data objects)2种方式，本文采用ADO的方式，该种方式在各种科技文献及网络教程中均有详细介绍，在此不再赘述。
为方便读取数据的管理，可以先自定义一个类，将从数据库中读取的井数据存储到类的实例中，包括井名称、坐标、地层情况等。从数据库读取数据之后，利用GOCAD开发套件提供的WellAPI生成Well对象，其接口函数为：
WellAPI::Well* create_well_from_vertical_path(const CString& name, const Point3d& loc, Coord max_depth, const HomogeneousCoordinateSystem& hcs, const Unit& zm_unit, bool specify_well_curve_intervals = false, Coord dzm = 0, EntityId::AccessLevel access_level = EntityId::PUBLIC_ENTITY )
其中，name、loc、max_depth分别代表井的名称、大地坐标以及井底高程。生成Well对象之后，还需将地层信息附加到Well对象上，通过WellAPI接口函数实现：
WellMarker* add_or_replace_marker(WellMarkerSet* marker_set, const CString& marker_name, DepthSystemType dtype, Coord z, Angle::AngleSystem atype, const CString& slope, bool replace, WellMarkerCollection* collection, const UniqueCStringList& forward_layer_name, bool set_gf, const CString& gf_name)
其中，marker_set、marker_name、z分别为Well对象、地层名称及其层底深度。同样的，测井曲线信息也可以通过对应的API接口函数附加到Well对象上。在创建Well对象时，Well对象就自动导入到了GOCAD中，通过GOCAD::store_created(const PtrList
图3 GOCAD中井数据的显示
需要注意的是，GOCAD仅支持UTF-8格式的中文字符，否则会出现乱码。而通过数据库读取出来的字符通常采用的是GB2312的编码格式，因此，在读取数据库的数据之后，需要对其中的中文字符进行处理，将其转为UTF-8编码格式。Qt库函数可以实现这一功能，具体做法是：利用Qt中的fromLocal8Bit()函数将中文字符从本地字符集GB转为Unicode，再利用toStdString()函数将其转为满足GOCAD的API接口函数的字符串类型。此外，GOCAD对于一些特殊字符，例如“&”、“*”等视为非法字符，应在钻孔编录中尽量避免使用这些符号或者在接口开发中对其加以替换。
3 地球物理勘探数据
由于施工环境或生产成本等因素的限制，部分工区钻探资料可能不足，地球物理勘探则成为了其重要的补充。地球物理勘探方法众多，对应的设备仪器多，数据文件类型也是各式各样，但是三维地质建模依赖的并不是这些原始的记录文件，而是经过处理或者反演过后的成果数据，这些数据一般可以以点数据的形式导出为文本文件。此外，地球物理勘探数据在进行成果解释时，通常会使用Golden Software公司的Surfer软件生成网格数据文件(*.grd)，而GOCAD中未提供这一数据格式的导入。
地球物理勘探可分为面积性勘探和剖面勘探，面积性勘探如重力勘探、磁法勘探等，其测点坐标一般为实测的大地坐标；剖面勘探如高密度电法勘探、电磁法勘探、微动勘探等通常沿测线进行，测点坐标为相对坐标。三维地质建模需要所有的数据源统一到相同的三维空间坐标系统下，因此，对于剖面勘探的点数据文件以及网格数据文件，还需进行额外的坐标匹配工作，将各测点的相对坐标转为大地坐标。以下介绍基于GOCAD开发套件开发数据接口，实现点数据文件以及网格数据文件的导入。
3.1 点数据
对于剖面勘探，点数据一般包括测线距离、深度、属性值(例如电阻率、波速等)3个物理量，为了匹配大地坐标，可在测线上实测所有测点或部分控制点的大地坐标，测量数据包括了测线距离、北坐标、东坐标、高程4个物理量。在接口开发中，需要同时读取点数据文件以及测量数据文件，通过简单的线性插值即可将测线上的坐标转换为三维空间上的大地坐标，再利用GOCAD开发套件提供的API函数进行点数据的导入。对于面积性勘探的点数据，其测点一般包含实测的大地坐标，无需进行坐标匹配，可直接利用API函数进行点数据的导入。
在GOCAD开发中，单个数据点对应于一个原子(Atom)对象，原子的空间坐标是其最基础的属性，也可以赋予其他属性值，若干个原子可以组成一个点子集(SubVSet)，若干个点子集又可以组成一个点集(VSet)。GOCAD中点数据集合即是以VSet对象的形式表现出来的[14]，具体实现过程如下：
VSet* vs = (VSet*)VSetAPI::create(name, *hcs, TypeId(“VSet”), false, coordinate_system);
//实例化一个VSet对象vs
PropertyAPI::create_property((AtomicGroup*)vs, property_name);//给vs添加一个属性
SubVSet* subset = vs->create_element();//在vs下创建一个点子集subset
for (int i = 0; i < N; i++) {
Point3d p(X[i], Y[i], Z[i]);//实例化空间坐标
VSetAPI::add_atom(subset, p, 1, V[i]); //添加原子对象
}//循环添加点数据到subset中
store_created(vs);//保存vs对象
首先，实例化一个VSet对象，并给该对象添加一个属性，然后在此对象上，创建一个SubVSet对象，之后，利用VSetAPI中的add_atom函数将原子对象连同其属性依次添加到SubVSet中，这样一个包含点数据的VSet对象就自动导入到了GOCAD中；通过GOCAD::store_created(const PtrList
图4 GOCAD中点数据的显示
3.2 网格数据
Surfer网格文件是地球物理勘探普遍采用的成果形式，GOCAD中尚未提供这一数据格式的直接导入[15-16]。在GOCAD中，面可以通过点集生成，因此，在GOCAD中导入网格文件可分为3步：①读取网格文件中的网格点数据；②导入点数据到GOCAD中；③利用点数据生成面。通过编程可以很容易地从网格文件中读取到点数据，点数据的导入可参照上节介绍进行。同样的，对于剖面勘探数据，还需先进行坐标匹配。这里重点介绍第三步的实现。
在GOCAD开发中，面为TSurf对象，GOCAD开发套件中的TSurfAPI提供了面的多种创建方法，例如，可以直接通过三维空间上的点集(VSet)来创建，也可以通过曲线(PLine)限定面的边界后，利用面内部的控制点来创建。GOCAD从点集生成面，采用的是Delaunay 三角剖分的方式，通过增加数据密度来增加正确生成更复杂表面的几率；为了找到确定点之间的连接性(三角剖分)的最佳投影平面，还可以人为指定投影平面。
在直接从点集生成面的过程中，GOCAD会自动建立点与点之间的拓扑关系，但只适用于简单规则的点集(没有扭曲、翻转、凹凸等)，对于不规则的点集则不太适用，必须利用曲线对点集的周围边界进行拟合，对曲线内部进行三角剖分，结合点集与曲线生成面。从Surfer的网格文件中提取的点数据是规则的网格点数据，但是如若经过坐标匹配，则可能出现扭曲变形，特别是存在测线弯曲以及地形起伏的情况下，因此，需要采用曲线与点集相结合的方式生成面。点集的生成可以参考上一节进行，由于Surfer网格文件中点数据的排列是规则的，可以很容易找到其边界点，通过这些边界点生成闭合的曲线。此外，Surfer成果剖面中通常存在“白化区域”，即空值区域，GOCAD中同样也有空值的概念，只需在读取白化文件时，将指定的白化区域的值设置为空值(一般为-99999)，即可达到“白化”的效果。生成点集和曲线后，利用TSurfAPI的create_from_points_and_curve函数生成TSurf对象。过程如下：
...//设置空值
...//生成VSet对象vs
PLine* curve = PLineAPI::create_from_points(pline_name, &points, true, *hcs, “”);
//生成闭合的边界曲线curve
TSurf* ts = TSurfAPI::create_from_points_and_curve(name, curve, vs, true, normal, true);
//根据点集及其边界曲线生成面
//GObjAPI::destroy(curve);删除curve对象
//GObjAPI::destroy(vs);删除vs对象
store_created(ts);//保存ts对象
由于点集和曲线在生成的过程中会自动导入到GOCAD中，如有必要，可以通过GObjAPI::destroy(GObj* obj)函数将其删除，仅保留面对象，通过GOCAD::store_created(const PtrList
至此，完成了网格数据到GOCAD平台的接口开发。图5为某测线微动谱比成果剖面(Surfer网格数据文件)导入到GOCAD中的实例。
图5 GOCAD中Surfer Grid的显示
4 结语
GOCAD具有强大的三维地质建模功能，数据来源是其建模的基础。本文基于GOCAD开发套件，对钻探数据、地球物理勘探数据接口开发做了简要介绍，实现了这些数据的快速高效导入，大大提高了三维地质建模的效率。以上的数据接口开发仅选取了部分典型的数据类别，实际上三维地质建模的数据来源多种多样，在GOCAD开发套件的基础上，可以根据实际生产需要开发多种数据的导入接口，在此不做过多介绍。
