CN109036567A - 一种基于子空间凝聚算法的软组织形变仿真方法 - Google Patents
一种基于子空间凝聚算法的软组织形变仿真方法 Download PDFInfo
- Publication number
- CN109036567A CN109036567A CN201810959632.9A CN201810959632A CN109036567A CN 109036567 A CN109036567 A CN 109036567A CN 201810959632 A CN201810959632 A CN 201810959632A CN 109036567 A CN109036567 A CN 109036567A
- Authority
- CN
- China
- Prior art keywords
- subspace
- soft tissue
- unit
- algorithm
- model
- 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.)
- Granted
Links
- 210000004872 soft tissue Anatomy 0.000 title claims abstract description 121
- 238000000034 method Methods 0.000 title claims abstract description 99
- 238000004088 simulation Methods 0.000 title claims abstract description 41
- 238000004422 calculation algorithm Methods 0.000 title claims abstract description 33
- 238000001514 detection method Methods 0.000 claims abstract description 45
- 230000008569 process Effects 0.000 claims abstract description 27
- 238000013334 tissue model Methods 0.000 claims abstract description 26
- 238000004364 calculation method Methods 0.000 claims abstract description 21
- 238000000926 separation method Methods 0.000 claims abstract description 12
- 230000002776 aggregation Effects 0.000 claims description 18
- 238000005054 agglomeration Methods 0.000 claims description 17
- 238000013507 mapping Methods 0.000 claims description 8
- 238000007781 pre-processing Methods 0.000 claims description 8
- 238000009877 rendering Methods 0.000 claims description 8
- 239000000203 mixture Substances 0.000 claims description 5
- 238000001914 filtration Methods 0.000 claims description 4
- 230000000694 effects Effects 0.000 abstract description 6
- 238000004458 analytical method Methods 0.000 abstract description 5
- 230000006872 improvement Effects 0.000 abstract description 3
- 241000208199 Buxus sempervirens Species 0.000 description 16
- 238000010276 construction Methods 0.000 description 8
- 238000009833 condensation Methods 0.000 description 7
- 230000005494 condensation Effects 0.000 description 7
- 238000009826 distribution Methods 0.000 description 5
- 238000001356 surgical procedure Methods 0.000 description 5
- 230000002452 interceptive effect Effects 0.000 description 4
- 210000000056 organ Anatomy 0.000 description 4
- 230000011218 segmentation Effects 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 230000006870 function Effects 0.000 description 3
- 238000005457 optimization Methods 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 239000003814 drug Substances 0.000 description 2
- 238000003780 insertion Methods 0.000 description 2
- 230000037431 insertion Effects 0.000 description 2
- 238000000968 medical method and process Methods 0.000 description 2
- 230000002980 postoperative effect Effects 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 238000011897 real-time detection Methods 0.000 description 2
- 238000012216 screening Methods 0.000 description 2
- 238000005481 NMR spectroscopy Methods 0.000 description 1
- 238000004220 aggregation Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 230000003750 conditioning effect Effects 0.000 description 1
- 238000013016 damping Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 238000009472 formulation Methods 0.000 description 1
- 230000036541 health Effects 0.000 description 1
- 238000003709 image segmentation Methods 0.000 description 1
- 238000007654 immersion Methods 0.000 description 1
- 230000000670 limiting effect Effects 0.000 description 1
- 239000002245 particle Substances 0.000 description 1
- 230000000704 physical effect Effects 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 230000000750 progressive effect Effects 0.000 description 1
- 230000008521 reorganization Effects 0.000 description 1
- 230000035807 sensation Effects 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H50/00—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics
- G16H50/50—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics for simulation or modelling of medical disorders
Landscapes
- Medical Informatics (AREA)
- Engineering & Computer Science (AREA)
- Public Health (AREA)
- Health & Medical Sciences (AREA)
- Pathology (AREA)
- Databases & Information Systems (AREA)
- Data Mining & Analysis (AREA)
- Biomedical Technology (AREA)
- Epidemiology (AREA)
- General Health & Medical Sciences (AREA)
- Primary Health Care (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Processing Or Creating Images (AREA)
Abstract
本发明涉及一种基于子空间凝聚算法的软组织形变仿真方法。该方法包括：a)数据预处理与软组织模型的三维重建；b)进行软组织的生物力学特性分析；c)采用AABB包围盒与三角形面片相交相结合实现碰撞检测；d)最后进行子空间凝聚算法完成软组织形变仿真计算，采用空间核方法实现空间映射，利用贪心算法及二维分离算法分离出最佳纯子空间。本发明结合了软组织的三维几何建模和实时的碰撞检测方法，最后采用改进子空间凝聚算法进行软组织形变仿真，能较好的调节软组织形变过程中真实性与实时性的均衡，防止模型形变过程出现畸变或停顿等现象，通过分析可具有良好的仿真效果且花费时间较短，对虚拟手术中的形变仿真有着重要的指导作用。
Description
技术领域
本发明涉及虚拟手术中软组织形变仿真领域，并结合子空间凝聚算法进行建模仿真，具体涉及一种基于子空间凝聚算法的软组织形变仿真方法；主要针对软组织形变过程中真实性与实时性的均衡问题，该方法在保证调节两者均衡情况下，同时能够实现模拟软组织深度形变的良好效果，保证虚拟手术的精确性。
背景技术
在医学领域，手术作为一种医疗手段，常用来维持患者的健康。随着医疗水平的提高，手术所涉及的范围不断地扩大，已经能够在人体的各个部位进行操作。在进行手术前，医生需根据病患的具体情况制定详细的手术方案，明确手术中的各个环节，预先进行手术的模拟。但由于受到了医生自身临床经验和复杂的人体构造的限制，传统手术模拟并不能完全根据医学影像来模拟真实的手术。因此，如何模拟出真实的手术成为了医学领域急需解决的问题。
虚拟手术是虚拟现实技术的重要应用之一，它是对现实手术中可能存在的各种过程的模拟，包括术前、术中、术后这一整个过程三个阶段。其综合了生物力学、运动力学、计算机图形学、传感技术等方面的知识，集视觉、触觉、听觉等各种器官感受于一身，能满足医生对病患真实情况的具体把握。利用虚拟手术进行手术前的模拟，能为医生呈现出一个具有沉浸感、联想性以及交互性的虚拟手术场景，使得医学上各种手术过程能很好地进行模拟和还原，为医学上的教学指导、术前方案制定、术中手术辅助以及术后身体调理提供技术上强大的支持，具有广泛的应用前景和深远的研究意义。
软组织形变在软组织虚拟手术系统中具有重要的地位。软组织普遍具有不均匀性、非线性和各向异性等一系列复杂的特性，所以软组织形变的研究是整个虚拟手术系统中的一个难点。根据软组织建模原理是否遵照物理原理，软组织形变模型可以分为：非物理模型(Non-physical DeformableModel)与物理模型(PhysicalDeformableModel)。非物理模型是通过几何方法来对形变过程进行仿真，这种模型的仿真速度快，但模型的还原度并不高，在计算机性能较差的年代得到普遍应用。随着不断提高的计算机性能和虚拟手术真实度的要求，由于物理模型还原出的模型更逼真，所以物理模型替代非物理模型得到了广泛的应用。非物理模型中具有代表性的类型主要是自由曲面形变模型和Chain-Mail模型这两类。而在物理模型中具有代表性的类型主要是弹簧质点模型、有限元模型以及质量张量模型。在虚拟手术系统中，当软组织受到外力后，要求虚拟的软组织形变仿真同时具有真实性和实时性的特性。一方面，要保证软组织形变的真实性，就需要保证模型具有生物力学特性，但如果具有生物力学特性将会伴随着计算量剧增的问题出现，同时仿真所需的时间也将大幅度的增加，这种情况下软组织形变的仿真过程将很难达到实时性的要求。另一方面，要保证软组织形变的实时性，就需要将计算量限定在一定范围内，这种情况下要求系统的计算步骤数要有一个限定值，计算量太小往往会伴随着很大的误差，这时仿真出来的软组织形变结果将很难满足真实性的要求。
弹簧质点模型是整个模型的质量离散到一个个的质点，不同的质点之间用弹簧连接来表示物体的弹力、阻尼力等表示物理特性的模型。该模型具有构建简单、方便、快捷等优点，但其仍有其很多的局限性，如：(1)仿真出来的模型具有一定的不真实性；(2)物体的仿真速度也存在较大的滞后性；(3)它的平衡受到破坏后总要再找寻新的平衡状态。
有限元模型是力学模型中对点进行离散取值的结果，是以数值计算堆积起来的建模方法。该模型一般较为准确，且单元的选取较为自由，但其还是有很大的局限性，如：(1)整个模型是相互关联的，任意一个网格的改变都会让受力分布重新改变调整；(2)该模型模拟软组织形变过程的迭代次数多、计算量大、花费时间多、实时性较差。
针对现有技术很难同时满足真实性和实时性的需求，因此，为了构造一个真实性和实时性均衡的软组织形变仿真，本发明提出了一种基于子空间凝聚的软组织形变仿真方法，其中碰撞检测采用AABB包围盒与三角形面片相交相结合法，同时采用空间核方法实现空间映射，利用贪心算法完成最佳子空间的产生并利用二维分离算法分离出最佳纯子空间以此来实现子空间凝聚算法的软组织形变仿真。为当下虚拟手术的应用提供广阔的前景。
发明内容
本发明的目的在于解决传统的软组织形变模型不能是实现真实性和实时性并存的问题，提供一种基于子空间凝聚算法的软组织形变仿真方法，针对软组织形变过程中真实性与实时性的均衡问题，该方法在保证调节两者均衡情况下，同时能够实现模拟软组织深度形变的良好效果，保证虚拟手术的精确性。
为实现上述目的，本发明的技术方案是：一种基于子空间凝聚算法的软组织形变仿真方法，包括如下步骤，
步骤S1、数据预处理与软组织模型的三维重建；
步骤S2、进行软组织的生物力学特性分析；
步骤S3、采用AABB包围盒与三角形面片相交相结合实现碰撞检测；
步骤S4、进行子空间凝聚算法完成软组织形变仿真计算，包括采用空间核方法实现空间映射，利用贪心算法完成最佳子空间的产生并利用二维分离算法分离出最佳纯子空间。
在本发明一实施例中，所述步骤S1的具体实现方式如下，
步骤S11、通过中值滤波对VHP图片中的软组织进行预处理，运用限定阈值的区域生长法提取出软组织的轮廓；
步骤S12、利用分区域绘制三角面片的方法对软组织模型进行三维重建；
步骤S13、利用OpenGL对重建后的软组织模型进行渲染；得出软组织模型重建后的效果图。
在本发明一实施例中，所述步骤S11中，运用限定阈值的区域生长法提取出软组织的轮廓，是利用软组织和周围环境间的灰度值分布的明显区别来实现的，由于软组织的像素值为周围环境最高像素值的8/20，因此，选择的限定阈值为图像中最高像素值7/20，其中限定阈值与相邻像素差之间的关系如式(1)所示：
|f(x,y)-f(x0,y0)|＜T0且f(x,y)＜T (1)
式中，(x0,y0)表示图片中软组织区域的生长点，T0为软组织所在区域的生长阈值，此外，(x,y)表示(x0,y0)邻域中的点，当(x,y)满足以上不等式关系时，(x0,y0)所在的区域M满足M＝M+(x,y)；其中的T满足等式(2)
T＝fmax(x,y)*7/20 (2)。
在本发明一实施例中，所述步骤S3中，采用AABB包围盒碰撞检测方法进行软组织模型内部区域处理，再用三角形面片相交的检测方法处理软组织模型边缘部分，从而使得检测过程中保证实时性的同时提高了碰撞检测的精确度。
在本发明一实施例中，所述三角形面片相交的检测方法的具体实现方式如下，
对于A、B两个碰撞模型中，两个三角形面片之间的相互关系只有两种可能：一种是相离，即这两个三角形间没有发生碰撞；一种是相交，即这两个三角形面片间发生了碰撞；
对于两个三角形面片位置判断公式如下：
其中，设I(x,y)是三角形N1N2N3里面的任何一点，t表示点N1到平面三角形的距离，d0＝M1-N1，d1＝N2-N1，d2＝N3-N1
判断原则如下：
(a)若求出x≥0且y≥0，则两个三角形相离，且距离t相对于向量为正方向；
(b)若求出的x≤0且y≤0，则两个三角形相离，且距离t相对于向量为负方向；
(c)非以上两种情况，则两个三角形相交，此时t的正负方向与x，y的绝对值的大小密切相关。
在本发明一实施例中，所述步骤S4的具体实现方式如下，
采用贪心法，产生候选子空间，再对候选子空间进行筛选；接着通过对先对几何距离进行判断来计算子空间凝聚力大小，根据子空间凝聚力大小为单元选择最佳子空间；而后采用二维空间计算分离最佳纯子空间；最后便可以进行软组织形变仿真；
其中子空间凝聚力计算如下：
子空间的组成与属性A和单元E相关，另外还有一个重要的元素就是凝聚力F；因此，子空间可以表示为S＝{A,E,F}，其中而整个子空间对某个单元的引力可以通过以下计算得到：
首先，可以计算单元M受到单元N的引力，其引力大小可用以下如公式(4)表示：
CM表示该子空间里单元M所包含的属性个数，CN表示该子空间里单元N所包含的属性个数，表示单元M和单元N的几何距离，满足公式(4)，G是万有引力中用到的引力常量；
接着，可以根据凝聚力的概念：一个包含n个单元的子空间，其中的任意一个单元会受到另外的n-1个单元的引力作用，该引力即为凝聚力，其大小可用公式(5)表示：
在本发明一实施例中，所述分离最佳子空间中分离时通过属性和单元两个方向将整个单元的集合分为四个部分，分别是：最佳纯子空间K、具有子空间K的属性而不具有子空间K的单元、具有子空间K的单元而不具有子空间K的属性、既不具有子空间K的属性也不具有它的单元。
相较于现有技术，本发明具有以下有益效果：本发明从虚拟手术中软组织数据图像的几何轮廓进行分割，三维重建，通过数字图像处理，完成碰撞检测，发明了基于子空间凝聚的算法来实现软组织形变仿真的方法，可以使得具有良好仿真效果且花费时间较短，同时能较好的调节软组织形变过程中真实性与实时性的均衡，对虚拟手术中的形变仿真有着重要的指导作用。
附图说明
图1为虚拟手术系统结构图。
图2为包围盒树的构建图。
图3为包围盒树的搜索方式图。
图4为包围盒树碰撞检测的流程图。
图5为两个三角形面片之间的相互关系。
图6为核方法的框架流程图。
图7为子空间凝聚优化方法的流程图。
具体实施方式
下面结合附图，对本发明的技术方案进行具体说明。
本发明的一种基于子空间凝聚算法的软组织形变仿真方法，包括如下步骤，
步骤S1、数据预处理与软组织模型的三维重建；
步骤S2、进行软组织的生物力学特性分析；
步骤S3、采用AABB包围盒与三角形面片相交相结合实现碰撞检测；
步骤S4、进行子空间凝聚算法完成软组织形变仿真计算，包括采用空间核方法实现空间映射，利用贪心算法完成最佳子空间的产生并利用二维分离算法分离出最佳纯子空间。
所述步骤S1的具体实现方式如下，
步骤S11、通过中值滤波对VHP图片中的软组织进行预处理，运用限定阈值的区域生长法提取出软组织的轮廓；
步骤S12、利用分区域绘制三角面片的方法对软组织模型进行三维重建；
步骤S13、利用OpenGL对重建后的软组织模型进行渲染；得出软组织模型重建后的效果图。
所述步骤S11中，运用限定阈值的区域生长法提取出软组织的轮廓，是利用软组织和周围环境间的灰度值分布的明显区别来实现的，由于软组织的像素值为周围环境最高像素值的8/20，因此，选择的限定阈值为图像中最高像素值7/20，其中限定阈值与相邻像素差之间的关系如式(1)所示：
|f(x,y)-f(x0,y0)|＜T0且f(x,y)＜T (1)
式中，(x0,y0)表示图片中软组织区域的生长点，T0为软组织所在区域的生长阈值，此外，(x,y)表示(x0,y0)邻域中的点，当(x,y)满足以上不等式关系时，(x0,y0)所在的区域M满足M＝M+(x,y)；其中的T满足等式(2)
T＝fmax(x,y)*7/20 (2)
所述步骤S3中，采用AABB包围盒碰撞检测方法进行软组织模型内部区域处理，再用三角形面片相交的检测方法处理软组织模型边缘部分，从而使得检测过程中保证实时性的同时提高了碰撞检测的精确度。
所述AABB包围盒碰撞检测方法的具体实现方式如下，
步骤S31、AABB的建立：选择一个包围盒树度数，包围盒树度数和树的深度成反比；在选择完树的度数、深度后，进行包围树的构建(包括自底向上、自顶向下或渐进插入法)；
步骤S32、构建包围盒：第一步，求出包围盒在三个坐标轴上投影的最大最小值，据此求出软组织模型外部的最大包围盒，显示构建软组织模型外部最大的包围盒的情况；第二步，根据均值法，将模型外部的最大包围盒一分为二，在分割面处选择两个子节点，并根据新的子节点重建新的包围盒；第三步，重复第一步和第二步，直到包围盒内元素的大小在误差范围内；
步骤S33、包围盒树的遍历：采用广度优先搜索BFS或深度优先搜索DFS进行包围盒树的遍历，其中，BFS搜索的顺序是按：根节点-所有子节点-所有叶节点来进行的，而DFS的搜索是按：根节点-子节点-叶节点-子节点-叶节点-……；
步骤S34、采用步骤S31至步骤S33对A，B两个碰撞模型分别构建出TreeA和TreeB两个包围盒树，并进行如下碰撞检测过程：
第一步，检测两个模型根节点是否发生碰撞，只有根节点发生碰撞后才能再执行下一步；
第二步，执行包围盒树的遍历，遍历过程从根节点到叶节点；直到所有的叶节点都没有检测到碰撞，结束本次碰撞检测，得出两模型没有发生碰撞的结论。
所述三角形面片相交的检测方法的具体实现方式如下，
对于A、B两个碰撞模型中，两个三角形面片之间的相互关系只有两种可能：一种是相离，即这两个三角形间没有发生碰撞；一种是相交，即这两个三角形面片间发生了碰撞；
对于两个三角形面片位置判断公式如下：
其中，设I(x,y)是三角形N1N2N3里面的任何一点，t表示点N1到平面三角形的距离，d0＝M1-N1，d1＝N2-N1，d2＝N3-N1
判断原则如下：
(a)若求出x≥0且y≥0，则两个三角形相离，且距离t相对于向量为正方向；
(b)若求出的x≤0且y≤0，则两个三角形相离，且距离t相对于向量为负方向；
(c)非以上两种情况，则两个三角形相交，此时t的正负方向与x，y的绝对值的大小密切相关。
所述步骤S4的具体实现方式如下，
采用贪心法，产生候选子空间，再对候选子空间进行筛选；接着通过对先对几何距离进行判断来计算子空间凝聚力大小，根据子空间凝聚力大小为单元选择最佳子空间；而后采用二维空间计算分离最佳纯子空间；最后便可以进行软组织形变仿真；
其中子空间凝聚力计算如下：
子空间的组成与属性A和单元E相关，另外还有一个重要的元素就是凝聚力F；因此，子空间可以表示为S＝{A,E,F}，其中而整个子空间对某个单元的引力可以通过以下计算得到：
首先，可以计算单元M受到单元N的引力，其引力大小可用以下如公式(4)表示：
CM表示该子空间里单元M所包含的属性个数，CN表示该子空间里单元N所包含的属性个数，表示单元M和单元N的几何距离，满足公式(4)，G是万有引力中用到的引力常量；
接着，可以根据凝聚力的概念：一个包含n个单元的子空间，其中的任意一个单元会受到另外的n-1个单元的引力作用，该引力即为凝聚力，其大小可用公式(5)表示：
所述分离最佳子空间中分离时通过属性和单元两个方向将整个单元的集合分为四个部分，分别是：最佳纯子空间K、具有子空间K的属性而不具有子空间K的单元、具有子空间K的单元而不具有子空间K的属性、既不具有子空间K的属性也不具有它的单元。
以下为本发明的具体实现过程。
完整的虚拟系统的结构，主要由以下几个部分组成：预处理模块、场景渲染模块、交互反馈模块，本发明的系统结构流程如图1所示。其中，预处理模块的作用是将研究部分与整个机体完整的分开并构建三维模型。而三维模型的构建是根据VHP或核磁共振医学方式获取某个人体器官的一组平面图片，通过图片中器官的分割与重组而构建某个人体器官的三维模型，将它作为后续步骤的操作对象；场景渲染模块的主要功能是完成操作系统中各个模型的可视化渲染、虚拟手术中场景的绘制以及软件平台的搭建；交互反馈模块是为了实现系统良好的交互性能也能同时提高整个系统的真实感和浸入感。具体实现过程如下：
(1)软组织的三维重建。采用限定阈值的区域生长法对VHP软组织数据图像的几何轮廓进行分割，并用分区域法绘制三角面片，实现二维图片到三维模型的重建，而后在OpenGL中对重建后的软组织面片模型进行渲染。
(2)子空间凝聚的软组织形变。根据形变的需求，对软组织几何模型进行了生物力学分析，使用AABB包围盒法进行实时碰撞检测。采用基于子空间凝聚的形变方法，将三维模型映射到二维平面圆形上，根据区域是否发生形变将其划分为全空间和子空间，分区域并行计算，以确保软组织形变仿真的真实性和实时性。
(3)子空间凝聚的软组织形变的优化。在软组织形变的过程中，针对子空间凝聚法实现过程中碰撞检测、子空间凝聚以及子空间映射存在的不足，分别发明一种优化改进方法：采用AABB包围盒与三角面片相交相结合的方法对碰撞检测进行优化，通过对空间采用核方法实现对空间映射的优化，通过最佳纯子空间的产生、筛选和分离对子空间凝聚进行优化。最后，可以通过实验仿真与分析，验证优化后的子空间凝聚模型具有更高的真实性。
上述步骤详述如下：
第一步，本发明需要对研究的数据进行预处理并完成几何建模。为了使软组织形变的效果更加真实，需建立一个真实度较高的几何模型。首先，通过中值滤波对VHP图片中的软组织进行预处理，本发明运用区域生长与阈值相结合的方法，提取出软组织的轮廓。接着利用分部分绘制三角面片的方法对软组织模型进行三维重建。然后利用OpenGL对重建后的软组织模型进行渲染。最后，得出软组织模型重建后的效果图。
本发明的分割方法为限定阈值的区域生长法，该方法运用在软组织图片分割中，能优化区域生长法的分割结果。运用该方法分割软组织区域的前提是设置的生长阈值要比某个区域内的相邻像素差大，才能保证生长区域内部的均匀性。具体实现过程如下：
由于软组织和周围环境间的灰度值分布有明显的区别，如果添加了限定效果较好的阈值条件，区域生长的范围只会在软组织的内部进行。根据医学上软组织的像素分布表可知，软组织的像素值一般为周围环境最高像素值的8/20，由此，本发明可以选择的阈值为图像中最高像素值的7/20。
限定的阈值与相邻像素差之间的关系如下所示：
|f(x,y)-f(x0,y0)|＜T0且f(x,y)＜T (1)
其中的(x0,y0)表示图片中某个软组织区域的生长点，T0为软组织所在区域的生长阈值。此外，(x,y)表示(x0,y0)邻域中的点，当(x,y)满足以上不等式关系时，(x0,y0)所在的区域M满足M＝M+(x,y)。其中的T满足等式(2)
T＝fmax(x,y)*7/20 (2)
第二步，完成软组织的生物力学特性分析后，进行软组织碰撞检测。实时检测碰撞包围盒检测算法是碰撞检测中常用的方法，其基本思想是用一个形状规则的包围盒近似地代替复杂的研究对象。常见方法为轴向包围盒(axis-aligned bounding box，AABB)和离散有向多面体(k-DOP)等。针对AABB包围盒的碰撞检测方法可能存在的碰撞误差问题，为了降低漏报和误报发生的概率，以减小碰撞检测的误差，在不引入过大计算量的前提下，本发明提出了一种AABB包围盒与三角形面片相交相结合的碰撞检测方法。
1、AABB的碰撞检测过程可分为三个步骤：AABB的建立、构建包围盒、包围盒树的遍历。
(1)AABB的建立：
首先，选择一个合适包围盒树度数，包围盒树度数和树的深度成反比。度数越大，深度就会越小，则树的遍历时间就会越少，但每个节点在访问时需要进行的计算就会增多。在一般情况下，包围盒树度数并没有明确规定，一般可以根据需要来选择。
在选择完树的度数、深度后，进行包围树的方法构建。包围树的构建方法一般有三种：自底向上、自顶向下以及渐进插入法，其构建过程如图2所示。
(2)构建包围盒，其过程又可分三步：第一步，求出包围盒在三个坐标轴上投影的最大最小值，据此求出软组织模型外部的最大包围盒，显示构建软组织模型外部最大的包围盒的情况。第二步，根据均值法，将模型外部的最大包围盒一分为二，在分割面处选择两个子节点，并根据新的子节点重建新的包围盒。第三步，重复第一和第二的步骤，直到包围盒内元素的大小在误差范围内。
(3)包围盒树的遍历。常用的包围盒树的遍历方法一般有两种，分别是：广度优先搜索(BFS)与深度优先搜索(DFS)，它们的搜索过程如图3所示。BFS搜索的顺序是按：根节点-所有子节点-所有叶节点来进行的，而DFS的搜索是按：根节点-子节点-叶节点-子节点-叶节点-……。
利用上述三个步骤对A，B两个碰撞模型分别构建出TreeA和TreeB两个包围盒树，它们间的检测过程可分为以下二步，其流程如图4所示：
第一，检测两个模型根节点是否发生碰撞，只有根节点发生碰撞后才能再执行下一步。
第二，执行包围盒树的遍历，遍历过程从根节点到叶节点。直到所有的叶节点都没有检测到碰撞，结束本次碰撞检测，得出两模型没有发生碰撞的结论。
2、三角网格面片间的碰撞检测方法的检测原理如下所述：
两个模型中，两个三角形网格面片之间的相互关系只有两种可能：一种是相离，即这两个三角形间没有发生碰撞，如图5(a)所示；一种是相交，即这两个三角形面片间发生了碰撞，如图5(b)所示。
用于两个三角形面片位置判断公式如下：
其中，设I(x,y)是三角形N1N2N3里面的任何一点，t表示点N1到平面三角形的距离，d0＝M1-N1，d1＝N2-N1，d2＝N3-N1
判断原则如下：
(a)若求出x≥0且y≥0，则两个三角形相离，且距离t相对于向量为正方向。
(b)若求出的x≤0且y≤0，则两个三角形相离，且距离t相对于向量为负方向。
(c)非以上两种情况，则两个三角形相交，此时t的正负方向与x，y的绝对值的大小密切相关。
三角形面片相交的检测方法能够提高结果的精确度，AABB包围盒碰撞检测方法又能保证检测的实时性，因此，本发明将这两种方法相结合，由AABB包围盒碰撞检测方法进行内部区域处理，再用三角形面片相交的检测方法处理边缘部分，从而使得检测过程中保证实时性的同时提高了碰撞检测的精确度。
第三步，进行空间映射。本发明采用核方法进行空间映射。核方法并非是特定的方法，而是一系列具有特定属性函数的处理方法的总称。数据通过核方法将数据以非线性的方式嵌入到某个具体的空间中去，并在空间中采用线性的运算方法求解非线性函数，其计算结果与非线性运算方法得出的结果是一致的。其实现流程如图6所示。
第四步，基于子空间凝聚力进行软组织形变仿真计算。子空间凝聚算法的流程如图7所示，分别是：产生候选子空间、筛选最佳纯子空间和分离最佳纯子空间。首先，本发明采用贪心法，将分类标准略微降低，并先将满足分类标准的单元都尽量地包括到集合中，产生候选子空间，再对候选子空间进行筛选；接着通过对先对几何距离进行判断来计算凝聚力大小，根据凝聚力大小为单元选择最佳的子空间。而后采用二维空间计算分离最佳纯子空间。最后便可以进行软组织形变仿真。
其中子空间凝聚力计算如下：
子空间的组成与属性(A)和单元(E)相关，另外还有一个重要的元素就是凝聚力(F)。因此，子空间可以表示为S＝{A,E,F}，其中而整个子空间对某个单元的引力可以通过以下计算得到：
首先，可以计算单元M受到单元N的引力，其引力大小可用以下如公式(4)表示：
CM表示该子空间里单元M所包含的属性个数，CN表示该子空间里单元N所包含的属性个数，表示单元M和单元N的几何距离，满足公式(4)，G是万有引力中用到的引力常量。
接着，可以根据凝聚力的概念：一个包含n个单元的子空间，其中的任意一个单元会受到另外的n-1个单元的引力作用，该引力即为凝聚力，其大小可用公式(5)表示：
较佳的，所述分离最佳子空间中分离时通过属性和单元两个方向将整个单元的集合分为四个部分，分别是：最佳纯子空间K、具有子空间K的属性而不具有子空间K的单元(搜索空间1)、具有子空间K的单元而不具有子空间K的属性(搜索空间2)、既不具有子空间K的属性也不具有它的单元。
以上是本发明的较佳实施例，凡依本发明技术方案所作的改变，所产生的功能作用未超出本发明技术方案的范围时，均属于本发明的保护范围。
Claims (7)
1.一种基于子空间凝聚算法的软组织形变仿真方法，其特征在于：包括如下步骤，
步骤S1、数据预处理与软组织模型的三维重建；
步骤S2、进行软组织的生物力学特性分析；
步骤S3、采用AABB包围盒与三角形面片相交相结合实现碰撞检测；
步骤S4、进行子空间凝聚算法完成软组织形变仿真计算，包括采用空间核方法实现空间映射，利用贪心算法完成最佳子空间的产生并利用二维分离算法分离出最佳纯子空间。
2.根据权利要求1所述的一种基于子空间凝聚算法的软组织形变仿真方法，其特征在于：所述步骤S1的具体实现方式如下，
步骤S11、通过中值滤波对VHP图片中的软组织进行预处理，运用限定阈值的区域生长法提取出软组织的轮廓；
步骤S12、利用分区域绘制三角面片的方法对软组织模型进行三维重建；
步骤S13、利用OpenGL对重建后的软组织模型进行渲染；得出软组织模型重建后的效果图。
3.根据权利要求2所述的一种基于子空间凝聚算法的软组织形变仿真方法，其特征在于：所述步骤S11中，运用限定阈值的区域生长法提取出软组织的轮廓，是利用软组织和周围环境间的灰度值分布的明显区别来实现的，由于软组织的像素值为周围环境最高像素值的8/20，因此，选择的限定阈值为图像中最高像素值7/20，其中限定阈值与相邻像素差之间的关系如式(1)所示：
|f(x,y)-f(x0,y0)|＜T0且f(x,y)＜T (1)
式中，(x0,y0)表示图片中软组织区域的生长点，T0为软组织所在区域的生长阈值，此外，(x,y)表示(x0,y0)邻域中的点，当(x,y)满足以上不等式关系时，(x0,y0)所在的区域M满足M＝M+(x,y)；其中的T满足等式(2)
T＝fmax(x,y)*7/20 (2)。
4.根据权利要求2所述的一种基于子空间凝聚算法的软组织形变仿真方法，其特征在于：所述步骤S3中，采用AABB包围盒碰撞检测方法进行软组织模型内部区域处理，再用三角形面片相交的检测方法处理软组织模型边缘部分，从而使得检测过程中保证实时性的同时提高了碰撞检测的精确度。
5.根据权利要求4所述的一种基于子空间凝聚算法的软组织形变仿真方法，其特征在于：所述三角形面片相交的检测方法的具体实现方式如下，
对于A、B两个碰撞模型中，两个三角形面片之间的相互关系只有两种可能：一种是相离，即这两个三角形间没有发生碰撞；一种是相交，即这两个三角形面片间发生了碰撞；
对于两个三角形面片位置判断公式如下：
其中，设I(x,y)是三角形N1N2N3里面的任何一点，t表示点N1到平面三角形的距离，d0＝M1-N1，d1＝N2-N1，d2＝N3-N1
判断原则如下：
(a)若求出x≥0且y≥0，则两个三角形相离，且距离t相对于向量为正方向；
(b)若求出的x≤0且y≤0，则两个三角形相离，且距离t相对于向量为负方向；
(c)非以上两种情况，则两个三角形相交，此时t的正负方向与x，y的绝对值的大小密切相关。
6.根据权利要求1所述的一种基于子空间凝聚算法的软组织形变仿真方法，其特征在于：所述步骤S4的具体实现方式如下，
采用贪心法，产生候选子空间，再对候选子空间进行筛选；接着通过对先对几何距离进行判断来计算子空间凝聚力大小，根据子空间凝聚力大小为单元选择最佳子空间；而后采用二维空间计算分离最佳纯子空间；最后便可以进行软组织形变仿真；
其中子空间凝聚力计算如下：
子空间的组成与属性A和单元E相关，另外还有一个重要的元素就是凝聚力F；因此，子空间可以表示为S＝{A,E,F}，其中而整个子空间对某个单元的引力可以通过以下计算得到：
首先，可以计算单元M受到单元N的引力，其引力大小可用以下如公式(4)表示：
CM表示该子空间里单元M所包含的属性个数，CN表示该子空间里单元N所包含的属性个数，表示单元M和单元N的几何距离，满足公式(4)，G是万有引力中用到的引力常量；
接着，可以根据凝聚力的概念：一个包含n个单元的子空间，其中的任意一个单元会受到另外的n-1个单元的引力作用，该引力即为凝聚力，其大小可用公式(5)表示：
7.根据权利要求6所述的一种基于子空间凝聚算法的软组织形变仿真方法，其特征在于：所述分离最佳子空间中分离时通过属性和单元两个方向将整个单元的集合分为四个部分，分别是：最佳纯子空间K、具有子空间K的属性而不具有子空间K的单元、具有子空间K的单元而不具有子空间K的属性、既不具有子空间K的属性也不具有它的单元。
Applications Claiming Priority (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810021793.3A CN108109699A (zh) | 2018-01-10 | 2018-01-10 | 一种基于子空间凝聚算法的软组织形变仿真方法 |
| CN2018100217933 | 2018-01-10 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109036567A true CN109036567A (zh) | 2018-12-18 |
| CN109036567B CN109036567B (zh) | 2021-11-09 |
Family
ID=62219836
Family Applications (2)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810021793.3A Withdrawn CN108109699A (zh) | 2018-01-10 | 2018-01-10 | 一种基于子空间凝聚算法的软组织形变仿真方法 |
| CN201810959632.9A Active CN109036567B (zh) | 2018-01-10 | 2018-08-22 | 一种基于子空间凝聚算法的软组织形变仿真方法 |
Family Applications Before (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810021793.3A Withdrawn CN108109699A (zh) | 2018-01-10 | 2018-01-10 | 一种基于子空间凝聚算法的软组织形变仿真方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (2) | CN108109699A (zh) |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102044086A (zh) * | 2010-11-30 | 2011-05-04 | 华北水利水电学院 | 一种软组织形变仿真方法 |
| CN103400023A (zh) * | 2013-06-28 | 2013-11-20 | 华北水利水电大学 | 软组织形变仿真方法 |
| CN105261069A (zh) * | 2015-11-05 | 2016-01-20 | 北京航空航天大学 | 基于gpu的软组织器官元球模型自动生成与碰撞检测的方法 |
| CN105302972A (zh) * | 2015-11-05 | 2016-02-03 | 北京航空航天大学 | 一种基于元球模型的软组织形变方法 |
| CN105513130A (zh) * | 2016-02-01 | 2016-04-20 | 福建师范大学福清分校 | 一种基于网格与无网格混合的软组织形变方法 |
| CN106504329A (zh) * | 2016-09-27 | 2017-03-15 | 西安科技大学 | 一种基于牙体长轴的质点弹簧模型的牙龈变形仿真方法 |
-
2018
- 2018-01-10 CN CN201810021793.3A patent/CN108109699A/zh not_active Withdrawn
- 2018-08-22 CN CN201810959632.9A patent/CN109036567B/zh active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102044086A (zh) * | 2010-11-30 | 2011-05-04 | 华北水利水电学院 | 一种软组织形变仿真方法 |
| CN103400023A (zh) * | 2013-06-28 | 2013-11-20 | 华北水利水电大学 | 软组织形变仿真方法 |
| CN105261069A (zh) * | 2015-11-05 | 2016-01-20 | 北京航空航天大学 | 基于gpu的软组织器官元球模型自动生成与碰撞检测的方法 |
| CN105302972A (zh) * | 2015-11-05 | 2016-02-03 | 北京航空航天大学 | 一种基于元球模型的软组织形变方法 |
| CN105513130A (zh) * | 2016-02-01 | 2016-04-20 | 福建师范大学福清分校 | 一种基于网格与无网格混合的软组织形变方法 |
| CN106504329A (zh) * | 2016-09-27 | 2017-03-15 | 西安科技大学 | 一种基于牙体长轴的质点弹簧模型的牙龈变形仿真方法 |
Non-Patent Citations (6)
| Title |
|---|
| 张辉荣 等: "面向分类数据的重叠子空间聚类算法SCCAT", 《西南大学学报（自然科学版）》 * |
| 徐海荣 等: "改进的区域生长算法在医学图像分割中的应用", 《生物医学工程研究》 * |
| 董红建 等: "生物软组织超弹性及夹持过程仿真分析", 《医用生物力学》 * |
| 陈寒青 等: "基于子空间凝聚的软组织形变仿真方法研究", 《佳木斯大学学报(自然科学版)》 * |
| 陈寒青: "基于子空间凝聚的软组织形变仿真方法研究", 《中国优秀硕士学位论文全文数据库 基础科学辑》 * |
| 陈晨 等: "基于包围盒和三角面片的碰撞检测优化算法", 《制造技术与机床》 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN108109699A (zh) | 2018-06-01 |
| CN109036567B (zh) | 2021-11-09 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106875462B (zh) | 一种基于元球模型和混合驱动方法的实时数字器官切割方法 | |
| CN108648548A (zh) | 一种脑神经外科虚拟手术训练系统 | |
| CN104821006B (zh) | 一种基于人体混合包围盒的动态服装仿真方法 | |
| Westwood | A GPU accelerated spring mass system for surgical simulation | |
| CN112102480B (zh) | 图像数据处理方法、装置、设备以及介质 | |
| Wu et al. | Interactive High-Resolution Boundary Surfaces for Deformable Bodies with Changing Topology. | |
| Fortmeier et al. | Image-based palpation simulation with soft tissue deformations using chainmail on the GPU | |
| Wang et al. | A review of collision detection for deformable objects | |
| Mosegaard et al. | GPU accelerated surgical simulators for complex morphology | |
| CN103699716B (zh) | 一种个性化三维医学图像驱动的器官虚拟显示方法 | |
| CN112785689A (zh) | 一种三维模型构建方法、装置、系统及存储介质 | |
| CN102682473B (zh) | 一种虚拟服装的实时物理建模方法 | |
| Peng et al. | Bi-potential and co-rotational formulations applied for real time simulation involving friction and large deformation | |
| Aras et al. | An analytic meshless enrichment function for handling discontinuities in interactive surgical simulation | |
| Wang et al. | A novel virtual cutting method for deformable objects using high‐order elements combined with mesh optimisation | |
| CN109036567A (zh) | 一种基于子空间凝聚算法的软组织形变仿真方法 | |
| CN112381825B (zh) | 用于病灶区图像几何特征提取的方法和相关产品 | |
| Tian et al. | A multi‐GPU finite element computation and hybrid collision handling process framework for brain deformation simulation | |
| Chen et al. | Hybrid optimization-based cutting simulation for soft objects | |
| CN118230967A (zh) | 一种带有力反馈设备的柔体组织仿真方法 | |
| Ye et al. | An improved collision detection and cutting algorithm of the soft tissue in virtual surgical simulation | |
| Yuan et al. | Real-time simulation of tissue cutting with CUDA based on GPGPU | |
| CN108511074B (zh) | 一种基于空间核映射和子空间聚集的软组织形变方法 | |
| JP3178808B2 (ja) | オブジェクト再構成方法及びオブジェクト近似方法及び空間描画方法 | |
| Hieber | Particle methods for flow-structure interactions |
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 |