CN109993830B - 一种软组织表面刺穿形变仿真方法及终端 - Google Patents
一种软组织表面刺穿形变仿真方法及终端 Download PDFInfo
- Publication number
- CN109993830B CN109993830B CN201910284423.3A CN201910284423A CN109993830B CN 109993830 B CN109993830 B CN 109993830B CN 201910284423 A CN201910284423 A CN 201910284423A CN 109993830 B CN109993830 B CN 109993830B
- Authority
- CN
- China
- Prior art keywords
- soft tissue
- tissue surface
- geometric model
- surface geometric
- position vector
- 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
- 210000004872 soft tissue Anatomy 0.000 title claims abstract description 284
- 238000000034 method Methods 0.000 title claims abstract description 52
- 238000004088 simulation Methods 0.000 title claims abstract description 49
- 239000013598 vector Substances 0.000 claims abstract description 133
- 239000002245 particle Substances 0.000 claims abstract description 116
- 230000014759 maintenance of location Effects 0.000 claims abstract description 47
- 238000004364 calculation method Methods 0.000 claims description 13
- 238000004590 computer program Methods 0.000 claims description 4
- 238000004321 preservation Methods 0.000 abstract description 13
- 238000010586 diagram Methods 0.000 description 6
- 238000004422 calculation algorithm Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 4
- 230000006870 function Effects 0.000 description 4
- 230000009286 beneficial effect Effects 0.000 description 3
- 230000008569 process Effects 0.000 description 3
- 239000000243 solution Substances 0.000 description 3
- 230000009471 action Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 239000011159 matrix material Substances 0.000 description 2
- 230000007704 transition Effects 0.000 description 2
- 238000010146 3D printing Methods 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 238000009472 formulation Methods 0.000 description 1
- 238000009432 framing Methods 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 238000002347 injection Methods 0.000 description 1
- 239000007924 injection Substances 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000012886 linear function Methods 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000009877 rendering Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000012549 training Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
- 238000000844 transformation Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30004—Biomedical image processing
Landscapes
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Computer Graphics (AREA)
- Geometry (AREA)
- Software Systems (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明提供的一种软组织表面刺穿形变仿真方法及终端，通过对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量，保证了在不丢失软组织表面受力形变的细节的同时实现软组织表面刺穿的实时仿真。
Description
技术领域
本发明涉及虚拟手术中软组织形变仿真领域，特别涉及一种软组织表面刺穿形变仿真方法及终端。
背景技术
虚拟手术在医学领域中具有广泛的应用，其将计算机图形学与虚拟技术结合，为医护人员提供手术方案制定与手术操作训练等功能。刺穿广泛应用于日常医学检查和诊断中，刺穿针刺入软组织表面到达目标位置，实现注射、取样等操作，因此软组织表面刺穿形变仿真是虚拟手术的关键技术之一。在传统虚拟手术系统中进行模拟时，软组织表面刺穿形变只简单的对模型进行拉伸形变，丢失了现实中软组织表面受力形变的一些细节效果，从而影响了模型的真实感。
申请号为CN201410723354.9的中国专利公开了一种基于子块的软组织表面变形追踪方法，该方法用于实现术中软组织变形矫正。
申请号为CN201810020385.6的中国专利公开了一种利用数字化空间重构及3D打印制备用于确定软组织表面替代皮瓣尺寸的切取导板的方法利用该发明制得的切取导板，能够精确地在供区切取替代皮瓣，从而实现以最小的皮瓣供区损伤获得最精确的受区适形功能重建。
申请号为CN201510020299.1的中国专利公开了一种基于凝胶拓印及图形像素技术获取软组织表面形变方法,涉及一种服装人体工效学研究软组织表面形变图与计算形变率获取软组织表面形变测量的方法,目的是解决获取软组织表面形变图可操作性差及动作或部位局限，利用此方法可以获取任何动作下的部位软组织表面形变图以及计算出网格软组织表面形变率。
尽管对于软组织表面形变方法的研究成果较为丰富，但是在某些方面依然存在一些问题，如低精度的软组织表面模型形变产生的细节不足，提高细节时需要产生额外的计算量，在保证局部形变的同时很难兼顾到全局的形变等问题。以上方法都只提到了软组织表面形变的仿真方法，但均未涉及软组织表面刺穿仿真方法。
本发明针对目前软组织表面形变存在的问题，结合虚拟手术中对软组织表面进行刺穿需要达到实时交互、真实感高等要求，以软组织表面为对象进行刺穿仿真研究，给了一种软组织表面刺穿形变仿真方法及终端。
发明内容
本发明所要解决的技术问题是：提供一种软组织表面刺穿形变仿真方法及终端，能够保证在不丢失软组织表面受力形变的细节的同时实现软组织表面刺穿的实时仿真。
为了解决上述技术问题，本发明采用的一种技术方案为：
一种软组织表面刺穿形变仿真方法，包括步骤：
S1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；
S2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；
S3、所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；
S4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量。
为了解决上述技术问题，本发明采用的另一种技术方案为：
一种软组织表面刺穿形变仿真终端，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现以下步骤：
S1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；
S2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；
S3、所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；
S4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量。
本发明的有益效果在于：通过对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量，保证了在不丢失软组织表面受力形变的细节的同时实现软组织表面刺穿的实时仿真。
附图说明
图1为本发明实施例的软组织表面刺穿形变仿真方法流程图；
图2为本发明实施例的软组织表面刺穿形变仿真终端的结构示意图；
图3为本发明实施例的软组织表面形变区域效果图；
图4为本发明实施例的曲面S所围成的区域；
图5为本发明实施例的软组织表面三角网格法向量；
图6为本发明实施例的余切Laplacian权值；
图7为本发明实施例的三角网格梯度示意图；
图8为本发明软组织表面刺穿形变仿真方法的效果图；
图9为本发明软组织表面刺穿仿真结果对比图；
标号说明：
1、软组织表面刺穿形变仿真终端；2、存储器；3、处理器。
具体实施方式
为详细说明本发明的技术内容、所实现目的及效果，以下结合实施方式并配合附图予以说明。
请参照图1，一种软组织表面刺穿形变仿真方法，包括步骤：
S1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；
S2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；
S3、所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量。
S4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量。
从上述描述可知，本发明的有益效果在于：通过对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量，保证了在不丢失软组织表面受力形变的细节的同时实现软组织表面刺穿的实时仿真。
进一步的，步骤S1包括：
S11、对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；
S12、根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程。
由上述描述可知，通过对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程，保证了软组织表面刺穿仿真的实时性。
进一步的，步骤S2还包括：
获取刺穿针与所述形变区域产生的碰撞信息。
由上述描述可知，通过获取刺穿针与所述形变区域产生的碰撞信息，便于后续对软组织表面几何模型进行求解。
进一步的，步骤S2具体为：
根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量。
由上述描述可知，通过根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量，保证了后续计算泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量的准确性。
进一步的，步骤S3包括：
S31、将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；
S32、将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量。
由上述描述可知，通过将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量，真实感强。
请参照图2，一种软组织表面刺穿形变仿真终端，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现以下步骤：
S1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；
S2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；
S3、所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量。
S4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量。
从上述描述可知，本发明的有益效果在于：通过对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量，保证了在不丢失软组织表面受力形变的细节的同时实现软组织表面刺穿的实时仿真。
进一步的，步骤S1包括：
S11、对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；
S12、根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程。
由上述描述可知，通过对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程，保证了软组织表面刺穿仿真的实时性。
进一步的，步骤S2还包括：
获取刺穿针与所述形变区域产生的碰撞信息。
由上述描述可知，通过获取刺穿针与所述形变区域产生的碰撞信息，便于后续对软组织表面几何模型进行求解。
进一步的，步骤S2具体为：
根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量。
由上述描述可知，通过根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量，保证了后续计算泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量的准确性。
进一步的，步骤S3包括：
S31、将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；
S32、将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量。
由上述描述可知，通过将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量，真实感强。
实施例一
请参照图1，一种软组织表面刺穿形变仿真方法，包括步骤：
S1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；
步骤S1包括：
S11、对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；
S12、根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程；
S2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；
步骤S2还包括：
获取刺穿针与所述形变区域产生的碰撞信息。
步骤S2具体为：
根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量；
S3、所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；
步骤S3包括：
S31、将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；
S32、将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；
S4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量。
实施例二
本实施例和实施例一的区别在于，本实施例将结合具体的应用场景进一步的说明本发明上述软组织表面刺穿形变仿真方法是如何实现的：
1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程
1.1、对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；
测地线是平面上直线在曲面上的拓展。在三角网格中，两个定点之间的测地线距离等于求解图的最短路径，使用精确的测地线算法求解得到测地线距离，D(s)表示起点到终点的最短路径长度，算法步骤如下：
步骤1.1.1：选择起点、终点。
步骤1.1.2：用dijkstra算法搜索边，计算经过搜索达到距离的上界。
步骤1.1.3：根据Euclidean disrance和曲面的距离函数得到距离的下界，且下界函数就满足Ls(p)+||p,et||≤Dst，其中||p,et||表示终点et到点p的距离，Dst表示测地线路径上一点到起点的距离，起点到当前窗口的一点p的曲面距离为Ls(p)。
步骤1.1.4：将窗口函数定义为W(b0,b1,d0,d1,K)，该窗口的值都为正数，K表示光源值，表示选择的起点传播与当前窗口光源之间的距离；由b0,b1确定半径方向上的窗口大小，如果当前窗口所在的边为||e||，当前窗口光源到b0,b1的距离分别为d0,d1，这时b0,b1与b0-b1可能成一条直线或为一个三角形，可以表示为|b0-b1|≤b0-b1≤|b0+b1|。
步骤1.1.5：如果窗口传播到的质点是衍生点或边界点这样的特殊点，这时需要将上一步的光源值与当前光源值进行比较，保证上一步光源值小于当前光源值，当窗口传播结束以后，使用衍生点或边界点进行窗口传播。
步骤1.1.6：结束窗口传播的条件为：窗口传播到终点所在平面，可以由公式L(s)＝min{Ls(p)+p,et}得到，如果L(s)pD(s)，那么L(s)＝D(s)。
步骤1.1.7：结束传播后，从终点开始进行遍历，最后得到遍历的路径，软组织表面形变区域如图3所示，(a)中阴影区域为刺穿形变区域，(b)为形变效果图。
1.2、根据散度定理得到软组织表面几何模型的体积V，并建立相应的体积保持约束方程；
图4为一个曲面S所围成的一个体积为V的区域且有对应关系
每个面都对应一个单位法向量。在软组织表面刺穿形变仿真中，软组织表面几何模型是一个封闭的结构，为曲面S包围的一个完整的边界。在仿真开始前先计算软组织表面几何模型的体积，然后在此后的刺穿形变仿真中的每一个时间步长执行恒定体积约束。软组织表面几何模型的体积和软组织表面几何模型的表面积分的三重积分关系用散度定理表示为：
∫∫Fv·NdS＝∫∫∫divFvdV (1)
其中，Fv表示软组织表面几何模型的连续可微向量场，N表示软组织表面几何模型表面三角网格的单位法向量，如图5所示，黑色线段表示软组织表面几何模型表面的单位法向量。
软组织表面几何模型的体积可以通过对其表面的三角网格求和的方式得到，假设向量Fv为位置向量r，因此3V＝∫∫r·NdS可以用如下公式计算：
其中，(x,y,z)表示三角网格的3个顶点坐标，Li表示区域坐标，A表示三角网格的面积。
在软组织表面几何模型中，软组织表面几何模型表面三角网格的单位法向量是一个常数，因此区域坐标Li可以用如下公式表示：
公式(3)中，a1,a2,a3表示非负整数，A表示软组织表面几何模型三角网格的面积。a1,a2,a3有三种情况：a1＝1,a2＝a3＝0，a2＝1,a1＝a3＝0，a3＝1，a1＝a2＝0。根据这三种情况，公式(3)可以表示为：
∫L1dxdy＝∫L2dxdy＝∫L3dxdy (4)
将公式(4)带入到公式(2)中，就能得到软组织表面几何模型的体积：
其中，nt表示软组织表面三角网格的网格数。根据公式(5)得到的软组织表面几何模型体积V在整个刺穿仿真过程中应保持不变。
2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量J；
具体的，获取刺穿针与所述形变区域产生的碰撞信息，并根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量J；
软组织表面几何模型的体积V，在每一刺穿仿真时间步长下体积V都将保持不变，因此就需要将体积V转化一个约束方程。本发明将体积V转化为约束方程的方法是使用拉格朗日乘数。对于软组织表面刺穿仿真来说，软组织表面几何模型的每一个质点都将受到体积V条件的约束。
首先，定义一个向量J表示软组织表面几何模型所有质点的三维坐标，如公式(6)所示，其中m表示软组织表面几何模型中质点的个数。
J＝[x1y1z1x2y2z2L xmymzm]T (6)
定义Φ(J,t)表示约束向量：
Φ(J,t)＝[Φ1(J,t)Φ2(J,t)LΦn(J,t)]T (7)
其中n表示软组织表面几何模型有n个约束，由于本文中只需要约束体积V，故n＝1。可将公式(6)改写为：
Φ(J,t)＝V-V0 (8)
其中V表示每个时间步长下软组织表面几何模型的体积，V0表示软组织表面几何模型原始体积。将公式(8)带入到牛顿第二定律中得到：
其中，F表示刺穿仿真中刺穿针对软组织表面施加的外力和虚拟弹簧的弹力作用在质点上的合力，M表示维度为3m×3m的对角质量矩阵，λ是包含拉格朗日乘数的维度为m×1的向量。将约束向量Φ(J,t)表示为：
Φ″+2αΦ′+β2Φ＝0 (10)
参数α和β确定Φ趋向于0的速度。对公式(4-15)的约束向量进行微分，可以得到：
整理后得到：
其中
在公式(13)中，下标J和t表示对应变量J和t的偏微分，将公式(9)代入到公式(13)中得到一个线性方程组，且需要求解得出拉格朗日乘数λ：
公式(14)可以使用显式数值计算的方法进行求解。但是显式欧拉法在大时间步长下有不稳定性的问题。出于稳定性考虑，使用隐式的方法进行求解，根据动力学公式可以得到：
J(t+Δt)＝J(t)+ΔtJ′(t+Δt) (16)
同样地，约束向量改写为隐式的表示方法并用泰勒级数展开：
将公式(16)代入到公式(15)得：
将公式(18)代入到公式(17)中消去J(t+Δt)得到：
在公式(19)中，λ为一个未知量。因为约束方程中只要约束软组织表面几何模型体积V，因此可以用简单的除法就能得到拉格朗日乘数λ的值。得到λ的值以后，将其带入公式(15)和公式(16)中，就能得到体积保持的下一时间步长的质点速度和位置。
3、所述软组织表面几何模型质点向量J通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；
3.1、将所述软组织表面几何模型质点位置向量J通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量Jb；
使用泊松形变算法，可以对软组织表面几何模型三角网格的三维坐标进行改变，使软组织表面在刺穿仿真中形变产生褶皱。把软组织表面几何模型的三角网格看作是一个标量场，形变指引场是软组织表面三角网格对应的梯度场，软组织表面几何模型的形变就可以通过改变三角网格所对应的梯度场产生形变，泊松方程公式表示为：
其中，w表示软组织表面几何模型三角网格的表面常矢量，矢量场的散度由▽·w表示，f表示待求解的形变后的软组织表面几何模型三角网格，f*提供了在边界
上可取的值，Δf为软组织表面三角网格的Laplacian算子，是该泊松方程的狄利克雷边界值。Laplacian坐标常用在三角网格的形变、平滑和去噪等方面。对于单个三角网格顶点，Laplacian算子定义为与其1-邻域顶点的质心偏差。将Laplacian坐标δi使用余切Laplacian权值表示为：
公式(21)中αij,βij表示边pij两个相邻三角网格的对角，如图6所示。
将软组织表面几何模型三角网格分别赋值坐标值，于是得到软组织表面几何模型三角网格上的梯度场，梯度场公式表示为：
如图7所示，公式(22)中bj,bk表示边pik,pij上的点。
软组织表面几何模型的三角网格矢量场的梯度的散度表示为：
其中，三角网格Tk的面积表示为Ak，▽Bik表示为点pi的分段线性函数在一个三角网格Tk上的梯度，▽Bik的值等于三角形Tk中pi对边的高的倒数。
对泊松方程细化得到：
对公式(24)表示的泊松方程进行求解的过程中，Laplacian矩阵只需要计算一次，之后的计算过程中作为一个定值。公式(4-5)等号左边第二项为泊松方程的未知数，也表示将软组织表面几何模型形变后的三维空间坐标。将形变后的三维坐标交给OpenGL图形库渲染就可以在屏幕上得到形变后的软组织表面几何模型。
软组织表面几何模型三角网格形变首先通过对三角网格的梯度场改变，然后能取得散度的改变，最后使用公式(24)得到得到泊松形变后的质点位置向量Jb。泊松网格编辑中最重要的是对梯度场的操作。
3.2、将所述泊松形变后软组织表面几何模型质点位置向量Jb带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量Jc；
4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量；
本发明所提出的一种软组织表面刺穿形变仿真方法及终端通过VS2015和OpenGL图形库实现，操作系统是Windows 7 64位。实验使用的硬件配置是Intel i7-i7-4790@3.60GHz四核处理器，8G内存。表1给出了实验C软组织表面几何模型刺穿形变数据对比，从表中可以看出：在相同的条件下，结合体积保持方法的形变计算量比单独使用泊松形变方法的计算量更大，但实时帧率两者接近，能达到实时性的要求。未使用泊松形变与体积保持方法的传统质点-弹簧方法计算量小，CPU占用率低，但是真实感不强。
表1软组织表面刺穿形变数据对比
图8中(a)表示第0帧，(b)表示第50帧，(c)表示第150帧，(d)表示第305帧，由图8可知，刺穿仿真到第50帧时，刺穿针与软组织表面发生碰撞，这时刺穿针对软组织表面几何模型施加外力使软组织表面开始形变；刺穿针继续对软组织表面施加外力，第150帧时软组织表面形变产生褶皱效果；当仿真到305帧时，软组织表面形变增大，产生更多的褶皱。使用基于泊松方程和体积保持约束方程的软组织表面刺穿形变仿真方法与未使用体积保持的泊松形变方法对软组织表面模型进行刺穿仿真对比，仿真时使用相同的软组织表面几何模型，相同的刺穿力和相同的虚拟弹簧弹性系数。
图9中(a)表示现有的未使用体积保持方程的泊松形变方法，(b)表示本申请基于泊松方程和体积保持约束方程，由图9可知，未使用体积保持的泊松形变方法在软组织表面形变时软组织表面不能产生平滑的过渡，如图9(a)所示。而基于泊松方程和体积保持约束方程的软组织表面刺穿形变仿真方法在软组织表面形变时能在形变区域的边缘产生平滑过渡，如图9(b)所示，在视觉上提升了刺穿仿真的真实感。
实施例三
请参照图2，一种软组织表面刺穿形变仿真终端，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现实施例一中的各个步骤：
综上所述，本发明提供的一种软组织表面刺穿形变仿真方法及终端，通过对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量，保证了在不丢失软组织表面受力形变的细节的同时实现软组织表面刺穿的实时仿真，通过对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程，保证了软组织表面刺穿仿真的实时性，通过将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量，真实感强。
以上所述仅为本发明的实施例，并非因此限制本发明的专利范围，凡是利用本发明说明书及附图内容所作的等同变换，或直接或间接运用在相关的技术领域，均同理包括在本发明的专利保护范围内。
Claims (4)
1.一种软组织表面刺穿形变仿真方法，其特征在于，包括步骤：
S1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；
S2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；
S3、所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；
S4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量；
步骤S1包括：
S11、对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；
S12、根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程；
步骤S2还包括：
获取刺穿针与所述形变区域产生的碰撞信息；
步骤S2具体为：
根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量；
使用拉格朗日乘数将软组织表面几何模型的体积V转化为体积保持约束方程；
其中，
nt表示软组织表面三角网格的网格数，Ai表示软组织表面几何模型第i个三角网格的面积，Nx表示软组织表面几何模型表面三角网格x轴的单位法向量，Ny表示软组织表面几何模型表面三角网格y轴的单位法向量，Nz表示软组织表面几何模型表面三角网格z轴的单位法向量，(x1,y1,z1)、(x2,y2,z2)和(x3,y3,z3)表示三角网格的3个顶点坐标。
2.根据权利要求1所述的软组织表面刺穿形变仿真方法，其特征在于，步骤S3包括：
S31、将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；
S32、将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量。
3.一种软组织表面刺穿形变仿真终端，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述程序时实现以下步骤：
S1、对软组织表面几何模型进行初始化，并建立相应的体积保持约束方程；
S2、对软组织表面几何模型进行数值求解，得到软组织表面几何模型质点位置向量；
S3、所述软组织表面几何模型质点向量通过泊松方程和体积保持约束方程进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量；
S4、根据所述泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量更新所述软组织表面几何模型质点位置向量；
步骤S1包括：
S11、对软组织表面几何模型进行初始化，并采用测地线的方法确定形变区域；
S12、根据散度定理得到软组织表面几何模型的体积，并建立相应的体积保持约束方程；
步骤S2还包括：
获取刺穿针与所述形变区域产生的碰撞信息；
步骤S2具体为：
根据所述碰撞信息采用质点-弹簧模型对软组织表面几何模型进行求解，得到软组织表面几何模型质点位置向量；
使用拉格朗日乘数将软组织表面几何模型的体积V转化为体积保持约束方程；
其中，
nt表示软组织表面三角网格的网格数，Ai表示软组织表面几何模型第i个三角网格的面积，Nx表示软组织表面几何模型表面三角网格x轴的单位法向量，Ny表示软组织表面几何模型表面三角网格y轴的单位法向量，Nz表示软组织表面几何模型表面三角网格z轴的单位法向量，(x1,y1,z1)、(x2,y2,z2)和(x3,y3,z3)表示三角网格的3个顶点坐标。
4.根据权利要求3所述的软组织表面刺穿形变仿真终端，其特征在于，步骤S3包括：
S31、将所述软组织表面几何模型质点位置向量通过泊松方程进行计算，得到泊松形变后软组织表面几何模型质点位置向量；
S32、将所述泊松形变后软组织表面几何模型质点位置向量带入体积保持约束方程中进行计算，得到泊松方程和体积保持约束方程条件下形变后的软组织表面几何模型质点位置向量。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910284423.3A CN109993830B (zh) | 2019-04-10 | 2019-04-10 | 一种软组织表面刺穿形变仿真方法及终端 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910284423.3A CN109993830B (zh) | 2019-04-10 | 2019-04-10 | 一种软组织表面刺穿形变仿真方法及终端 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109993830A CN109993830A (zh) | 2019-07-09 |
| CN109993830B true CN109993830B (zh) | 2022-05-03 |
Family
ID=67132869
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910284423.3A Active CN109993830B (zh) | 2019-04-10 | 2019-04-10 | 一种软组织表面刺穿形变仿真方法及终端 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109993830B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113610875A (zh) * | 2021-09-03 | 2021-11-05 | 成都天地罔极科技有限公司 | 软组织填充效果的仿真方法及系统 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103400023A (zh) * | 2013-06-28 | 2013-11-20 | 华北水利水电大学 | 软组织形变仿真方法 |
| CN106934189A (zh) * | 2015-12-29 | 2017-07-07 | 中国科学院深圳先进技术研究院 | 外科手术软组织变形的仿真方法及装置 |
| CN107330972A (zh) * | 2017-06-28 | 2017-11-07 | 华中科技大学鄂州工业技术研究院 | 模拟生物力学特性的实时软组织形变方法和系统 |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US8140304B2 (en) * | 2007-07-13 | 2012-03-20 | Hyeong-Seok Ko | Method of cloth simulation using linear stretch/shear model |
-
2019
- 2019-04-10 CN CN201910284423.3A patent/CN109993830B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103400023A (zh) * | 2013-06-28 | 2013-11-20 | 华北水利水电大学 | 软组织形变仿真方法 |
| CN106934189A (zh) * | 2015-12-29 | 2017-07-07 | 中国科学院深圳先进技术研究院 | 外科手术软组织变形的仿真方法及装置 |
| CN107330972A (zh) * | 2017-06-28 | 2017-11-07 | 华中科技大学鄂州工业技术研究院 | 模拟生物力学特性的实时软组织形变方法和系统 |
Non-Patent Citations (1)
| Title |
|---|
| "穿刺针尖光纤力传感信号的小波变换分析";杨唐文，陈盼飞，韩建达，徐卫良;《光学 精密工程》;20180831;全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109993830A (zh) | 2019-07-09 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Bender et al. | Position-based Methods for the Simulation of Solid Objects in Computer Graphics. | |
| Weber et al. | Complex barycentric coordinates with applications to planar shape deformation | |
| Sobiecki et al. | Comparison of curve and surface skeletonization methods for voxel shapes | |
| US20180122125A1 (en) | Animating a virtual object in a virtual world | |
| CN108648548A (zh) | 一种脑神经外科虚拟手术训练系统 | |
| Xu et al. | Dynamic harmonic fields for surface processing | |
| Jin et al. | General constrained deformations based on generalized metaballs | |
| González et al. | Computational vademecums for the real-time simulation of haptic collision between nonlinear solids | |
| Cai et al. | Graphical Simulation of Deformable Models | |
| Peng et al. | Higher Dimensional Vector Field Visualization: A Survey. | |
| Li et al. | Fast simulation of deformable characters with articulated skeletons in projective dynamics | |
| CN109993830B (zh) | 一种软组织表面刺穿形变仿真方法及终端 | |
| Celikcan et al. | Example‐Based Retargeting of Human Motion to Arbitrary Mesh Models | |
| CN116580164B (zh) | 一种面向单视角三维人体重建的着装特征学习方法 | |
| Segall et al. | Hele-shaw flow simulation with interactive control using complex barycentric coordinates. | |
| Volino et al. | Stop-and-go cloth draping | |
| Tian et al. | A multi‐GPU finite element computation and hybrid collision handling process framework for brain deformation simulation | |
| Li et al. | Multi-resolution modeling of shapes in contact | |
| Takayama et al. | A sketch-based interface for modeling myocardial fiber orientation that considers the layered structure of the ventricles | |
| Huang et al. | A survey on fast simulation of elastic objects | |
| CN110060779B (zh) | 一种软组织表面刺穿仿真方法及装置 | |
| CN107292942B (zh) | 一种权值c2连续的线性混合形状编辑方法 | |
| Kil et al. | 3D warp brush modeling | |
| Gois et al. | Curvature-driven modeling and rendering of point-based surfaces | |
| US20240221287A1 (en) | Method and apparatus for rendering virtual garment |
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 | ||
| CP03 | Change of name, title or address | ||
| CP03 | Change of name, title or address |
Address after: No.1 campus new village, Rongcheng Town, Fuqing City, Fuzhou City, Fujian Province Patentee after: Fujian Normal University of Technology Country or region after: China Address before: No.1 campus new village, Rongcheng Town, Fuqing City, Fuzhou City, Fujian Province Patentee before: FUQING BRANCH OF FUJIAN NORMAL University Country or region before: China |