CN110610478A - 一种基于邻域拓扑的医学图像三维重建方法 - Google Patents
一种基于邻域拓扑的医学图像三维重建方法 Download PDFInfo
- Publication number
- CN110610478A CN110610478A CN201910671883.1A CN201910671883A CN110610478A CN 110610478 A CN110610478 A CN 110610478A CN 201910671883 A CN201910671883 A CN 201910671883A CN 110610478 A CN110610478 A CN 110610478A
- Authority
- CN
- China
- Prior art keywords
- cube
- topology
- cubes
- algorithm
- medical image
- 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
- 238000000034 method Methods 0.000 title claims abstract description 30
- 238000012545 processing Methods 0.000 claims abstract description 35
- 230000008859 change Effects 0.000 claims abstract description 21
- 238000009826 distribution Methods 0.000 claims abstract description 13
- 230000007704 transition Effects 0.000 claims abstract description 11
- 238000003708 edge detection Methods 0.000 claims description 3
- IYLGZMTXKJYONK-ACLXAEORSA-N (12s,15r)-15-hydroxy-11,16-dioxo-15,20-dihydrosenecionan-12-yl acetate Chemical compound O1C(=O)[C@](CC)(O)C[C@@H](C)[C@](C)(OC(C)=O)C(=O)OCC2=CCN3[C@H]2[C@H]1CC3 IYLGZMTXKJYONK-ACLXAEORSA-N 0.000 abstract 2
- IYLGZMTXKJYONK-UHFFFAOYSA-N ruwenine Natural products O1C(=O)C(CC)(O)CC(C)C(C)(OC(C)=O)C(=O)OCC2=CCN3C2C1CC3 IYLGZMTXKJYONK-UHFFFAOYSA-N 0.000 abstract 2
- 238000002591 computed tomography Methods 0.000 description 11
- 238000010586 diagram Methods 0.000 description 7
- 230000008569 process Effects 0.000 description 5
- 230000000694 effects Effects 0.000 description 4
- 230000008439 repair process Effects 0.000 description 4
- 238000004026 adhesive bonding Methods 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 2
- 210000000988 bone and bone Anatomy 0.000 description 2
- 238000003745 diagnosis Methods 0.000 description 2
- 238000003384 imaging method Methods 0.000 description 2
- 230000006872 improvement Effects 0.000 description 2
- 238000004898 kneading Methods 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 210000001519 tissue Anatomy 0.000 description 2
- 238000012285 ultrasound imaging Methods 0.000 description 2
- 208000006735 Periostitis Diseases 0.000 description 1
- 239000000654 additive Substances 0.000 description 1
- 230000000996 additive effect Effects 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 238000005520 cutting process Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000002059 diagnostic imaging Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 238000002595 magnetic resonance imaging Methods 0.000 description 1
- 238000010603 microCT Methods 0.000 description 1
- 238000009206 nuclear medicine Methods 0.000 description 1
- 210000000056 organ Anatomy 0.000 description 1
- 210000003460 periosteum Anatomy 0.000 description 1
- 238000002600 positron emission tomography Methods 0.000 description 1
- 238000007670 refining Methods 0.000 description 1
- 238000009877 rendering Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000002603 single-photon emission computed tomography Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
- 238000000844 transformation Methods 0.000 description 1
- 210000000623 ulna Anatomy 0.000 description 1
- 239000013598 vector Substances 0.000 description 1
- 238000012800 visualization Methods 0.000 description 1
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
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
- G06T7/0012—Biomedical image inspection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/13—Edge detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/136—Segmentation; Edge detection involving thresholding
-
- 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
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Medical Informatics (AREA)
- Quality & Reliability (AREA)
- Radiology & Medical Imaging (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Health & Medical Sciences (AREA)
- General Health & Medical Sciences (AREA)
- Computer Graphics (AREA)
- Geometry (AREA)
- Software Systems (AREA)
- Image Processing (AREA)
Abstract
本发明公开了一种基于邻域拓扑的医学图像三维重建方法。载入原始医学图像序列，获得待标定图像等值线及其等值线图像，构建四个概率参数进而获得目标阈值。读取像素信息构建立方体集，判断立方体顶点的灰度值与目标阈值的大小关系，按顶点分布特征预分类为拓扑相似和拓扑变化两类；拓扑相似类，其邻域内若不存在拓扑变化则采用移动立方体算法，若存在拓扑变化则采用共有面拓扑一致处理的过渡算法；拓扑变化类，则采用按拓扑指定方式和拓扑相异方式划分的移动四面体算法。本发明使移动立方体算法和移动四面体算法融合起来，既避免了等值线拓扑连接的二义性，保证了精度，又显著减少计算量，提升医学图像序列的在线三维重建速度。
Description
技术领域
本发明属于医学与工程科学结合领域，涉及医学图像处理以及三维重建过程中几何结构模型的处理方法。
背景技术
医学图像序列通常是指描述医学对象的切片化的相互平行的图像序列，获得方法可以为电子计算机断层扫描(Computed Tomography,CT)、超声波成像(UltrasonicImaging,UI)、核磁共振成像(Magnetic Resonance Imaging,MRI)、核医学成像(NuclearMedicine Imaging,NMI)、多普勒血流成像(Doppler Flow Imaging,DFI)、正电子辐射断层摄影(Positron Emission Tomography,PET)、单光子发射计算机断层成像(Single PhotonEmission Computed Tomography,SPECT)等。从二维医学图像序列重建三维医学组织，其重建精度与准确性直接决定了医学分析诊断的效果。而重建的三维几何模型作为数据驱动源头，可进一步进行后续个性化定制治疗。当前医学图像三维重建的核心研究聚焦在提高重建精度与准确性，同时提高重建效率，实现远程在线快速诊断。
针对此问题，有学者提出了用鞍点保证拓扑正确的快速等值面提取技术，它利用线性插值函数求取鞍点，并利用鞍点的特征解决拓扑结构不相同的问题，但它将每一个立方体都划分为12个以上的四面体，增加了计算量；还有学者提出了基于移动四面体的曲面重建方法，它基于四面体内等值线的唯一性，使曲面的细节得到加强，但在处理大数据点云时要耗费更多的时间；还有一种基于统一计算设备架构(Compute Unified DeviceArchitecture,CUDA)的改进移动立方体法(Marching Cube，MC)算法，但无法直接读入太大的数据量。美国通用电气公司的Lorenson等(Lorensen W，Cline H A.High resolution 3Dsurface on construction algorithm.ACM Computer Graphics，1987，21(4):163-170)提出了面绘制的代表算法MC，它把相邻切片的8个像素点构成一个体素，在每个体素内用三角面片拟合等值面，虽然它存在二义性等问题，但是它可以较为精确地拟合目标等值面。Doi在《IEICE Transactions》(1991,E74(1):214-224.)发表论文“An Efficient Method ofTriangulating Equi-Valued Surfaces by Using Tetrahedral Cells”，提出移动四面体法，该算法首先将立方体剖分成四面体，然后在其中构造等值面进而重建三维表面。S.P.Liou等人在《Computer Vision》(1997,46(3):237-252.)发表论文“An Approach toThree Dimensional Image Segmentation”，提出一种使用三维四次多项式函数来估计图像的局部三维灰度分布，并以此为基础进行三维分割的方法。美国IBM公司的Gyezuec等(Gyezuec A,Taubin G Lazarus F,et a1.Converting sets of polygons to manifoldsurfaces by cutting and stitching//Proceedings of the Conference onVisualization'98,1998:383-390.)提出了将顶点和边分割再缝合的修复算法，通过复制非流形边和顶点，再将这些复制按某种策略粘合的方式来修复非流形，作为最具代表性的修复算法之一，该算法在简单的网格与接近原模型之间做出了很好的权衡，并且根据粘合策略的不同可以产生具有不同几何观感的网格，但是对于某些特例的修复效果并不好。Sidky等人在《Journal of X-ray science and technology》(2006,14(2):119-139.)发表论文“Accurate image reconstruction from few-views and limited-angle data indivergent-beam CT”，提出了可用于少投影和有限角投影CT重建的凸集投影自适应最速下降(ASD-POCS)法。Johnston等人在《SPIE Medical Imaging》(2010:3684-3689)发表论文“GPU-based iterative reconstruction with total variation minimization formicro-CT”，提出采用图形处理单元(GPU)进行辅助计算。Hu等人在《Biomedicalengineering online》(2014,13(1):70.)发表论文“Accurate image reconstructionfrom few-views and limited-angle data in divergent-beam CT”，提出采用预扫描CT图像或解析重建图像作为迭代初始条件，提高了约束重建的速度。
发明内容
为了解决背景技术中存在的问题，本发明提出了一种基于邻域拓扑的医学图像三维重建方法，该方法利用共有面拓扑一致处理的过渡算法将效率较高、精度较低的移动立方体算法与精度较高、效率较低的移动四面体算法融合起来，从而在保证较高精度的情况下，减少医学图像三维重建所需时间，提高三维重建的效率。
本发明所采用的技术方案是：
第一步：采集并载入三维重建医学对象的原始医学图像帧序列，根据三维重建医学对象属性，初始设定与三维重建医学对象灰度值接近的若干不同相近的阈值；
第二步：选取原始医学图像帧序列中待标定CT图像，CT图像包含有完整的三维重建医学对象的轮廓，根据步骤一设定的阈值利用移动正方形算法(Marching Squares,MS)对CT图像处理获得等值线及其等值线图像，然后对比不同阈值下的等值线，得到最优阈值作为目标阈值δ；
所述的三维重建医学对象为医学图像中的医学上的生物组织或者器官。所述的医学图像为骨骼的序列化医学图像。
第三步：读取原始医学图像帧序列中每层图像的像素信息，以原始医学图像帧序列中最底层图像的左下角为原点O建立直角坐标系OXYZ，即全局坐标系GCS(GlobalCoordinate System)，X轴正方向沿图像右方向建立，Y轴正方向沿图像上方向建立，Z轴遵循右手坐标系建立，且像素间隔设定为1，由原始医学图像帧序列相邻层图像中紧邻排列且位于正方体顶点上的八个像素点构成一个立方体，由所有的立方体构建立方体集C，且所有立方体均建立在X、Y、Z三轴的正方向，初始化k＝0；
第四步：对立方体集C中的立方体采用以下方式进行处理：以立方体中X，Y，Z三轴坐标均为最小值的点为标记点，遍历标记点的Z轴坐标值为k的立方体，判断立方体8个顶点的灰度值与目标阈值δ的大小关系，按照顶点分布特征将所有立方体划分到拓扑相似的立方体集C1和拓扑变化的立方体集C2；
第五步：对于两种立方体集中的每个立方体，采用以下方式进行处理：
对于立方体集C2中每个立方体，对标记点的X、Y、Z坐标值相加为奇数的立方体使用按拓扑指定方式划分的移动四面体算法(MT，Marching Tetrahedral)进行处理；
对于立方体集C2中每个立方体，对标记点X，Y，Z坐标值相加为偶数的立方体使用按与拓扑指定方式相异的拓扑相异方式划分的移动四面体算法进行处理；
对立方体集C1中不与立方体集C2中任何一个立方体相邻的立方体，使用移动立方体算法(MC，Marching Cubes)进行处理；
对于立方体集C1中与立方体集C2中立方体相邻的立方体，使用共有面拓扑一致处理的过渡算法进行处理；
经过处理后的立方体内部的三角面片形成了网络模型，将k加1；
第六步：重复第四步至第五步，直到k＝M-1，使得原始医学图像帧序列中的每层图像均得到处理，M表示原始医学图像帧序列中的图像总数量，最后通过各立方体之间邻域拓扑的一致性消除立方体集等值面的二义性，得到完整的医学图像三维重建模型。
所述第二步中，对比不同阈值下的等值线，其具体步骤如下：
第2.1步：利用Canny边缘检测算子提取目标边缘，将其作为参考边缘图像；
第2.2步：将使用移动正方形算法所获得的等值线图像与参考边缘图像作对比，得到不同阈值下的精确率(P)和召回率(R)：
其中，为第一概率参数，为第二概率参数，为第三概率参数，为第四概率参数；上述四个概率参数计算如下：
其中，W和H分别为图像的宽和高，R1和R0分别表示参考边缘图像中的边缘点和非边缘点，T1和T0分别表示等值线图像中的边缘点和非边缘点；
第2.3步：计算出不同阈值对应的F-Score：
其中，F-Score表示为等值线相似度参数；
第2.4步：选取等值线相似度参数F-Score最大时所对应的阈值作为目标阈值δ。
所述第四步中，按照顶点分布特征判断的分类方法为：
不存在任何一个顶点灰度值或仅存在一个顶点灰度值小于等于目标阈值δ的立方体归属于拓扑相似立方体集C1；
对于存在有且仅有两个顶点P1、P2灰度值小于等于目标阈值δ的立方体，两个顶点P1、P2位于立方体同一面上的对角处，则立方体归属于拓扑变化立方体集C2，否则归属于拓扑相似立方体集C1；
对于存在有且仅有三个顶点P1、P2、P3灰度值小于等于目标阈值δ的立方体，三个顶点P1、P2、P3位于立方体的同一面上，则立方体归属于拓扑相似立方体集C1，否则归属于拓扑变化立方体集C2；
对于存在有且仅有四个顶点P1、P2、P3、P4灰度值小于等于目标阈值δ的立方体，四个顶点P1、P2、P3、P4中的三个点位于立方体的同一面上且另一个点和前面三个点其中任意一点共边，或者四个顶点P1、P2、P3、P4位于立方体的同一面上，则立方体归属于拓扑相似立方体集C1，否则归属于拓扑变化立方体集C2。
若立方体中存在四个以上顶点灰度值小于等于目标阈值δ，按照上述方式反推判断存在几个顶点灰度值小于等于目标阈值δ，即8个顶点状态置反，即将灰度值大于目标阈值δ的顶点标记为灰度值小于等于目标阈值δ的顶点，将灰度值小于等于目标阈值δ的顶点标记为灰度值大于目标阈值δ的顶点，然后按前4种情况处理。
所述的移动四面体算法将立方体划分为5个四面体，且具有两种划分方式。
所述的拓扑指定方式具体为：以与标记点共边的三个顶点以及与标记点立体对角的顶点的四个点为顶点划分为中心四面体，以中心四面体将立方体中除了中心四面体以外的其余部分划分获得其他四个四面体，从而共划分获得五个四面体；
所述的拓扑相异方式具体为：以与标记点共面且处于对角的三个顶点以及标记点自身的四个点为顶点划分为中心四面体，以中心四面体将立方体中除了中心四面体以外的其余部分划分获得其他四个四面体，从而共划分获得五个四面体；
所述的共有面拓扑一致处理的过渡算法是在移动立方体算法基础上进行改变，具体是对于立方体集C1中立方体和立方体集C2中立方体的共有面进行特殊处理，在共有面上额外添加一条与立方体集C2中该立方体的共有面上的划分对角线相同的对角线，将额外添加的该对角线上的等值点加入到移动立方体算法中与其他等值点一起进行处理，获得立方体中的三角面片。
本发明的有益效果是：
1、本发明判断立方体顶点的灰度值与目标阈值的大小关系，按顶点分布特征预分类为拓扑相似立方体集和拓扑变化立方体集，提升了医学图像序列的在线三维重建速度，有助于提高医学定制骨骼增材制造的生产效率。
2、本发明针对了医学图像三维重构阈值选取的难点，构建四个概率参数，将F-score最大时的阈值作为目标阈值，提升了医学图像三维重建对于不同医学对象的适配性。
3、本发明使移动立方体算法和移动四面体算法融合起来，既避免了移动立方体算法的等值面拓扑连接的二义性，同等于移动四面体算法的精度，又显著减少计算量，有助于实现高精度的远程在线快速医学分析。
附图说明
图1是本发明方法流程总图。
图2是不同阈值下的F-Score和精确率P示意图。
图3是本发明的最底层图像上坐标系的示意图。
图4是本发明中拓扑相似立方体集C1中立方体的顶点分布特征。
图5是本发明中拓扑变化立方体集C2中立方体的顶点分布特征。
图6是本发明中移动四面体算法所采用的拓扑指定方式示意图。
图7是本发明中移动四面体算法所采用的拓扑相异方式示意图。
图8是本发明中拓扑相似立方体集C1中一个立方体使用按拓扑指定方式划分的移动四面体算法时的示意图。
图9是本发明中拓扑相似立方体集C1中一个立方体使用按拓扑相异方式划分的移动四面体算法时的示意图。
图10是本发明中拓扑变化立方体集C2中一个立方体使用按拓扑指定方式划分的移动四面体算法时的示意图。
图11是本发明中拓扑变化立方体集C2中一个立方体使用按拓扑相异方式划分的移动四面体算法时的示意图。
图12是本发明中共有面拓扑一致处理的过渡算法的示意图。
图13是本发明的采用F-score选取的阈值重建的三维桡骨茎突和尺骨茎突。
图14是本发明的未采用F-score选取的阈值重建的三维桡骨茎突和尺骨茎突。
图15是本发明与移动四面体算法所耗时间对比图。
具体实施方式
下面结合附图和实施例对本发明作进一步说明。
本发明具体实施例如下：
第一步：以人体桡骨茎突和尺骨茎突作为三维重建医学对象，初步设定阈值为骨膜的灰度值及与该灰度值相近的几个阈值，选取130HU、135HU、140HU、145HU和150HU作为初步设定的阈值；
第二步：采集并载入三维重建医学对象的原始骨骼医学图像帧序列，选取121层CT图像，选取原始骨骼医学图像帧序列中待标定CT图像，CT图像包含有完整的三维重建医学对象的轮廓，根据步骤一设定的阈值利用移动正方形算法(Marching Squares,MS)对CT图像处理获得等值线及其等值线图像，然后对比不同阈值下的等值线，具体如下：
第2.1步：利用Canny边缘检测算子提取桡骨和尺骨的边缘，将其作为参考边缘图像；
第2.2步：将使用移动正方形算法所获得的等值线图像与参考边缘图像作对比，得到不同阈值下的精确率(P)和召回率(R)：
其中，为第一概率参数，为第二概率参数，为第三概率参数，为第四概率参数；上述四个概率参数计算如下：
其中，W和H分别为图像的宽和高，R1和R0分别表示参考边缘图像中的边缘点和非边缘点，T1和T0分别表示等值线图像中的边缘点和非边缘点；
第2.3步：计算出不同阈值对应的F-Score：
其中，F-Score表示为等值线相似度参数；
经计算得出如下数据：
| 阈值(HU) | 130 | 135 | 140 | 145 | 150 |
| 精确率P | 0.727 | 0.742 | 0.757 | 0.759 | 0.755 |
| F-Score(％) | 2.66 | 2.68 | 2.72 | 2.75 | 2.71 |
第2.4步：选取等值线相似度参数F-Score最大时所对应的阈值145HU作为目标阈值δ。
第三步：读取原始骨骼医学图像帧序列中每层图像的像素信息，以原始骨骼医学图像帧序列中最底层图像的左下角为原点O建立直角坐标系OXYZ，即全局坐标系GCS(Global Coordinate System)，X、Y轴方向如图3所示，Z轴遵循右手坐标系建立，且像素间隔设定为1，即相邻像素之间的坐标值之差为1，每个像素视为点无坐标长度，由原始骨骼医学图像帧序列相邻层图像中紧邻排列且位于正方体顶点上的八个像素点构成一个立方体，由所有的立方体构建立方体集C，且所有立方体均建立在X、Y、Z三轴的正方向，初始化k＝0；
第四步：对立方体集C中的立方体采用以下方式进行处理：以立方体中X，Y，Z三轴坐标均为最小值的点为标记点，遍历标记点的Z轴坐标值为k的立方体，判断立方体8个顶点的灰度值与目标阈值δ145HU的大小关系，按照顶点分布特征将所有立方体划分到拓扑相似的立方体集C1和拓扑变化的立方体集C2，具体如下：
不存在任何一个顶点灰度值或仅存在一个顶点灰度值小于等于目标阈值δ的立方体归属于拓扑相似立方体集C1；
对于存在有且仅有两个顶点P1、P2灰度值小于等于目标阈值δ的立方体，两个顶点P1、P2位于立方体同一面上的对角处，则立方体归属于拓扑变化立方体集C2，否则归属于拓扑相似立方体集C1；
对于存在有且仅有三个顶点P1、P2、P3灰度值小于等于目标阈值δ的立方体，三个顶点P1、P2、P3位于立方体的同一面上，则立方体归属于拓扑相似立方体集C1，否则归属于拓扑变化立方体集C2；
对于存在有且仅有四个顶点P1、P2、P3、P4灰度值小于等于目标阈值δ的立方体，四个顶点P1、P2、P3、P4中的三个点位于立方体的同一面上且另一个点和前面三个点其中任意一点共边，或者四个顶点P1、P2、P3、P4位于立方体的同一面上，则立方体归属于拓扑相似立方体集C1，否则归属于拓扑变化立方体集C2。
若立方体中存在四个以上顶点灰度值小于等于目标阈值δ，按照上述方式反推判断存在几个顶点灰度值小于等于目标阈值δ，即8个顶点状态置反，即将灰度值大于目标阈值δ的顶点标记为灰度值小于等于目标阈值δ的顶点，将灰度值小于等于目标阈值δ的顶点标记为灰度值大于目标阈值δ的顶点，然后按前4种情况处理。
如图4和图5所示，分别为本发明中拓扑相似立方体集C1中立方体的顶点分布特征和拓扑变化立方体集C2中立方体的顶点分布特征。其中球状顶点处灰度值小于等于目标阈值δ，其余顶点处灰度值大于目标阈值δ。若立方体经过旋转和互补处理(所有顶点状态置反)后符合图4所示顶点分布特征则为立方体集C1内立方体；若立方体经过旋转和互补处理后符合图5所示顶点分布特征则为立方体集C2内立方体。
第五步：对于两种立方体集中的每个立方体，采用以下方式进行处理：
对于立方体集C2中每个立方体，对标记点的X、Y、Z坐标值相加为奇数的立方体使用按拓扑指定方式划分的移动四面体算法(MT，Marching Tetrahedral)进行处理；
拓扑指定方式具体为：以与标记点共边的三个顶点以及与标记点立体对角的顶点的四个点为顶点划分为中心四面体，以中心四面体将立方体中除了中心四面体以外的其余部分划分获得其他四个四面体，从而共划分获得五个四面体。如图6所示，实线为可见线，虚线为不可见线，图中的P6为标记点，划分的中心四面体如图6所示，图6中按拓扑指定方式将立方体P1P5P7P3-P2P6P8P4划分为四面体P1-P5P2P3、四面体P7-P5P3P8、四面体P6-P5P2P8、四面体P4-P2P3P8和四面体P5-P2P3P8。
对于立方体集C2中每个立方体，对标记点X，Y，Z坐标值相加为偶数的立方体使用按与拓扑指定方式相异的拓扑相异方式划分的移动四面体算法进行处理；
拓扑相异方式具体为：以与标记点共面且处于对角的三个顶点以及标记点自身的四个点为顶点划分为中心四面体，以中心四面体将立方体中除了中心四面体以外的其余部分划分获得其他四个四面体，从而共划分获得五个四面体。如图7所示，图中的P6为标记点，划分的中心四面体如图7所示，图7中按拓扑相异方式将立方体P1P5P7P3-P2P6P8P4划分为四面体P5-P1P6P7、四面体P3-P1P4P7、四面体P2-P1P4P6、四面体P8-P4P6P7和四面体P1-P4P6P7。
通过按照上述两种划分方式划分的移动四面体算法处理后，可得出拓扑相似立方体集C1和拓扑变化立方体集C2之间的重要区别：
对立方体集C1中的立方体使用按不同划分方式划分的移动四面体算法后立方体的6个面上的等值线不变、变为折线或变为线段。
如图8和图9所示，该立方体为拓扑相似立方体集C1中一个立方体，其顶点P3、P4和P6处灰度值小于等于目标阈值δ，其余顶点处灰度值大于目标阈值δ。对图8中立方体使用按拓扑指定方式划分的移动四面体算法后面P5P6P2P1和面P7P8P4P3内等值线分别与图9中使用按拓扑相异方式划分的移动四面体算法后面P5P6P2P1和面P7P8P4P3内等值线完全相同，分别为线段P9P16和线段P12P13。图8中面P1P5P7P3内等值线的折线P13P14P15变为图9中面P1P5P7P3内等值线的线段P13P14，图8中面P2P6P8P4内等值线的折线P10P11P12变为图9中面P2P6P8P4内等值线的线段P11P12。图8中面P5P6P8P7内等值线的线段P9P10变为图9中面P5P6P8P7内等值线的折线P9P10P11，图8中面P1P2P4P3内等值线的线段P15P16变为图9中面P1P2P4P3内等值线的折线P14P15P16。
对立方体集C2中的立方体使用按不同划分方式划分的移动四面体算法后立方体的一些面上的等值线完全不同。
如图10和图11所示，该立方体为拓扑变化立方体集C2中一个立方体，其顶点P3、P6和P8处灰度值小于等于目标阈值δ，其余顶点处灰度值大于目标阈值δ。
对图10中立方体使用按拓扑指定方式划分的移动四面体算法后面P7P8P4P3内等值线的线段P9P10和线段P11P12与图11中使用按拓扑相异方式划分的移动四面体算法后面P7P8P4P3内等值线的线段P10P12和线段P9P11完全不同。
对立方体集C1中不与立方体集C2中任何一个立方体相邻的立方体，使用移动立方体算法(MC，Marching Cubes)进行处理；
对于立方体集C1中与立方体集C2中立方体相邻的立方体，使用共有面拓扑一致处理的过渡算法划分为多个多面体；立方体相邻是指两个立方体具有共面，也认为其中一个立方体位于另一立方体的邻域；
共有面拓扑一致处理的过渡算法是在移动立方体算法基础上进行改变，具体是对于立方体集C1中立方体和立方体集C2中立方体的共有面进行特殊处理，在共有面上额外添加一条与立方体集C2中该立方体的共有面上的划分对角线相同的对角线，将额外添加的该对角线上的等值点加入到移动立方体算法中与其他等值点一起进行处理，获得立方体中的三角面片。
如图12中立方体P5P7P3P1-P6P8P4P2为立方体集C2中使用按拓扑指定方式划分的移动四面体算法处理的立方体，立方体P7P11P9P3-P8P12P10P4为立方体集C1中的立方体，面P7P8P4P3为两个立方体的共有面，折线P13P14P15为该共有面上拓扑一致的等值线。对立方体P7P11P9P3-P8P12P10P4使用共有面拓扑一致处理的过渡算法，在共有面上添加与立方体P5P7P3P1-P6P8P4P2相同的划分对角线P3P8，在移动立方体算法的基础上额外计算线段P3P8上的等值点P14。
经过处理后的立方体内部的三角面片形成了网络模型，由此完成实现了最底层立方体内部的三维重建，将k加1；
第六步：重复第四步至第五步，直到k＝M-1，使得原始骨骼医学图像帧序列中的每层图像均得到处理，M表示原始骨骼医学图像帧序列中的图像总数量，最后通过各立方体之间邻域拓扑的一致性消除立方体集等值面的二义性，得到完整的医学图像三维重建模型。实施例结果如图13所示，为实施例经过六个步骤重建的三维桡骨茎突和尺骨茎突，选取的摄像机视角的单位方向向量为(0.43，0.75,-0.5)。若选择150HU为目标阈值进行重建，则实施例结果如图14所示，该三维重建模型的下方茎突有一些缺陷，出现了明显的孔洞，说明根据相似度参数F-score最大时选取的阈值更优，重建效果更好。
经统计，立方体集C中共有31334520个立方体，其中使用移动立方体算法的立方体共31331700个，占99.9910％，使用共有面拓扑一致处理的过渡算法的立方体共2313个，占0.0074％，使用移动四面体算法的立方体共507个，占0.0016％。其中使用移动四面体算法的少量立方体单位耗时较多，其余立方体则单位耗时较少。
在相同的硬件计算平台上，采用本发明方法，耗时为2461.35毫秒(ms)，采用移动四面体方法，耗时为3963.33毫秒(ms)。经比较，本发明中的医学图像三维重建方法比同等精度的移动四面体算法节约了38％的计算时间，如图15所示，效率有明显提升。随着重建精度要求的提升，及重建医学对象复杂程度的提升，需要载入巨大数目的医学图像数据，且为了达到最好效果需多次重建，采用本发明中的医学图像三维重建方法可以显著节约重建时间。
Claims (5)
1.一种基于邻域拓扑的医学图像三维重建方法，其特征在于：该方法包括以下步骤：
第一步：采集并载入三维重建医学对象的原始医学图像帧序列，根据三维重建医学对象属性，初始设定若干不同阈值；
第二步：选取原始医学图像帧序列中待标定CT图像，CT图像包含有完整的三维重建医学对象的轮廓，根据步骤一设定的阈值利用移动正方形算法(Marching Squares,MS)对CT图像处理获得等值线及其等值线图像，然后对比不同阈值下的等值线，得到最优阈值作为目标阈值δ；
第三步：读取原始医学图像帧序列中每层图像的像素信息，以原始医学图像帧序列中最底层图像的左下角为原点O建立直角坐标系OXYZ，即全局坐标系GCS(Global CoordinateSystem)，X轴正方向沿图像右方向建立，Y轴正方向沿图像上方向建立，Z轴遵循右手坐标系建立，且像素间隔设定为1，由原始医学图像帧序列相邻层图像中紧邻排列且位于正方体顶点上的八个像素点构成一个立方体，由所有的立方体构建立方体集C，且所有立方体均建立在X、Y、Z三轴的正方向，初始化k＝0；
第四步：对立方体集C中的立方体采用以下方式进行处理：以立方体中X，Y，Z三轴坐标均为最小值的点为标记点，遍历标记点的Z轴坐标值为k的立方体，判断立方体8个顶点的灰度值与目标阈值δ的大小关系，按照顶点分布特征将所有立方体划分到拓扑相似的立方体集C1和拓扑变化的立方体集C2；
第五步：对于两种立方体集中的每个立方体，采用以下方式进行处理：
对于立方体集C2中每个立方体，对标记点的X、Y、Z坐标值相加为奇数的立方体使用按拓扑指定方式划分的移动四面体算法(MT，Marching Tetrahedral)进行处理；
对于立方体集C2中每个立方体，对标记点X，Y，Z坐标值相加为偶数的立方体使用按与拓扑指定方式相异的拓扑相异方式划分的移动四面体算法进行处理；
对立方体集C1中不与立方体集C2中任何一个立方体相邻的立方体，使用移动立方体算法(MC，Marching Cubes)进行处理；
对于立方体集C1中与立方体集C2中立方体相邻的立方体，使用共有面拓扑一致处理的过渡算法进行处理；
经过处理后的立方体内部的三角面片形成了网络模型，将k加1；
第六步：重复第四步至第五步，直到k＝M-1，使得原始医学图像帧序列中的每层图像均得到处理，M表示原始医学图像帧序列中的图像总数量，最后通过各立方体之间邻域拓扑的一致性消除立方体集等值面的二义性，得到完整的医学图像三维重建模型。
2.根据权利要求1所述的一种基于邻域拓扑的医学图像三维重建方法，其特征在于：所述第二步中，对比不同阈值下的等值线，其具体步骤如下：
第2.1步：利用Canny边缘检测算子提取目标边缘，将其作为参考边缘图像；
第2.2步：将使用移动正方形算法所获得的等值线图像与参考边缘图像作对比，得到不同阈值下的精确率(P)和召回率(R)：
其中，为第一概率参数，为第二概率参数，为第三概率参数，为第四概率参数；上述四个概率参数计算如下：
其中，W和H分别为图像的宽和高，R1和R0分别表示参考边缘图像中的边缘点和非边缘点，T1和T0分别表示等值线图像中的边缘点和非边缘点；
第2.3步：计算出不同阈值对应的F-Score：
其中，F-Score表示为等值线相似度参数；
第2.4步：选取等值线相似度参数F-Score最大时所对应的阈值作为目标阈值δ。
3.根据权利要求1所述的一种基于邻域拓扑的医学图像三维重建方法，其特征在于：所述第四步中，按照顶点分布特征判断的分类方法为：
不存在任何一个顶点灰度值或仅存在一个顶点灰度值小于等于目标阈值δ的立方体归属于拓扑相似立方体集C1；
对于存在有且仅有两个顶点P1、P2灰度值小于等于目标阈值δ的立方体，两个顶点P1、P2位于立方体同一面上的对角处，则立方体归属于拓扑变化立方体集C2，否则归属于拓扑相似立方体集C1；
对于存在有且仅有三个顶点P1、P2、P3灰度值小于等于目标阈值δ的立方体，三个顶点P1、P2、P3位于立方体的同一面上，则立方体归属于拓扑相似立方体集C1，否则归属于拓扑变化立方体集C2；
对于存在有且仅有四个顶点P1、P2、P3、P4灰度值小于等于目标阈值δ的立方体，四个顶点P1、P2、P3、P4中的三个点位于立方体的同一面上且另一个点和前面三个点其中任意一点共边，或者四个顶点P1、P2、P3、P4位于立方体的同一面上，则立方体归属于拓扑相似立方体集C1，否则归属于拓扑变化立方体集C2。
4.根据权利要求1所述的一种基于邻域拓扑的医学图像三维重建方法，其特征在于：所述的拓扑指定方式具体为：以与标记点共边的三个顶点以及与标记点立体对角的顶点的四个点为顶点划分为中心四面体，以中心四面体将立方体中除了中心四面体以外的其余部分划分获得其他四个四面体，从而共划分获得五个四面体；
所述的拓扑相异方式具体为：以与标记点共面且处于对角的三个顶点以及标记点自身的四个点为顶点划分为中心四面体，以中心四面体将立方体中除了中心四面体以外的其余部分划分获得其他四个四面体，从而共划分获得五个四面体。
5.根据权利要求1所述的一种基于邻域拓扑的医学图像三维重建方法，其特征在于：所述的共有面拓扑一致处理的过渡算法是在移动立方体算法基础上进行改变，具体是对于立方体集C1中立方体和立方体集C2中立方体的共有面进行特殊处理，在共有面上额外添加一条与立方体集C2中该立方体的共有面上的划分对角线相同的对角线，将额外添加的该对角线上的等值点加入到移动立方体算法中与其他等值点一起进行处理，获得立方体中的三角面片。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910671883.1A CN110610478B (zh) | 2019-07-24 | 2019-07-24 | 一种基于邻域拓扑的医学图像三维重建方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910671883.1A CN110610478B (zh) | 2019-07-24 | 2019-07-24 | 一种基于邻域拓扑的医学图像三维重建方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110610478A true CN110610478A (zh) | 2019-12-24 |
| CN110610478B CN110610478B (zh) | 2021-05-25 |
Family
ID=68890223
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910671883.1A Active CN110610478B (zh) | 2019-07-24 | 2019-07-24 | 一种基于邻域拓扑的医学图像三维重建方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110610478B (zh) |
Cited By (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111310736A (zh) * | 2020-03-26 | 2020-06-19 | 上海同岩土木工程科技股份有限公司 | 一种保护区内车辆卸载堆土的快速识别方法 |
| CN111476881A (zh) * | 2020-03-18 | 2020-07-31 | 创业慧康科技股份有限公司 | 基于结构相似性水平集算法的人体组织器官三维重建方法 |
| CN113113090A (zh) * | 2021-03-10 | 2021-07-13 | 长安大学 | 一种dna多面体拓扑结构的构筑方法 |
| CN113838557A (zh) * | 2021-09-15 | 2021-12-24 | 王其景 | 一种医学影像三维重建模拟方法及系统 |
| CN115359049A (zh) * | 2022-10-19 | 2022-11-18 | 首都师范大学 | 基于非线性扩散模型的有限角ct图像重建方法及装置 |
| CN116912427A (zh) * | 2023-09-12 | 2023-10-20 | 武汉工程大学 | 基于标记点三角特征聚类的三维扫描重建方法及系统 |
| CN117333615A (zh) * | 2023-10-09 | 2024-01-02 | 上海交通大学 | 分层扫描数据重构四面体网格模型的方法及系统 |
| CN117934717A (zh) * | 2024-01-25 | 2024-04-26 | 哈尔滨理工大学 | 一种基于二维超声切片影像的胎儿头部三维影像重建系统 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20130163836A1 (en) * | 2011-12-23 | 2013-06-27 | Stmicroelectronics S.R.L. | Computing the mass of an object |
| CN103871107A (zh) * | 2014-02-26 | 2014-06-18 | 浙江大学 | 一种等值面提取精度的定量测量方法 |
| CN107895364A (zh) * | 2017-10-31 | 2018-04-10 | 哈尔滨理工大学 | 一种用于虚拟手术术前规划的三维重建系统 |
| CN109636915A (zh) * | 2018-12-19 | 2019-04-16 | 哈尔滨理工大学 | 一种基于边增长的医学图像三维重建等值面网格提取方法 |
| CN109934922A (zh) * | 2019-03-14 | 2019-06-25 | 哈尔滨理工大学 | 一种基于改进mc算法的三维重建方法 |
-
2019
- 2019-07-24 CN CN201910671883.1A patent/CN110610478B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20130163836A1 (en) * | 2011-12-23 | 2013-06-27 | Stmicroelectronics S.R.L. | Computing the mass of an object |
| CN103871107A (zh) * | 2014-02-26 | 2014-06-18 | 浙江大学 | 一种等值面提取精度的定量测量方法 |
| CN107895364A (zh) * | 2017-10-31 | 2018-04-10 | 哈尔滨理工大学 | 一种用于虚拟手术术前规划的三维重建系统 |
| CN109636915A (zh) * | 2018-12-19 | 2019-04-16 | 哈尔滨理工大学 | 一种基于边增长的医学图像三维重建等值面网格提取方法 |
| CN109934922A (zh) * | 2019-03-14 | 2019-06-25 | 哈尔滨理工大学 | 一种基于改进mc算法的三维重建方法 |
Non-Patent Citations (2)
| Title |
|---|
| 张晔: ""基于改进的MC与MT方法的三维重建的研究与实现"", 《万方学位论文》 * |
| 徐敬华等: ""基于邻接拓扑的流形网格模型层切多连通域构建方法"", 《计算机辅助设计与图形学学报》 * |
Cited By (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111476881A (zh) * | 2020-03-18 | 2020-07-31 | 创业慧康科技股份有限公司 | 基于结构相似性水平集算法的人体组织器官三维重建方法 |
| CN111310736A (zh) * | 2020-03-26 | 2020-06-19 | 上海同岩土木工程科技股份有限公司 | 一种保护区内车辆卸载堆土的快速识别方法 |
| CN113113090A (zh) * | 2021-03-10 | 2021-07-13 | 长安大学 | 一种dna多面体拓扑结构的构筑方法 |
| CN113113090B (zh) * | 2021-03-10 | 2024-01-30 | 长安大学 | 一种dna多面体拓扑结构的构筑方法 |
| CN113838557A (zh) * | 2021-09-15 | 2021-12-24 | 王其景 | 一种医学影像三维重建模拟方法及系统 |
| CN115359049A (zh) * | 2022-10-19 | 2022-11-18 | 首都师范大学 | 基于非线性扩散模型的有限角ct图像重建方法及装置 |
| CN116912427A (zh) * | 2023-09-12 | 2023-10-20 | 武汉工程大学 | 基于标记点三角特征聚类的三维扫描重建方法及系统 |
| CN116912427B (zh) * | 2023-09-12 | 2023-11-24 | 武汉工程大学 | 基于标记点三角特征聚类的三维扫描重建方法及系统 |
| CN117333615A (zh) * | 2023-10-09 | 2024-01-02 | 上海交通大学 | 分层扫描数据重构四面体网格模型的方法及系统 |
| CN117934717A (zh) * | 2024-01-25 | 2024-04-26 | 哈尔滨理工大学 | 一种基于二维超声切片影像的胎儿头部三维影像重建系统 |
| CN117934717B (zh) * | 2024-01-25 | 2024-08-23 | 哈尔滨理工大学 | 一种基于二维超声切片影像的胎儿头部三维影像重建系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110610478B (zh) | 2021-05-25 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110610478B (zh) | 一种基于邻域拓扑的医学图像三维重建方法 | |
| Jones et al. | 3D distance fields: A survey of techniques and applications | |
| US9147239B2 (en) | Computing the mass of an object | |
| US9367958B2 (en) | Method and apparatus for correction of errors in surfaces | |
| CN111696210A (zh) | 一种基于三维点云数据特征轻量化的点云重构方法及系统 | |
| CN104966317B (zh) | 一种基于矿体轮廓线的三维自动建模方法 | |
| CN109191510B (zh) | 一种病理切片的3d重建方法及其装置 | |
| d'Otreppe et al. | Generating smooth surface meshes from multi‐region medical images | |
| CN104318554B (zh) | 基于三角剖分优化匹配的医学影像刚性配准方法 | |
| Liu et al. | High-quality textured 3D shape reconstruction with cascaded fully convolutional networks | |
| CN105825471A (zh) | 一种基于Unity 3D的三维体表面重构与渲染方法 | |
| CN111105476B (zh) | 一种基于Marching Cubes的CT图像三维重建方法 | |
| CN114612318A (zh) | 基于文物ct图像轮廓线的三维建模方法、系统及设备 | |
| Chen et al. | Image vectorization with real-time thin-plate spline | |
| Van Nguyen et al. | Geometric modeling: background for processing the 3d objects | |
| CN109509259B (zh) | 一种医学图像三维重建等值面网格提取方法 | |
| Moriconi et al. | High quality surface reconstruction in radiotherapy: cross-sectional contours to 3D mesh using wavelets | |
| CN108597589B (zh) | 模型生成方法、目标检测方法及医学成像系统 | |
| CN111369662A (zh) | Ct图像中血管的三维模型重建方法及系统 | |
| CN111833448A (zh) | 源自局部体积的腔室重建 | |
| Bertram et al. | Non-manifold mesh extraction from time-varying segmented volumes used for modeling a human heart. | |
| Rahmawati et al. | Modification Rules for Improving Marching Cubes Algorithm to Represent 3D Point Cloud Curve Images. | |
| CN114758073A (zh) | 一种基于rgbd输入和柔性配准口腔数字化系统 | |
| Guo et al. | Three dimension reconstruction of medical images based on an improved marching cubes algorithm | |
| Jensen et al. | Deep Active Latent Surfaces for Medical Geometries |
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 |