CN118097055A - 一种三维地质结构模型的网格化表示方法及系统 - Google Patents

Abstract

本发明公开一种三维地质结构模型的网格化表示方法及系统，充分考虑地质结构和属性的地下分布具有极强的非均匀性、不连续性以及空间多尺度性且表现出纵向分层特性和多属性场耦合的相关性，构建统一描述三维地质结构的第一数组、第二数组和第三数组，实现复杂地质结构的矢量栅格一体化表达，能够减少模型占用的内存空间和减少计算量，能够更准确地刻画三维复杂地质，并且易于组织管理大规模数据，整体上能够为大规模三维复杂地质数据的空间分析及可视化提供支撑。

Description

一种三维地质结构模型的网格化表示方法及系统

技术领域

本发明涉及三维地质结构表达技术领域，更具体地，涉及一种三维地质结构模型的网格化表示方法及系统。

背景技术

前人对三维地质模型开展过许多研究，它是地质空间内地质构造以及地质体对象边界的数字化表征模型，具有确定性、可视性和可修改性等特点。目前主要有结构模型和属性模型两种方法实现地质模型的构建和表达，其中结构模型侧重于地质对象几何信息构建，主要内容是利用表达地质体的相关原始数据，在计算机中建立描述地质构造形态、地质体对象之间关系的模型。

三维地质结构模型目前主要有两种方法进行表达，一种是使用三角网表达边界来划分地下空间；一种是网格化方法，强调将地质空间离散为网格化单元，通常为常规六面体元素。其中三角网（TIN）最为流行，主要原因是TIN对复杂表面的刻画能力以及保持对象间边界拓扑一致性方面具有很好的优势，但是三角形数据量规模大影响可视化以及分析的效率。因此，复杂地质模型的网格化表达在科学研究和工程应用中发挥着越来越重要的作用，广泛用于自然资源的评估和预测、岩石物理性质的估计、地质统计学研究。

地质模型的网格化表示方法使用不同的网格填充单元来划分空间。规则网格通过标准长方体或者平行六面体来划分空间，结构简单，在大规模数据需求下内存过大；PillarGrid的结构基于一组从地质模型的顶部到底部延伸的垂直或倾斜的柱子。柱子中的单元是由8个相邻顶点定义的六面体，顶点之间相互独立，单元可以发生形变，能够表达常见的断层，但是需要记录8个顶点，结构复杂；四面体网格自适应表达能力强，但需要记录复杂的网格信息，包括每个网格的顶点信息、连接方式，内存占用量比较大且空间查询困难；Stack-Based Representation of Terrains（SBRT)水平方向上仍是规则网格，在垂直方向上记录岩层边界位置和对应的岩性；Vertically Projected Triangulated Network（VPTN）在水平方向上做三角网剖分，垂直方向上根据地层界面切割出三棱柱，记录切割点高程和地层岩性；SBRT和VPTN的内存更为紧凑，结构简单，但是无法模拟复杂的地形现象，如Y型断层和分支断层。

随着多源异构地质数据矢量栅格一体化集成以及地质数据空间查询分析需求增多，基于不规则三角网的方法不能满足要求，地质模型的网格化表达方法成为当前的研究热点。数据采集技术的进步使地球数字化变得越来越准确。由此产生的大量数据显著影响计算资源和性能：存储和可视化所需的内存大小，以及传输和处理。这极大地影响了整体模型的表达。同时，由于地质信息具有极强的空间非均匀性、不连续性、多尺度和纵向分层的多属性场耦合等特性，目前的网格化方法比如规则网格、PillarGrid和SBRT模型等，仍然无法实现复杂地质结构的精细表达。构建统一表达地质结构的一体化空间数据模型是三维地质信息科学领域需要解决的基础性问题之一。

发明内容

本发明的首要目的是提供一种三维地质结构模型的网格化表示方法，解决现有技术中大量数据影响计算资源和性能，并极大地影响了整体模型的表达的技术问题；本发明的次要目的是提供一种三维地质结构模型的网格化表示系统。

本发明的第一发明提供一种三维地质结构模型的网格化表示方法，包括以下步骤：

将预设区域的三维地质属性场离散化为网格单元；

根据所述预设区域的地质勘察数据，构建每一网格单元的不同地质界面片段的空间关系二叉树、第一数组和第二数组，根据所述预设区域的地质勘察数据和每一网格单元的不同地质界面片段的空间关系二叉树，构建第三数组，其中，所述第一数组存储不同网格单元的数据在所述第二数组和第三数组的存储位置以及网格单元中地质界面片段的数量，所述第二数组存储网格单元中地质界面片段的高程数据，所述第三数组存储地质界面片段的属性信息和地质界面片段在空间关系二叉树的右子树信息；

根据所述预设区域的地质勘察数据，将所述地质界面片段网格化；

根据所述空间关系二叉树描述的不同地质界面片段空间关系，将网格化后的不同地质界面片段投影至所述网格单元中，每一网格单元纵向上记录该网格单元内的不同地质界面片段的空间关系，并根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息。

在上述技术手段中，充分考虑地质结构和属性的地下分布具有极强的非均匀性、不连续性以及空间多尺度性且表现出纵向分层特性和多属性场耦合的相关性。基于InterfaceGrid设计理论模型和数据结构，构建统一描述三维地质结构的InterfaceGrid数据模型，实现复杂地质结构的矢量栅格一体化表达。

进一步的，所述网格单元表达为单元柱。

进一步的，所述第一数组包括：

令预设区域的三维地质属性场离散化为Mx×My个网格单元的规模，则所述第一数组为一个Mx×My×2的三维数组，以行优化顺序进行线性化，在索引（x，y）位置处有：

式中，表示三维的第一数组，/>表示在(x,y)位置处的网格单元在所述第二数组和第三数组的存储位置，/>表示在(x,y)位置处的网格单元中地质界面片段的数量，/>表示一维的第一数组。

进一步的，所述第二数组中同一网格单元中不同地质界面片段的顺序对应于该网格单元的空间关系二叉树的先序序列。

进一步的，所述第二数组存储网格单元中地质界面片段的高程数据时，将浮点型高程值在整数范围内离散化，包括：

式中，表示离散化后的高程值，/>表示浮点型高程值，/>表示网格单元底层的高程值，/>表示网格单元顶层的高程值。

进一步的，根据所述预设区域的地质勘察数据，将所述地质界面片段网格化，包括：

所述预设区域的地质勘察数据包括某个地质界面片段上的采样点集合，根据所述采样点集合采用薄板样条插值算法网格化的拟合地质界面片段。

进一步的，根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息，包括：

确定待查询点所在的网格单元，通过所述第二数组、第三数组获取待查询点所在的网格单元的地质界面片段数据；

令待查询点所在的网格单元待搜索的地质界面片段为N，index表示当前处理的地质界面片段的索引值，初始值为0，最大值为N-1；invalidindex表示无效查询区域，即在该网格单元的空间关系二叉树的先序序列中，invalidIndex以外的地质界面片段数据属于已被排除的范围；top和bottom分别记录待查询点的地质界面片段索引位置；

进行递归处理：对于当前处理的地质界面片段数据，比较待查询点与地质界面片段的空间位置关系，若待查询点位于当前处理的地质界面片段的上方，则更新bottom为当前处理的地质界面片段的索引位置，index的值加1，更新invalidIndex为位于当前处理的地质界面片段数据下方的地质界面片段索引；若待查询点不位于当前处理的地质界面片段的上方，则更新top为当前处理的地质界面片段数据的索引位置，index更新为位于当前处理的地质界面片段的右子树索引；

当满足以下任一条件时结束递归：

1）index与invalidIndex的值相同；

2）index为无效值。

进一步的，在添加和删除地质界面片段时，执行以下步骤：

根据预设的先序序列和第三数组中对应的右子树信息构造链式二叉树；

根据添加和删除地质界面片段的信息，更新所述链式二叉树；

根据添加和删除地质界面片段后的地质勘察数据和更新后的链式二叉树，构建更新后的第一数组、第二数组和第三数组。

本发明的第二方面提供一种三维地质结构模型的网格化表示系统，包括：

离散化模块，所述离散化模块将预设区域的三维地质属性场离散化为网格单元；

数组构建模块，所述数组构建模块根据所述预设区域的地质勘察数据，构建每一网格单元的不同地质界面片段的空间关系二叉树、第一数组和第二数组，根据所述预设区域的地质勘察数据和每一网格单元的不同地质界面片段的空间关系二叉树，构建第三数组，其中，所述第一数组存储不同网格单元的数据在所述第二数组和第三数组的存储位置以及网格单元中地质界面片段的数量，所述第二数组存储网格单元中地质界面片段的高程数据，所述第三数组存储地质界面片段的属性信息和地质界面片段在空间关系二叉树的右子树信息；

网格化模块，所述网格化模块根据所述预设区域的地质勘察数据，将所述地质界面片段网格化；

建模模块，所述建模模块根据所述空间关系二叉树描述的不同地质界面片段空间关系，将网格化后的不同地质界面片段投影至所述网格单元中，每一网格单元纵向上记录该网格单元内的不同地质界面片段的空间关系，并根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息。

进一步的，根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息，包括：

确定待查询点所在的网格单元，通过所述第二数组、第三数组获取待查询点所在的网格单元的地质界面片段数据；

令待查询点所在的网格单元待搜索的地质界面片段为N，index表示当前处理的地质界面片段的索引值，初始值为0，最大值为N-1；invalidindex表示无效查询区域，即在该网格单元的空间关系二叉树的先序序列中，invalidIndex以外的地质界面片段数据属于已被排除的范围；top和bottom分别记录待查询点的地质界面片段索引位置；

进行递归处理：对于当前处理的地质界面片段数据，比较待查询点与地质界面片段的空间位置关系，若待查询点位于当前处理的地质界面片段的上方，则更新bottom为当前处理的地质界面片段的索引位置，index的值加1，更新invalidIndex为位于当前处理的地质界面片段数据下方的地质界面片段索引；若待查询点不位于当前处理的地质界面片段的上方，则更新top为当前处理的地质界面片段数据的索引位置，index更新为位于当前处理的地质界面片段的右子树索引；

当满足以下任一条件时结束递归：

1）index与invalidIndex的值相同；

2）index为无效值。

与现有技术相比，本发明技术方案的有益效果是：

本发明提出的与目前提出的结构模型不仅具有同等的能力，并且从理论上能够对整个复杂地质模型的表达起到有力的补充作用。本发明通过使用第一数组、第二数组和第三数组组成的数据模型能够减少模型占用的内存空间和减少计算量，能够更准确地刻画三维复杂地质，并且易于组织管理大规模数据，整体上能够为大规模三维复杂地质数据的空间分析及可视化提供支撑。

附图说明

图1为本发明实施例提供的一种三维地质结构模型的网格化表示方法的流程示意图；

图2为本发明实施例提供的InterfaceGrid数据模型示意图；

图3为本发明实施例提供的InterfaceGrid数据模型构建流程图；

图4为本发明实施例提供的同一视角下SBRT模型和InterfaceGrid数据模型的对比示意图，其中，图4的（a）和图的4（c）为SBRT模型体绘制效果图，图4的（b）和图4的（d）为InterfaceGrid数据模型体绘制效果图；

图5为本发明实施例提供的单元结构示意图，其中，图5的（a）为PillarGrid模型的单元结构示意图，图5的（b）为InterfaceGrid数据模型的单元结构示意图；

图6为本发明实施例提供的一种三维地质结构模型的网格化表示系统的示意图。

具体实施方式

附图仅用于示例性说明，不能理解为对本专利的限制；

为了更好说明本实施例，附图某些部件会有省略、放大或缩小，并不代表实际产品的尺寸；

对于本领域技术人员来说，附图中某些公知结构及其说明可能省略是可以理解的。

下面结合附图和实施例对本发明的技术方案做进一步的说明。

实施例

如图1所示为一种三维地质结构模型的网格化表示方法的流程示意图，包括以下步骤：

将预设区域的三维地质属性场离散化为网格单元；

根据所述预设区域的地质勘察数据，构建每一网格单元的不同地质界面片段的空间关系二叉树、第一数组和第二数组，根据所述预设区域的地质勘察数据和每一网格单元的不同地质界面片段的空间关系二叉树，构建第三数组，其中，所述第一数组存储不同网格单元的数据在所述第二数组和第三数组的存储位置以及网格单元中地质界面片段的数量，所述第二数组存储网格单元中地质界面片段的高程数据，所述第三数组存储地质界面片段的属性信息和地质界面片段在空间关系二叉树的右子树信息；

根据所述预设区域的地质勘察数据，将所述地质界面片段网格化；

根据所述空间关系二叉树描述的不同地质界面片段空间关系，将网格化后的不同地质界面片段投影至所述网格单元中，每一网格单元纵向上记录该网格单元内的不同地质界面片段的空间关系，并根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息。

本发明实施例提出一种InterfaceGrid数据模型，该InterfaceGrid数据模型结构设计的基本思路是根据地质数据确定定义域，在定义域内构建全局模型，模型单元表达为单元柱（Pillar）。水平方向将网格抽象为XYGrid数组，Pillar则表达为地质界面（Interface）和空间关系（二叉树）的组合，如图2所示。

对于InterfaceGrid数据模型需要表达的两部分信息，第一部分是高程和属性数据，可以由静态数组表达；第二部分地质界面片段的关系信息，直观上表现为链式二叉树。由于链式存储的内存开销较大，本发明实施例使用静态数组表达链式二叉树，同时设计静态数组与链式二叉树之间相互转换的算法以支持数据更新。

进一步的实施例中，所述第一数组（XYGrid数组）包括：

令预设区域的三维地质属性场离散化为Mx×My个网格单元的规模，则所述第一数组为一个Mx×My×2的三维数组，以行优化顺序进行线性化，在索引（x，y）位置处有：

式中，表示三维的第一数组，/>表示在(x,y)位置处的网格单元在所述第二数组和第三数组的存储位置，/>表示在(x,y)位置处的网格单元中地质界面片段（Interface）的数量，/>表示一维的第一数组。

进一步的实施例中，所述第二数组（Interfaces数组）存储网格单元中地质界面片段的高程数据时，将浮点型高程值在整数范围内离散化，使每个高程值节省4字节的内存消耗，包括：

式中，表示离散化后的高程值，/>表示浮点型高程值，/>表示网格单元底层的高程值，/>表示网格单元顶层的高程值。

进一步的实施例中，所述Interfaces数组中同一网格单元中不同地质界面片段的顺序对应于该网格单元的空间关系二叉树的先序序列。

进一步的实施例中，Interfaces数组的内存大小由所有网格单元中Interface的数量决定。

具体的，在Interfaces数组和InterfaceProp数组（第三数组）中同一地质界面片段的数据分组存储；在Interfaces数组中，一个地质界面片段的数据对应一组数据，包括地质界面片段的四个顶点的高程值；在InterfaceProp数组中，一个地质界面片段的数据对应一组数据，包括Interface属性编码及右子树信息；XYGrid数组中指示存储位置的索引为组数序号。

本实施例提供在InterfaceGrid数据模型中，利用基于曲面的边界表示方法，在特定空间关系约束下使用统一网格化流程和组织方式来构造、存储以及管理地下曲面对象的拓扑结构和几何结构，基于InterfaceGrid理论表达将模型的构建划分为以下具体步骤，构建流程图如图3所示：

（1）制定建模参数。计算层面控制点数据的AABB包围框作为定义域，根据区域地质背景以及建模源数据分析结果设置水平方向分辨率Rxy，三维地质属性场离散化为网格单元Grid，称为基网格。

（2）空间关系定义。使用已有的钻孔、地调、物探资料以及经验推断出地层层序关系及断层的接触关系，将层面对象抽象为节点并以二叉树森林定义空间关系。

（3）地质界面网格化。在地质勘查数据中获得地下地质构造的采样点信息，可以利用采样点信息可以了解地质构造以及地质对象的几何形态，但是地质构造中其他点的信息是未知的。因此，需要利用地质勘察获取的采样点信息拟合网格化曲面来预测未知点信息。本发明中采用薄板样条插值算法利用控制点数据拟合地质界面，已知某个地质界面控制点集合为，则薄板样条插值函数为：

其中，为系数，/>，/>为待插值点/>与地质界面控制点中第/>个点的距离,插值点的采样密度与InterfaceGrid数据模型的水平方向分辨率一致。

（4）一体化建模。步骤（3）中生成的网格化地质界面投影到地质网格中形成Interface，空间关系包括地层之间、断层之间、地层与断层之间的关系由（2）定义的空间关系森林逐个网格点嵌入基网格生成Interface，建立整体模型的过程。整体模型在水平方向上为均匀网格，纵向上记录建模工区内的地质界面，地层之间的层序关系、尖灭关系，地层与断层以及断层之间的接触关系由（2）中的定义的二叉树表示。

在进一步的实施例中，在可视化过程中，需要快速有效地获取到模型空间任意位置或者任意区域的属性信息。根据InterfaceGrid数据模型的结构特点，可以设计一种稳健高效的空间数据查询方法- QueryPoint。

QueryPoint总体上是一个递归函数，首先确定待查询点point所在单元柱pillar。根据索引数据获得地质界面数据，即Interfaces数组和InterfaceProp数组中对应的数据。假设待搜索单元柱的地质界面（Interface）数量为N， index表示正在处理的interface索引值，初始值为0，最大值为N-1；invalidindex记录无效查询区域，即在先序序列中，invalidIndex索引位置以外的Interface数据属于已被排除的范围；top和bottom分别记录待查询点的地质界面索引位置。

对于正在处理的地质界面（curInterface）数据，使用双线性插值计算查询点位置在Interface上对应的高度值（z）。比较待查询点与interface的空间位置关系，若待查询点在较上面的位置，更新bottom为curInterface，index向右滑动即数值加1。同时由于当前处理的Interface已经作为底部，该Interface空间位置以下的所有Interface将无需判断，因此更新invalidIndex为位于curInterface下面的地质界面索引（即curInterface的右子树位置）；否则更新top为curInterface，index更新为位于curInterface的右子树索引。

递归结束条件有两个，满足其一则程序结束：1、index与invalidIndex指示位置重合；2、index为无效值（空间位置超出模型数据底部边界）。最终，待查询点位于top和bottom对应的两个地质界面组成的地层中，该地层的属性信息即查询点的属性信息。算法伪代码如下：

function QueryPoint(point, index, invalidindex, top, bottom)

1. Read Interfaces and InterfaceProp from curInterfaceGrid

2. Set pos = GetPosFromGrid(curInterfaces) ,Set index = 0

3. Set invalidIndex = GetNumOflnerfaces(curInterfaces)

4. Set cuInterface = getlnterfaceFromGrid(pos)

5. Set zOnInterface = ComputeZoflnterface(curInterface)

6. If pointz>zOnInterface then

7. If curInterfaceProp[index×2 +1] is not -1

8. Set inalidIndex = curInterfaceProp[index × 2 + 1]

9. Set bottom = index

10. index ++

11. If index = invalidindex

12. return

13. QueryPoint(point, index, invalidindex, top, bottom)

14. else

15. Set top = index

16. Set index = curInterfaceProp[index × 2 + 1]

17. If index = -1

18. return

19. QueryPoint(point,index,invalidindex,top,bottom)

单元柱空间可以看做是二叉树空间分区，二叉树空间分区是一种常见的空间划分方法，单元柱空间可以看做是二叉树空间分区，二叉树空间分区是一种常见的空间划分方法，它递归地将n维空间划分为关于划分超平面的子空间。Interface是单元柱中用于划分空间的平面。通常，在二叉树空间分区中需要建立链式二叉树才能进行查询，从数据构建一棵平衡二叉树的时间复杂度通常是O(nlogn)。本研究提出的查询算法的显著优势在于无需构建二叉树，利用两个序列数组进行递归查询。QueryPoint函数的空间复杂度主要由递归调用所引起的系统栈空间的消耗。每次递归调用都会将当前的函数栈帧压入系统栈中，并在递归结束后将该函数栈帧弹出，因此系统栈的最大深度取决于平衡二叉树的高度h，时间复杂度为O(log(n))，其中n是平衡二叉树中的节点数。表1所示为本发明实施例方法与常规构建二叉树查询方法的复杂度对比分析。表1所示为本发明实施例方法与常规构建二叉树查询方法的复杂度对比分析。

表1 查询算法复杂度分析

在进一步的实施例中，考虑界面之间复杂的交切关系，在添加和删除地质界面数据时，需要由静态数组构建二叉树，执行常规二叉树数据更新算法实现。本研究设计基于InterfaceGrid数据结构的二叉树创建算法（BuildBiTree）。

BuildBiTree函数的功能是使用给定的先序遍历序列和各先序结点对应的右子树序列构造链式二叉树。其中，先序遍历序列按照根节点、左子树、右子树的顺序排列，右子树序列存储节点右子树在先序遍历序列中的位置。函数伪代码如下：

function BuildBiTree(index, node, treeOrder, rightTreeOrder, flag)

1. If index>treeOrder.size()-1 or flag[index]

2.return

3. Create a tree-node interfaceNode of binary tree

4. Set interfaceNode ->prop = treeOrder[index]

5. If rightTreeOrder[index] is not -1

6. BuildBiTree(rightTreeOrder[index],node->downInterfaceNode,treeOrder, rightTreeOrder, flag)

7. If flag[index + 1] is false

8. BuildBiTree(index+1,node->upInterfaceNode,TreeOrder,rightTreeOrder, flag)

BuildBiTree函数的复杂度可以分析为O(n)，其中n为树中结点的数量。该函数通过递归的方式实现了对树的构造。在每一次递归中，函数都会创建一个新的结点，并记录节点属性。接着，函数会根据rightTreeOrder数组的信息递归地构造结点的右子树，再递归地构造结点的左子树。由于在每一次递归中都会对Flag数组进行标记，保证了每个结点只会被构造一次。因此，对于n个结点的树，该函数的时间复杂度为O(n)。因此，该函数的时间复杂度与二叉树的节点个数成线性关系。表2所示为本研究二叉树创建方法与先序和中序序列构建二叉树方法的复杂度对比分析。

表2构建树算法复杂度分析

根据地质空间数据矢量栅格一体化集成以及对数据查询、分析和计算的需求，采用网格化的空间数据表达方法已成为地学数据集成的主要研究趋势。针对规则网格、角点网格、广义三棱柱等无法实现复杂地质结构的表达问题，本发明实施例考虑地质结构和属性的地下分布具有极强的非均匀性、不连续性以及空间多尺度性且表现出纵向分层特性和多属性场耦合的相关性，提出了界面网格数据模型，应用它统一描述地质结构和多属性场的地下分布实现地质模型矢量栅格的一体化表达，从而实现地质空间数据查询、分析、可视化及应用。

作为矢量栅格一体化模型，界面网格具备栅格模型和矢量模型的特点和优势。一方面，界面网格结构简单，空间坐标隐含在网格索引中，利用本发明实施例设计的查询算法纵向上实现了二分查找，查询效率高，利于进行空间分析计算和在线可视化；另一方面，很大程度上减少了由于网格尺寸导致的空间误差，建模精度高，有利于大尺度精细建模。本研究测试验证表明，与SBRT模型相比，界面网格模型的可视化结果边界无锯齿问题，能够更加真实地描述地质边界情况，提高了三维地质结构模型的准确性和可靠性。

本发明实施例使用Citcoms全球球面网格剖分成功构建地球表面到莫霍面9个岩石圈分层的界面网格模型数据并实现可视化。实验显示，界面网格模型可以实现全球化圈层结构的网格化描述，同时界面网格模型相比较PillarGrid数据量更小，更加适合基于网络的大数据地学应用服务。

在具体的模型对比测试中，三维地质构造的表达能力是数据模型可视化效果的重要评价标准，SBRT模型由于其数据量和空间查询的优势被广泛应用于在线可视化领域。本研究应用测试数据三维网格模型对比验证SBRT模型和界面网格模型的体绘制效果。它共有6个岩层，岩性属性在三维网格内呈现层状分布，存在岩层尖灭现象。水平方向上划分网格规模为，网格间距为/>，测试环境为CPU：Intel Core i7，显卡：NVidiaGeofore MX150，采用谷歌浏览器，体绘制效果如图4所示。

图4可以看出，SBRT模型作为规则网格的纵向压缩，继承了规则网格的特点，体绘制时在岩层边界处会产生明显的锯齿或台阶问题，影响可视化效果，不能准确描述地质情况；界面网格模型在SBRT模型的基础上增加了岩层连接约束，岩性边界无锯齿现象，更符合真实数据效果，提高了地质模型的准确性。

在绘制效率方面，SBRT模型的绘制效率是40帧率，界面网格数据模型的绘制效率是25帧率。由于界面网格模型的数据量比SBRT模型大，界面网格模型体绘制的显示帧率比SBRT模型低。随着GPU硬件的升级，计算机显存提升明显，因此内存导致的一定的帧率降低是可接受的。

全球球面网格剖分测试：

测试数据LITHO1.0模型的球面是一个正二十面体球面三角网，在每个三角形单元的径向上描述了岩石圈的9个亚层(水、冰、上层沉积物、中层沉积物、下层沉积物、上地壳、中部地壳、下地壳和岩石圈地幔)。针对全球岩石圈模型多源异构地质空间信息整合以及地质大数据一体化组织需求，本发明实施例来自LITHO1.0模型的9个亚层数据，利用相同分辨率的Citcoms网格建立了PillarGrid模型和界面网格模型，验证界面网格模型。PillarGrid和界面网格均为柱状网格，基于一组从地质模型的顶部到底部延伸的垂直或倾斜的柱子。PillarGrid的一个单元由其 8 个相邻顶点定义，与界面网格相比，它在纵向上重复记录了相邻单元的共用顶点，两种不同结构的单元对比如图5所示。两种模型的单元数量和内存空间占用量如表3所示。

界面网格模型中的单元数量比PillarGrid模型的单元数量增加了约8%，这是由于界面网格的单元为四边形面，PillarGrid模型的单元为六面体，同一位置纵向上面单元与体单元表达范围的差异决定。界面网格模型中的内存占用量比PillarGrid模型减少了约36.8%，对于不同区域的模型，这一比例可能会有所不同，与建模区域的地层条件有关。一般来说，在相同的运行环境中，平均层数越大，界面网格模型与PillarGrid模型的比值越小，总体上界面网格模型中的内存占用量约是PillarGrid模型的2/3。

表3 模型数据量对比

本发明实施例还提供一种三维地质结构模型的网格化表示系统，如图6所示，包括：

离散化模块，所述离散化模块将预设区域的三维地质属性场离散化为网格单元；

数组构建模块，所述数组构建模块根据所述预设区域的地质勘察数据，构建每一网格单元的不同地质界面片段的空间关系二叉树、第一数组和第二数组，根据所述预设区域的地质勘察数据和每一网格单元的不同地质界面片段的空间关系二叉树，构建第三数组，其中，所述第一数组存储不同网格单元的数据在所述第二数组和第三数组的存储位置以及网格单元中地质界面片段的数量，所述第二数组存储网格单元中地质界面片段的高程数据，所述第三数组存储地质界面片段的属性信息和地质界面片段在空间关系二叉树的右子树信息；

网格化模块，所述网格化模块根据所述预设区域的地质勘察数据，将所述地质界面片段网格化；

建模模块，所述建模模块根据所述空间关系二叉树描述的不同地质界面片段空间关系，将网格化后的不同地质界面片段投影至所述网格单元中，每一网格单元纵向上记录该网格单元内的不同地质界面片段的空间关系，并根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息。

进一步的实施例中，根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息，包括：

确定待查询点所在的网格单元，通过所述第二数组、第三数组获取待查询点所在的网格单元的地质界面片段数据；

令待查询点所在的网格单元待搜索的地质界面片段为N，index表示当前处理的地质界面片段的索引值，初始值为0，最大值为N-1；invalidindex表示无效查询区域，即在该网格单元的空间关系二叉树的先序序列中，invalidIndex以外的地质界面片段数据属于已被排除的范围；top和bottom分别记录待查询点的地质界面片段索引位置；

进行递归处理：对于当前处理的地质界面片段数据，比较待查询点与地质界面片段的空间位置关系，若待查询点位于当前处理的地质界面片段的上方，则更新bottom为当前处理的地质界面片段的索引位置，index的值加1，更新invalidIndex为位于当前处理的地质界面片段数据下方的地质界面片段索引；若待查询点不位于当前处理的地质界面片段的上方，则更新top为当前处理的地质界面片段数据的索引位置，index更新为位于当前处理的地质界面片段的右子树索引；

当满足以下任一条件时结束递归：

1）index与invalidIndex的值相同；

2）index为无效值。

相同或相似的标号对应相同或相似的部件；

附图中描述位置关系的用语仅用于示例性说明，不能理解为对本专利的限制；

显然，本发明的上述实施例仅仅是为清楚地说明本发明所作的举例，而并非是对本发明的实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或变动。这里无需也无法对所有的实施方式予以穷举。凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明权利要求的保护范围之内。

Claims (10)

1.一种三维地质结构模型的网格化表示方法，其特征在于，包括以下步骤：

将预设区域的三维地质属性场离散化为网格单元；

根据所述预设区域的地质勘察数据，构建每一网格单元的不同地质界面片段的空间关系二叉树、第一数组和第二数组，根据所述预设区域的地质勘察数据和每一网格单元的不同地质界面片段的空间关系二叉树，构建第三数组，其中，所述第一数组存储不同网格单元的数据在所述第二数组和第三数组的存储位置以及网格单元中地质界面片段的数量，所述第二数组存储网格单元中地质界面片段的高程数据，所述第三数组存储地质界面片段的属性信息和地质界面片段在空间关系二叉树的右子树信息；

根据所述预设区域的地质勘察数据，将所述地质界面片段网格化；

根据所述空间关系二叉树描述的不同地质界面片段空间关系，将网格化后的不同地质界面片段投影至所述网格单元中，每一网格单元纵向上记录该网格单元内的不同地质界面片段的空间关系，并根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息。

2.根据权利要求1所述的三维地质结构模型的网格化表示方法，其特征在于，所述网格单元表达为单元柱。

3.根据权利要求1所述的三维地质结构模型的网格化表示方法，其特征在于，所述第一数组包括：

令预设区域的三维地质属性场离散化为Mx×My个网格单元的规模，则所述第一数组为一个Mx×My×2的三维数组，以行优化顺序进行线性化，在索引（x，y）位置处有：

式中，表示三维的第一数组，/>表示在(x,y)位置处的网格单元在所述第二数组和第三数组的存储位置，/>表示在(x,y)位置处的网格单元中地质界面片段的数量，/>表示一维的第一数组。

4.根据权利要求3所述的三维地质结构模型的网格化表示方法，其特征在于，所述第二数组中同一网格单元中不同地质界面片段的顺序对应于该网格单元的空间关系二叉树的先序序列。

5.根据权利要求4所述的三维地质结构模型的网格化表示方法，其特征在于，所述第二数组存储网格单元中地质界面片段的高程数据时，将浮点型高程值在整数范围内离散化，包括：

式中，表示离散化后的高程值，/>表示浮点型高程值，/>表示网格单元底层的高程值，/>表示网格单元顶层的高程值。

6.根据权利要求5所述的三维地质结构模型的网格化表示方法，其特征在于，根据所述预设区域的地质勘察数据，将所述地质界面片段网格化，包括：

所述预设区域的地质勘察数据包括某个地质界面片段上的采样点集合，根据所述采样点集合采用薄板样条插值算法网格化的拟合地质界面片段。

7.根据权利要求6所述的三维地质结构模型的网格化表示方法，其特征在于，根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息，包括：

确定待查询点所在的网格单元，通过所述第二数组、第三数组获取待查询点所在的网格单元的地质界面片段数据；

令待查询点所在的网格单元待搜索的地质界面片段为N，index表示当前处理的地质界面片段的索引值，初始值为0，最大值为N-1；invalidindex表示无效查询区域，即在待查询点所在的网格单元的空间关系二叉树的先序序列中，invalidIndex以外的地质界面片段数据属于已被排除的范围；top和bottom分别记录待查询点的地质界面片段索引位置；

进行递归处理：对于当前处理的地质界面片段数据，比较待查询点与地质界面片段的空间位置关系，若待查询点位于当前处理的地质界面片段的上方，则更新bottom为当前处理的地质界面片段的索引位置，index的值加1，更新invalidIndex为位于当前处理的地质界面片段数据下方的地质界面片段索引；若待查询点不位于当前处理的地质界面片段的上方，则更新top为当前处理的地质界面片段数据的索引位置，index更新为位于当前处理的地质界面片段的右子树索引；

当满足以下任一条件时结束递归：

1）index与invalidIndex的值相同；

2）index为无效值。

8.根据权利要求7所述的三维地质结构模型的网格化表示方法，其特征在于，在添加和删除地质界面片段时，执行以下步骤：

根据预设的先序序列和第三数组中对应的右子树信息构造链式二叉树；

根据添加和删除地质界面片段的信息，更新所述链式二叉树；

根据添加和删除地质界面片段后的地质勘察数据和更新后的链式二叉树，构建更新后的第一数组、第二数组和第三数组。

9.一种三维地质结构模型的网格化表示系统，其特征在于，包括：

离散化模块，所述离散化模块将预设区域的三维地质属性场离散化为网格单元；

数组构建模块，所述数组构建模块根据所述预设区域的地质勘察数据，构建每一网格单元的不同地质界面片段的空间关系二叉树、第一数组和第二数组，根据所述预设区域的地质勘察数据和每一网格单元的不同地质界面片段的空间关系二叉树，构建第三数组，其中，所述第一数组存储不同网格单元的数据在所述第二数组和第三数组的存储位置以及网格单元中地质界面片段的数量，所述第二数组存储网格单元中地质界面片段的高程数据，所述第三数组存储地质界面片段的属性信息和地质界面片段在空间关系二叉树的右子树信息；

网格化模块，所述网格化模块根据所述预设区域的地质勘察数据，将所述地质界面片段网格化；

建模模块，所述建模模块根据所述空间关系二叉树描述的不同地质界面片段空间关系，将网格化后的不同地质界面片段投影至所述网格单元中，每一网格单元纵向上记录该网格单元内的不同地质界面片段的空间关系，并根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息。

10.根据权利要求9所述的三维地质结构模型的网格化表示系统，其特征在于，根据所述第一数组、第二数组和第三数组查询网格单元中任意点的属性信息，包括：

确定待查询点所在的网格单元，通过所述第二数组、第三数组获取待查询点所在的网格单元的地质界面片段数据；

令待查询点所在的网格单元待搜索的地质界面片段为N，index表示当前处理的地质界面片段的索引值，初始值为0，最大值为N-1；invalidindex表示无效查询区域，即在待查询点所在的网格单元的空间关系二叉树的先序序列中，invalidIndex以外的地质界面片段数据属于已被排除的范围；top和bottom分别记录待查询点的地质界面片段索引位置；

进行递归处理：对于当前处理的地质界面片段数据，比较待查询点与地质界面片段的空间位置关系，若待查询点位于当前处理的地质界面片段的上方，则更新bottom为当前处理的地质界面片段的索引位置，index的值加1，更新invalidIndex为位于当前处理的地质界面片段数据下方的地质界面片段索引；若待查询点不位于当前处理的地质界面片段的上方，则更新top为当前处理的地质界面片段数据的索引位置，index更新为位于当前处理的地质界面片段的右子树索引；

当满足以下任一条件时结束递归：

1）index与invalidIndex的值相同；

2）index为无效值。