CN108776956B - 一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法 - Google Patents
一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法 Download PDFInfo
- Publication number
- CN108776956B CN108776956B CN201810339104.3A CN201810339104A CN108776956B CN 108776956 B CN108776956 B CN 108776956B CN 201810339104 A CN201810339104 A CN 201810339104A CN 108776956 B CN108776956 B CN 108776956B
- Authority
- CN
- China
- Prior art keywords
- vertex
- model
- spectrogram
- grid
- dimensional
- 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
- 238000000034 method Methods 0.000 title claims abstract description 80
- 238000001514 detection method Methods 0.000 claims abstract description 15
- 238000004458 analytical method Methods 0.000 claims abstract description 7
- 230000008447 perception Effects 0.000 claims abstract description 5
- 239000013598 vector Substances 0.000 claims description 28
- 239000011159 matrix material Substances 0.000 claims description 10
- 238000000354 decomposition reaction Methods 0.000 claims description 3
- 238000012886 linear function Methods 0.000 claims description 3
- 210000000056 organ Anatomy 0.000 abstract description 3
- 238000002474 experimental method Methods 0.000 description 13
- 238000012545 processing Methods 0.000 description 12
- 238000001914 filtration Methods 0.000 description 11
- 238000011282 treatment Methods 0.000 description 8
- 238000012800 visualization Methods 0.000 description 8
- 210000000988 bone and bone Anatomy 0.000 description 6
- 230000004069 differentiation Effects 0.000 description 6
- 230000008569 process Effects 0.000 description 6
- 238000003745 diagnosis Methods 0.000 description 5
- 238000005516 engineering process Methods 0.000 description 5
- 230000014759 maintenance of location Effects 0.000 description 5
- 230000002146 bilateral effect Effects 0.000 description 4
- 230000000052 comparative effect Effects 0.000 description 4
- 238000002591 computed tomography Methods 0.000 description 4
- 238000010586 diagram Methods 0.000 description 4
- 238000009499 grossing Methods 0.000 description 4
- 230000008859 change Effects 0.000 description 3
- 238000013461 design Methods 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 210000003692 ilium Anatomy 0.000 description 3
- 238000005457 optimization Methods 0.000 description 3
- 238000004321 preservation Methods 0.000 description 3
- 210000002303 tibia Anatomy 0.000 description 3
- 230000002411 adverse Effects 0.000 description 2
- 238000013079 data visualisation Methods 0.000 description 2
- 238000012938 design process Methods 0.000 description 2
- 230000006870 function Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000009877 rendering Methods 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 230000003595 spectral effect Effects 0.000 description 2
- 241000276398 Opsanus tau Species 0.000 description 1
- 230000002159 abnormal effect Effects 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 238000004195 computer-aided diagnosis Methods 0.000 description 1
- 230000008602 contraction Effects 0.000 description 1
- 201000010099 disease Diseases 0.000 description 1
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 1
- 230000009977 dual effect Effects 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 210000002391 femur head Anatomy 0.000 description 1
- 238000005286 illumination Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 238000001959 radiotherapy Methods 0.000 description 1
- 238000005096 rolling process Methods 0.000 description 1
- 238000001356 surgical procedure Methods 0.000 description 1
- 210000001519 tissue Anatomy 0.000 description 1
- 210000000689 upper leg Anatomy 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- 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
- G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/20—Image enhancement or restoration using local operators
-
- 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/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10081—Computed x-ray tomography [CT]
-
- 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
- G06T2207/30008—Bone
Landscapes
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Graphics (AREA)
- Geometry (AREA)
- Software Systems (AREA)
- Apparatus For Radiation Diagnosis (AREA)
- Measuring And Recording Apparatus For Diagnosis (AREA)
Abstract
本发明公开了一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法,包括如下步骤:(1)构建特定方向感知的特征检测方法,准确识别梯田型噪声;(2)建立三维医学网格模型离散拉普拉斯‑贝尔特拉米算子,执行谱图分析,构建谱图空间;(3)构建基于非均匀谱图编码的三维低通滤波器去除高频随机噪声以及梯田型噪声。本发明不仅能去除三维医学模型高频噪声,还能有效去除梯田型噪声,同时本发明能很好地保持模型体积,光顺的结果能最大限度逼近于人体真实器官。
Description
技术领域
本发明涉及医学数据可视化领域,具体涉及一种三维医学模型光顺的方法。
背景技术
通过一组二维计算机断层扫描数据(CT)重建的三维医学模型目前已广泛地应用于解剖三维重建、立体仿真手术与导航、放射治疗规划、计算机辅助诊断以及个性化假肢再造等领域。然而由于机器设备、扫描分辨率以及环境噪声的影响,重建的三维医学模型不仅含有高频噪声,普遍还会含有一种特殊的梯田型噪声。这种梯田型噪声不仅会使组织器官呈现出非正常结构,同时还会干扰医生对疾病的诊断,对后续的临床治疗也会产生及其不利的影响。因此,采用网格光顺方法去除高频噪声,尤其是梯田型噪声是三维医学数据可视化必不可少的步骤。
目前已有的网格光顺方法中,通用的网格光顺方法包括各向同性和各向异性网格光顺方法。各向同性方法例如Laplace光顺[D.A.Field,Laplacian smoothing anddelaunay triangulations,Communications in applied numerical methods 4(6)1988,709–712]虽然能有效去除噪声,但是模型会产生体积收缩并且无法保持特征。为了能在网格光顺的同时保持特征防止变形,国内外学者在各向异性网格光顺方面展开了相关工作。Fleishman等人[S.Fleishman,I.Drori,D.Cohen-Or,Bilateral mesh denoising,in:ACMTrans actions on Graphics(TOG),Vol.22,ACM,2003,pp.950–953.]将二维图像处理中双边滤波器的思想拓展到三维,提出基于双边滤波去噪算子的特征保持的网格模型光顺方法。Jones等人[T.R.Jones,F.Durand,M.Desbrun,Non-iterative,feature-preservingmesh smoothing,in:ACM Transactions on Graphics(TOG),Vol.22,ACM,2003,943–949]提出了非迭代的网格光顺方法,该方法通过控制顶点邻域的大小来保持网格模型的特征。Hildebrandt等人[K.Hildebrandt,K.Polthier,Anisotropic filtering of non-linearsurface features,in:Computer Graphics Forum,Vol.23,Wiley Online Library,2004,391–400]提出了各向异性设定平均曲率流的方法(Prescribed Mean Curvature,PMC),该方法可以在有效去除噪声的同时保证光顺后的网格表面能收敛到预测的曲率分布网格表面,网格模型体积以及特征都得到保持。Zheng等人[Y.Zheng,H.Fu,O.-C.Au,C.-L.Tai,Bilateral normal filtering for mesh denoising,Visualization and ComputerGraphics,IEEE Transactions on 17(10)(2011)1521–1530.]提出了基于双面法向量滤波的网格光顺方法。He等人[L.He,S.Schaefer,Mesh denoising via l 0 minimization,ACMTransactions on Graphics(TOG)32(4),2013,64]提出了基于L0最小化的光顺方法。Zhu等人[L.Zhu,M.Wei,J.Yu,W.Wang,J.Qin,P.-A.Heng,Coarse-to-fine normal filteringfor feature-preserving mesh denoising based on isotropic subneighborhoods,in:Computer Graphics Forum,Vol.32,Wiley Online Library,2013,371–380]提出了多尺度的可以同时保持尖锐以及平缓特征的网格光顺方法。Cheng等人[X.Cheng,M.Zeng,X.Liu,Feature-preserving filtering with l 0 gradient minimization,Computers&Graphics 38,2014,150–157]提出基于特征保持的近似L0梯度最小化光顺方法。Wang etal.[P.S.Wang,X.M.Fu,Y.Liu,X.Tong,S.L.Liu,B.Guo,Rolling guidance normal filterfor geometric processing,Acm Transactions on Graphics 34(6)(2015)173.]提出RGN滤波器有效消除小尺度几何特征保存网格模型的大尺度特征。Lu等人[X.Lu,Z.Deng,W.Chen,A robust scheme for feature-preserving mesh denoising,IEEETransactions on Visualization&Computer Graphics 22(3)(2016)1181.]采用基于二次优化的方法进行特征检测,并基于特征边迭代更新顶点位置,获得光滑表面,同时保持特征。Yadav等人[S.K.Yadav,U.Reitebuch,K.Polthier,Mesh denoising based on normalvoting tensor and binary optimization,IEEE Transactions on Visualization&Computer Graphics PP(99)(2017)]提出了一种基于法向张量投票和二进制优化的网格去噪方法。Wei等人[M.Wei,L.Liang,W.M.Pang,J.Wang,W.Li,H.Wu,Tensor voting guidedmesh denoising,IEEE Transactions on Automation Science&Engineering 14(2)(2017)931–945.]提出了一种基于张量投票的网格去噪方法,该方法结合面法线和二次曲面拟合出分段光滑曲面。
然而,无论是各向同性的网格光顺方法还是各向异性的网格光顺方法均不能很好地处理三维医学模型的网格光顺问题。原因是应用于数字化医学诊疗的三维医学模型对网格模型光顺技术提出了较高的要求:首先网格光顺技术不仅需去除三维医学模型高频噪声,还应去除梯田型噪声,因为真实人体骨骼器官不存在锐利的边以及尖锐的点;同时,经光顺处理的三维医学模型不能产生体积收缩,应最大限度逼近于人体真实骨骼。而各向同性的网格光顺方法虽然能去除噪声,但是会产生过度光顺的现象,模型体积收缩严重,不能满足医学诊疗的需求;各种各向异性的网格光顺方法在网格光顺的过程中会将三维医学模型梯田型噪声识别为模型的特征加以保持而非去除,同样不能满足医学诊疗的需求。
直接针对三维医学网格模型光顺处理的方法中,Whitaker等人[R.T.Whitaker,Reducing aliasing artifacts in iso-surfaces of binary volumes,in:Proceedingsof the 2000 IEEE symposium on Volume visualization,ACM,2000,pp.23–32.]、Nielson等人[G.M.Nielson,Dual marching cubes,in:Proceedings of the conferenceon Visualization’04,IEEE Computer Society,2004,pp.489–496.]和Gibson等人[S.F.Gibson,Constrained elastic surface nets:Generating smooth surfaces frombinary segmented data,in:Medical Image Computing and Computer-AssistedInterventation MICCAI’98,Springer,1998,pp.888–898.]的方法虽然可以获得令人满意的光顺效果,但是遗憾的是这些方法在执行过程中都离不开原始的二维医学切片图像数据的处理,并非可以独立执行的数字几何处理方法。Bade等人[R.Bade,J.Haase,B.Preim,Comparison of fundamental mesh smoothing algorithms for medical surfacemodels.,in:SimVis,Vol.6,Citeseer,2006,pp.289–304.]对三维医学网格模型的光顺问题进行了较为深入的讨论。其中提到,三维低通滤波器能在一定程度上解决三维医学模型的噪声去除的问题。三维低通滤波器中,Taubin[G.Taubin,A signal processingapproach to fair surface design,in:Proceedings of the 22nd annual conferenceon Computer graphics and interactive techniques,ACM,1995,pp.351–358.]的方法能获得较好的结果,但是该方法运行的过程中要进行繁琐的参数设置,而这些参数又不具备直观意义,因此给操作者带来极大不便。同时,该方法仅能去除较平缓的梯田型噪声,对于尖锐的梯田型噪声依然无法有效去除。Levy等人[B.Vallet,B.L′evy,Spectral geometryprocessing with manifold harmonics,in:Computer Graphics Forum,Vol.27,WileyOnline Library,2008,pp.251–260.]的三维低通滤波器可以很好地去除模型的高频噪声并保持模型体积,但是该低通滤波器会将梯田型噪声识别为模型几何结构加以重建,从而导致梯田型噪声无法被有效清除。
发明内容
本发明的目的在于提供一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法。该方法可以有效去除三维医学模型高频噪声以及梯田型噪声,同时有效保持模型体积,获得的三维医学模型能最大程度逼近于人体真实骨骼器官,满足数字化临床医学诊疗的需求。
为达上述目的,本发明采用的技术方案如下:
一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法,其特征在于,包括如下步骤:(1)构建特定方向感知的特征检测方法,准确识别梯田型噪声;(2)建立三维医学网格模型离散拉普拉斯-贝尔特拉米算子,执行谱图分析,构建谱图空间;(3)构建基于非均匀谱图编码的三维低通滤波器去除高频随机噪声以及梯田型噪声,获得光顺的三维医学模型。
本发明所述步骤(1)为构建特定方向感知的特征检测方法以准确识别沿z轴方向分布的梯田型噪声:
首先,构建笛卡尔坐标系;
接着,定义方向感知的特征检测方法,公式如下:
其中,k为网格顶点i一邻域邻接三角形,为一邻域邻接三角形面法向量;计算顶点i一邻域邻接三角形面法向量与x轴之间的角度差,记为找出网格顶点i一邻域三角形沿x轴方向的最大偏离以及最小偏离计算沿x轴方向的最大偏离以及最小偏离之间的差,记为δx,i;同理,在y轴以及z轴方向也执行上述相同的操作,得到δy,i与δz,i。
本发明所述步骤(1)设置一个合适的阈值τ,将δx,i>τ的顶点标识为x轴方向特征点,将δy,i>τ的顶点标识为y轴方向特征点,将δz,i>τ的顶点标识为z轴方向特征点;由于三维医学模型的梯田型噪声仅沿z轴方向分布,因此,若三维医学模型顶点i的δz,i>τ,则顶点i为梯田型噪声点;对于梯田型噪声点,为其设置标签ωi=1,非梯田型噪声的网格顶点,为其设置标签ωi=0。
本发明所述步骤(2)为构建离散拉普拉斯-贝尔特拉米算子,执行三维医学模型谱图分析,构建谱图空间:
首先,构建离散拉普拉斯-贝尔特拉米操作算子Δ(Laplace-Beltrami operatorΔ,以下简称LBO操作算子)。LBO操作算子是定义在黎曼流形上的作用于标量函数f的二阶微分操作算子,在网格曲面上,LBO操作算子可以定义为:
wij=(cotαij+cotβij)/2Ai
其中,i与j为网格顶点;设V为顶点集合,E为边集,N(i)={j∈V|(i,j)∈E}是网格顶点i的邻接点;wij=(cotαij+cotβij)/2Ai是边(i,j)的权重,αij与βij是共享同一条边(i,j)的两个邻接三角形的对角,Ai为网格顶点i一邻域Vorionor面积;将网格模型离散拉普拉斯-贝尔特拉米算子写成矩阵L的形式如下:
其中,i与j为网格顶点;设V为顶点集合,E为边集,N(i)={j∈V|(i,j)∈E}是顶点i的邻接点;wij=(cotαij+cotβij)/2Ai为边(i,j)的权重;矩阵L特征分解的结果是产生一系列特征值λi与特征向量fi(1≤i≤n),这里n为网格顶点数,特征值与特征向量成对出现(λi,fi),特征向量两两正交;将模型顶点几何信息看作是信号,将其投射到正交的特征向量构造谱图空间,首先需要标准化特征向量:
Φi=fi/‖fi‖
其中,标准化的特征向量Φi构建一个矩阵,矩阵的第i行提供了网格顶点i空间几何坐标(xi,yi,zi)的一个嵌入,Φi的第k个元素Φi(k)为网格顶点i的分段线性函数;
接着,利用所得到的标准化的特征向量Φi构造网格模型顶点的谱图空间,公式如下:
本发明所述步骤(2)中的k最大值被设置为1000。基于傅里叶变换的原理,低频因子对应于模型的基本几何外形,高频因子对应于模型的细节特征。由于本发明设计的是一个三维低通滤波器,因此需要确定所采用的低频因子数量值。注意到,若k值取值过小,则不能很好地描述模型的几何外形,且重建的模型体积收缩严重;若k值取值过大,运算效率会降低。本发明设计过程中,记录了从采用前100个频率因子到前1500个频率因子进行重建时网格模型的体积变化。结果显示,当k小于900时,重建的模型体积丢失较多,而当k超过1000时,大部分模型体积收缩率都低于3%。因此,本发明k最大值被设置为1000。
本发明所述步骤(3)为构建基于非均匀谱图编码的三维低通滤波器,有效去除高频噪声以及梯田型噪声获得光顺的三维医学模型:
m=m*(1-0.8ωi),
根据在步骤(1)中的定义,网格顶点若为梯田型噪声顶点,顶点标签ωi被设置为1(ωi=1),则在本步骤中仅有前200个(0.2m,m=1000)频率因子被采用重建网格梯田型噪声顶点,梯田型噪声被有效抑制;而在网格非梯田型噪声顶点,由于在步骤(1)中顶点标签ωi被设置为0(ωi=0),因此在本步骤中前1000个(m=1000)频率因子被采用重建网格顶点,高频噪声被有效消除,同时模型体积得以有效保持。
相比于现有的技术,本发明具有以下优点:
通过一组二维计算机断层扫描数据(CT)重建的三维医学模型在医学数据可视化方面有广泛的应用。但是,原始的三维医学模型不仅包含高频噪声还含有一种特异的梯田型噪声。使用现有的网格光顺技术均无法有效清除该梯田型噪声。原因是现有的通用的网格光顺技术会将梯田型噪声识别为特征加以保持;现有的三维低通过滤器会将梯田型噪声识别为模型固有结构加以重建而无法清除。本发明提供了一种基于非均匀谱图编码的三维低通滤波方法,采取在梯田型网格顶点以及非梯田型网格顶点差异化几何信息重建的方法,本发明不仅能有效去除三维医学模型高频噪声,还能有效去除现有网格光顺技术无法去除的梯田型噪声,同时还能很好地保持三维医学模型体积,获得高度光顺的能满足医学临床数字化诊疗需求的三维医学模型。
附图说明
图1为本发明实施例提供的一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法的流程示意图;
图2为原始三维医学模型以及本发明实施例光顺效果示意图,(a)为原始三维医学模型,(b)本发明光顺后的三维医学模型,(图2(a)与(b)中左图均为前视图,右图均为后视图);
图3为三维网格模型方向感知的特征检测方法示意图,(a)三维网格模型及笛卡尔坐标系(Cartesian coordinate system),(b)x轴导向的特征点,(c)y轴导向的特征点,(d)z轴导向的特征点,(阈值:30°);
图4为三维医学网格模型梯田型噪声检测示意图;
图5为本发明实施例边权wij设置示意图,(a)网格顶点i及其一邻域顶点,(b)αij与βij为共享同一条边(i,j)两个三角形的对角;
图6为采用不同频率因子进行重建时三维医学网格模型体积的变化,水平坐标轴:参与网格重建的频率因子的个数,垂直坐标轴:重建的模型相比于原始模型体积收缩率;
图7为基于均匀谱图编码的低通滤波器与本发明实施例基于非均匀谱图编码的低通滤波器对钩骨模型光顺实验结果线框模型对比图,其中(a)原始网格模型及局部放大图,(b)基于均匀谱图编码的低通滤波器实验结果线框模型,(c)本发明实施例基于非均匀谱图编码的低通滤波器实验结果线框模型;
图8为基于均匀谱图编码的低通滤波器与本发明实施例基于非均匀谱图编码的低通滤波器对三维医学模型光顺实验结果对比图,其中(a)原始三维医学网格模型,(b)基于均匀谱图编码的低通滤波方法实验结果,(c)本发明实施例基于非均匀谱图编码的低通滤波方法实验结果;
图9为采用APSS方法、Zheng等人的方法、Taubin等人的方法以及本发明实施例方法对三维医学模型网格光顺对比实验结果,其中(a)原始三维医学网格模型,(b)APSS方法光顺的结果,(c)Zheng等人方法光顺的结果,(d)Taubin等人方法光顺的结果,(e)本发明实施例方法光顺的结果;
图10为Guennebaud等人的APSS方法、Zheng等人的方法、Taubin等人的方法及本发明实施例方法的体积保持对比实验结果。
具体实施方式
下面结合附图及具体实施方式对本发明的技术方案进行详细说明。
本发明实施例提供的基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法,包括如下步骤:
步骤1、首先构建特定方向感知的特征检测方法以准确识别沿z轴方向分布的梯田型噪声:
受激光扫描环境、机器物理特性、扫描物体表面、光照以及扫描分辨率等因素的影响,通常原始的三维医学模型会含有一种典型的梯田型噪声,如图2(a)所示。梯田型噪声的存在会对后续的医学诊疗产生极其不利的影响。考虑到梯田型噪声均沿z轴方向分布,本发明提供以下特定方向感知的特征检测方法以准确识别沿z轴方向分布的梯田型噪声。
首先,构建笛卡尔坐标系,如图3(a)所示。
接着,定义方向感知的特征检测方法,公式如下:
其中,k为网格顶点i一邻域邻接三角形,为一邻域邻接三角形面法向量。计算顶点i一邻域邻接三角形面法向量与x轴之间的角度差,记为找出网格顶点i一邻域三角形沿x轴方向的最大偏离以及最小偏离计算沿x轴方向的最大偏离以及最小偏离之间的差,记为δx,i。同理,在y轴以及z轴方向也执行上述相同的操作,得到δy,i与δz,i。
设置一个合适的阈值τ(图3中,模型阈值τ设为30°),将δx,i>τ的顶点标识为x轴方向特征点(如图3(b)所示)、将δy,i>τ的顶点标识为y轴方向特征点(如图3(c)所示),将δz,i>τ的顶点标识为z轴方向特征点(如图3(d)所示)。
注意到,三维医学模型的梯田型噪声仅沿z轴方向分布。因此,设置一个阈值τ,若三维医学模型顶点i的δz,i>τ,则顶点i为梯田型噪声点,如图4所示。对于梯田型噪声点,为其设置标签ωi=1,非梯田型噪声的网格顶点,为其设置标签ωi=0。本发明在六个三维医学模型实施例上进行梯田型噪声检测。六个三维医学模型实施例包括:钩骨模型,第三近端指骨模型,髂骨模型,胫骨模型,股骨头模型,股骨髁状突模型。每一个三维医学模型梯田型噪声检测阈值的设置见表1。
步骤2、构建离散拉普拉斯-贝尔特拉米算子,执行三维医学模型谱图分析,构建谱图空间:
首先,构建离散拉普拉斯-贝尔特拉米操作算子Δ(Laplace-Beltrami operatorΔ,以下简称LBO操作算子)。LBO操作算子是定义在黎曼流形上的作用于标量函数f的二阶微分操作算子。在网格曲面上,LBO操作算子可以定义为:
wij=(cotαij+cotβij)/2Ai
其中,i与j为网格顶点。设V为顶点集合,E为边集,N(i)={j∈V|(i,j)∈E}是网格顶点i的邻接点。wij=(cotαij+cotβij)/2Ai是边(i,j)的权重,αij与βij是共享同一条边(i,j)的两个邻接三角形的对角,Ai为网格顶点i一邻域Vorionor面积,边权wij的定义如图5所示。将网格模型离散拉普拉斯-贝尔特拉米算子写成矩阵L的形式如下:
其中,i与j为网格顶点。设V为顶点集合,E为边集,N(i)={j∈V|(i,j)∈E}是顶点i的邻接点。wij=(cotαij+cotβij)/2Ai为边(i,j)的权重。矩阵L特征分解的结果是产生一系列特征值λi与特征向量fi(1≤i≤n),这里n为网格顶点数。特征值与特征向量成对出现(λi,fi),特征向量两两正交。本发明将模型顶点几何信息看作是信号,将其投射到正交的特征向量构造谱图空间。首先需要标准化特征向量:
Φi=fi/‖fi‖
其中,标准化的特征向量Φi构建一个矩阵,矩阵的第i行提供了网格顶点i空间几何坐标(xi,yi,zi)的一个嵌入,Φi的第k个元素Φi(k)为网格顶点i的分段线性函数。
接着利用所得到的标准化的特征向量构造网格模型顶点的谱图空间,公式如下:
其中,xi=(xi,yi,zi)(1≤i≤n)为网格顶点空间几何坐标,这里n为网格顶点数, 为谱图空间几何频率因子。基于傅里叶变换的原理,低频因子对应于模型的基本几何外形,高频因子对应于模型的细节特征。由于本发明设计的是一个三维低通滤波器,因此需要确定所采用的低频因子数量值。注意到,若k值取值过小,则不能很好地描述模型的几何外形,且重建的模型体积收缩严重;若k值取值过大,运算效率会降低。本发明设计过程中,记录了从采用前100个频率因子到前1500个频率因子进行重建时网格模型的体积变化,如图6所示。结果显示,当k小于900时,重建的模型体积丢失较多,而当k超过1000时,大部分模型体积收缩率都低于3%。因此,本发明k最大值被设置为1000。本发明采用Arnolidi方法[G.H.Golub,C.F.Van Loan,Matrix computations,Vol.3,JHU Press,2012]计算每一个医学模型实施例的前1000个特征值-特征向量对为每一个模型构建谱图空间,模型谱图分析所用的计算时间见表1。
表1实施例谱图空间计算时间以及参数设置
步骤3、构建基于非均匀谱图编码的三维低通滤波器,有效去除高频噪声以及梯田型噪声获得光顺的三维医学模型:
使用在步骤1中获取的ωi以及在步骤2中获取的频率因子构建基于非均匀谱图编码的三维低通滤波器对原始医学网格模型进行光顺去噪。本发明与现有的基于均匀谱图编码的三维低通滤波器[B.Vallet,B.L′evy,Spectral geometry processing with manifold harmonics,in:ComputerGraphics Forum,Vol.27,Wiley Online Library,2008,pp.251–260.]作对比实验以验证本发明的有效性。采用现有的基于均匀谱图编码的三维低通滤波器对钩骨模型进行网格光顺的结果如图7(b)所示。现有的三维低通滤波器虽然可以去除医学模型高频噪声,但是无法去除梯田型噪声。原因是该低通滤波器在全部网格顶点均采用相同频率构造几何信息,梯田型噪声不可避免地会被当成几何模型固有结构加以重建而无法被去除。为了解决这个问题,本发明构建了一种新的基于非均匀谱图编码的三维低通滤波器,公式如下:
m=m*(1-0.8ωi),
其中,m为参与医学模型三维重建的频率因子的个数,m值为1000, 为在步骤2中获取的频率因子,Φi(k)为标准化的特征向量,xi=(xi,yi,zi)(1≤i≤n)为重建后的顶点几何坐标,这里n为网格顶点数。根据在步骤1中的定义,网格顶点若为梯田型噪声顶点,顶点标签ωi被设置为1(ωi=1),则在本步骤中仅有前200(0.2m,m=1000)个频率因子被采用重建网格梯田型噪声顶点。采用前200个频率因子重建梯田型噪声网格顶点的原因是,在实验的过程中发现,在本发明所采用的所有实施例上,当采用前200个频率因子重建网格时都会得到一个高度光顺的不包括任何高频噪声以及梯田型噪声的网格模型,梯田型噪声被完全抑制。而在网格非梯田型噪声顶点,由于在步骤1中顶点标签ωi被设置为0(ωi=0),因此在本步骤中前1000(m=1000)个频率因子被采用重建网格顶点,高频噪声被有效消除的同时模型体积得以有效保持。采用本发明的基于非均匀谱图编码的三维低通滤波器对钩骨模型进行网格光顺的结果如图7(c)所示。正是由于采用了这种在梯田型顶点以及非梯田型顶点差异化重建的方法,本发明不仅能消除三维医学模型高频噪声还能有效清除梯田型噪声,同时还能很好地保持模型体积。
除钩骨模型外,本发明在更多的三维医学模型实施例上,包括第三近端指骨模型,髂骨模型,胫骨模型,股骨头模型以及股骨髁状突模型与基于均匀谱图编码的三维低通滤波器作对比实验以验证本方法的有效性。运算时间以及参数设置见表1。实验结果如图8所示。虽然两种方法都能去除高频随机噪声,但是基于均匀谱图编码的三维低通滤波器不可避免地会将梯田型噪声视为模型固有的结构加以重建,如图8(b)所示;而本发明提供的基于非均匀谱图编码的三维低通滤波器不仅能去除高频噪声,还能针对性地去除梯田型噪声,得到高度光顺的医学模型,如图8(c)所示。
本发明在上述六个三维医学模型实施例上与三种体积保持的网格光顺方法APSS方法[G.Guennebaud,M.Gross,Algebraic point set surfaces,in:ACM Transactions onGraphics(TOG),Vol.26,ACM,2007,pp.23.]、Zheng等人的方法[Y.Zheng,H.Fu,O.C.Au,C.L.Tai,Bilateral normal filtering for mesh denoising,Visualization andComputer Graphics,IEEE Transactions on 17(10),2011,pp.1521–1530.]以及Taubin的方法[G.Taubin,A signal processing approach to fair surface design,in:Proceedings of the 22nd annual conference on Computer graphics andinteractive techniques,ACM,1995,pp.351–358.]作对比实验以验证本发明在进行网格光顺处理时体积保持的有效性。实验结果如图9所示。模型从上到下,采用APSS方法做光顺处理时,参数设置分别为5,4,5,5,5和9,采用Taubin的方法做光顺处理时,迭代次数设置分别为100,16,200,200,120和100。针对每一个模型,本发明记录模型体积保持率,实验结果见表2以及图10。
表2本发明与其他三种网格光顺方法体积保持对比实验结果
APSS网格光顺方法可以很好地保持模型体积,但是在光顺髂骨模型时,模型被轻微放大,Zheng等人的方法亦可以很好地保持体积,如表2中的数据所示。注意到,作为特征保持的网格光顺算法,APSS网格光顺方法(如图9(b)所示)以及Zheng等人的方法(如图9(c)所示)都仅能去除模型高频噪声,而无法消除三维医学模型梯田型噪声。原因是两个方法都会将原始三维医学网格模型的梯田型噪声识别为模型的特征加以保持甚至增强而非去除。Taubin的方法在光顺同时亦具有保持模型体积的特性,但是其仅能去除较为缓和的梯田型噪声,对于显著的尖锐的梯田型噪声不能很好地去除,如图9(d)所示。而本发明提供的基于非均匀谱图编码的三维低通滤波器能够在梯田型顶点以及非梯田型顶点进行差异性几何信息重建。结果是,本发明不仅能去除高频噪声,还能针对性地去除梯田型噪声,得到高度光顺的三维医学模型,如图9(e)所示。同时,表2中的数据表明,在对钩骨模型,第三近端指骨模型,髂骨模型,胫骨模型,股骨髁状突模型等做光顺处理后,模型的体积保持率均在98%以上,股骨头模型的体积保持率为97%以上。本发明提供的基于非均匀谱图编码的三维低通滤波方法在能有效去除三维医学模型高频噪声以及梯田型噪声的同时,还能很好地保持模型的体积。
最后应说明的是:以上所述仅为本发明的优选实施例而己,并不用于限制本发明,尽管参照上述实施例对本发明进行了详细的说明,对于本领域的技术人员来说,其依然可以对实施例所记载的技术方案进行修改,或者对其中部分技术特征进行等同替换。凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。
Claims (1)
1.一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法,其特征在于,包括如下步骤:(1)构建特定方向感知的特征检测方法,准确识别梯田型噪声;(2)建立三维医学网格模型离散拉普拉斯-贝尔特拉米算子,执行谱图分析,构建谱图空间;(3)构建基于非均匀谱图编码的三维低通滤波器去除高频随机噪声以及梯田型噪声,获得光顺的三维医学模型;其中
所述步骤(1)为构建特定方向感知的特征检测方法以准确识别沿z轴方向分布的梯田型噪声:首先,构建笛卡尔坐标系;接着,定义方向感知的特征检测方法,公式如下:
其中,k为网格顶点i一邻域邻接三角形,为一邻域邻接三角形面法向量;计算顶点i一邻域邻接三角形面法向量与x轴之间的角度差,记为找出网格顶点i一邻域三角形沿x轴方向的最大偏离以及最小偏离计算沿x轴方向的最大偏离以及最小偏离之间的差,记为δx,i;同理,在y轴以及z轴方向也执行上述相同的操作,得到δy,i与δz,i;
所述步骤(1)设置一个合理的阈值τ,将δx,i>τ的顶点标识为x轴方向特征点,将δy,i>τ的顶点标识为y轴方向特征点,将δz,i>τ的顶点标识为z轴方向特征点;由于三维医学模型的梯田型噪声仅沿z轴方向分布,因此,若三维医学模型网格顶点i的δz,i>τ,则顶点i为梯田型噪声点;对于梯田型噪声点,为其设置标签ωi=1,非梯田型噪声点,为其设置标签ωi=0;
所述步骤(2)为构建离散拉普拉斯-贝尔特拉米算子,执行三维医学模型谱图分析,构建谱图空间:首先,构建离散拉普拉斯-贝尔特拉米操作算子Δ,在网格曲面上,离散拉普拉斯-贝尔特拉米操作算子Δ定义为:
wij=(cotαij+cotβij)/2Ai
其中,i与j为网格顶点;设V为顶点集合,E为边集,N(i)={j∈V|(i,j)∈E}是顶点i的邻接顶点;wij=(cotαij+cotβij)/2Ai是边(i,j)的权重,αij与βij是共享同一条边(i,j)的两个邻接三角形的对角,Ai为顶点i一邻域Vorionor面积;将网格模型离散拉普拉斯-贝尔特拉米算子写成矩阵L的形式如下:
其中,i与j为网格顶点;设V为顶点集合,E为边集,N(i)={j∈V|(i,j)∈E}是顶点i的邻接点;wij=(cotαij+cotβij)/2Ai为边(i,j)的权重;矩阵L特征分解的结果是产生一系列特征值λi与特征向量fi,这里1≤i≤n,n为网格顶点数,特征值与特征向量成对出现(λi,fi),特征向量两两正交;将模型顶点几何信息看作是信号,将其投射到正交的特征向量构造谱图空间,首先需要标准化特征向量:
Φi=fi/‖fi‖
其中,标准化的特征向量Φi会构建一个矩阵,矩阵的第i行提供了网格顶点i空间几何坐标(xi,yi,zi)的一个嵌入,Φi的第k个元素Φi(k)为网格顶点i的分段线性函数;
接着,利用所得到的标准化的特征向量Φi构造网格模型顶点的谱图空间,公式如下:
所述步骤(3)为构建基于非均匀谱图编码的三维低通滤波器,有效去除高频噪声以及梯田型噪声获得光顺的三维医学模型:
m=m*(1-0.8ωi),
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201810339104.3A CN108776956B (zh) | 2018-04-16 | 2018-04-16 | 一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201810339104.3A CN108776956B (zh) | 2018-04-16 | 2018-04-16 | 一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN108776956A CN108776956A (zh) | 2018-11-09 |
CN108776956B true CN108776956B (zh) | 2023-02-28 |
Family
ID=64033750
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201810339104.3A Active CN108776956B (zh) | 2018-04-16 | 2018-04-16 | 一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN108776956B (zh) |
Citations (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2010073251A2 (en) * | 2008-12-25 | 2010-07-01 | Medic Vision - Brain Technologies Ltd. | Denoising medical images |
Family Cites Families (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US7502499B2 (en) * | 2004-11-02 | 2009-03-10 | Siemens Medical Solutions Usa, Inc. | System and method for filtering noise from a medical image |
US9891072B2 (en) * | 2014-12-08 | 2018-02-13 | Here Global B.V. | Method and apparatus for providing a map display based on velocity information |
-
2018
- 2018-04-16 CN CN201810339104.3A patent/CN108776956B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2010073251A2 (en) * | 2008-12-25 | 2010-07-01 | Medic Vision - Brain Technologies Ltd. | Denoising medical images |
CN102203826A (zh) * | 2008-12-25 | 2011-09-28 | 梅迪奇视觉成像解决方案有限公司 | 医学图像的降噪 |
Non-Patent Citations (2)
Title |
---|
Comparison of fundamental mesh smoothing algorithms for medical surface models;Bade R,et al.;《In: SimVis》;20061231;全文 * |
基于局部特征概率密度估计的三维模型特征提取方法;孙挺 等;《计算机科学》;20151231;全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN108776956A (zh) | 2018-11-09 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Han et al. | A review of algorithms for filtering the 3D point cloud | |
US6813373B1 (en) | Image segmentation of embedded shapes using constrained morphing | |
CN109285225B (zh) | 一种基于医学影像的虚拟现实辅助手术的建立方法 | |
Cerrolaza et al. | Automatic multi-resolution shape modeling of multi-organ structures | |
CN110956698B (zh) | 一种在血管数字模型中生成肿瘤的方法、系统及电子设备 | |
Ha et al. | Low dose CT image restoration using a database of image patches | |
US20230169732A1 (en) | A Method and System for Enforcing Smoothness Constraints on Surface Meshes from a Graph Convolutional Neural Network | |
Zhong et al. | 4D cone-beam CT reconstruction using multi-organ meshes for sliding motion modeling | |
Wei et al. | Morphology-preserving smoothing on polygonized isosurfaces of inhomogeneous binary volumes | |
Duan et al. | A subdivision-based deformable model for surface reconstruction of unknown topology | |
CN108776956B (zh) | 一种基于非均匀谱图编码三维低通滤波器的医学模型光顺的方法 | |
Muraki et al. | A survey of medical applications of 3D image analysis and computer graphics | |
Zifan et al. | The use of the Kalman filter in the automated segmentation of EIT lung images | |
Zhang et al. | Skeleton Based As-Rigid-As-Possible Volume Modeling. | |
CN112862975B (zh) | 骨骼数据处理方法、系统、可读存储介质和设备 | |
Leonardi et al. | 3D reconstruction from CT-scan volume dataset application to kidney modeling | |
CN114119928A (zh) | 一种基于网格操作的肺内器官三维模型优化方法及系统 | |
Roerdink | Mathematical morphology in computer graphics, scientific visualization and visual exploration | |
Hadwiger et al. | State of the art report 2004 on GPU-based segmentation | |
Feng et al. | Guided normal filter for 3D point clouds | |
Audette et al. | A review of mesh generation for medical simulators | |
Hamzah et al. | 3D model visualization for brain tumour using MRI images | |
Boissonnat et al. | From segmented images to good quality meshes using delaunay refinement | |
Guo et al. | 3D medical model low-pass filtering based on non-uniform spectral synthesis | |
Wang et al. | Three-dimensional surface mesh optimization and centerline extraction of vasculatures for endovascular intervention simulation |
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 |