CN114359500B - 一种面向滑坡危险范围预测的三维建模与可视化方法 - Google Patents
一种面向滑坡危险范围预测的三维建模与可视化方法 Download PDFInfo
- Publication number
- CN114359500B CN114359500B CN202210229728.6A CN202210229728A CN114359500B CN 114359500 B CN114359500 B CN 114359500B CN 202210229728 A CN202210229728 A CN 202210229728A CN 114359500 B CN114359500 B CN 114359500B
- Authority
- CN
- China
- Prior art keywords
- landslide
- terrain
- node
- sph
- particle
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000007794 visualization technique Methods 0.000 title claims abstract description 10
- 238000004088 simulation Methods 0.000 claims abstract description 84
- 238000000034 method Methods 0.000 claims abstract description 76
- 238000009877 rendering Methods 0.000 claims abstract description 52
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 33
- 230000000007 visual effect Effects 0.000 claims abstract description 33
- 238000009825 accumulation Methods 0.000 claims abstract description 30
- 230000000694 effects Effects 0.000 claims abstract description 22
- 230000004927 fusion Effects 0.000 claims abstract description 22
- 238000012800 visualization Methods 0.000 claims abstract description 21
- 238000005457 optimization Methods 0.000 claims abstract description 10
- 239000002245 particle Substances 0.000 claims description 77
- 238000004364 calculation method Methods 0.000 claims description 17
- 238000011156 evaluation Methods 0.000 claims description 9
- 238000012545 processing Methods 0.000 claims description 8
- 239000012530 fluid Substances 0.000 claims description 7
- 230000001133 acceleration Effects 0.000 claims description 6
- 238000013507 mapping Methods 0.000 claims description 5
- 238000013459 approach Methods 0.000 claims description 4
- 230000009471 action Effects 0.000 claims description 3
- 230000008859 change Effects 0.000 claims description 3
- 230000010354 integration Effects 0.000 claims description 3
- 238000010276 construction Methods 0.000 claims description 2
- 238000005243 fluidization Methods 0.000 abstract description 2
- 230000008569 process Effects 0.000 description 12
- 238000005516 engineering process Methods 0.000 description 9
- 238000010586 diagram Methods 0.000 description 8
- 230000002829 reductive effect Effects 0.000 description 7
- 230000006870 function Effects 0.000 description 5
- 238000004458 analytical method Methods 0.000 description 3
- 230000005484 gravity Effects 0.000 description 3
- 239000000463 material Substances 0.000 description 3
- 241000581364 Clinitrachus argentatus Species 0.000 description 2
- 230000003044 adaptive effect Effects 0.000 description 2
- 230000003628 erosive effect Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 150000003254 radicals Chemical class 0.000 description 2
- 238000010187 selection method Methods 0.000 description 2
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 2
- 241001441571 Hiodontidae Species 0.000 description 1
- 239000011324 bead Substances 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000000903 blocking effect Effects 0.000 description 1
- 239000003086 colorant Substances 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 230000008021 deposition Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000004880 explosion Methods 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 230000001788 irregular Effects 0.000 description 1
- 230000009191 jumping Effects 0.000 description 1
- 238000002372 labelling Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 230000008520 organization Effects 0.000 description 1
- 230000036961 partial effect Effects 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 230000000750 progressive effect Effects 0.000 description 1
- 238000007670 refining Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 239000011435 rock Substances 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
- 239000002689 soil Substances 0.000 description 1
- 239000000126 substance Substances 0.000 description 1
- 239000013589 supplement Substances 0.000 description 1
- 238000010200 validation analysis Methods 0.000 description 1
- 230000003313 weakening effect Effects 0.000 description 1
Images
Classifications
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Processing Or Creating Images (AREA)
Abstract
本发明公开了一种面向滑坡危险范围预测的三维建模与可视化方法,将滑坡概化为流化模型,采用SPH法构建滑坡危险范围数值模型;提出顾及滑坡运动堆积特征的滑坡可视化优化算法,该算法采用Metaball方法渲染SPH数值模型,LOD方法优化滑坡地形建模,以解决大型滑坡计算资源不足的问题。本发明方法提出的三维SPH与Metaball融合建模方法在取得较好的数值模拟结果的同时,可有效保证滑坡可视化模拟的效果,能够提供具有较高准确性的滑坡危险范围预测结果,提出的四叉树多尺度地形渲染算法能有效简化地形渲染,提升程序运行速度。
Description
技术领域
本发明涉及滑坡预测模拟技术领域,具体涉及一种面向滑坡危险范围预测的三维建模与可视化方法。
背景技术
高速远程滑坡冲淤距离及危险范围预测对于滑坡灾害防治至关重要,具有很高的研究价值。高速远程滑坡具有滑动时间短、速度快、滑动距离长等特点,对下游社区和财产造成重大损害。SPH模型可应用于高速远程滑坡冲淤距离及危险范围的预测。然而,现有的SPH模型在实际仿真应用中存在诸多不足。
从高速远程滑坡的运动堆积特征来看,滑体从启动到静止既有崩滑后破碎的过程,也有受到雨水冲刷后崩滑体粘合的作用。尽管SPH在数值模拟方面具有一定的优势,但SPH滑坡数值模型是由数万个SPH离散粒子组成的滑体,因此完全基于SPH方法的高速远程滑坡堆积过程三维可视化,无法真实地模拟滑坡运动堆积特征。
利用滑坡模拟预测滑坡危险范围是滑坡灾害防治、减少人员和财产损失的重要手段。尽管滑坡危险范围数值模拟预测的准确性需要重点考虑,但能够反映滑坡危险性关键参数及其运动堆积过程的可视化模拟也不容忽视。此外,由于当前的滑坡模拟通常将大规模地形数据和对应的影像数据一次性调入内存并全部渲染,这是导致滑坡模拟程序运行缓慢且仅能模拟小尺度滑坡区域的重要原因。
通过以上分析不难看出,当前缺乏一种既能保证滑坡危险范围预测准确性及可靠的运行速度,又能通过三维可视化的方法真实再现滑坡堆积过程。
发明内容
为解决现有技术中存在的问题,本发明提供了一种面向滑坡危险范围预测的三维建模与可视化方法,将滑坡概化为流化模型,采用SPH法构建滑坡危险范围数值模型;提出顾及滑坡运动堆积特征的滑坡可视化优化算法,该算法采用Metaball方法渲染SPH数值模型,LOD方法优化滑坡地形建模,以解决大型滑坡计算资源不足的问题,解决了上述背景技术中提到的问题。
为实现上述目的,本发明提供如下技术方案:一种面向滑坡危险范围预测的三维建模与可视化方法,包括如下步骤:
S1、采用SPH法构建滑坡危险范围数值模型模拟滑坡运动;
S2、采用Metaball与SPH融合可视化建模,为SPH粒子生成多边形包围体,再现滑坡运动堆积特征;
S3、使用四叉树多尺度地形渲染算法进行地形渲染,采用不同尺度的地形模型表示地形,利用塔形四叉树结构管理地形数据,通过动态调度地形数据块以及不同尺度地形数据块的组合渲染。
优选的,所述步骤S1中的SPH法具体是指:通过连续性密度法对粒子密度进行更新,公式如下:
粒子的运动方程如下式所示:
由此,得到粒子的加速度,加速度改变粒子速度,速度改变粒子位置,通过使用蛙跳法进行时间积分,即可完成对滑坡运动模拟。
优选的,所述步骤S2中Metaball与SPH融合可视化建模具体包括:设整个三维空间中存在一个能量场,能量场中任意一个空间点具有的能量值等于所有SPH粒子对该点造成的影响的总和,通过使用与SPH法中的密度求和法相同的式子进行计算,能量值计算公式如下:
优选的,所述的等值面的构建采用MC算法,所述MC算法包括如下步骤:
S22、遍历所有正方体单元,对于每一个正方体单元,计算其八个顶点各自的能量
值,并将这些能量值分别与预先设定的能量值相比较,若大于或等于,则说明该顶点位
于等值面的内部,标记为0;若小于,则说明该顶点位于等值面的外部,标记为1;
将八个顶点的标记情况按顶点编号顺序存入一个八位二进制数中,该二进制数的每一位分别对应一个顶点的标记情况,
建立边搜索表EdgeTable和三角形搜索表TriTable,EdgeTable保存正方体单元的边与等值面相交的情况,TriTable保存在正方体单元内部生成三角形面片的情况的三角形顶点连接顺序;
S23、将步骤S22中得到的八位二进制数转为十进制数并代入EdgeTable中得到一个十二位二进制数,新的十二位二进制数的每一位分别对应正方体单元的一条边,表示该正方体单元的哪些边与等值面相交,然后进行线性插值,计算出这些边各自与等值面的交点坐标;
S24、查找TriTable,依据TriTable提供的连接顺序,将步骤S23中所得的交点连接成一个个三角形面片,用这些三角形面片组成的多边形面来逼近等值面。
优选的,所述线性插值的公式如下:
优选的,所述步骤S3中的使用四叉树多尺度地形渲染算法进行地形渲染具体包括如下步骤:
S31、创建队列1和队列2, 队列1用于存放正在处理的当前层次的地形节点(每一个地形数据块可看作一个四叉树节点,即地形节点),队列2用于存放待处理的下一个层次的地形节点;
S32、采用广度优先遍历四叉树的原则从根节点开始对四叉树进行遍历;在遍历四叉树时,一个四叉树节点即为一个地形节点,该节点所对应的高程点集直接从滑坡区域DEM数据中选取;
S33、将遍历到的节点送入队列1,判断该节点在视景体内是否可见,若是不可见的,则直接出队列,不予处理;若是可见的,则使用节点细分评价体系判断该节点是否需要细分,若不需要细分,则直接渲染该节点;若需要细分,则将该节点出列,进行细分,并将该节点细分产生的下一个层次的四个子节点送入队列2,等待下一轮的处理;
S34、当程序运行到队列1为空时,则交换队列1和队列2中的节点,并重复步骤S31至步骤S33,直到所有需要处理的节点均已处理完毕为止;对于那些在视景体内可见的叶子节点,则直接进行渲染。
优选的,所述步骤S32中选取的具体步骤如下:
将DEM数据存入一个二维数组中,每个地形数据块可由9个高程
特征点构成,用一个大小的二维数组进行存储,通过行索引值和列索引值进
行访问,则任意地形数据块点集中某一高程点的高程值与原始地形数据块
(未经分层分块的地形数据块)的高程值相对应,因此可通过映射关系计
算出某一地形数据块行索引值和列索引值在原始地形数据块中对应的行索引值row和
列索引值col,计算公式如下:
优选的,所述步骤S33的节点细分评价体系具体公式如下:
本发明的有益效果是:本发明方法提出的三维SPH与Metaball融合建模方法,有机结合metaball技术、SPH建模以及LOD地形优化算法,充分发挥三项技术的各自优势,使得滑坡模拟不会因为SPH与Metaball融合建模技术由于模拟计算量过大,加之地形渲染耗费大量的计算资源,导致滑坡模拟无法实现,更为重要的是本发明在取得较好的数值模拟结果的同时,有效保证了滑坡可视化模拟的效果。
附图说明
图1为本发明方法步骤流程图;
图2为实施例中SPH粒子互相靠近时外部等值面的变化情况图;
图3为实施例中SPH滑坡粒子群外部生成矩形包围盒并分割包围盒示意图;
图4为实施例中滑坡运动堆积特征的可视化模拟图,图4(a)为滑坡碎屑化过程,图4(b)为滑坡碎屑物质重新粘合过程;
图5为实施例中地形节点点集示意图;
图6为实施例中滑坡地形简化渲染效果图;
图7(a)为0s时运动堆积状态厚度截图,图7(b)为10s时运动堆积状态厚度截图,图7(c)为20s时运动堆积状态厚度截图,图7(d)为40s时运动堆积状态厚度截图,图7(e)为60s时运动堆积状态厚度截图,图7(f)为100s时运动堆积状态厚度截图;
图8 为实施例中模拟滑体运动足迹与实际滑坡堆积范围对比图;
图9为实施例中茂县滑坡灾害场景模拟与实际情况对比图,图9(a)是当t=0s时对茂县滑坡灾害场景的可视化模拟效果,图9(b)是当t=20s时对茂县滑坡灾害场景的可视化模拟效果,图9(c)是当t=30s时对茂县滑坡灾害场景的可视化模拟效果,图9(d)是当t=40s时对茂县滑坡灾害场景的可视化模拟效果,图9(e)是当t=100s时对茂县滑坡灾害场景的可视化模拟效果,图9(f)滑坡实际情况;
图10为实施例中滑坡模拟程序渲染帧率随运行时间变化图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
与三维数值模拟相比,二维数值模拟具有便捷和易于实现的优势,因此目前滑坡数值模拟仍以二维数值模拟为主,可以更有效地研究滑坡灾害机理,但在滑坡危险范围的预测和可视化表达方面仍然存在局限性,特别是二维数值模拟展现的二维可视化信息难以展现滑坡至灾的三维细节信息。三维SPH数值模拟结合三维可视化模拟是滑坡模拟和预测向着多维度、多源信息表达方向发展可行有效的方案,具有很好的发展前景。
在滑坡模拟三维可视化方面,Metaball技术是一种适用于建立可变形表面的技术,与SPH融合建模时常被用于构建三维流体表面,本发明人尝试将metaball技术与SPH融合建模,从三维可视化模拟效果上看,metaball技术与SPH融合建模可为SPH粒子生成多边形包围体,当SPH粒子相互靠近时,相邻SPH粒子的包围体会互相融合生成新的包围体,由于包围体在可视化时能够模拟滑坡碎屑物质的质感且便于进行纹理映射,因此可弥补SPH粒子在滑坡可视化表达方面的不足,是对SPH滑坡数值模拟的重要补充,但是该项融合建模技术由于模拟计算量过大,且地形渲染将耗费大量的计算资源,因此简单的融合metaball技术与SPH建模方法无法在大规模滑坡模拟中得以实际应用。
另一方面,大规模滑坡滑程远、堆积范围广,而地形是滑坡模拟的边界,因此滑坡模拟存在的另一个重要问题是滑坡地形数据的调度和可视化问题。当前,复杂地形边界条件下的滑坡模拟常常将大规模数字高程模型数据一次性加载到计算机内存中并全部渲染,再加上滑坡模型复杂度高、科学计算复杂,因此常造成滑坡模拟程序运行速度缓慢甚至无法运行的后果。基于以上问题,发明人通过长时间的实验与论证,克服了上面的问题,具体方案如下:
本发明提供一种技术方案:一种面向滑坡危险范围预测的三维建模与可视化方法,如图1所示,包括如下步骤:
采用SPH法构建滑坡危险范围数值模型模拟滑坡运动;
采用Metaball与SPH融合可视化建模,为SPH粒子生成多边形包围体,再现滑坡运动堆积特征;
使用四叉树多尺度地形渲染算法进行地形渲染,采用不同尺度的地形模型表示地形,利用塔形四叉树结构管理地形数据,通过动态调度地形数据块以及不同尺度地形数据块的组合渲染。
SPH法滑坡建模
在SPH方法中有两种密度近似法:一种是对密度直接应用SPH近似法,称为密度求和法;另一种方法是连续性密度法。由于连续性密度法对于具有强间断问题的模拟(如爆炸、高速冲击等)可以得到较好的结果。
因此,本发明选择连续性密度法(Monaghan 1994)对粒子密度进行更新:
粒子的运动方程(Becker and Teschner 2007)如下式所示:
由式(3)和式(4)可得到粒子的加速度,加速度改变粒子速度,速度改变粒子位置,如此,通过使用简单易用的蛙跳法(Leap-Frog,LF)(Desbrun and Cani1996)进行时间积分,即可模拟滑坡运动。
可视化渲染
从高速远程滑坡的运动堆积特征来看,滑体从启动到静止既有崩滑后破碎的过程,也有受到雨水冲刷后崩滑体粘合的作用。尽管SPH在数值模拟方面具有一定的优势,但SPH滑坡数值模型是由数万个SPH离散粒子组成的滑体,因此完全基于SPH方法的高速远程滑坡堆积过程三维可视化,无法真实地模拟滑坡运动堆积特征。尽管Metaball方法无法提升滑坡危险范围数值模拟的正确性,但通过使用Metaball技术能够从可视化角度真实再现滑坡运动堆积特征。
Metaball与SPH融合可视化建模方法如下:
设整个三维空间中存在一个能量场,能量场中任意一个空间点具有的能量值等于所有SPH粒子对该点造成的影响的总和,通过使用与SPH法中的密度求和法相同的式子进行计算,能量值计算公式如下:
等值面的构建采用Marching Cubes(MC)算法,该算法最早由Lorensen和Cline于1987年提出(Lorensen and Cline 1987),是从三维离散数据场中提取等值面的经典算法,具有生成等值面质量好、实现简单等优点。
第二步,遍历所有正方体单元,对于每一个正方体单元,计算其八个顶点各自的能
量值,并将这些能量值分别与预先设定的能量值相比较,若大于或等于,则说明该顶点
位于等值面的内部,标记为0;若小于,则说明该顶点位于等值面的外部,标记为1。接着将
八个顶点的标记情况按一定的顶点编号顺序存入一个八位二进制数中,该二进制数的每一
位分别对应一个顶点的标记情况。
由于每个顶点存在两种标记情况,因此八个顶点共存在有256种情况,基于这256种情况建立一个边搜索表EdgeTable和一个三角形搜索表TriTable,EdgeTable保存正方体单元的边与等值面相交的256种情况,TriTable保存在正方体单元内部生成三角形面片的256种情况的三角形顶点连接顺序。
第三步,将第二步中得到的八位二进制数转为十进制数代入EdgeTable中得到一个十二位二进制数,新的二进制数的每一位分别对应正方体单元的一条边,表示该正方体单元的哪些边与等值面相交,然后通过式(7)进行线性插值,计算出这些边各自与等值面的交点坐标:
最后一步,通过查找TriTable,依据TriTable提供的连接顺序,将第三步中所得的
交点连接成一个个三角形面片,用由这些三角形面片组成的多边形面来逼近等值面。理论
上来说,单元宽度设置得越小,则分割出来的正方体单元越多,最后生成的多边形
面越接近等值面,但也意味着计算量的急剧上升,因此选择一个合适的也很重要。
由于等值面完全包裹SPH粒子,因此在三维空间中可称之为包围体。对于单个SPH
粒子,Marching Cubes算法中的单元宽度设置得越小,则生成的多边形包围体越近
似于球体。在模拟水流时,通常会取较小的以模拟光滑的水珠。但对于滑坡来说,其
碎屑物质通常是不规则的且棱角分明的。因此在模拟滑坡时, 可以取较大的值,如
此既可生成粗糙的包围体,模拟碎屑物质的质感,又可减少Marching Cubes算法的计算量,
提升程序运行速度。
由于Metaball的特性,当SPH粒子相互远离时其对应的多边形包围体会互相分离,基于这个特性,可模拟高速远程滑坡崩滑后由于滑体内部粘聚力减弱以及重力作用下不断破碎、碎屑化的可视化过程,如图4中图4(a)所示。与之相反的,当SPH粒子相互靠近时其对应的多边形包围体会互相融合,依据这个特性,可模拟滑体受到雨水冲刷作用后滑体碎屑化以后重新粘合的可视化过程,如图4(b)所示。
相较于粒子,用多边形包围体对滑体进行渲染更符合实际情况,既可利用相邻包围体互相融合的特性表达滑坡整体运动趋势,又可利用相互远离的包围体互相分离的特性模拟部分滑坡滑动后碎屑化的过程,从可视化角度反映滑坡运动堆积特征。
滑坡地形建模
地形是滑坡模拟的边界,是滑坡场景的重要组成部分。本发明在渲染场景时,利用滑坡区域数字高程模型DEM数据构建地形,将滑坡区域数字正射影像图(DigitalOrthophoto Map,DOM)作为纹理映射至地形上,生成滑坡地形。
为了保证滑坡模拟的顺利进行,需要将DEM数据一次性全部读入内存中,如果在模拟程序运行时将DEM数据全部渲染将会大大降低程序运行速度。使用LOD技术是解决上述问题的有效途径。常用的LOD算法有:基于四叉树的LOD算法、实时优化自适应网格算法、递进网格算法以及状态无关单通过自适应细化算法。其中,由于四叉树结构方便对地形数据进行分层分块和管理,因此基于四叉树的LOD算法便成了首选。
为此,本发明以基于四叉树的LOD算法思想为框架,提出了一种四叉树多尺度地形渲染算法,采用不同尺度的地形模型表示地形,依据人类视觉原理,在靠近观察者的区域用大量小尺度地形模型渲染高分辨率地形,在远离观察者的区域用少量大尺度地形模型渲染低分辨率地形,运用塔形四叉树结构管理地形数据,通过动态调度地形数据块以及不同尺度地形数据块的组合渲染,达到简化地形渲染的目的,在保证可视化效果的前提下,提高地形渲染效率,提升程序运行速度。
滑坡地形数据的组织与调度
本发明采用塔形四叉树结构对DEM和DOM数据进行组织与管理。将整个滑坡地形区域按四叉树结构划分成一个个小的地形数据块,每一个地形数据块可看成一个四叉树节点,称为地形节点。从四叉树根节点开始,除了叶子节点外,每个节点均可细分产生下一个层次的四个子节点。
每个地形节点所对应的高程点集直接从滑坡区域DEM数据中选取,如图5所示,DOM数据同理。
选取方法如下:将DEM数据存入一个二维数组中,数组中的每个元素用于存
储其对应点位的高程值,数组的行数和列数与DEM数据的行数和列数相同。本发明采用三角
形扇的形式绘制每个地形节点,每个地形节点均由9个高程点构成,可用一个大小的
二维数组进行存储,通过行索引值和列索引值进行访问。则任意地形数据块点集中某一高程点的高程值与原始地形数据块(未经分层分块的地形数据块)
的高程值相对应,因此可通过映射关系计算出某一地形数据块行索引
值和列索引值在原始地形数据块中对应的行索引值row和列索引值col,计算公式如下:
地形节点的调度主要由两个方面的内容来决定:①地形节点的细分程度;②地形节点的取舍。
针对第一点,本发明构建了一个地形节点细分评价体系(Liu et al. 2010),用于对地形节点细分与否做出判断,具体公式如下:
针对第二点,本发明在每次绘制地形节点前,都对节点进行了视景体裁剪,仅保留处于视景体内部的节点,而舍弃位于视景体外部的节点,并由此判断需要调度的节点对象。另外,由于地形节点本身的不规则性,本发明使用节点包围盒(Decoret et al.1999)来替代原有节点进行视景体裁剪判断。
四叉树多尺度地形渲染算法
为了满足滑坡数值模拟过程中边界处理的需要,在程序运行过程中,需要将滑坡区域DEM数据完整读入内存中,但在渲染时没有必要渲染全部数据,仅在保证可视化效果的前提下渲染部分数据即可。
为此,本发明提出了一种四叉树多尺度地形渲染算法控制地形渲染:
1)创建两条队列(队列1用于存放正在处理的当前层次的地形节点,队列2用于存放待处理的下一个层次的地形节点)。
2)采用广度优先遍历四叉树的原则从根节点开始对四叉树进行遍历。在遍历四叉树时,一个四叉树节点即为一个地形节点,该节点所对应的高程点集直接从滑坡区域DEM数据中选取,DOM数据同理,选取方法如上所述。
3)将遍历到的节点送入队列1,判断该节点在视景体内是否可见。若是不可见的,则直接出队列,不予处理;若是可见的,则使用节点细分评价体系(见上面的内容)判断该节点是否需要细分。若不需要细分,则直接渲染该节点;若需要细分,则将该节点出列,进行细分,并将该节点细分产生的下一个层次的四个子节点送入队列2,等待下一轮的处理。
4)当程序运行到队列1为空时,说明某一个层次的所有节点均已处理完毕,则交换队列1和队列2中的节点,并重复如上处理步骤,直到所有需要处理的节点均已处理完毕为止。对于那些在视景体内可见的叶子节点,则直接进行渲染。
通过利用以上算法控制地形渲染,可有效减少不必要的渲染资源开销,提升程序运行速度,保证滑坡场景的渲染效率。对于不同层次地形节点拼接处产生的裂缝(Zhang etal. 2005),采用跳点法(Decoret et al.1999)进行消除。如图6所示,图6用格网模式展示了滑坡地形简化渲染效果。
实验验证分析
本发明以2017年四川茂县滑坡为例,搭建了滑坡灾害场景模拟系统,对滑坡危险范围进行预测。实验平台硬件参数如下:处理器:英特尔 Core i5-9400F,显卡:NvidiaGeForce GTX 1050 Ti。编程环境为C++和OpenGL。
1)滑坡危险范围数值模拟与分析
茂县滑坡从小规模垮塌到主滑体启动、失稳运动、到最终停积整个过程用时约
100s,其中主滑坡仅持续60s,滑坡启动体积约为,停止时平均堆积厚度大于
10m,最大堆积厚度大于32m,图7中黑色轮廓线圈出了滑坡堆积范围(Ouyang et al.
2017b; Xu et al. 2017)。
本发明对茂县滑坡的模拟结果如图7(a)-(f)所示:图中展示了模拟滑体在0s、10s、20s、40s、60s、100s这些不同时刻的运动堆积状态(时间序列按照从左到右从上到下的次序),滑体的不同颜色代表不同的堆积厚度。0s至10s期间失稳的山体在重力作用下开始逐渐垮塌,10s至40s期间滑体沿着斜坡坡面高速运动并不断解体破碎,40s左右滑体前部到达坡脚平坦的地面,由于受到阻挡,滑体运动速度开始减缓,40s至60s期间滑体在平坦的地面上不断堆积,60s后整个滑坡趋于停止,60s至100s期间是滑坡的最后调整阶段,仍有一些碎石掉落,100s后滑坡几乎停止。模拟滑体在t=100s时的平均堆积厚度为7.3m。模拟滑体的运动堆积过程与滑坡实际情况吻合较好。
为了准确度量模拟滑坡堆积范围与实际滑坡堆积范围的吻合程度,本发明在滑坡模拟过程中每间隔1s便对模拟滑体的运动状态进行截图(参照图7(a)-(f),在模拟结束后将这些截图叠加合成一张模拟滑体运动足迹图,该图中记录有每个SPH粒子的运动足迹。接着,这里引入一个模拟滑坡拟合度的概念,用以描述模拟滑坡与实际滑坡的相似程度,如图8所示,其具体数值通过对比模拟滑体运动足迹图与实际滑坡堆积范围图的相似程度来确定,计算公式如下:
式中,为实际滑坡在图像中所占的像素数量(对应图8中的黑色区域加灰色区
域),为模拟滑坡在图像中重合于实际滑坡的像素数量(对应图8中的灰色区域),为
模拟滑坡在图像中不重合于实际滑坡的像素数量(对应图8中的白色区域),为模拟滑坡
拟合度值,该值越接近1说明模拟滑坡堆积范围与实际滑坡堆积范围吻合程度越高,模拟滑
坡与实际滑坡相似程度越高,滑坡模拟效果越好,反之则越差。
通过使用式(11)进行计算可得本发明模拟滑坡拟合度值为0.8542,模拟滑坡堆积范围与实际滑坡堆积范围的吻合程度较高,滑坡模拟效果较好。
2)可视化渲染与分析
图9为实施例中茂县滑坡灾害场景模拟与实际情况对比图,图9(a)-(f)中黑色轮廓线圈出了滑坡堆积范围。图9(a)至图9(e)展示了滑坡启动后的0s至100s期间对滑坡灾害场景的可视化模拟效果。图9(a)中用裂缝线圈出了模拟滑体初始化轮廓。图9(a)至图9(e)中使用了SPH与Metaball融合建模方法模拟滑体,相较于图7(a)-(f)中使用纯SPH方法模拟滑体,在进行纹理映射后,前者的可视化效果明显优于后者。图9(e)展示了当t=100s时对茂县滑坡灾害场景的可视化模拟效果,与实际情况(图9(f))相比,二者相似程度较高,滑坡堆积范围吻合较好。滑坡模拟所用SPH粒子数量为45630个,模拟场景中滑坡地形三角面的渲染数量是随着视点移动而动态变化的,但总体控制在10万个以下。滑坡模拟程序渲染帧率变化情况如图10所示。
通过对比未做简化地形渲染操作和使用本发明方法简化地形渲染这两种情况下程序渲染帧率随运行时间变化的情况,可以看出:在使用本发明方法简化地形渲染后,整个程序的运行速度得到明显提升。另外,当程序运行时间超过100s后,SPH数值计算程序停止运行,仅剩可视化渲染程序在持续运行,此时,未做简化地形渲染操作的情况下程序渲染帧率仅维持在18帧/秒左右,而使用本发明方法简化地形渲染的情况下程序渲染帧率能维持在60帧/秒以上,且两者间的差距会随地形数据量的增加而不断扩大,该现象对于程序运行时间未超过100s时亦是如此。
本发明主要由三类建模方法构成,分别为SPH法滑坡建模、Metaball可视化渲染以及LOD滑坡地形建模方法构成。滑坡模拟主要通过数值模拟预测滑坡危险范围,可视化模拟反映灾情信息,数值模拟的准确性是滑坡危险范围预测最为核心的问题,而模拟的准确性和可视化灾情信息都需要模拟发挥模拟计算性能。由于地形(由DEM地形数据和DOM遥感数据构成)是滑坡模拟的边界,因此SPH与地形的交互耦合必不可少,由于大规模滑坡如茂县滑坡模拟时SPH数值模拟的粒子数必须达到一定数量(本次模拟达到了45630个),才能保证模型达到85.42%的危险范围预测的正确率,如果没有LOD地形优化方法,地形渲染将会占用大量的计算资源,导致分配给SPH数值模拟的计算资源大幅减少,势必造成模拟所用的SPH粒子数的大幅减少,滑坡危险范围预测将因为无法发挥模拟计算性能,导致正确率大大降低。同理,SPH/Metaball融合可视化建模时,LOD地形优化方法也是发挥模拟计算性能,有效发挥可视化表达灾情信息最可靠的保证。没有SPH/Metaball与LOD三类建模方法的深度融合,加载的SPH粒子数大幅减少,势必造成可视化模拟滑坡堆积时滑坡岩土的分离和融合的模拟非常粗糙,完全不能模拟出滑坡的灾情信息和运动堆积效果。数值模拟和可视化模拟需要将LOD优化方法、SPH以及Metaball的深度融合,是准确模拟滑坡危险范围、高保真表达滑坡灾情信息最为有效的手段,可创造性的解决了因模拟性能不佳所造成的模拟准确性问题,以及滑坡灾情信息可视化表达力不强的问题。
综上,本发明提出的三维SPH与Metaball融合建模方法在取得较好的数值模拟结果的同时,可有效保证滑坡可视化模拟的效果,能够提供具有较高准确性的滑坡危险范围预测结果,提出的四叉树多尺度地形渲染算法能有效简化地形渲染,提升程序运行速度,本发明方法是行之有效的。通过2017年四川茂县高速远程滑坡危险范围模拟验证了模型在数值模拟和可视化模拟方面的卓越特性,为滑坡模拟及危险范围预测向着多维度、多源信息表达方向发展提供了行之有效的方案。
尽管参照前述实施例对本发明进行了详细的说明,对于本领域的技术人员来说,其依然可以对前述各实施例所记载的技术方案进行修改,或者对其中部分技术特征进行等同替换,凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。
Claims (2)
1.一种面向滑坡危险范围预测的三维建模与可视化方法,其特征在于,包括如下步骤:
S1、采用SPH法构建滑坡危险范围数值模型模拟滑坡运动;
S2、采用Metaball与SPH融合可视化建模,为SPH粒子生成多边形包围体,再现滑坡运动堆积特征;
S3、使用四叉树多尺度地形渲染算法进行地形渲染,采用不同尺度的地形模型表示地形,利用塔形四叉树结构管理地形数据,通过动态调度地形数据块以及不同尺度地形数据块的组合渲染;
所述步骤S1中的SPH法具体是指:通过连续性密度法对粒子密度进行更新,公式如下:
式中,ρi为粒子i的密度,mj为粒子j的质量,vi和vj分别代表粒子i和粒子j的速度,W(rij,h)为搜索半径为h的光滑核函数,rij=ri-rj,ri和rj分别代表粒子i和粒子j的位置坐标;
粒子的运动方程如下式所示:
式中,Pi和Pj分别代表粒子i和粒子j的压强,g代表外力作用;
将人工粘度Πij引入到计算中,保证流体模拟的数值稳定性,公式如下:
其中,Πij的定义如下:
由此,得到粒子的加速度,加速度改变粒子速度,速度改变粒子位置,通过使用蛙跳法进行时间积分,即可完成对滑坡运动模拟;
所述步骤S2中Metaball与SPH融合可视化建模具体包括:设整个三维空间中存在一个能量场,能量场中任意一个空间点具有的能量值等于所有SPH粒子对该点造成的影响的总和,能量值计算公式如下:
E(ri)=∑jmjw(ri-rj,h);
设定一个能量值E′,则可以在SPH粒子外部生成一个由具有E′能量值的空间点构建的等值面,当SPH粒子互相靠近时,空间点的能量值随之发生变化,进而影响生成的等值面的形状,形成融合的效果;
所述的等值面的构建采用MC算法,所述MC算法包括如下步骤:
S21、为SPH粒子群生成一个矩形包围盒,并按设置的单元宽度width将该包围盒分割成一个个规则的正方体单元;
S22、遍历所有正方体单元,对于每一个正方体单元,计算其八个顶点各自的能量值,并将这些能量值分别与预先设定的能量值E′相比较,若大于或等于E′,则说明该顶点位于等值面的内部,标记为0;若小于E′,则说明该顶点位于等值面的外部,标记为1;
将八个顶点的标记情况按顶点编号顺序存入一个八位二进制数中,该二进制数的每一位分别对应一个顶点的标记情况,
建立边搜索表EdgeTable和三角形搜索表TriTable,EdgeTable保存正方体单元的边与等值面相交的情况,TriTable保存在正方体单元内部生成三角形面片的情况的三角形顶点连接顺序;
S23、将步骤S22中得到的八位二进制数转为十进制数并代入EdgeTable中得到一个十二位二进制数,新的十二位二进制数的每一位分别对应正方体单元的一条边,表示该正方体单元的哪些边与等值面相交,然后进行线性插值,计算出这些边各自与等值面的交点坐标;
S24、查找TriTable,依据TriTable提供的连接顺序,将步骤S23中所得的交点连接成一个个三角形面片,用这些三角形面片组成的多边形面来逼近等值面;
所述线性插值的公式如下:
其中,r为正方体单元的边与等值面的交点三维坐标,r1和r2分别为边的两个端点的三维坐标,E1和E2分别为边的两个端点的能量值;
所述步骤S3中的使用四叉树多尺度地形渲染算法进行地形渲染具体包括如下步骤:
S31、创建队列1和队列2,队列1用于存放正在处理的当前层次的地形节点,队列2用于存放待处理的下一个层次的地形节点;
S32、采用广度优先遍历四叉树的原则从根节点开始对四叉树进行遍历, 在遍历四叉树时,一个四叉树节点即为一个地形节点,该节点所对应的高程点集直接从滑坡区域DEM数据中选取;
S33、将遍历到的节点送入队列1,判断该节点在视景体内是否可见,若是不可见的,则直接出队列,不予处理;若是可见的,则使用节点细分评价体系判断该节点是否需要细分,若不需要细分,则直接渲染该节点;若需要细分,则将该节点出列,进行细分,并将该节点细分产生的下一个层次的四个子节点送入队列2,等待下一轮的处理;
S34、当程序运行到队列1为空时,则交换队列1和队列2中的节点,并重复步骤S31至步骤S33,直到所有需要处理的节点均已处理完毕为止;对于那些在视景体内可见的叶子节点,则直接进行渲染;
所述步骤S32中节点点集选取的具体方法如下:将DEM数据存入一个二维数组dem[row][col]中,每个地形数据块可由9个高程特征点构成,用一个3×3大小的二维数组cell进行存储,通过行索引值i和列索引值j进行访问,则任意地形数据块点集cell中某一高程点的高程值cell[i][j]与原始地形数据块的高程值dem[row][col]相对应,因此,可通过映射关系计算出某一地形数据块行索引值i和列索引值j在原始地形数据块中对应的行索引值row和列索引值col,计算公式如下:
式中,N为当前地形节点的索引值,L为当前地形节点的层次值。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210229728.6A CN114359500B (zh) | 2022-03-10 | 2022-03-10 | 一种面向滑坡危险范围预测的三维建模与可视化方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210229728.6A CN114359500B (zh) | 2022-03-10 | 2022-03-10 | 一种面向滑坡危险范围预测的三维建模与可视化方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114359500A CN114359500A (zh) | 2022-04-15 |
CN114359500B true CN114359500B (zh) | 2022-05-24 |
Family
ID=81095264
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210229728.6A Active CN114359500B (zh) | 2022-03-10 | 2022-03-10 | 一种面向滑坡危险范围预测的三维建模与可视化方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114359500B (zh) |
Families Citing this family (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115577571B (zh) * | 2022-11-25 | 2023-03-10 | 海南浙江大学研究院 | 一种控制滑坡体体积量的海底滑坡海啸的建模方法 |
CN116504032B (zh) * | 2023-06-28 | 2023-09-22 | 湖南科技大学 | 一种基于实景三维的滑坡危险性监测预警方法及系统 |
CN117688784B (zh) * | 2024-02-02 | 2024-05-17 | 中国科学院地质与地球物理研究所 | 基于等效摩擦系数和趋势分析的滑坡威胁范围成图方法 |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110738721A (zh) * | 2019-10-12 | 2020-01-31 | 四川航天神坤科技有限公司 | 基于视频几何分析的三维场景渲染加速方法及系统 |
CN111475937A (zh) * | 2020-04-03 | 2020-07-31 | 中国地质科学院地质力学研究所 | 一种流固二相流流化滑坡的模拟仿真方法 |
CN113436316A (zh) * | 2021-06-28 | 2021-09-24 | 哈尔滨理工大学 | 一种基于空间投影的适用于sph算法的邻近粒子搜索方法 |
Family Cites Families (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106875462B (zh) * | 2017-01-13 | 2020-07-07 | 北京航空航天大学 | 一种基于元球模型和混合驱动方法的实时数字器官切割方法 |
CN110930509A (zh) * | 2019-11-27 | 2020-03-27 | 王程 | 线性四元树多层模型驱动的dem即时可视化方法 |
CN113762090B (zh) * | 2021-08-16 | 2024-03-08 | 国网浙江省电力有限公司湖州供电公司 | 一种特高压密集输电通道灾害监测预警方法 |
-
2022
- 2022-03-10 CN CN202210229728.6A patent/CN114359500B/zh active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110738721A (zh) * | 2019-10-12 | 2020-01-31 | 四川航天神坤科技有限公司 | 基于视频几何分析的三维场景渲染加速方法及系统 |
CN111475937A (zh) * | 2020-04-03 | 2020-07-31 | 中国地质科学院地质力学研究所 | 一种流固二相流流化滑坡的模拟仿真方法 |
CN113436316A (zh) * | 2021-06-28 | 2021-09-24 | 哈尔滨理工大学 | 一种基于空间投影的适用于sph算法的邻近粒子搜索方法 |
Non-Patent Citations (3)
Title |
---|
Research on Three-Dimensional Visualization System for Landslide and Mud-rock Flows based on Direct3D and SPH Methods;Jian Gao 等;《Applied Mechanics and Materials》;20141008;第675-677卷;1179-1183 * |
面向大规模滑坡灾害模拟的地形建模与三维可视化;吕奕杰 等;《武汉大学学报(信息科学版)》;20190722;第45卷(第3期);467-474 * |
高速远程滑坡颗粒流研究进展;李坤 等;《地球科学》;20211102;第47卷(第3期);893-912 * |
Also Published As
Publication number | Publication date |
---|---|
CN114359500A (zh) | 2022-04-15 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN114359500B (zh) | 一种面向滑坡危险范围预测的三维建模与可视化方法 | |
CN110738721B (zh) | 基于视频几何分析的三维场景渲染加速方法及系统 | |
CN112347546A (zh) | 基于轻量级设备的bim渲染方法、设备和计算机可读存储介质 | |
Gao et al. | Feature suppression based CAD mesh model simplification | |
US8570322B2 (en) | Method, system, and computer program product for efficient ray tracing of micropolygon geometry | |
Peytavie et al. | Arches: a framework for modeling complex terrains | |
CN108537869B (zh) | 一种基于级联纹理的圆锥追踪动态全局光照方法 | |
EP1754196B1 (en) | Resource management for rule-based procedural terrain generation | |
Ernst et al. | Early split clipping for bounding volume hierarchies | |
CN108919954B (zh) | 一种动态变化场景虚实物体碰撞交互方法 | |
CN109118588B (zh) | 一种基于块分解的彩色lod模型自动生成方法 | |
CN103714575A (zh) | 一种sph与动态表面网格相结合的流体仿真方法 | |
US7439970B1 (en) | Computer graphics | |
WO2024037116A1 (zh) | 三维模型的渲染方法、装置、电子设备及存储介质 | |
CN112906125A (zh) | 铁路固定设施bim模型轻量化加载方法 | |
Onoue et al. | An interactive deformation system for granular material | |
He et al. | All range and heterogeneous multi-scale 3D city models | |
Wiemann et al. | Automatic Map Creation For Environment Modelling In Robotic Simulators. | |
CN114170394B (zh) | 一种海量倾斜数据在Web端的展示优化方法及装置 | |
CN110930509A (zh) | 线性四元树多层模型驱动的dem即时可视化方法 | |
CN107767458B (zh) | 不规则三角网曲面几何拓扑一致分析方法及系统 | |
Glander et al. | Techniques for generalizing building geometry of complex virtual 3D city models | |
CN115033972A (zh) | 一种建筑主体结构批量单体化方法、系统及可读存储介质 | |
CN114419256A (zh) | 基于多级抽壳算法的城市级bim数据轻量化方法及系统 | |
Chang et al. | Hierarchical simplification of city models to maintain urban legibility. |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |