含透镜体的地质体真三维建模算法

李昌领，张 虹，朱良峰

（1.中国矿业大学 环境与测绘学院，江苏 徐州221008；2.华东师范大学 地理信息科学教育部重点实验室，上海200062）

0 引 言

目前，针对地层、地层中透镜体的三维建模研究不太多见，针对二者整体统一建模的研究更加少见。Lamon［1］等采用B－reps模型对包含透镜体、断层等地质体进行了表面建模；程丹等［2］基于水文剖面图设计了一种 “顶面—尖灭线—透镜体”的三步建模方法构造透镜体模型，对复杂地质体中透镜体的建模进行了非常有益的探索。该算法侧重于表面建模，还没有考虑真三维实体的本质特征；其次把透镜体从复杂地质体中 “分离出来”，单独建模，乎略了周围的地质环境，没有考虑与地层或煤层的相依相伴关系。也有学者同时考虑了实体建模和地质体关系的问题，提出了虚拟钻孔链的方法［3］，但人工刻意干预成份多，自动化程度不高。

本文在以上算法的基础上，以表现地质体的三维实体性为目标，强调透镜体和地层的整体统一建模，提出了“虚拟广义三棱柱”的新模型，设计了一种基于该模型的真三维建模算法。

1 基本概念

定义虚拟广义三棱柱。虚拟广义三棱柱。使用GTP体元［4］对含透镜体等特殊地质体建模时，在不增加虚拟钻孔前提下，一个GTP体元中可能会出现二种或二种以上的不同地质体，若这些不同地质体的几何形态、拓扑关系均能被该GTP体元准确的识别、表达，称具有这些特性的GTP体元为虚拟广义三棱柱体元 （virtual generalized triprism，VGTP），典型的VGTP体元如图1所示。

图1 典型的虚拟广义三棱柱VGTP体元

图1（a）的VGTP中只含单一地层，它和普通的GTP没有差别；图1（b）的VGTP包含两种地层，这是由于地层1的出露线位于VGTP顶面三角面内造成的；图1（c）表示VGTP内部包含透镜体 （或透镜体的一部分），因而由地层和透镜体两种地质体构成；图1（d）则是图1（b）、（c）两种情况的混合，同一个VGTP中包含了3种地质体。

由此可知，每个VGTP中可能包含二种或二种以上的不同地质体，具体到实体表示时，每个VGTP一次只能表达其中一种地质体，如图1（c）中的VGTP，虽然其中包含两种地质体，但只能二选一：用于表示普通地层时，VGTP代表了地层，用于表示透镜体时，VGTP代表了透镜体 （为简单起见，用VGTP中包含透镜体的多面体来表示透镜体），图1（b）中含有两种地层，VGTP要么表示地层1，要么表示地层2。所以生成VGTP时应同时生成该VGTP的地质描述，即表示的是其中的哪一种地质体。

2 含透镜体的地质体三维模型

包含或部分包含透镜体的复杂地层三维模型如图2所示。

图2 含透镜体的地质体三维模型

从图2中看出，含透镜体的地质体作为一个整体，由地层和其中的透镜体构成。普通地层由多个VGTP表达，每个VGTP几何空间可能只表示单一地层，也可能因包含了透镜体而只有部分空间表示地层。无论是VGTP的全部空间 （含单一地层）还是部分空间 （剔除了透镜体后的剩余空间）最终会剖分成多个四面体 （tetrahedron，简记为TEN），用多个TEN的整体来实现一个VGTP中的普通地层。透镜体使用多面体模型表示，通过一定算法，也剖分成四面体TEN。TEN严格地表达了地质实体空间，非常适合于三维可视化、空间分析、空间统计、空间计算等。

VGTP和TEN是并存的，二者缺一不可，有VGTP就必然存在了相应的TEN，离开了TEN，VGTP便只有了概念上的意义，没有实用性可言。透镜体多面体和TEN也有同样的关系。

TEN在模型中不是临时存在的，而是整个模型不可分割的一部分。

3 建模算法

3.1 基本思想

总的来说，根据地层与透镜体相依关系，考虑二者的整体建模。使用VGTP来描述普通地层，使用多面体模型来描述透镜体，然后分别剖分成能够精确表示这二种地质体的四面体。

具体来说，设定建模范围，取得建模区域边界面；生成透镜体的三角网表面模型 （从实体角度来看，是一种多面体模型）；建立地层层序表，地层编号；加密钻孔点；根据钻孔数据生成地表Delaunay三角网 （D－TIN）；将D－TIN沿钻孔迹线向下扩展生成地层的VGTP；若VGTP中包含或部分包含透镜体，把 “含透镜体”的信息写入VGTP地质描述中。

根据VGTP地质描述，若VGTP中包含透镜体，则按一定算法把透镜体空间从VGTP空间中剔除，再把代表地层的剩余不规则空间剖分成四面体；若VGTP只包含单一地层，把整个VGTP空间剖分成四面体。

对于表示透镜体的多面体，按一定算法也剖分成结构一致的四面体。

3.2 建模流程

建模流程如图3所示。

图3 建模流程

3.3 透镜体的多面体建模

根据程丹、杨钦等著文献［2］中提出的方法，使用三角面构造透镜体的表面，得到透镜体的多面体模型。具体实现算法详见相关文献，这里不再赘述。

3.4 地层的逻辑层建模

地层层面编号［5－7］，钻孔加密［8］，地层层面相交处理后，首先在最上层对钻孔点进行Delaunay三角化［9］，形成Delaunay三角网，然后沿钻孔迹线方向向下生成各地层的VGTP［10］。

生成VGTP过程中，若VGTP中包含或部分包含透镜体，要把 “含透镜体”信息、相关的透镜体模型编号信息写入VGTP地质描述中。

接下来还要做更精确的建模，即把VGTP中地层所占空间、透镜体多面体分别剖分成四面体。

3.5 VGTP的四面体剖分算法

VGTP的四面体剖分算法是要剔除其中的透镜体空间，再把剩下空间剖分成四面体。本文算法两步合并，一次性完成四面体的剖分。

3.5.1 基本约定

四面体的每个三角面的法线方向为垂直于三角面指向四面体内部的方向；沿三角面法线方向看过去 （向四面体的内部看），每个三角面的顶点按顺时针排列。

3.5.2 构成四面体的第四点的选点原则

已知一个三角面及其法向量，和位于法线一侧的若干候选空间点，需从中选择一点作为构成四面体的第四点，相对于三角面而言，选择第四点的原则是：计算过空间点和三角面的外接球球心 （四面体外心），与三角面距离最小的外心对应的空间点就是最佳第四点。

3.5.3 剖分算法

对于地层模型中的每个VGTP，在不增加额外点或增加少数中间点 （由于出现了Schnhardt Polyhedron［11］情况，导致该类多面体不能直接被四面体化）的情况下，把其内部表示的地层空间剖分成四面体的算法如下：

若VGTP中只包含单一地层，则把整个VGTP看作一个多面体，使用Shewchuk的gift－wrapping算法［12］或切角算法［13］剖分成四面体。

若VGTP包含透镜体或与透镜体有相交关系，用VGTP的边界面裁减掉位于VGTP外面的透镜体。此时，地层模型抽象成这样的几何模型：VGTP面构成外环面、内部透镜体表面构成若干内环面的含内孔的多面体，。

在文献 Ghadyani的 LAST RESORT算法［14］和Shewchuk的gift－wrapping算法的基础上，提出了一种空间约束域直接Delaunay四面体生长算法，在不增加额外顶点（或增加少量额外顶点）的前提下，以外环面和内环面为约束边界面条件，一次性地把其中的不规则地层空间生长为四面体。剖分算法描述如下：

（1）初始化一个用于存放三角面信息的队列；

（2）从外环面中任选一个三角面；

（3）根据3.5.2节选点原则，从该三角面处直接可见且位于该三角面法线一侧的内、外环面上的所有顶点中，选择一个最优的点，组成一个四面体，加入结果集；新生成的3个三角面的反面加入队列；

（4）若队列为空，则转 （7），否则从队列中移出一个三角面；

（5）对于新移出的三角面，从两方面找出构成四面体的第四点的所有候选邻接点。①若该三角面有一条边位于外环面 （或内环面）上，则该三角面所在的外环面 （或内环面）上的候选点只有一个，即含该边的三角形中与该边相对的顶点 （还有一个顶点，但已构成了四面体）；②若该三角面只有一点位于外环面 （或内环面）上，则找出该点的所有邻接点，再把位于该三角面的法线反向一侧的邻接点去除，剩下的邻接点就是候选空间点；

（6）根据3.5.2节选点原则，从以上候选空间点中找出一个最优的点与该三角面形成一个新的四面体，加入结果集；新生成的三角面的反面进入队列；转第 （4）步；

（7）四面体剖分完毕。

算法执行完毕，若还有少量空间没有剖分成四面体，应加入少量中间点、中间三角面，重新执行 （4）～ （7），直至多面体空间全部剖分成四面体。

3.6 透镜体模型的四面体剖分算法

透镜体模型是多面体模型，可以任选以下两种算法剖分成四面体：Shewchuk的gift－wrapping算法、切角算法。

4 实验结果

参照了开滦精煤股份有限公司某矿业分公司的地质数据，设计了一个地层模型，内部包含了两个不规则的透镜体，顶点数为518，初始边界面数为1024。使用本文算法对其进行了约束Delaunay四面体剖分，结果如图4所示。

图4 含两个透镜体地层约束Delaunay四面体剖分结果

5 结束语

以往算法把透镜体从复杂地质体中 “分离出来”单独建模，忽略了其所在地层。本文充分考虑二者的相依关系，进行整体全面建模。其次，在以往面模型的基础上，考虑地质体真三维实体的本质，以体元模型进行三维实体建模。在使用VGTP建模时，克服了GTP无法同时跨越地层和透镜体的问题。最后，算法不需交互加入位于地层与透镜体交界上的虚拟钻孔点链，从而减少手工干预，增强了算法的自动化程度。本算法初步应用于研究区含透镜体的复杂地质建模，取得了令人满意的结果。

［1］Lamon A M，Jones NL.Building solid models from boreholes and user－defined cross－section ［J］.Computer and Geosciences，2003，29 （3）：547－555.

［2］CHENG Dan， YANG Qin. Three－dimensional geological modeling based on hydro sections ［J］.Journal of Beijing University of Aeronautics and Astronautics，2007，33 （11）：1362－1366 （in Chinese）. ［程丹，杨钦.基于水文剖面的三维地质建模方法 ［J］.北京航空天大学学报，2007，33 （11）：1362－1366.］

［3］LUO Zhiyong.Research and development of 3Dvisualization system for geological surveys ［D］.Chengdu：Chengdu University of Technology，2008 （in Chinese）. ［罗智勇.面向地质勘查的三维可视化系统研制与开发 ［D］.成都：成都理工大学，2008.］

［4］WU Lixin.Topological relations embodied in a generalized triprism （GTP）model for a 3Dgeoscience modeling system ［J］.Computers and Geosciences，2004，30 （4）：405－418.

［5］LUO Zhiyong，YANG Wunian.3Dmodeling and visualization of geological bodies based on borehole dada ［J］.Science of Surveying and Mapping，2008，33 （2）：130－132 （in Chinese）.［罗智勇，杨武年.基于钻孔数据的三维地质建模与可视化研究 ［J］.测绘科学，2008，33 （2）：130－132.］

［6］XIA Yanhua，BAI Shiwei.Study of preprocessing of drilling data for building 3Dstrata model ［J］.Rock and Soil Mechanics，2012，33 （4）：1223－1226 （in Chinese）. ［夏艳华，白世伟.三维地层建模钻孔数据预处理研究 ［J］.岩土力学，2012，33 （4）：1223－1226.］

［7］MING Jing.Quick construction and update of three－dimensional geological models based on boreholes ［J］.Geography and Geo－Information Science，2012，28 （5）：55－60 （in Chinese）. ［明镜.基于钻孔的三维地质模型快速构建及更新 ［J］.地理与地理信息科学，2012，28 （5）：55－60.］

［8］ZHU Liangfeng，PAN Xin，WU Xincai，et al.Construction method and actualizing techniques of 3Dvisual model for geological faults ［J］.Journal of Software，2008，19 （8）：2004－2017（in Chinese）.［朱良峰，潘信，吴信才，等.地质断层三维可视化模型的构建方法与实现技术 ［J］.软件学报，2008，19 （8）：2004－2017.］

［9］SHI Wenzhong，WU Lixin，LI Qingquan，et al.Models and algorithms for three dimensional spatial information system ［M］.Beijing：Publishing House of Electronics Industry，2007：50－58 （in Chinese）.［史文中，吴立新，李清泉，等.三维空间信息系统模型与算法 ［M］.北京：电子工业出版社，2007：50－58.］

［10］HU Gang.Research of 3Dvisualization of geologic body based on borehole data ［J］.Geomatics World，2012 （2）：64－69（in Chinese）.［胡刚.基于钻井数据的地质体3维可视化研究 ［J］.地理信息世界，2012 （2）：64－69.］

［11］Edward Verbree.Piecewise linear complex representation through conforming delaunay tetrahedronization ［C］／／Proceedings of 4th International Conference on Geographic Information Science.New York：Springer，2006：385－387.

［12］Shewchuk J R.Constrained delaunay tetrahedralizations and provably good boundary recovery ［C］／／Proceedings of 11th International Meshing Roundtable.New York：Springer，2002：193－204.

［13］LI Changling，ZHANG Hong，ZHU Liangfeng.Improved algorithm for dividing an arbitrary polyhedron into tetrahedrons ［J］.Computer Engineering and Applications，2012，48（25）：20－23 （in Chinese）.［李昌领，张虹，朱良峰.一种任意多面体剖分成四面体的改进算法 ［J］.计算机工程与应用，2012，48 （25）：20－23.］

［14］Ghadyani H，Sullivan J，WU Ziji.Boundary recovery for delaunay tetrahedral meshes using local topological transformations ［J］.Finite Elements in Analysis and Design，2010，46（1－2）：74－83.