CN101436303B - 从物体三维图像中获取四面体网格的方法 - Google Patents

从物体三维图像中获取四面体网格的方法 Download PDF

Info

Publication number
CN101436303B
CN101436303B CN2008102077799A CN200810207779A CN101436303B CN 101436303 B CN101436303 B CN 101436303B CN 2008102077799 A CN2008102077799 A CN 2008102077799A CN 200810207779 A CN200810207779 A CN 200810207779A CN 101436303 B CN101436303 B CN 101436303B
Authority
CN
China
Prior art keywords
image
value
sampled point
sampling
voronoi
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.)
Expired - Fee Related
Application number
CN2008102077799A
Other languages
English (en)
Other versions
CN101436303A (zh
Inventor
顾力栩
吕思哲
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shanghai Jiaotong University
Original Assignee
Shanghai Jiaotong University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Shanghai Jiaotong University filed Critical Shanghai Jiaotong University
Priority to CN2008102077799A priority Critical patent/CN101436303B/zh
Publication of CN101436303A publication Critical patent/CN101436303A/zh
Application granted granted Critical
Publication of CN101436303B publication Critical patent/CN101436303B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Image Processing (AREA)

Abstract

本发明涉及一种图像处理技术领域的从物体三维图像中获取四面体网格的方法。本发明直接在运行的计算机上进行,首先读入物体三维图像,然后依次通过图像前处理模块、采样模块以及网格重建模块处理,得到四面体网格。其中图像前处理模块对输入的二值图像反复进行二值腐蚀操作获得图像的距离变换结果,在距离变换的结果上计算出密度函数,并将其输出到采样模块。采样模块根据密度函数对图像进行采样,得到所需要的空间采样点,并将其输出到网格重建模块,网格重建模块由上述采样点构建最终的Delaunay网格。本发明易于实现,在网格质量上相比于现有方法,拥有良好的表现。

Description

从物体三维图像中获取四面体网格的方法
技术领域
本发明涉及的是一种图像处理技术领域的方法,具体是一种从物体三维图像中获取四面体网格的方法。
背景技术
随着科学技术的高速发展,高科技医疗设备的不断涌现为医疗的现代化提供了越来越多的帮助。在众多的计算机辅助医疗系统中,获取建立系统分析数据所需要的四面体网格,是系统必须具备的功能模块。该模块需要运用计算机图形学、计算几何、数字图像处理相关专业知识,将以像素为基本表示单位的图像,转换生成为以四面体为基本单位的网格数据。在计算机辅助医疗系统中,需要将真实病人的图像数据(如CT、MRI等),转换为计算机图形学范畴内的四面体网格,再进行后继的分析处理。生成的网格的质量将严重影响后继处理的效果及正确性。因此,国际上众多研究机构在该技术上倾注了大量的资源和精力。
当前医学数据四面体网格重建,一般分为两个阶段。第一阶段为数据的表面网格重建,即抽取出图像的边缘信息,而后通过一系列处理,将该边缘信息转换为三角形或其他形式的表面网格。代表方法有Marching Cube、α-Shape、CRUST等。第二阶段为四面体网格重建阶段,即在前一阶段的基础上,用四面体网格填充表面网格,得到最终的四面体网格。代表方法很多,一般说来,可以分为基于栅格的方法、基于Delaunay网格的方法、前进波方法等。一般的医学仿真软件,通过以上两个阶段的处理,得到最终所需要的四面体网格。
然而,上述处理方法存在以下几个问题:首先,两个阶段的衔接存在问题。一般情况下,第二阶段对第一阶段的处理结果有一定的要求,如必须生成封闭的表面网格数据、质量要达到某些指标等等。在实际应用场景中,这些要求第一阶段的技术未必能够百分之一百的达到,可能导致最终无法正确地生成网格或者生成低质量的网格。其次,绝大多数的处理技术,输出的网格的正确性及质量,严重依赖于用户手工设定的参数值。这样一来,在实际应用中,用户往往需要凭经验、多次试验才能得到所需要的一定质量的网格,带来了相当的麻烦。再者,上述方法,基本上来源于计算机辅助设计制造领域,并未针对医学图像数据的特点,进行相应的优化和处理,导致难以生成高质量的网格。最后,由于实现上述技术的复杂性较高,在应用到实际的工程领域时,亦会存在不小的困难。
经对现有技术的文献检索,N.Amenta等在Proceedings of the 25th annualconference on Computer graphics and interactive techniques(98年出版的第25届计算机图形与交互技术年会论文集),1998.pages 415-421,上发表的“A new Voronoi-based surface reconstruction algorithm”(《一种基于Voronoi图的表面重建方法》)中提出了一种利用插入某些特殊的空间点,通过两次构建Voronoi图的方式建立表面网格数据。该方法存在的问题是,最终输出的表面三角形网格,不能保证为流型网格,且可能出现空洞,不能保证为封闭网格。这对于后继的造型技术,如四面体划分等,存在很大的障碍。再如,Qiang Du等人在International Journal for Numerical Methods in Engineering杂志Volume 56 Issue 9,Pages 1355-1373中发表的“Tetrahedral mesh generationand optimization based on Centroidal Voronoi tessellations”(基于质心Voronoi图的四面体网格划分及优化)。该文主要针对从表面网格到四面体网格的划分与网格质量优化,其主要问题在于该方法单纯的追求四面体的质量,虽然对四面体的质量进行了局部优化,但是,在全局上并未能保证四面体的采样率与物体几何形态相关。此外,该文仅给出了简单数据的测试结果,并未给出对于复杂形态物的建模效果,其在复杂场景下的实用性值得怀疑。
发明内容
本发明的目的在于克服现有方法中的不足,提供一种从物体三维图像中获取四面体网格的方法,直接从图像数据中进行网格重建,同时,针对医学图像数据的特点,考虑复杂图像的优化,生成Delaunay网格,方便后继处理,保证网格质量。
本发明是通过以下技术方案实现的。本发明直接在运行的计算机上进行,首先读入物体三维图像,然后依次通过图像前处理模块、采样模块以及网格重建模块处理,得到物体三维图像四面体网格。其中图像前处理模块对输入的二值图像反复进行二值腐蚀操作获得图像的距离变换结果,在距离变换的结果上计算出密度函数,并将其输出到采样模块。采样模块根据密度函数对图像进行采样,得到所需要的空间采样点,并将其输出到网格重建模块,网格重建模块由上述采样点构建最终的Delaunay网格。该方案稳定可靠,自动化程度高,且实现方便,非常容易在实际的项目中加以应用。
其中:
所述的图像距离变换,指的是求取图像每个像素点到物体边界的最短距离。该距离可以直观的理解为,以相应像素位置为圆心,能够在物体内部放置的最大的球体的半径的大小。距离变换所得到的数值,反映了在不同位置,物体的几何形态特征。半径越大,说明物体在该像素位置的几何形态变化越平缓,反之,说明该像素位置的物体几何形态变化越剧烈。
所述的质心Voronoi图,是一种特殊的Voronoi图。输入一组采样点,每个采样点可以得到一个对应的Voronoi单元,单元内的所有空间点,到该采样点的距离,小于到其他采样点的距离。所有的Voronoi单元,形成了对空间的一个划分,该划分被成为对空间的Voronoi图。质心Voronoi图是Voronoi图的一个特例,即,给定空间内的密度函数,用于划分Voronoi图的所有采样点,恰好位于其对应的Voronoi单元的质心位置。
所述的Delaunay网格,是一种特殊的网格结构。三维Delaunay网格由四面体组成。每个四面体的外接球内,不包含除该四面体的四个顶点之外的其他采样点。在数学上,Delaunay网格与Voronoi图呈对偶关系。用于生成Voronoi图的采样点,即为Delaunay网格中四面体的顶点。当两个Voronoi单元相接壤,则连接该区域对应的采样点,将得到对应的Delaunay网格内的一条边。
本发明包括如下步骤:
第一步,将物体的三维图像数据,以二值的形式读入系统;
第二步,图像前处理模块接收输入的二值图像,并对其反复进行二值腐蚀操作,将第k次图像二值腐蚀中被腐蚀的像素置为k,直至所有有值的像素点都被赋值为止,得到图像距离变换的结果,k为自然数;
第三步,图像前处理模块在获得图像距离变换的结果基础上,遍历所有像素点,获得每个像素点的密度值,并将其输出到采样模块;
第四步,采样模块接收图像前处理模块输出的像素点的密度值,在图像内随机选取初始采样点,然后进行迭代操作;
第五步,采样模块每次迭代获得当前采样点集的Voronoi图,从而得到图像中每个Voronoi单元的质心及质量,将获得的Voronoi单元中对应的采样点位移到该质心位置;
第六步,检测Voronoi单元的质量,根据其质量控制在图像不同区域的采样率,如果质量小于用户设定的参数值,则删除该采样点;如果质量大于用户设定的参数值,则在该采样点附近随机插入新的样点,每次迭代得到的采样点趋近于Voronoi单元的质心位置,采样点的变化将会逐步减小,直至最终几乎没有发生位移。如果插入或者删除的采样点数量小于给定的某个值,且采样点的位移和小于指定值,则采样模块结束处理,将处理的结果输出到网格重建模块;
第七步,网格重建模块接收采样模块输出的采样点,由边界采样点构造对应的Voronoi图,以此得到所有的Voronoi单元,并找到每个Voronoi单元的两个极点,将找到的所有的位于物体外部的极点,记为集合P,然后根据所有的采样点以及极点构建Delaunay网格;
第八步,检查构造出来的Delaunay网格内的所有四面体,如果四面体至少有一个顶点在集合P中,则该四面体在物体的外部,否则就在物体内部,删除外部四面体,保留内部四面体,得到最终需要的四面体网格数据并输出。
所述极点,即该边界采样点所在的Voronoi单元中(呈细长状),在长轴方向上,距离最远的两个顶点。
所述每个像素点的密度值通过物体密度函数得到,当当前像素点位于物体外部时,密度函数值为0。否则,根据下述公式进行计算:
&rho; ( x , y , z ) = tan [ &alpha; ( 1.0 - d ( x , y , z ) - 1 max { d ( x , y , z ) - 1 } ) ] tan &alpha; , ( 1.0 &le; &alpha; < &Pi; 2 )
其中,ρ(x,y,z)为密度函数,d(x,y,z)为前述的距离变换函数,α为用户输入的参数。
本发明利用本发明进行物体三维图像数据的四面体网格重建,具备如下优点:(1)在工程应用中便于操作和实现。在实际应用中,可以在非常短的时间内对该技术进行实现并集成至系统中。(2)网格质量良好,非常适合于进行几何造型及其后继处理。该方法生成的网格,是标准的Delaunay网格,具备Delaunay网格应有的所有的数学特性。且,网格的分布合理均匀,在物体中心附近,网格体积较大,越接近边界,网格体积越小。(3)一步构建网格,省去了用户需要多步进行构建、反复调节参数的麻烦,避免了原有表面重建、四面体网格划分两步走的数据无法交接的问题。(4)该方法针对医学图像数据的特点,非常适合于复杂形态的物体的网格获取,尤其是人体各个组织器官。
附图说明
图1本发明方法流程图
具体实施方式
下面结合附图对本发明的实施例作详细说明:本实施例在以本发明技术方案为前提下进行实施,给出了详细的实施方式和具体的操作过程,但本发明的保护范围不限于下述的实施例。
本实施例在CPU为Pentuim M 1.8GHz,显卡为NVIDIA Geforce 9800GTS,内存为2.0GB的计算机中实现,编程语言为C++。
本实施例的实施流程如图1所示。
测试用两个图像具体参数情况如表1所示:
表1实施例中各种参数情况
情况  图像大小  四面体个数  顶点数
图像A  168x142x221  137025  25905
图像B  168x142x221  384348  69050
第一步,将图像A和图像B的上述三维图像数据,以二值的形式读入系统;
第二步,图像前处理模块接收输入的二值图像,并对其反复进行二值腐蚀操作,将第k次图像二值腐蚀中被腐蚀的像素置为k,直至所有有值的像素点都被赋值为止,得到图像距离变换的结果,k为自然数;
第三步,图像前处理模块在获得图像距离变换的结果基础上,遍历所有像素点,获得每个像素点的密度值,并将其输出到采样模块;
第四步,采样模块接收图像前处理模块输出的像素点的密度值,在图像内随机选取初始采样点,然后进行迭代操作;
第五步,采样模块每次迭代获得当前采样点集的Voronoi图,从而得到图像中每个Voronoi单元的质心及质量,将获得的Voronoi单元中对应的采样点位移到该质心位置;
第六步,检测Voronoi单元的质量,根据其质量控制在图像不同区域的采样率,如果质量小于用户设定的参数值,则删除该采样点;如果质量大于用户设定的参数值,则在该采样点附近随机插入新的样点,每次迭代得到的采样点趋近于Voronoi单元的质心位置,采样点的变化将会逐步减小,直至最终几乎没有发生位移。如果插入或者删除的采样点数量小于给定的某个值,且采样点的位移和小于指定值,则采样模块结束处理,将处理的结果输出到网格重建模块;
第七步,网格重建模块接收采样模块输出的采样点,由边界采样点构造对应的Voronoi图,以此得到所有的Voronoi单元,并找到每个Voronoi单元的两个极点,将找到的所有的物体位于外部的极点,记为集合P,然后根据所有的采样点以及极点构建Delaunay网格;
第八步,检查构造出来的Delaunay网格内的所有四面体,如果四面体至少有一个顶点在集合P中,则该四面体在物体的外部,否则就在物体内部,删除外部四面体,保留内部四面体,得到最终需要的四面体网格数据并输出。
所述极点,即该边界采样点所在的Voronoi单元中(呈细长状),在长轴方向上,距离最远的两个顶点。
所述每个像素点的密度值通过物体密度函数得到,当当前像素点位于物体外部时,密度函数值为0。否则,根据下述公式进行计算:
&rho; ( x , y , z ) = tan [ &alpha; ( 1.0 - d ( x , y , z ) - 1 max { d ( x , y , z ) - 1 } ) ] tan &alpha; , ( 1.0 &le; &alpha; < &Pi; 2 )
其中,ρ(x,y,z)为密度函数,d(x,y,z)为前述的距离变换函数,α为用户输入的参数。
图像A的实例与已有方法的比较。方法一为栅格化网格方法,即先将物体所在包围盒划分为相同大小的立方体,而后再将每个小立方体划分为四面体。方法二使用原有的两步走的策略。使用TightCocone方法进行表面重建后,再使用NETGEN方法进行四面体网格划分。
表2中,计算、统计了以下几个质量衡量指标。这些指标,数值越高,表明网格质量越高。
参数1:
Figure G2008102077799D00071
衡量基于该网格进行线性插值后的函数,与真实函数的误差限大小。其中,V为四面体的体积,rmc为包围该四面体的最小包围球的半径大小。
参数2:
Figure G2008102077799D00072
衡量及于该网格进行线性插值后的函数的一阶导,与真实函数的一阶导之间的误差限。其中Am为四面体的三角形面的面积,lij为四面体的边长。
参数3:
Figure G2008102077799D00073
衡量该网格对有限元中刚度矩阵的条件数的影响。
表2图像A实例实施结果及与传统方法的比较情况
Figure G2008102077799D00081
从表2中可以知道,不论从网格的何种指标上来看,本实施例所得的网格,都具备更加合理优异的全局分布。即,高质量的网格所占的比例,远远高于其他方法。如,参数3中,本实施例高质量网格的比例为37.033%,大于其他两方法的23.732%和14.230%。这表明,用此方法得到的网格,在后继的其他模块(如数值计算)的相关任务中,将比其他方法更具优势。

Claims (2)

1.一种从物体三维图像中获取四面体网格的方法,其特征在于包括如下步骤:
第一步,将物体的三维图像数据,以二值的形式读入系统;
第二步,图像前处理模块接收输入的二值图像,并对其反复进行二值腐蚀操作,将第k次图像二值腐蚀中被腐蚀的像素置为k,直至所有有值的像素点都被赋值为止,得到图像距离变换的结果,k为自然数;
第三步,图像前处理模块在获得图像距离变换的结果基础上,遍历所有像素点,获得每个像素点的密度值,并将其输出到采样模块;
所述每个像素点的密度值通过物体密度函数得到,当当前像素点位于物体外部时,密度函数值为0,否则,根据下述公式进行计算:
&rho; ( x , y , z ) = tan [ &alpha; ( 1.0 - d ( x , y , z ) - 1 max { d ( x , y , z ) - 1 ) ] tan &alpha; ( 1.0 &le; &alpha; < &Pi; 2 )
其中,ρ(x,y,z)为密度函数,d(x,y,z)为前述的距离变换函数,α为用户输入的参数;
第四步,采样模块接收图像前处理模块输出的像素点的密度值,在图像内随机选取初始采样点,然后进行迭代操作;
第五步,采样模块每次迭代获得当前采样点集的Voronoi图,从而得到图像中每个Voronoi单元的质心及质量,将获得的Voronoi单元中对应的采样点位移到该质心位置;
第六步,检测Voronoi单元的质量,根据其质量控制在图像不同区域的采样率,如果质量小于用户设定的参数值,则删除该采样点,如果质量大于用户设定的参数值,则在该采样点附近随机插入新的样点,每次迭代得到的采样点趋近于Voronoi单元的质心位置,采样点的变化将会逐步减小,直至最终几乎没有发生位移;如果插入或者删除的采样点数量小于给定的某个值,且采样点的位移和小于指定值,则采样模块结束处理,将处理的结果输出到网格重建模块;
第七步,网格重建模块接收采样模块输出的采样点,由边界采样点构造对应的Voronoi图,以此得到所有的Voronoi单元,并找到每个Voronoi单元的两个极点,将找到的所有的位于物体外部的极点,记为集合P,然后根据所有的采样点以及极点构建Delaunay网格;
第八步,检查构造出来的Delaunay网格内的所有四面体,如果四面体至少有一个顶点在集合P中,则该四面体在物体的外部,否则就在物体内部,删除外部四面体,保留内部四面体,得到最终需要的四面体网格数据并输出。
2.根据权利要求1所述的从物体三维图像中获取四面体网格的方法,其特征是,第七步中,所述极点,是指边界采样点所在的Voronoi单元中,在长轴方向上,距离最远的两个顶点。
CN2008102077799A 2008-12-25 2008-12-25 从物体三维图像中获取四面体网格的方法 Expired - Fee Related CN101436303B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN2008102077799A CN101436303B (zh) 2008-12-25 2008-12-25 从物体三维图像中获取四面体网格的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN2008102077799A CN101436303B (zh) 2008-12-25 2008-12-25 从物体三维图像中获取四面体网格的方法

Publications (2)

Publication Number Publication Date
CN101436303A CN101436303A (zh) 2009-05-20
CN101436303B true CN101436303B (zh) 2011-07-20

Family

ID=40710733

Family Applications (1)

Application Number Title Priority Date Filing Date
CN2008102077799A Expired - Fee Related CN101436303B (zh) 2008-12-25 2008-12-25 从物体三维图像中获取四面体网格的方法

Country Status (1)

Country Link
CN (1) CN101436303B (zh)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101706966B (zh) * 2009-11-06 2011-11-16 上海第二工业大学 基于二维图像和多点统计方法的多孔介质三维重构方法
CN104036552B (zh) * 2014-06-23 2017-05-10 中国科学院自动化研究所 基于最远点优化的蓝噪声网格生成方法
CN104063903B (zh) * 2014-07-08 2016-09-14 清华大学 三维实体模型的四面体网格生成方法和装置
CN104504757A (zh) * 2014-12-15 2015-04-08 北京工业大学 一种基于限制点集方法的血管壁网格划分方法
US10783708B2 (en) 2017-03-14 2020-09-22 Siemens Industry Software Inc. Systems and methods for determining mass properties of a modeled object
CN107967397B (zh) * 2017-12-12 2020-12-29 北京航空航天大学 一种基于有限元分析的飞行器结构质心漂移量高精度设计方法
CN108510505B (zh) * 2018-03-30 2022-04-01 南京工业大学 一种基于双边格的高分辨率图像的图切分图像分割方法

Also Published As

Publication number Publication date
CN101436303A (zh) 2009-05-20

Similar Documents

Publication Publication Date Title
CN101436303B (zh) 从物体三维图像中获取四面体网格的方法
CN110288695B (zh) 基于深度学习的单帧图像三维模型表面重建方法
Ye et al. Reverse innovative design—an integrated product design methodology
KR100900824B1 (ko) 스케치 기반 3차원 모델 생성 장치 및 방법
Al Akhras et al. Isogeometric analysis-suitable trivariate NURBS models from standard B-Rep models
US20090295803A1 (en) Image processing method
Maadasamy et al. A hybrid parallel algorithm for computing and tracking level set topology
CN106548484A (zh) 基于二维凸包的产品模型散乱点云边界特征提取方法
CN102999937A (zh) 心脏散乱点云数据曲面重建的方法
Duczek et al. The finite cell method for tetrahedral meshes
Rossignac et al. HelSweeper: Screw-sweeps of canal surfaces
Huang et al. Automatic CAD model reconstruction from multiple point clouds for reverse engineering
Cao et al. Computation of the medial axis of planar domains based on saddle point programming
Fu et al. An algorithm for finding intersection between ball B-spline curves
Jahanshahloo et al. Reconstruction of 3D shapes with B-spline surface using diagonal approximation BFGS methods
Raseli et al. The Construction of Cubic Bezier Curve
Mascarenhas et al. Isocontour based visualization of time-varying scalar fields
Lizier et al. Comparing techniques for tetrahedral mesh generation
Batuhan Arisoy et al. Free form surface skinning of 3d curve clouds for conceptual shape design
CN112614046B (zh) 一种二维平面上绘制三维模型的方法及装置
Liu et al. Reconstruction of surgical instruments in virtual surgery system
Ma et al. 3D Facial Similarity Measure Based on Deformation Field
Liu et al. End-to-End Surface Reconstruction for Touching Trajectories
Buczek Area collapse algorithm computing new curve of 2D geometric objects
CN116486035A (zh) 一种基于特征的cad体网格模型编辑优化系统

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
C17 Cessation of patent right
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20110720

Termination date: 20131225