CN111612743B - 一种基于ct图像的冠状动脉中心线提取方法 - Google Patents
一种基于ct图像的冠状动脉中心线提取方法 Download PDFInfo
- Publication number
- CN111612743B CN111612743B CN202010332170.5A CN202010332170A CN111612743B CN 111612743 B CN111612743 B CN 111612743B CN 202010332170 A CN202010332170 A CN 202010332170A CN 111612743 B CN111612743 B CN 111612743B
- Authority
- CN
- China
- Prior art keywords
- point
- points
- marked
- template
- voxel
- 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
- 210000004351 coronary vessel Anatomy 0.000 title claims abstract description 45
- 238000000605 extraction Methods 0.000 title claims abstract description 27
- 210000004204 blood vessel Anatomy 0.000 claims abstract description 54
- 238000000034 method Methods 0.000 claims abstract description 33
- 230000002792 vascular Effects 0.000 claims abstract description 19
- 238000007670 refining Methods 0.000 claims abstract description 7
- 238000009792 diffusion process Methods 0.000 claims description 16
- 239000011159 matrix material Substances 0.000 claims description 14
- 230000008569 process Effects 0.000 claims description 13
- 239000013598 vector Substances 0.000 claims description 9
- 238000001914 filtration Methods 0.000 claims description 6
- 238000012545 processing Methods 0.000 claims description 6
- 230000004907 flux Effects 0.000 claims description 4
- 230000001629 suppression Effects 0.000 claims description 3
- 238000012216 screening Methods 0.000 claims description 2
- 239000000284 extract Substances 0.000 abstract description 5
- 238000005457 optimization Methods 0.000 abstract description 3
- 230000000694 effects Effects 0.000 abstract description 2
- 238000012217 deletion Methods 0.000 description 6
- 230000037430 deletion Effects 0.000 description 6
- 208000029078 coronary artery disease Diseases 0.000 description 5
- 238000012795 verification Methods 0.000 description 5
- 208000024172 Cardiovascular disease Diseases 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 238000004445 quantitative analysis Methods 0.000 description 3
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 2
- 210000001367 artery Anatomy 0.000 description 2
- 230000017531 blood circulation Effects 0.000 description 2
- 238000002591 computed tomography Methods 0.000 description 2
- 238000004195 computer-aided diagnosis Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 238000012804 iterative process Methods 0.000 description 2
- 210000000056 organ Anatomy 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 206010003211 Arteriosclerosis coronary artery Diseases 0.000 description 1
- 208000005475 Vascular calcification Diseases 0.000 description 1
- 206010057469 Vascular stenosis Diseases 0.000 description 1
- 210000000709 aorta Anatomy 0.000 description 1
- 210000001992 atrioventricular node Anatomy 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008901 benefit Effects 0.000 description 1
- 239000008280 blood Substances 0.000 description 1
- 210000004369 blood Anatomy 0.000 description 1
- 230000036770 blood supply Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 208000026758 coronary atherosclerosis Diseases 0.000 description 1
- 230000009977 dual effect Effects 0.000 description 1
- 230000002526 effect on cardiovascular system Effects 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 230000005484 gravity Effects 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 230000036541 health Effects 0.000 description 1
- 210000005240 left ventricle Anatomy 0.000 description 1
- 230000000670 limiting effect Effects 0.000 description 1
- 230000002503 metabolic effect Effects 0.000 description 1
- 230000004060 metabolic process Effects 0.000 description 1
- 230000002107 myocardial effect Effects 0.000 description 1
- 235000015097 nutrients Nutrition 0.000 description 1
- 230000035764 nutrition Effects 0.000 description 1
- 235000016709 nutrition Nutrition 0.000 description 1
- 230000008520 organization Effects 0.000 description 1
- 231100000915 pathological change Toxicity 0.000 description 1
- 230000036285 pathological change Effects 0.000 description 1
- 238000013146 percutaneous coronary intervention Methods 0.000 description 1
- 230000035790 physiological processes and functions Effects 0.000 description 1
- 238000013139 quantization Methods 0.000 description 1
- 210000005241 right ventricle Anatomy 0.000 description 1
- 230000000630 rising effect Effects 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 210000001013 sinoatrial node Anatomy 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Images
Classifications
-
- 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
- G06T3/00—Geometric image transformations in the plane of the image
- G06T3/60—Rotation of whole images or parts thereof
-
- 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
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- 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/187—Segmentation; Edge detection involving region growing; involving region merging; involving connected component labelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/60—Analysis of geometric attributes
- G06T7/68—Analysis of geometric attributes of symmetry
-
- 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/20—Special algorithmic details
- G06T2207/20024—Filtering details
-
- 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/30101—Blood vessel; Artery; Vein; Vascular
-
- 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/30172—Centreline of tubular or elongated structure
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)
- Geometry (AREA)
- Image Analysis (AREA)
- Apparatus For Radiation Diagnosis (AREA)
- Image Processing (AREA)
Abstract
本发明公开了一种基于CT图像的冠状动脉中心线提取方法。现有技术由于多种条件限制，不能较完整的、自动的提取出CT图像中的冠状动脉中心线；本发明如下：一、消除被处理的CT图像的散粒噪声。二、对CT图像内的所有像素点进行血管区域增强。三、在CT图像的血管区域内使用区域生长法，对血管进行分割，使得CT图像变为二值图像。四、将CT图像的根据原始数据中的空间分辨率进行三维重构，得到三维空间中的血管数据。五、对血管数据进行细化，得到血管的初步的中心线。六、将中心线转换为有向图，删除初步提取的中心线的细小分支。本发明提高了基于CT图像的冠状动脉中心线的提取精度与效率，优化效果较为明显。
Description
技术领域
本发明属于医学图像处理领域，具体涉及到一种利用计算机断层扫描(ComputedTomography，CT)图像的基于十二方向拓扑细化算法和Dijkstra算法的对冠状动脉中心线进行自动提取的方法。
背景技术
心脏作为人体中，最重要的器官之一，负责向循环系统提供动力，以保证身体各处正常的新陈代谢，而冠状动脉负责为心脏提供所需养料，带走代谢废物，是保证心脏正常生理功能重要器官。近年来，由于人们日益重视维护心脑血管健康，冠心病(又称冠状动脉粥样硬化)作为一种高发的心血管疾病，得到国际社会的广泛关注。《中，国心血管病报告2018》中，指出经推算我国现有冠心病患者1100万人，且死亡率自2012年有明显的上升趋势。
冠状动脉附着于心脏表面，为心脏提供足量的营养与能量，它源起升主动脉根部，向左、向右各发散出左冠状动脉(Left Coronary Artery，LCA)和右冠状动脉(RightCoronary Artery，RCA)，两支动脉紧紧包裹着心脏，形似一顶皇冠，倒置与心脏之上。LCA保证心脏前壁和侧壁的心肌供血，在其分支处10-25mm后分出左前降支(LeftAnteriorDescending Artery，LAD)和回旋支(Left CircumflexArtery，LCX)，前者在左右心室间延伸至心尖，LCA分支之前的部分被称为左主干(Left MainArtery，LM)；RCA负责向右室、右房、左室部分、后壁、窦房结、房室结供血，从中，又发散出锐缘支(Right MarginalArtery，MA)和后降支(Posterior DescendingArtery，PDA)，其中，LCA与RCA为冠状动脉的主要血管，其余的分支都是由它们发散出去。
CTA图像是将血管特征利用计算机重构出来，可以在血管的横截面上提供三维血管的原始层切图像，可根据此图像分析冠状动脉的病变情况。现有技术中，不能自动提取出CT层切图像中的三维血管；近年来基于冠心病的计算机辅助诊断系统(Computer AidedDiagnosis，CAD)成为研究的热点，而冠状动脉中心线的提取又是CAD研究中的一个难点。在CAD中冠状动脉中心线是一个定量分析冠心病等心血管疾病重要参数，是对血管狭窄分析、血管钙化的定量分析、血管斑块类型的判定等多种心血管疾病定量分析的起点。但是，目前现有的冠状动脉中心线提取方法不能很好的适应冠状动脉树上所有分支的形态及走向。因此，获得一种能够在各种情况下精确提取冠状动脉中心线的方法对于CAD来说尤为重要。此外，冠状动脉中心线对多种冠状动脉的可视化重建也有着不可替代的作用，尤其是在曲面重建技术和多平面重建技术中，发挥了巨大的作用。同时，冠状动脉的中心线还可以为冠状动脉手术治疗提供术前的路径规划和术中，的路径引导，如经皮冠脉介入治疗和冠脉旁路移植术。因此，获取精准的冠状动脉中心线对于临床治疗有极大的意义。
发明内容
本发明的目的在于提供一种基于CT图像的冠状动脉中心线提取方法。
本发明的具体步骤如下：
步骤一、消除被处理的CT层切图像的散粒噪声。
步骤二、对CT层切图像内的所有体素点进行血管区域增强。
步骤三、经过血管增强的图像血管区域的轮廓被增强，在血管区域内使用区域生长法，对CT层切图像中的每一层图像中的血管区域进行分割；使得CT层切图像变为二值图像。
步骤四、对经步骤三分割后的CT层切图像的每层图像分别根据原始数据中的空间分辨率进行三维重构，得到三维空间中的血管数据。
步骤五、对血管数据进行细化，得到血管的初步的中心线；血管数据上各个体素点的值为1，作为目标体素点；对所有目标体素点进行如下操作：提取出目标体素点的六个6-邻接点、十二个18-邻接点和八个26-邻接点；对步骤四重构得到的血管数据采用十二方向的迭代细化处理；每个迭代步骤均对血管数据的所有边界点依次进行十二个方向的细化处理；十二个方向分别为目标体素点指向十二个18-邻接点的方向。
步骤六、将步骤五中提取出的中心线转换为有向图。在有向图中，仅有一个邻接体素的体素点视为端点，有且仅有两个邻接体素的体素点视为普通连接点，存在多于两个邻接体素的体素点视为节点；相邻的两个节点或端点之间的各普通连接点形成一条有向边。根据各有向边上的节点数量筛选出中心线和分支。
作为优选，步骤五中，若一轮细化处理中删除一个或多个体素点，则再进行一轮迭代，直至迭代到一轮细化处理中未有能够被删除的体素点。
作为优选，步骤五中，用U、D、W、E、S、N表示点目标体素点p的6-邻接点；用US、NE、WD、ES、UW、ND、SW、UN、ED、NW、UE、SD分别表示十二个18-邻接点；组成18-邻接点的两个字母表示该18-邻接点与那两个6-邻接点相邻。用点A、B、C、F、G、H、I、J分别表示八个26-邻接点；点A为点U、S、W所在平面的交点；点B为点U、W、N所在平面的交点；点C为点U、N、E所在平面的交点；点F为点U、S、E所在平面的交点；点G为点D、S、W所在平面的交点；点H为点D、W、N所在平面的交点；点I为点D、N、E所在平面的交点；点J为点D、S、E所在平面的交点；十二个方向的细化处理顺序为：US、NE、WD、ES、UW、ND、SW、UN、ED、NW、UE、SD。十二方向的细化处理各自对应一组模板；一组模板内有14个模板。在一个方向的细化处理中，若目标体素点的3×3×3邻域与对应的14个模板中任意一个或多个相符合，则将目标体素点删除。
US方向的十四个模板如下：
第1个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点均标记为“○”，D点标记为“●”，其余各点均标记为“X”。
第2个模板中，目标体素点p的3×3×3邻域中的点A、US、F、SW、S、ES、G、SD、J均标记为“○”，N点标记为“●”，其余各点均标记为“X”。
第3个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点和点SW、S、SE、G、SD、J均标记为“○”，N点和D点标记为“·”，ND点标记为“●”，其余各点均标记为“X”。
第4个模板中，目标体素点p的3×3×3邻域中点UW和SW标记为“V”，UE点和ES点标记为“W”，N点和D点标记为“●”，点U、A、US、F、S标记为“○”，其余各点均标记为“·”。
第5个模板中，目标体素点p的3×3×3邻域中点UW和SW标记为“V”，点UE和ES标记为“Z”，点F、N、D标记为“●”，点U、A、US、S标记为“○”，其余各点均标记为“·”。
第6个模板中，目标体素点p的3×3×3邻域中UW点和SW点标记为“Z”，UE点和ES点标记为“V”，点A、N、D标记为“●”，点U、F、US、S标记为“○”，其余各点均标记为“·”。
第7个模板中，目标体素点p的3×3×3邻域中点UW和SW标记为“V”，点E、N、D标记为“●”，点U、A、US、S标记为“○”，其余各点均标记为“·”。
第8个模板中，目标体素点p的3×3×3邻域中UE点和ES点标记为“V”，点W、N、D标记为“●”，点U、F、US、S标记为“○”，其余各点均标记为“·”。
第9个模板中，目标体素点p的3×3×3邻域中UW点和SW点标记为“Z”，点A、N、E、D标记为“●”，点U、US、S标记为“○”，其余各点均标记为“·”。
第10个模板中，目标体素点p的3×3×3邻域中点UE和ES标记为“Z”，点F、N、W、D标记为“●”，点U、US、S标记为“○”，其余各点均标记为“·”。
第11个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点和点SW、S、G、SD均标记为“○”，ND点和ED点标记为“●”，其余各点均标记为“·”。
第12个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点和点ES、S、SD、J均标记为“○”，ND点和WD点标记为“●”，其余各点均标记为“·”。
第13个模板中，目标体素点p的3×3×3邻域中前侧的九个体素点和点B、UN、UW、U、A、US均标记为“○”，NE点和ND点标记为“●”，其余各点均标记为“·”。
第14个模板中，目标体素点p的3×3×3邻域中前侧的九个体素点和点UN、C、U、UE、US、F均标记为“○”，，NW点和ND点标记为“●”，其余各点均标记为“·”。
符号“○”表示该位置为背景点；符号“●”表示该位置为目标点；符号“·”表示该位置上既可以为目标点也可以为背景点；同一个模板中的字母“X”表示在这些位置上至少有一个为目标点；同一个模板中的字母“V”和“W”各自表示这些位置上至少有一个为背景点；同一个模板中有两个不同的位置标注了字母“Z”，则表示这两个位置上必须一个为目标点，另一个为背景点。
剩余十一个方向对应的模板均通过US方向对应的十四个模板通过旋转、镜像对称得到；具体如下：以目标体素点p为坐标原点，点p到点E为X轴正方向，p到点N为Y轴正方向，p到点U为Z轴正方向建立空间直角坐标系；
UN方向的模板由US方向的模板根据AGJF平面镜像对称得到；UW方向的模板由US方向的模板根据AGIC平面镜像对称得到；UE方向的模板由US方向的模板根据BHJF平面镜像对称得到；SW方向的模板由US方向的模板根据AJIB平面镜像对称得到；ES方向的模板由US方向的模板根据FGHC平面镜像对称得到；NW方向的模板由UN方向的模板根据AJIB平面镜像对称得到；NE方向的模板由UN方向的模板根据FGHC平面镜像对称得到；WD方向的模板由SW方向的模板根据BGJC平面镜像对称得到；SD方向的模板由US方向的模板根据HGJI平面镜像对称得到；ND方向的模板由UN方向的模板根据HGJI平面镜像对称得到；ED方向的模板由UE方向的模板根据HGJI平面镜像对称得到。
作为优选，步骤六中，设置一个长度阈值，将小于长度阈值的分支作为需要被删除的细小分支，大于长度阈值的分支为需要被保留的中心线起点和终点所在的有向边。
作为优选，步骤一的具体过程如下：采用异性扩散滤波器进行散粒噪声的消除。CT层切图像上的任意点x的各向异性扩散方程的表达式如式(1)所示：
式(1)中，
表示图像在点x的灰度，
表示CT层切图像中，点x的坐标所对应的向量，t表示迭代次数；
表示梯度算子，
表示
的梯度，div(·)表示散度算子；
表示扩散函数，其表达式如式(2)所示：
式(2)中，k为与噪声级和边界强度相关的常数；求解式(2)的偏微分方程，得到消除散粒噪声后的滤波结果
其表达式如式(3)所示：
u(t+Δt)＝u(t)+Δt(Φeast-Φwest+Φnorth-Φsouth)#(3)
式(3)中，Φeast、Φwest、Φnorth、Φsouth分别为点x的东侧、西侧、上侧、下侧邻接像素点的扩散通量。Δt表示迭代步长。
作为优选，步骤二中，以CT层切图像内所有体素点分别作为被分析点A进行如下操作：
2-1、本发明利用泰勒级数的展开式分析体素点的局部特征；被分析点A在所在图像I中的泰勒级数展开式I(A+ΔA)如式(4)所示：
式(4)中，I(A)表示图像I中的被分析点A；
表示图像I在被分析点A处的梯度，Δ(·)表示拉普拉斯算子；AT表示被分析点A的转置；H(A)表示被分析点A的Hessian矩阵，其表达式如式(5)所示。
式(5)中，矩阵H(A)内各个元素均为图像I在被分析点A处的二阶偏导；矩阵H(A)存在三个特征值λ1、λ2和λ3，以及三个特征值各自对应的特征向量为
和
三个特征值满足|λ1|≤|λ2|≤|λ3|的关系。
2-2、利用基于Hessian矩阵的Frangi血管增强函数对血管区域进行增强，对非血管区域进行抑制，血管相似性函数vo(λ)如式(6)所示：
式(6)中，RA、RB为椭球体的两个定义参数；S为对比度系数；α、β和γ为三个调参系数。RA、RB和S分别用式(7)、式(8)和式(9)表示：
本发明具有的有益效果是：
本发明提出了一种在三维空间中，的针对冠状动脉的拓扑细化算法，该算法在提取中心线的同时尽可能的保持了冠状动脉原有的拓扑结构；同时，针对细化过程中，所产生的细小分支，本发明采用Dijkstra算法予以消除。本发明提高了基于CT图像的冠状动脉中心线的提取精度与效率，优化效果较为明显。
附图说明
图1为本发明验证时使用的原始CTA图像；
图2a为本发明验证过程中，对原始CTA图像进行各项异性扩散滤波后的结果；
图2b为本发明验证过程中，进行血管增强之后的结果；
图2c为本发明验证过程中，对分割出的冠状动脉血管进行三维重构后的结果；
图3a为本发明提出的十二方向细化算法的十二个方向的示意图；
图3b为本发明提出的十二方向细化算法的18-邻接点和26-邻接点的示意图；
图4a-4n为本发明所用进行十二方向拓扑细化算法时所需的细化模板；
图5为本发明验证过程中，对冠状动脉中心线提取的初步结果；
图6为本发明进行细小分支删除后的冠状动脉中心线提取的最终结果。
具体实施方式
以下结合附图对本发明作进一步说明。
本发明的基础思想为：对采集的原始CT心脏图像中，的冠状动脉进行分割，之后利用十二方向上的拓扑细化方法对分割后的冠状动脉中心线进行提取，并对提取过程中，产生的细小分支结合Dijkstra算法进行消除。
一种基于CT图像的冠状动脉中心线提取方法的具体步骤如下：
步骤一、对原始的CT层切图像利用各向异性扩散滤波器消除散粒噪声；CT层切图像上的任意点x的各向异性扩散方程的表达式如式(1)所示：
式(1)中，
表示图像在点x的灰度，
表示CT层切图像中，点x的坐标所对应的向量，表示迭代次数；本实施例中，为10，
表示梯度算子，
表示
的梯度，div(·)表示散度算子；
表示扩散函数，用来控制扩散速率，其表达式如式(2)所示：
式(2)中，为与噪声级和边界强度相关的常数，本实施例中选取
表示点x的梯度图像中的最大值，经计算k＝5。在4-邻域图像中，求解式(2)的偏微分方程，得到本发明在消除散粒噪声后的滤波结果
其表达式如式(3)所示：
步骤二、对CT层切图像内的所有体素点进行血管区域增强。以CT层切图像内所有体素点分别作为被分析点A进行如下操作：
2-1、本发明利用泰勒级数的展开式分析体素点的局部特征；被分析点A在所在图像I(CT层切图像的其中一层图像)中的泰勒级数展开式I(A+ΔA)如式(4)所示：
式(4)中，被分析点A的坐标为A＝(x，y，z)；
表示图像I在被分析点A处的梯度，Δ(·)表示拉普拉斯算子，AT表示被分析点A的转置，H(A)表示被分析点A的Hessian矩阵，被分析点A的Hessian矩阵表达式如式(5)所示。
式(5)中，Iij(A)表示图像I在点A处的二阶偏导，其求偏导顺序为先对i方向求偏导，再对j方向求偏导。i,j分别取x、y、z；矩阵H(A)存在三个特征值λ1、λ2和λ3，以及三个特征值各自对应的特征向量为
和
且这三个特征向量在三维空间中，相互正交且垂直。在三个特征值中满足|λ1|≤|λ2|≤|λ3|的关系；其中，绝对值最大的特征值λ3所对应的特征向量
表示在三维曲面中，该点的最大曲率方向；而绝对值最小的特征值λ1所对应的特征向量
表示了该点的血管走向。
2-2、利用基于Hessian矩阵的Frangi血管增强函数对血管区域进行增强，对非血管区域进行抑制，血管相似性函数vo(λ)如式(6)所示：
式(6)中，RA、RB为椭球体的两个定义参数，它们的大小可以反应出图中，组织的结构；RA表示垂直于血流方向的横截面与椭球体半长轴之间的关系；RB则表示椭球体体积与垂直于血流方向的横截面之间的关系；S为对比度系数，表示了图像中的感兴趣区域(RegionofInterest,ROI)与背景间的关系，感兴趣区域与背景对比度越高，S其值就越大；α、β和γ为三个调参系数，通过影响RA、RB和S在血管相似性函数中所占的比重来决定式(6)的灵敏度；作为本实施例中的优选，α和β的取值为0.5，γ的取值为Hessian矩阵行列式的值的一半。RA、RB和S分别用式(7)、式(8)和式(9)表示：
步骤三、经过血管增强的图像血管区域的轮廓被增强，在血管区域内使用区域生长法，对CT层切图像中的每一层图像中的血管区域进行分割，使得CT层切图像变为二值图像。本发明中，区域生长的种子点即为人工指定的血管区域内的像素点；与该种子点的灰度值相近的灰度值即为区域生长的条件；目标区域和背景在相交的部分会产生灰度值的变化，从而产生边界，此时发生变化的灰度值就是区域生长的终止条件。
步骤四、对经步骤三分割后的CT层切图像的每层图像分别根据原始数据中的空间分辨率进行三维重构，得到三维空间中的血管数据。
步骤五、通过拓扑细化算法对CT层切图像中的血管数据进行细化，得到血管的初步的中心线；拓扑细化算法通常在二值化图片上进行操作，二值化图像只存在数值为0和1的两种像素，1代表的是白色，0代表的是黑色。通常将图中，值为1的像素定义为目标点，值为0的像素定义为背景点即在图像中，由目标点组成的白色部分为目标区域，由背景点组成的黑色部分为背景区域。用细化操作提取中心线的基本思想就是将无关图形拓扑结构的目标点删除(将该点的值由1置0)，从而得到一条位于为物体中心的单体素宽的中心线。所以这些能够被删除的点被称为简化点。血管数据上各个体素点的值为1，作为目标体素点；对所有目标体素点进行如下操作。
5-1.对计算目标体素点与CT层切图像中其他体素点的欧氏距离d如式(10)所示。
式(10)中，x、y、z分别表示两个体素的三个坐标之间的差值。
若这两个体素点相邻接，则根据距离的不同，这两个体素的关系分为6-邻接、18-邻接和26-邻接，欧式距离分别为1,
和
因此，根据计算的欧氏距离能够找出目标体素点的六个6-邻接点、十二个18-邻接点和八个26-邻接点如图3a和3b所示；图3a用U、D、W、E、S、N表示点目标体素点p的6-邻接点，其分别位于目标体素点的上、下、西、东、南、北方向；该6个方向即为三维立体网格中的6个主要方向；图3a中用12个圆圈标识表示目标体素点p的十二个18-邻接点；用五角星标识表示目标体素点p的八个26-邻接点。图3b中，用US、NE、WD、ES、UW、ND、SW、UN、ED、NW、UE、SD分别表示十二个18-邻接点；组成18-邻接点的两个字母表示该18-邻接点与那两个6-邻接点相邻。图3b中用点A、B、C、F、G、H、I、J分别表示八个26-邻接点；点A为点U、S、W所在平面的交点；点B为点U、W、N所在平面的交点；点C为点U、N、E所在平面的交点；点F为点U、S、E所在平面的交点；点G为点D、S、W所在平面的交点；点H为点D、W、N所在平面的交点；点I为点D、N、E所在平面的交点；点J为点D、S、E所在平面的交点；
本发明对三维重构后的血管采用十二方向细化处理，即每个迭代步骤根据图3b所示的十二个方向进行并行的细化过程，图中，所示的十二个方向即为一个立方体的十二条边(也即中心体素点到十二个18-邻接体素点的方向)，每个方向都是由前述的六个主要方向两两配对组成的，八个26-邻接体素点(即立方体的八个顶点)由A、B、C、F、G、H、I、J分别表示。本发明中，选择的删除方向顺序为：US、NE、WD、ES、UW、ND、SW、UN、ED、NW、UE、SD。根据这种删除顺序又可以将这十二个方向平均分为四组，每组三个方向，所以第一组为US、NE、WD；第二组为ES、UW、ND；第三组为SW、UN、ED；第四组为NW、UE、SD。以这种顺序删除简化点的最大优势在于每个小组中，又包含了6个主要方向，细化的过程可以尽可能保证目标区域会在各个方向上均匀的收缩。
以第一组中，的US方向为例，在该组的迭代过程中，会率先在US方向上进行血管数据的所有边界点的匹配，即分别在U和S方向上删除某些边界点，当完成US方向上的删除之后，原始的血管区域的边界已经发生了变化，产生了新的边界，NE方向的删除操作将会在新的边界上继续进行。余下的各迭代过程以此类推。需要被细化的冠状动脉会在每个方向上从外层开始，向内被层层腐蚀，直至十二个子迭代过程中，不再有能够被删除的目标点时，主迭代过程结束，中心线提取完成。
十二方向的细化处理各自对应一组模板；一组模板内有14个模板。每轮细化处理中，均在用一个方向的模板对血管数据上所有的边界点分别进行处理后，再用下一个方向的模板对血管数据上所有的边界点分别进行处理，直到十二个方向处理完；若一轮细化处理中删除(即置为0)一个或多个体素点，则再进行一轮细化处理，直至迭代到一轮细化处理中未有能够被删除的体素点。
在一个方向对一个模板体素点的处理中，若目标体素点的3×3×3邻域(26个各邻接点)与该方向对应的14个模板中任意一个或多个相符合，则将目标体素点p是可以在细化过程中删除的简单点，将其删除。
US方向的十四个模板如图4a-4n所示；点p为图中各模板的中心点。图中，符号“○”表示该位置为背景点；符号“●”表示目标点，因为每个模板的中心点都是点p，而点p也是目标点所以用相同的符号表示；符号“·”表示该位置上可以为目标点也可以为背景点，如图4c和图4d中所示；图中，同一个模板中的字母“X”表示在这些位置上至少有一个为目标点，如图4a和图4b中所示；同一个模板中的字母“V”和“W”各自表示这些位置上至少有一个为背景点，如图4d、图4g、图4h所示；其中，在模板图4d中，既包含两个标记为“V”的点也包含两个标记为“W”的点，则表明两个标记为“V”的点中，至少有一个为背景点，两个标记为“W”的点中，也至少有一个为背景点，而不是这四个点中，任意两个为背景点；当有两个不同的位置标注了字母“Z”，则表示这两个位置上必须为属性不同的点(如一个为目标点，另一个为背景点)，但是这两个点的位置不是固定的，如图4f、图4j等。
模板图4a中，目标体素点p的3×3×3邻域中顶层的九个体素点均标记为“○”，D点标记为“●”，其余各点均标记为“X”。
模板图4b中，目标体素点p的3×3×3邻域中的点A、US、F、SW、S、ES、G、SD、J均标记为“○”，N点标记为“●”，其余各点均标记为“X”。
模板图4c中目标体素点p的3×3×3邻域中顶层的九个体素点和点SW、S、SE、G、SD、J均标记为“○”，N点和D点标记为“·”，ND点标记为“●”，其余各点均标记为“X”。
模板图4d中目标体素点p的3×3×3邻域中点UW和SW标记为“V”，UE点和ES点标记为“W”，N点和D点标记为“●”，点U、A、US、F、S标记为“○”，其余各点均标记为“·”。
模板图4e中目标体素点p的3×3×3邻域中点UW和SW标记为“V”，点UE和ES标记为“Z”，点F、N、D标记为“●”，点U、A、US、S标记为“○”，其余各点均标记为“·”。
模板图4f中目标体素点p的3×3×3邻域中UW点和SW点标记为“Z”，UE点和ES点标记为“V”，点A、N、D标记为“●”，点U、F、US、S标记为“○”，其余各点均标记为“·”。
模板图4g中目标体素点p的3×3×3邻域中点UW和SW标记为“V”，点E、N、D标记为“●”，点U、A、US、S标记为“○”，其余各点均标记为“·”。
模板图4h中目标体素点p的3×3×3邻域中UE点和ES点标记为“V”，点W、N、D标记为“●”，点U、F、US、S标记为“○”，其余各点均标记为“·”。
模板图4i中目标体素点p的3×3×3邻域中UW点和SW点标记为“Z”，点A、N、E、D标记为“●”，点U、US、S标记为“○”，其余各点均标记为“·”。
模板图4j中目标体素点p的3×3×3邻域中点UE和ES标记为“Z”，点F、N、W、D标记为“●”，点U、US、S标记为“○”，其余各点均标记为“·”。
模板图4k中，目标体素点p的3×3×3邻域中顶层的九个体素点和点SW、S、G、SD均标记为“○”，ND点和ED点标记为“●”，其余各点均标记为“·”。
模板图4l中，目标体素点p的3×3×3邻域中顶层的九个体素点和点ES、S、SD、J均标记为“○”，ND点和WD点标记为“●”，其余各点均标记为“·”。
模板图4m中目标体素点p的3×3×3邻域中前侧的九个体素点和点B、UN、UW、U、A、US均标记为“○”，NE点和ND点标记为“●”，其余各点均标记为“·”。
模板图4n中目标体素点p的3×3×3邻域中前侧的九个体素点和点UN、C、U、UE、US、F均标记为“○”，，NW点和ND点标记为“●”，其余各点均标记为“·”。
剩余十一个方向对应的模板均通过US方向对应的十四个模板通过旋转、镜像对称得到；具体如下：以目标体素点p为坐标原点，点p到点E为X轴正方向，p到点N为Y轴正方向，p到点U为Z轴正方向建立空间直角坐标系；
UN方向的模板由US方向的模板根据AGJF平面镜像对称得到；UW方向的模板由US方向的模板根据AGIC平面镜像对称得到；UE方向的模板由US方向的模板根据BHJF平面镜像对称得到；SW方向的模板由US方向的模板根据AJIB平面镜像对称得到；ES方向的模板由US方向的模板根据FGHC平面镜像对称得到；NW方向的模板由UN方向的模板根据AJIB平面镜像对称得到；NE方向的模板由UN方向的模板根据FGHC平面镜像对称得到；WD方向的模板由SW方向的模板根据BGJC平面镜像对称得到；SD方向的模板由US方向的模板根据HGJI平面镜像对称得到；ND方向的模板由UN方向的模板根据HGJI平面镜像对称得到；ED方向的模板由UE方向的模板根据HGJI平面镜像对称得到。
步骤六、为了消除细化过程中产生的细小分支，本发明采用Dijkstra算法对细小分支进行限制和消除。Dijkstra是贪婪算法的一种实现形式，其特点是迭代从起点向外层层扩张，直至到达终点。
首先，将步骤五中提取出的中心线转换为由节点和有向边组成的有向图。在有向图中，仅有一个邻接体素的体素点视为端点，有且仅有两个邻接体素的体素点视为普通连接点，存在多于两个邻接体素(仅考虑值为1的目标点)的体素点视为节点；相邻的两个节点或端点之间的各普通连接点形成一条有向边。
根据有向图和节点的数量N建立维度为n×n的邻接矩阵
zij为第i个节点到第j个节点的路径长度，i＝1,2,…,n，j＝1,2,…,n。若第i个节点与第j个节点不在同一条中心线上，则zij的值为无穷大。
若一个节点直接到另一个节点的路径长度小于所有其他经过第三个节点的路径长度，则该两个节点为相邻节点；当所有节点的相邻节点均被确定后，将相邻节点依次连接，即获得不包含细小分支的中心线，从而区分出细小分支和中心线；之后，将区分出的细小分支删除。由于起点和终点的位置也只与一个节点相连，本发明为了较好的将起点与终点所在的有向边保存在中心线中，会根据每根血管的实际情况确定一个固定的长度阈值，将小于长度阈值的分支作为需要被删除的细小分支，长于该长度的分支就是需要被保留的中心线起点和终点所在的有向边。
举例来说，假设从节点B到节点C存在两条路径，一条是经过aBC到达节点C，另一条是先通过aBD到达D节点，再通过aDC到达节点C，由于B、C节点为相邻节点，所以选择第一条路径的距离会小于第二条路径，所以节点B、C之间的最短路径为aBC。依此方法逐一找出两节点之间的最短路径即可完成细小分支的删除。
为证明上述方法的可行性，下面采用具体的冠状动脉CT图像进行举例说明。
本发明验证所采用的CT数据均来自鹿特丹冠状动脉算法评价框架(RotterdamCoronary ArteryAlgorithmEvaluationFramework，RCAAEF)提供的公共数据集。该数据集包含32例数据，这些数据于2005年6月至2006年6月采自荷兰鹿特丹的伊拉斯谟医学中，心(Erasmus MC)，其中，20例采自西门子64排螺旋CT，12例采自西门子双源CT。两种CT设备扫描时所用的球管电压均为120KV；采集方式为心电图脉冲采集；扫描仪宽度为0.6mm；三维数据重建所用卷积核为B46f；所获得的图像中，每一个体素对应实际的大小的平均值为0.32×0.32×0.4mm3。在每一个CT层切图像中，包含4根冠状动脉的主要分支，其中，前三根血管的选取是固定的，分别为右冠状动脉、左前降支和左回旋支，第四根血管根据心脏的特点可能选取第一对角支、第二对角支、选择性斜行冠状动脉、第一钝缘支、后降支和锐缘支动脉中，的一支。专家在这些血管中，对中心线进行了人工标记，并提供了对中心线上的每一点给出了允许的误差范围，可将专家标记的坐标作为中心线的实际位置。在标记过程中，为了保证标记的准确性，每一根血管都有两位专家进行标注，标注信息如表1所示。表中，记录了专家标注的五个病灶点的位置信息，前三列代表的是病灶在世界系坐标中，的位置信息，第四列所对应的是提取出的中心线的点在图像中，对应的允许误差范围，也就是说，如果算法提取出的中心线上的某一点落在以表中，前三列为圆心，第四列为半径构建出来的圆形区域内，那么就认为算法对于中心线的提取在该点的位置是准确的。
表1标注数据展示
实验过程如下：先对原始的CT图像(如图1所示)，进行各项异性扩散滤波，滤波结果如图2a所示；对图2a进行血管增强处理，血管增强的结果如图2b所示；在经过血管增强处理的各层切图像中，利用区域生长法分割出血管区域，再根据原始数据的空间分辨率重构出三维的血管，重构结果如图2c所示；而后用本发明步骤五中，提出的十二方向的拓扑细化算法，对三维血管进行中心线提取操作，中心线提取结果如图4所示，图中，灰色部分为重构出的冠状动脉区域，中间的黑色曲线为提取出的中心线，通过观察发现图中提取出的中心线上存在着明显的分支，这些分支将会干扰中心线提取的精度。为了量化中心线提取的精度，本发明分别对数据集中，所有血管数据进行中心线提取，然后对提取结果进行量化分析，并随机选取了6根血管的分析结果进行展示，如表2所示。在表2中，所涉及的量化参数包括中心线提取的重叠率(OV)、提取的中心线与参考标准间的平均距离(AD)和在提取过程中，第一次出现错误前的重叠率(OF)。
表2未删除细小分支的冠状动脉中心线提取结果
为了消除提取过程中，的细小分支，本发明在图5的基础上用Dijkstra算法对毛刺进行消除，所得的结果如图6所示。通过图5与图6的对比不难发现，本发明的算法对提取过程中，产生的细小分支有很好的限制作用，大大减少了中，心上的毛刺，使提取出的中心线更加的平滑，更加符合中心线的实际形状。为了更加直观的展示该算法对提取精度的优化，本发明对上述相同6根血管中心线进行细小分支删除的操作，对操作后的结果进行了定量分析，分析结果如表3所示。
根据表2和表3中，的数据，经统计后发现，与未进行细小分支移除的拓扑细化算法相比，运用了细小分支删除算法后，中心线的OV提高了2％，AD减少了38.2％，OF提高了16％。对数据集中，的所有血管的OV、AD、OF进行计算，得到平均重叠率、平均距离和平均第一次出错前的重叠率分别为0.985、1.974和0.832。
本发明的所用的十二方向的拓扑细化方法可以有效的提取出冠状动脉的中心线；再者，结合Dijkstra算法使中心线的提取结果进一步提高。另外，实验中，只给出了一个冠状动脉CT图像的数据集的结果，该算法也同样适用于其他冠状动脉CT图像数据集。
Claims (5)
1.一种基于CT图像的冠状动脉中心线提取方法，其特征在于：步骤一、消除被处理的CT层切图像的散粒噪声；
步骤二、对CT层切图像内的所有体素点进行血管区域增强；
步骤三、在经过步骤二增强的血管区域内使用区域生长法，对CT层切图像中的每一层图像中的血管区域进行分割，使得CT层切图像变为二值图像；
步骤四、对经步骤三分割后的CT层切图像的每层图像分别根据原始数据中的空间分辨率进行三维重构，得到三维空间中的血管数据；
步骤五、对血管数据进行细化，得到血管的初步的中心线；血管数据上各个体素点的值为1，作为目标体素点；对所有目标体素点进行如下操作：提取出目标体素点的六个6-邻接点、十二个18-邻接点和八个26-邻接点；对步骤四重构得到的血管数据采用十二方向的迭代细化处理；每个迭代步骤均对血管数据的所有边界点依次进行十二个方向的细化处理；十二个方向分别为目标体素点指向十二个18-邻接点的方向；
用U、D、W、E、S、N表示点目标体素点p的6-邻接点；用US、NE、WD、ES、UW、ND、SW、UN、ED、NW、UE、SD分别表示十二个18-邻接点；组成18-邻接点的两个字母表示该18-邻接点与那两个6-邻接点相邻；用点A、B、C、F、G、H、I、J分别表示八个26-邻接点；点A为点U、S、W所在平面的交点；点B为点U、W、N所在平面的交点；点C为点U、N、E所在平面的交点；点F为点U、S、E所在平面的交点；点G为点D、S、W所在平面的交点；点H为点D、W、N所在平面的交点；点I为点D、N、E所在平面的交点；点J为点D、S、E所在平面的交点；十二个方向的细化处理顺序为：US、NE、WD、ES、UW、ND、SW、UN、ED、NW、UE、SD；十二方向的细化处理各自对应一组模板；一组模板内有14个模板；在一个方向的细化处理中，若目标体素点的3×3×3邻域与对应的14个模板中任意一个或多个相符合，则将目标体素点删除；
US方向的十四个模板如下：
第1个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点均标记为“○”，D点标记为“●”，其余各点均标记为“X”；
第2个模板中，目标体素点p的3×3×3邻域中的点A、US、F、SW、S、ES、G、SD、J均标记为“○”，N点标记为“●”，其余各点均标记为“X”；
第3个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点和点SW、S、SE、G、SD、J均标记为“○”，N点和D点标记为“·”，ND点标记为“●”，其余各点均标记为“X”；
第4个模板中，目标体素点p的3×3×3邻域中点UW和SW标记为“V”，UE点和ES点标记为“W”，N点和D点标记为“●”，点U、A、US、F、S标记为“○”，其余各点均标记为“·”；
第5个模板中，目标体素点p的3×3×3邻域中点UW和SW标记为“V”，点UE和ES标记为“Z”，点F、N、D标记为“●”，点U、A、US、S标记为“○”，其余各点均标记为“·”；
第6个模板中，目标体素点p的3×3×3邻域中UW点和SW点标记为“Z”，UE点和ES点标记为“V”，点A、N、D标记为“●”，点U、F、US、S标记为“○”，其余各点均标记为“·”；
第7个模板中，目标体素点p的3×3×3邻域中点UW和SW标记为“V”，点E、N、D标记为“●”，点U、A、US、S标记为“○”，其余各点均标记为“·”；
第8个模板中，目标体素点p的3×3×3邻域中UE点和ES点标记为“V”，点W、N、D标记为“●”，点U、F、US、S标记为“○”，其余各点均标记为“·”；
第9个模板中，目标体素点p的3×3×3邻域中UW点和SW点标记为“Z”，点A、N、E、D标记为“●”，点U、US、S标记为“○”，其余各点均标记为“·”；
第10个模板中，目标体素点p的3×3×3邻域中点UE和ES标记为“Z”，点F、N、W、D标记为“●”，点U、US、S标记为“○”，其余各点均标记为“·”；
第11个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点和点SW、S、G、SD均标记为“○”，ND点和ED点标记为“●”，其余各点均标记为“·”；
第12个模板中，目标体素点p的3×3×3邻域中顶层的九个体素点和点ES、S、SD、J均标记为“○”，ND点和WD点标记为“●”，其余各点均标记为“·”；
第13个模板中，目标体素点p的3×3×3邻域中前侧的九个体素点和点B、UN、UW、U、A、US均标记为“○”，NE点和ND点标记为“●”，其余各点均标记为“·”；
第14个模板中，目标体素点p的3×3×3邻域中前侧的九个体素点和点UN、C、U、UE、US、F均标记为“○”，NW点和ND点标记为“●”，其余各点均标记为“·”；
符号“○”表示这个位置为背景点；符号“●”表示这个位置为目标点；符号“·”表示这个位置上既可以为目标点也可以为背景点；同一个模板中的字母“X”表示在这些位置上至少有一个为目标点；同一个模板中的字母“V”和“W”各自表示这些位置上至少有一个为背景点；同一个模板中有两个不同的位置标注了字母“Z”，则表示这两个位置上必须一个为目标点，另一个为背景点；
剩余十一个方向对应的模板均通过US方向对应的十四个模板通过旋转、镜像对称得到；具体如下：以目标体素点p为坐标原点，点p到点E为X轴正方向，p到点N为Y轴正方向，p到点U为Z轴正方向建立空间直角坐标系；
UN方向的模板由US方向的模板根据AGJF平面镜像对称得到；UW方向的模板由US方向的模板根据AGIC平面镜像对称得到；UE方向的模板由US方向的模板根据BHJF平面镜像对称得到；SW方向的模板由US方向的模板根据AJIB平面镜像对称得到；ES方向的模板由US方向的模板根据FGHC平面镜像对称得到；NW方向的模板由UN方向的模板根据AJIB平面镜像对称得到；NE方向的模板由UN方向的模板根据FGHC平面镜像对称得到；WD方向的模板由SW方向的模板根据BGJC平面镜像对称得到；SD方向的模板由US方向的模板根据HGJI平面镜像对称得到；ND方向的模板由UN方向的模板根据HGJI平面镜像对称得到；ED方向的模板由UE方向的模板根据HGJI平面镜像对称得到；
步骤六、将步骤五中提取出的中心线转换为有向图；在有向图中，仅有一个邻接体素的体素点视为端点，有且仅有两个邻接体素的体素点视为普通连接点，存在多于两个邻接体素的体素点视为节点；相邻的两个节点或端点之间的各普通连接点形成一条有向边；根据各有向边上的节点数量筛选出中心线和分支。
2.根据权利要求1所述的一种基于CT图像的冠状动脉中心线提取方法，其特征在于：步骤五中，若一轮细化处理中删除一个或多个体素点，则再进行一轮迭代，直至迭代到一轮细化处理中未有能够被删除的体素点。
3.根据权利要求1所述的一种基于CT图像的冠状动脉中心线提取方法，其特征在于：步骤六中，设置一个长度阈值，将小于长度阈值的分支作为需要被删除的细小分支，大于长度阈值的分支为需要被保留的中心线起点和终点所在的有向边。
4.根据权利要求1所述的一种基于CT图像的冠状动脉中心线提取方法，其特征在于：步骤一的具体过程如下：采用异性扩散滤波器进行散粒噪声的消除；CT层切图像上的任意点x的各向异性扩散方程的表达式如式(1)所示：
式(1)中，
表示图像在点x的灰度，
表示CT层切图像中，点x的坐标所对应的向量，t表示迭代次数；
表示梯度算子，
表示
的梯度，div(·)表示散度算子；
表示扩散函数，其表达式如式(2)所示：
式(2)中，k为与噪声级和边界强度相关的常数；求解式(2)的偏微分方程，得到消除散粒噪声后的滤波结果
其表达式如式(3)所示：
u(t+Δt)＝u(t)+Δt(Φeast-Φwest+Φnorth-Φsouth) (3)
式(3)中，Φeast、Φwest、Φnorth、Φsouth分别为点x的东侧、西侧、上侧、下侧邻接像素点的扩散通量；Δt表示迭代步长。
5.根据权利要求1所述的一种基于CT图像的冠状动脉中心线提取方法，其特征在于：步骤二中，以CT层切图像内所有体素点分别作为被分析点A进行如下操作：
2-1、利用泰勒级数的展开式分析体素点的局部特征；被分析点A在所在图像I中的泰勒级数展开式I(A+ΔA)如式(4)所示：
式(4)中，I(A)表示图像I中的被分析点A；
表示图像I在被分析点A处的梯度，Δ(·)表示拉普拉斯算子；AT表示被分析点A的转置；H(A)表示被分析点A的Hessian矩阵，其表达式如式(5)所示；
式(5)中，矩阵H(A)内各个元素均为图像I在被分析点A处的二阶偏导；矩阵H(A)存在三个特征值λ1、λ2和λ3，以及三个特征值各自对应的特征向量为
和
三个特征值满足|λ1|≤|λ2|≤|λ3|的关系；
2-2、利用基于Hessian矩阵的Frangi血管增强函数对血管区域进行增强，对非血管区域进行抑制，血管相似性函数vo(λ)如式(6)所示：
式(6)中，RA、RB为椭球体的两个定义参数；S为对比度系数；α、β和γ为三个调参系数；RA、RB和S分别用式(7)、式(8)和式(9)表示：
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010332170.5A CN111612743B (zh) | 2020-04-24 | 2020-04-24 | 一种基于ct图像的冠状动脉中心线提取方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010332170.5A CN111612743B (zh) | 2020-04-24 | 2020-04-24 | 一种基于ct图像的冠状动脉中心线提取方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111612743A CN111612743A (zh) | 2020-09-01 |
| CN111612743B true CN111612743B (zh) | 2023-05-02 |
Family
ID=72201498
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010332170.5A Active CN111612743B (zh) | 2020-04-24 | 2020-04-24 | 一种基于ct图像的冠状动脉中心线提取方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111612743B (zh) |
Families Citing this family (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112215941B (zh) * | 2020-09-04 | 2021-06-25 | 北京科技大学 | 一种基于微分几何的主动脉扭曲的评估方法及系统 |
| CN112465772B (zh) * | 2020-11-25 | 2023-09-26 | 平安科技（深圳）有限公司 | 眼底彩照图像血管评估方法、装置、计算机设备和介质 |
| CN112598639A (zh) * | 2020-12-23 | 2021-04-02 | 济南国科医工科技发展有限公司 | 一种虚拟支气管镜的拼接点多分支路径规划方法及系统 |
| CN112308857B (zh) * | 2020-12-25 | 2021-05-11 | 数坤（北京）网络科技有限公司 | 一种血管中心线的确定方法、装置及可读存储介质 |
| CN113222964B (zh) * | 2021-05-27 | 2021-11-12 | 推想医疗科技股份有限公司 | 一种冠脉中心线提取模型的生成方法及装置 |
| CN113450364B (zh) * | 2021-06-29 | 2022-05-13 | 湖南大学 | 一种基于三维通量模型的树状结构中心线提取方法 |
| CN113628346B (zh) * | 2021-06-30 | 2023-08-08 | 上海市胸科医院 | 一种自由浏览vb的方法以及进行标记的方法和系统 |
| CN113689388B (zh) * | 2021-08-03 | 2024-02-06 | 慧影医疗科技(北京)股份有限公司 | 用于主动脉曲面重建的三维中心线起始点定位方法及装置 |
| CN114972859A (zh) * | 2022-05-19 | 2022-08-30 | 推想医疗科技股份有限公司 | 像元的分类方法、模型训练方法、装置、设备和介质 |
| CN114677374B (zh) * | 2022-05-27 | 2022-09-02 | 杭州键嘉机器人有限公司 | 一种下颌神经管中心线提取与半径计算的方法 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101393644A (zh) * | 2008-08-15 | 2009-03-25 | 华中科技大学 | 一种肝门静脉血管树建模方法及其系统 |
| CN106447645A (zh) * | 2016-04-05 | 2017-02-22 | 天津大学 | 增强ct图像中冠脉钙化检测及量化装置和方法 |
Family Cites Families (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US9129417B2 (en) * | 2012-02-21 | 2015-09-08 | Siemens Aktiengesellschaft | Method and system for coronary artery centerline extraction |
| CN106296664B (zh) * | 2016-07-30 | 2019-10-08 | 上海联影医疗科技有限公司 | 血管提取方法 |
| CN106815853B (zh) * | 2016-12-14 | 2020-06-05 | 海纳医信（北京）软件科技有限责任公司 | 对眼底图像中视网膜血管的分割方法和装置 |
| EP3352135B1 (en) * | 2017-01-19 | 2019-09-11 | INESC TEC - Instituto de Engenharia de Sistemas e Computadores, Tecnologia e Ciência | Method and apparatus for segmentation of blood vessels |
| CN107248155A (zh) * | 2017-06-08 | 2017-10-13 | 东北大学 | 一种基于swi图像的脑静脉血管分割方法 |
| CN108765385B (zh) * | 2018-05-14 | 2021-12-07 | 广东药科大学 | 一种双源ct冠状动脉自动提取方法 |
| CN109461138B (zh) * | 2018-09-29 | 2020-10-27 | 深圳睿心智能医疗科技有限公司 | 参数计算方法、系统、可读存储介质及计算机设备 |
| CN110706246B (zh) * | 2019-10-15 | 2022-08-02 | 上海微创卜算子医疗科技有限公司 | 一种血管图像分割方法、装置、电子设备和存储介质 |
-
2020
- 2020-04-24 CN CN202010332170.5A patent/CN111612743B/zh active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101393644A (zh) * | 2008-08-15 | 2009-03-25 | 华中科技大学 | 一种肝门静脉血管树建模方法及其系统 |
| CN106447645A (zh) * | 2016-04-05 | 2017-02-22 | 天津大学 | 增强ct图像中冠脉钙化检测及量化装置和方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN111612743A (zh) | 2020-09-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111612743B (zh) | 一种基于ct图像的冠状动脉中心线提取方法 | |
| WO2020001217A1 (zh) | 一种基于卷积神经网络的ct图像中带夹层主动脉分割方法 | |
| CN106097298B (zh) | 基于球空间划分的冠状动脉自动分段与解剖学标记方法 | |
| Sun et al. | Graph-based IVUS segmentation with efficient computer-aided refinement | |
| WO2018120644A1 (zh) | 血管提取方法及系统 | |
| WO2018000564A1 (zh) | 血管中心线提取方法及系统 | |
| JP6542129B2 (ja) | リアルタイムの診断上有用な結果 | |
| US8358819B2 (en) | System and methods for image segmentation in N-dimensional space | |
| US8571278B2 (en) | System and methods for multi-object multi-surface segmentation | |
| US7995810B2 (en) | System and methods for image segmentation in n-dimensional space | |
| CN109478327B (zh) | 用于在任意视场计算断层摄影血管造影（cta）中自动检测全身动脉的方法 | |
| JP2007289704A (ja) | 大動脈瘤の半自動分析のためのシステム及び方法 | |
| CN111815663A (zh) | 一种基于Hessian矩阵和灰度法的肝血管分割系统 | |
| CN112862833A (zh) | 血管分段方法、电子装置和存储介质 | |
| Tan et al. | An approach for pulmonary vascular extraction from chest CT images | |
| JP5105997B2 (ja) | 医用画像処理装置、医用画像診断装置、及びプログラム | |
| Xiong et al. | Comprehensive modeling and visualization of cardiac anatomy and physiology from CT imaging and computer simulations | |
| JP4807819B2 (ja) | 画像処理装置 | |
| KR20140120236A (ko) | 심근 및 심혈관 정보의 통합 분석 방법 | |
| Feuerstein et al. | Adaptive model based pulmonary artery segmentation in 3D chest CT | |
| Roy et al. | MDL-IWS: multi-view deep learning with iterative watershed for pulmonary fissure segmentation | |
| KR20230092947A (ko) | 의료용 이미지 내에서 적어도 하나의 관상 구조체를 세그멘트화하고 식별하기 위한 방법 및 시스템 | |
| Novikov et al. | Automated anatomy-based tracking of systemic arteries in arbitrary field-of-view CTA scans | |
| Rossignac et al. | Pearling: 3D interactive extraction of tubular structures from volumetric images | |
| Whited et al. | Pearling: 3D interactive extraction of tubular structures from volumetric images |
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 |