【word】 面向医学假体CAD的Loop细分曲面拟合系统

【word】 面向医学假体CAD的Loop细分曲面拟合系统 面向医学假体CAD的Loop细分曲面拟合 系统 中国机械工程第18卷第17期2007年9月上半月 面向医学假体CAD的Loop细分曲面拟合系统 刘泗岩廖文和刘浩 南京航空航天大学,南京,210016 摘要:为实现基于Loop细分曲面表示的医学假体CAD造型,提出了一个完整的Loop细分曲面拟 合系统.输入为CT图像经三维重建并经滤噪与平滑处理所得的密集三角网格,再用QEM简化算法 与提出的网格优化算法进行简化优化得到基网格,将提出的优化算法与QEM算法无缝集成,最后用顶 点迭代调整细分拟合算法输出细分控制网格.该控制网格可作为后续设计与加工模块的输入,对其作 若干次Loop细分后即得到逼近于原始网格的光滑细分曲面. 关键词:医学假体CAD;Loop细分;细分曲面拟合;逆向工程;曲面重建 中图分类号:TP391.7文章编号:1004—132X(2007)17—2066,05 ALoopSubdivisionSurfaceFittingSystemforProsthesisCAD LiuSiyanLiaoWenheLiuHao NanjingUniversityofAeronauticsandAstronautics,Nanjing,20016 Abstract:TowardstheimplementationofprosthesisCADsystembasedonsubdivisionsurface representation,acompleteLoopsubdivisionsurfacefittingsystemwasputforward.Theinputdense surfacingfromCTimages,whichw triangularmeshwastheresultofMCiso— asthensmoothed,sim— plifiedandoptimizedtogetthebasemeshforsubdivisionsurfacefitting.QEMsimplificationalgo— rithmwasused,andanewoptimizationalgorithmwasputforward,whichwasseamlesslyintegrated withQEM.Then,thevertexiterativeadjustingsubdivisionsurfacefittingalgorithmwasimplemented tooutputthesubdivisioncontro1mesh,whichcanbeusedastheinputoftheafterwardsdesignand fabricationmodules.Fromthecontro1mesh,smoothsubdivisionsurfacewhichapproximatetheorigi- na1meshcanbeobtainedaftersevera1timesofLoopsubdivision. Keywords:medicalprosthesisCAD;Loopsubdivision;subdivisionsurfacefitting;reverseengi— neering;surfacereconstruction 0引言 随着生产,生活与医疗水平的提高,应用 CAD/CAM技术精确制作各种假体用于临床治 疗的需求在不断增长,这一需求促进了医学假体 CAD/CAM技术u的发展.近年来,细分曲面已 成为CG/CAGD领域的研究热点.,细分曲面 用于医学假体CAD前景广阔. Loop细分是基于三角网格的逼近型细分,该 细分模式因其简单性而特别适合工程应用j. Loop细分曲面拟合可分为数学优化法l3与顶点 迭代调整法l7两类.在数学优化法中,Hoppe 等首先扩展了Loop细分规则用于处理尖锐特 征,并提出了基于密集点云的细分曲面拟合方法, 但计算过程非常耗时.Ma等],李桂清等用 Garland等_1的QEM算法进行网格简化与优化 得到基网格,然后用最小二乘拟合得到控制网格. 最近李登高等又针对方程组的共轭梯度求解提 收稿日期:2OO6,O619 基金项目:国家863高技术研究发展 计划 项目进度计划表范例计划下载计划下载计划下载课程教学计划下载 资助项目 (2005AA420240) ? 2066? 出一种加速近似算法. Suzuki等采用顶点迭代调整法通过不断 调整控制网格顶点,使细分曲面逐渐逼近原始网 格,避免了方程组的求解.此法收敛速度很快,只 是拟合精度受基网格影响比较大,基网格中的不 良三角形容易造成拟合细分曲面的自交与扭 曲J.Kobbeh等[81提出了”变形薄膜(shrink wrapping)”算法,该算法除了对控制网格顶点施 加向原始网格的拉伸作用外,还用Laplacian算 子对网格顶点进行”松弛”,使拉伸后的网格顶点 分布趋于均匀,以避免或减轻自交扭曲现象,此法 引入了额外的平滑算子,通过使平滑算子系数逐 渐趋于零来保证算法的收敛性.周海等通过引 入Laplacian算子的切向分量实现了一定精度下 的收敛. 考虑到基网格质量对细分曲面拟合结果具有 重要影响,本文提出了一种与QEM简化算法集 成的网格优化算法来获取简化优化基网格的方 法,用迭代法拟合得到了光顺平滑的细分曲面,为 基于细分曲面表示的假体CAD/CAM系统开发 面向医学假体CAD的Loop细分曲面拟合系统——刘泗岩廖文和刘 浩 打下了初步基础. 1系统总 流程 快递问题件怎么处理流程河南自建厂房流程下载关于规范招聘需求审批流程制作流程表下载邮件下载流程设计 细分曲面拟合总流程如图1所示,分为三维 重建,网格处理与细分曲面拟合3个模块. 细分曲面拟合模块网格处理模块} @合化 图1细分曲面拟合系统流程 三维重建模块输入的是DICOM格式的CT 图像,用MC算法抽取等值面得到三角网格模 型.对CT图像进行滤噪分割等处理能够减少三 维重建的游离三角片噪声与不需要的其他组织数 据,处理方法一般采用手工勾勒与阈值分割相结 合的方法.本文先用文献[11]的方法进行图像预 处理,然后用一种”种子增长”三维分割方法来滤 除游离的三角片噪声,用三维重建网格后处理方 法弥补MC算法的三维重建质量的不足.实现 过程是:在目标组织的三维模型上拾取一个种子 三角片,对该三角片作标记,并向周围相邻三角片 扩散,直到所有相邻的三角片标记完毕,把未标记 的三角片作为噪声或背景删除,从而得到连通的 三角网格. 用网格处理模块对密集网格进行平滑,简化 与优化处理,得到细分拟合所需的基网格. 2基网格的生成 细分曲面拟合的基本思路是,首先对原始网 格进行简化得到拓扑相同,几何近似的稀疏网格, 称为基网格,然后通过细分重新加密网格,并按一 定几何精度逼近原网格.基网格的质量直接影响 细分曲面拟合的质量与效率,理想的基网格应能 在保持特征的同时尽可能地简化,且所含三角形 应尽量接近等边三角形. 目前用于网格简化的算法虽有很多,但同时 考虑网格简化及三角形形状优化的成功算例并不 多见.Hoppe等口2_的网格优化算法能够得到高 质量的网格,但计算比较耗时.Garland等口..的 QEM算法在计算高效率的同时保持了比较好的 质量,成为广泛应用的简化算法,但QEM算法未 考虑形状优化问题.李桂清等?5对QEM算法进 行扩展,加上了形状优化因子,但形状因子的加入 却有损于QEM算法固有的数学意义.本文提出 了两步法的简化与优化算法,即先用QEM算法 进行简化,然后对结果进行优化.优化算法与 QEM算法相集成,直接使用QEM算法的数据结 构,在保持QEM算法简捷高效优点的同时,获得 了较好的优化结果. 2.1OEM网格简化算法 Garland等..的QEM网格简化算法是一种 基于二次误差度量(quadricerrormetric,QEM) 的边折叠简化方法.所谓边折叠,是指将边 的两端点合并到新顶点的过程,此操作记为 一, 如图2a所示. 设为其中某三角形F所在平面-厂的单位 法向量,则任意顶点到平面-厂距离的平方可表 示为 ()一(n7)+2()+(1) 式中,d为原点到平面的有向距离. 定义边的二次误差算子Q对于任意 顶点,有 . ()一(+)()(2) FkC-{F7,FF:1 上式的计算结果实际上是新顶点到顶点 和v,相关平面集中所有平面的距离的平方和. Q(v)为二次曲面,对其求极小值可确定新顶点 v,若极小值不存在则取边的中点为新顶点v. QEM网格简化算法的基本原理是,通过计算 每条边的折叠代价,按从大到小的顺序压人堆栈, 每次从栈顶弹出折叠代价最小的边来执行折叠, 同时对其邻域内的边进行更新.对边折叠一 v,QEM误差按式(2)进行计算.若考虑网格中三 角形大小的影响,还需要对式(2)进行面积加权 处理,详见文献[10]. 根据QEM的几何意义,该算法具有自动保 持特征的特点,但对尖锐棱边与边界边还需要添 加约束平面F,以将折叠后的新顶点限制在棱边 或边界边上,约束平面的二次元算子为 QFc—KQF 式中,K为惩罚系数. Garland等口阳简单地取K一1000进行计算, 我们发现,当K>lO时即能比较严格地保持边 界,但会使边界处含有很多狭长的三角形,当 0.05<K<2时即可获得比较好的边界保持效 果,同时能够减轻边界处三角形的狭长化. 2.2网格优化技术 2.2.1网格质量评定准则 本文采用三角形正则度作为 评价 LEC评价法下载LEC评价法下载评价量规免费下载学院评价表文档下载学院评价表文档下载 网格质量的 ? 2067? 中国机械工程第18卷第17期2007年9月上半月 依据,三角形正则度是指三角形接近正三角形的 程度.正三角形的正则度为1,三角形越狭长,正 则度就越接近零. 定义若三角形F的三内角为a,,y,则定 义三角形F的正则度为R(F)一2(cos口+cos+ cos7—1). 2.2.2拓扑调整操作 Hoppe等用边折叠,边分裂与边交换3种 局部拓扑调整操作进行网格优化,操作过程如图 2所示. (b) (c) 图2局部拓扑调整操作 边折叠操作能有效删除网格中的针形三角 形,边交换操作更适合优化扁平三角形.边分裂操 作一般会带来新的狭长三角形,尽管新产生的这 些狭长三角形可用再次边折叠而删除,但不如直 接运用一次边交换更有效,边分裂后再进行边折 叠操作相当于一次边交换操作,如图2b,图2c 所示. Hoppe的网格优化法由于计算太复杂而限制 了其实用性,本文采用边折叠与边交换两种操作 设计了一种简捷而实用的网格优化算法. 2.2.3网格优化算法总流程 优化算法的总流程为:初始化时按正则度标 记出网格中的不良三角形集合{Fi,指定正则度 阈值r,如果R(F)<r则将三角形F添加到 {F;然后对{F中的每个三角形计算边折叠 或边交换操作的执行优先度p(F),按从大到小的 顺序依次执行.设三角形F进行边折叠的优先度 为P(F),进行边交换的优先度为P(F),若 P(F)>声(F),则需对F进行边折叠操作,并令 (F)一P(F),反之亦然.三角形F的边折叠操 作是指对F的最小内角所对的边进行折叠;F的 边交换操作是指对F的最大内角所对的边进行交 ? 2068? 换(边交换只限于内部光滑边)操作.若F的最大 内角所对的边恰为棱边或边界边,则直接令 P(F):0. 2.2.4拓扑调蔓操作的优先度计算 拓扑调整操作的优先度计算有两个原则,一 是要提高正则度,AR应尽可能地大,二是几何误 差/xG应尽可能地小.本文拓扑调整操作的优先 度用下式计算: pc(F)一7凸:} (3) P(F)一JL,sJ 式中,为比值系数,>0. 正则度变化AR用拓扑调整前后邻域结构内 的三角形最小正则度的变化量来表示,进行调整 的前提是要保证AR>0,否则不作调整. 2.2.5.&,lW’kU-差的计算 2.2.5.1边折叠的几何误差计算 优化过程中的边折叠以三角形正则度提高为 目的,同时限制最大几何误差不超过 规定 关于下班后关闭电源的规定党章中关于入党时间的规定公务员考核规定下载规定办法文件下载宁波关于闷顶的规定 的值. Remi等_1采用边折叠后的新顶点到相关平面集 的最大距离平方来表示边折叠误差.QEM实际上 是此方法的变形与扩展,只是用距离的平方和代 替了最大距离平方,并提供了一种基于二次元 (Quadric)的简捷实现方法.本文参考Remi等与 Garland等的方法,用新顶点到相关平面集的最 大距离来计算边折叠误差,相当于Garland等所 用二次元的最大平方根,因而称为max~/Q误差 度量,即 AG:一max(max,_(,max,,F?sar()F?( 2.2.5.2边交换的几何误差计算 由于边交换操作与边分裂加边折叠操作等 效,而进行边分裂操作时,边的几何关系是不变 的,所以可以把边交换操作的几何误差计算转化 为边折叠操作的几何误差计算. 对图2b中边分裂操作的结果再进行边折叠 一’,或一操作,即可得到图2c所示的边 交换结果.取两种能进行边折叠操作,且几何误差 较小的作为边交换的几何误差,即 /xGs—min(~/”.)(),~/.()) 3细分曲面拟合 3.1细分极限位置 网格顶点在细分极限曲面上的位置计算是细 分拟合算法的基础,Loop等常用的细分模式已有 成熟的计算公式.].带特征的Loop细分将顶点 面向医学假体CAD的Loop细分曲面拟合系统——刘泗岩廖文和刘 浩 分为光滑顶点,刺点,棱边与边界顶点,角点等类 型.其中角点在细分中的位置不变,在细分拟合中 也不必再作调整. 对光滑顶点与刺点,极限位置按下式计算: 一n+?(一Pi)(4) ? 一 ‘+) :1E5一(3十1c.s).] 设顶点’,在棱边或边界上的两相邻顶点为 ‘,H,Vi+1,则 =V+[(V一,V)+(V一V)](5) 3.2顶点迭代调整细分拟合算法 3.2.1逼近调整 设原始密集网格为M,基网格为S.,细分k次 后的网格为S,细分极限曲面为S..,细分拟合就 是实现映射Ss,使得S..按一定精度逼 近M.设S上任意顶点’,在S上的相应位置为 ‘,,’,在M上的最近点为m;M上任意顶点P在 S..上的最近点为.注意最近点不一定是网格的 顶点,也可能在边或者三角形面内. 要使S..理想地逼近网格M,显然应有m一 ‘,,代入式(4),得 +(一Pi)一ml一0(6) ?} 顶点迭代法的思路是不断调整顶点’,的位 置,从而使’,逼近到m.取每次调整量?’,为 ?v一1(m一)0<?1(7) 每次调整’,到新位置’,一’,十?’,,代入式 (4)得到逼近调整操作为 +(J,l一)一(一p)(8) ?{} 用来控制逼近速度,当一1时逼近速度最 快,但会造成局部扭曲,我们取一0.65. 3.2.2优化调整 上述顶点调整会使网格中的三角形质量下 降,因此用切向Laplacian平滑算子进行纠正,从 而使顶点分布重新趋于均匀.注意此调整操作只 针对内部光滑点,对边界,棱边与刺点均不作此 操作. 若顶点v的相邻顶点个数为,v,为任意邻 点,Laplacian算子为 L(v)一‚?(vj—v)” ?】 L(v)可分解为法向分量Lt(‘,)与切向分量 L(‘,),L(‘,)使顶点’,向网格内部移动,引起网 格收缩;L(‘,)使顶点’,沿切向滑动,使顶点分布 趋于均匀,起到网格优化的作用.若顶点v的法向 量为,l,L(‘,)计算公式为 L(v)一L(v)一[L(v)nIn 顶点’,的法向量n的估算采用’,相邻三角 形法向量的面积加权平均. 从而得到优化调整操作为 一 +正(p)(9) 其中,/1?E0,1]是优化调整因子,/1取值太小起 不到优化效果,如取大值又会使拟合误差很小(例 如取e.<10),算法不收敛.本文初始设置/1= 0.6,以后每经过一次迭代,/1值缩减75,经过5 次迭代后,/1—0.6x0.75?0.19,此后设/1—0, 即不再进行优化调整,因为此时细分网格顶点已 充分接近原始网格. 3.2.3误差评估 S..与M的逼近精度用双向L误差表示,其 中S..到M的L.误差为 E】一Ilmi一ll. ‘. 式中,ff为网格&所含顶点的数量. M到S..的L误差为 E2一一Pll.(10) .p?M 由于在计算E前S已经作了若干次细分, 计算E时,再对S细分2,3次,网格就已非常 密集,故用S抖(一2或h一3)代替细分极限曲 面S.当k<3且基网格S.比较稀疏时取h一3, 其他情况均取h一2. E-误差控制的是S网格顶点的细分极限位 置点对原始网格M的插值精度,它并不能保证S 细分后添加的新顶点也贴近网格M.通过对S进 一 步细分,当网格足够密集后,用E来验证细分 曲面对原始网格M的逼近精度. 3.2.4算法流程 细分拟合模块的输人为初始网格M及其简 化与优化后所得基网格S.,拟合过程包括外部的 细分循环与内部的逼近循环两个子过程.在内部 的逼近循环中,对网格S的每个顶点先后按式 (8),式(9)进行迭代逼近与优化调整,使S网格 顶点的细分极限点到网格M的逼近误差E.?e; 在外部的细分循环中,对S细分一次得到S计., 再对S抖-进行逼近循环,细分.次后验证E,根 据基网格S.的疏密程度,一般取为2或.为 3,S.较稀疏时取一3,S.较密集时取.一2.如 果E2>e,则需要对S继续细分与迭代调整,直 到E?ez.如果S细分次以后仍然E>e,可 ? 2O69? 中国机械工程第18卷第17期2007年9月上半月 能继续细分到程序崩溃也不能收敛,程序给出提 示信息”不能拟合到此精度”而返回.的取值根 据基网格S.的规模与密度而定,一般可取4,6, 对稀疏的基网格,可取z一6,对规模比较大或者 比较密集的基网格,取为4或5. 4结果与总结 图3为股骨头细分曲面拟合结果.图3a为CT 图像三维重建并经滤噪平滑等处理后所得原始网 格M,含有20015个三角形;图3b为M经过本文 简化优化算法所得基网格s.,含45个三角形;图 3c为s.经过两次细分与三次迭代调整所得控制 网格S,取e=10,,e一10,此精度已能满足 假体制作的需要,考虑到CT测量与三维重建中 的误差,取更高的拟合精度是没必要的;图3d为 S细分两次的结果,它对M的实际逼近精度 E,=0.0034. (c)控制网格(d)细分曲面 图3股骨头曲面拟合结果 拟合结果所得控制网格可以作为后续设计与 加工的输入,设计模块通过对控制网格的编辑来 得到所需医学假体.假体的控制网格可以进一步 输入到加工模块,再经过若干次细分后,进行切片 处理与RP制造. 本文从CT图像三维重建开始,完整地实现 了一个面向医学假体CAD应用的Loop细分曲 面拟合系统.提出了与QEM网格简化算法相集 成的网格优化算法,以此得到基网格,然后用顶点 迭代调整细分拟合法得到细分控制网格,为基于 Loop细分曲面的假体CAD/CAM系统的开发打 下了初步基础. ? 2070? 参考文献: [1]SunW,StarlyB,NamJ,eta1.Bio,CADModel [2] [3] [4] [5] [6] [7] [8] [9] ingandItsApplicationsinComputer--aidedTissue Engineering[J]Computer--aidedDesign,2005(7): 1—18. LoopCSmoothSubdivisionSurfacesBasedonTri— angles[D].Utah:UniversityofUtah,1987. HoppeH,DeRoseT,DuchampT,etalPiece— wiseSmoothSurfaceReconstruction[C]//ACM SIGGRAPH.NewYork:ACMPress,1994:295— 302 MaWY,MaXH,TSOSK,eta1.ADirectAp— proachforSubdivisionSurfaceFittingfromaDense TriangleMesh[J]Computer--aidedDesign,2004, 36(6):525—536 李桂清,马维银,鲍虎军带尖锐特征的Loop细分 曲面拟合系统[J].计算机辅助设计与图形学, 2005,17(6):1179—1185. 李登高,秦开怀Loop细分曲面的优化拟合算法 [J].计算机辅助设计与图形学,2006,18(6): 755—759 SuzukiH.TakeuchoiS.KanalT.SubdivisionSur— faceFittingtOaRangeofPoints[c]//Proceedings ofthe7thPacificConferenceonComputerGraphics andApplicationsWashington:IEEEComputerSo— ciety,1999:158—167. KobbeltLP,VorsatzJ,LabsikU,eta1.AShrink WrappingApproachtORemeshingPolygonalSur— faces[J].ComputerGraphicsForum,1999,18(1): 209—237. 周海,周来水,王占东,等.散乱数据点的细分曲面 重构算法及实现[J].计算机辅助设计与图形学学 报,2003,15(10):1287—1292. [10]GarlandM,HeckbertPS.SurfaceSimplification UsingQuadricErrorMetrics[C]//ACMSIG— GRAPH.NewYork:ACMPress,1997:209—216. [11]秦绪佳,欧宗瑛,纪凤欣.医学图像的交互分割及 三维表面重建[J].工程图学,2001,22(2):94 - 10l_ [12]HoppeH,DeRoseT,DuchampT,eta1.Mesh Optimization[C]//ACMSIGGRAPH.Anaheim: ACMPress,1993:19-26. [13]RemiR,JarekR.Full—rangeApproximationof TriangulatedPolyhedra[J].ComputerGraphics Forum,1996,15(3):67—76.(编辑何成根) 作者简介:刘泗岩,男,1975年生.南京航空航天大学机电学院 博士研究生.主要研究方向为CAD/CAM及其在生物医学工程 中的应用.发表 论文 政研论文下载论文大学下载论文大学下载关于长拳的论文浙大论文封面下载 2篇.廖文和,男,1965年生.南京航空航 天大学机电学院教授,博士研究生导师.刘浩,男,1975年生. 南京航空航天大学机电学院讲师. ?一?