血管提取方法与流程

本发明涉及医学图像处理领域，特别是涉及一种血管提取方法。

背景技术：

血管成像技术包括CT血管造影(CTA，CT angiography)，磁共振血管造影(Magnetic Resonance Angiography，MRA)，以及MR非造影成像等，可以帮助医生诊断血管的各种疾病，如钙化、狭窄、动脉瘤、硬脑膜等。通过血管成像技术获取的血管图像，主要是三维图像，并不能给医生直观的感受。因此需要从医学图像中提取血管并以三维显示技术展示血管的形态，以更好地辅助医生对病灶进行分析，提高医疗诊断的准确性和科学性，有利于制定最优的治疗方案及手术规划，对医学研究具有重要的意义。

但是在实际血管提取中，对于紧贴周围骨骼区域，例如穿颅的颈内动脉，紧贴髂骨的髂动脉，由于空间位置错综复杂，同时增强后的血管CT值和骨的CT值重叠，通常的血管提取方法，例如区域生长，水平集等很容易混淆血管和骨骼，造成血管分割失败。

为解决上述技术问题，现有技术中采用的血管提取方法有：

1.基于血管模型计算血管的中心线，根据此中心线，用单纯网格拟合实际图像中的血管结构，将血管分割出来。此方法虽然可以准确提取血管，但由于网格计算非常复杂，导致分割速度很慢，阻碍了实际应用的推广。

2.基于非模型的方法，例如采用Graph Cut的算法实现了颈内动脉的提取。该算法经过手动确定起始点后，计算起始点之间的Dijkstra距离。以该距离为中心，标记血管感兴趣区域。最后通过Graph Cut算法，优化感兴趣区域的边界实现血管提取。由于该算法计算复杂，导致血管提取速度较慢。还有采用基于Ray-Casting的动态追踪实现了穿颅血管的提取。该算法经过手动确定起始点后，通过Ray-Casting向四周放射线，获得血管的边界，再通过椭圆拟合实现血管提取。该方法鲁棒性差，不同数据的血管提取需要不同的参数，同时该方法无法提取旋转跨度很大的椎动脉。

3.血管剪影方法，病人需要同一部位扫描2次(不注射造影剂扫描和注射造影剂扫描)，得到非CTA图像和CTA图像。非CTA图像中，血管CT值低，CTA图像中血管CT值升高。这样两种图像配准后相减就可以将血管提取出来。但是该方法需要对病人扫描2次，比较耗时，同时给病人带来更多的扫描辐射量。

因此，有必要改进现有的血管提取方法。

技术实现要素：

本发明所要解决的技术问题提出一种方法，其能有效改进现有血管提取方法，以有效获取与紧贴骨骼的血管组织。

本发明为解决上述技术问题所采用的技术方案为一种血管提取方法，包括如下步骤：

获取由若干层切片图像构成的三维图像，选取一层切片图像，基于所述切片图像中像素点的灰度信息获取若干个连通域；

基于所述连通域的圆度确定种子点；

以所述种子点为起点，基于阈值a对所述三维图像进行区域生长，获取第一图像；

以所述种子点为起点，基于第一阈值对所述第一图像进行区域生长，获取血管掩膜；

将所述第一图像与血管掩膜相减获取第一减影图像；

粗分割所述第一减影图像，获取初始掩膜；

形态学处理所述初始掩膜，获取骨骼掩膜；

将所述第一图像与所述骨骼掩膜相减获取第二减影图像，所述第二减影图像包括若干个连通域；

基于血管生长条件连接所述若干个连通域，获取血管提取结果。

优选的，所述选取一层切片图像为沿Z轴方向上的终止层切片图像。

优选的，所述确定种子点包括：

计算所述连通域中所有像素点的边界距离场值；

选取最大边界距离场值为该连通域的半径，根据圆面积公式计算该连通域的圆形面积，并计算所述圆形面积与所述连通域实际面积的比值，所述比值为所述连通域的圆度；

选取所述切片图像中圆度与1差值的绝对值的的最小值所在的连通域，确定该连通域中最大边界距离场值所在的像素点为所述种子点。

优选的，所述获取血管掩膜包括：

计算所述第一图像中各体素点的边界距离场值；

根据所述种子点基于所述第一阈值进行区域生长，获取所述血管掩膜；

其中，所述第一阈值为大于3.0个体素点的宽度。

优选的，获取所述血管掩膜还包括：对所述血管掩膜作膨胀操作。

优选的，所述粗分割包括：

计算所述第一减影图像上所有体素点的边界距离场值；

选取边界距离场值大于阈值b的体素点为骨骼种子点，根据所述第一减影图像中体素点的灰度信息或边界距离场值进行区域生长，获取所述初始掩膜；

其中，所述阈值b为大于3.0个体素点的宽度；所述体素点的灰度信息为灰度值大于500的常数，所述体素点的边界距离场值为大于1.5个体素点的宽度。

优选的，形态学处理包括：

以阈值c为宽度腐蚀所述初始掩膜，获取若干个连通域；

选取最大体积的连通域为骨骼区域；

以阈值d为宽度膨胀所述骨骼区域，获取第一骨骼掩膜。

优选的，所述获取骨骼掩膜还包括填实操作，包括如下步骤：

根据所述种子点基于阈值d对所述三维图像进行区域生长，获取第二图像；

逐层填实所述第二图像中各连通域；

将填实后的第二图像与所述第一骨骼掩膜进行图像叠加，获取所述骨骼掩膜。

优选的，所述填实操作还包括：

提取叠加后的骨骼掩膜中最大体积的连通域，对该连通域经闭操作获取所述骨骼掩膜。

优选的，所述血管生长条件为：

计算所述第二减影图像包括若干个连通域的边界距离场值；

根据所述源距离场，并且排除所述骨骼掩膜中各个连通域所在区域提取血管中心线，获取血管提取结果。

与现有技术相比，本发明的优点在于：通过基于所述连通域的圆度确定种子点，根据所述种子点基于阈值进行三维区域生长，获取包含血管和骨骼的第一图像；将所述第一图像与所述骨骼掩膜相减，获取包含有主血管和少部分碎骨的连通域的第二减影图像；基于血管生长条件，排除所述骨骼区域，通过多区域区域生长的方式连接所述血管连通域，有效获取与骨骼粘连的血管。

【附图说明】

图1为本发明一实施例中血管提取方法的流程图；

图2为本发明一实施例中确定种子点的方法流程图；

图3为本发明一实施例中获取血管掩膜的方法流程图；

图4为本发明一实施例中获取血管掩膜的结果示意图；

图5为本发明一实施例中获取初始掩膜的方法流程图；

图6为本发明一实施例中获取骨骼掩膜的方法流程图；

图7a～7c为本发明一实施例中填实操作的各阶段结果示意图；

图8为本发明一实施例中提取血管的方法流程示意图；

图9为本发明一实施例中血管提取的结果示意图。

【具体实施方式】

在下面的描述中阐述了很多具体细节以便于充分理解本发明。但是本发明能够以很多不同于在此描述的其它方式来实施，本领域技术人员可以在不违背本发明内涵的情况下做类似推广，因此本发明不受下面公开的具体实施的限制。

其次，本发明利用示意图进行详细描述，在详述本发明实施例时，为便于说明，所述示意图只是实施例，其在此不应限制本发明保护的范围。使本发明的上述目的、特征和优点能够更为明显易懂，下面结合附图和实施例对本发明的具体实施方式做详细的说明。

本实施例采用一种的血管提取，如图1所示，包括如下步骤：

S10：选取三维图像中的切片图像，基于灰度信息获取若干个连通域。本实施例中，所述三维图像包含若干个二维切片图像。所述三维图像中的最小处理单元为体素点，切片图像中最小处理单元为像素点。所述三维图像可以通过各类模态的成像系统扫描采集获得三维或二维图像，也可以通过诸如存储系影像归档和通信系统(Picture Archiving and Communication Systems，PACS)等内部或外部存储系统传输获得。所述模态包括但不限于磁共振成像(MRI)、磁共振血管造影(MRA)、计算机断层扫描(CT)、正电子发射断层扫描(Positron Emission Tomography，PET)等一种或多种的组合。

所述三维图像包含各种组织结构，例如包含各类血管区域，例如头颈部血管，胸腹部血管和下肢血管，所述血管区域包含主血管和微血管，主血管可以是动脉血管、静脉血管等血管，或者上述的任意组合。动脉血管可以是头动脉、颈动脉、胸动脉、腹动脉、腋动脉、肱动脉等，或者上述的任意组合。还包括各类骨骼组织，例如髂骨、小腿骨、胫骨等下肢骨骼组织。部分血管组织紧贴骨骼组织(例如胫动脉)，或者穿过骨骼组织(例如椎骨中穿行的左右椎动脉)，并且由于造影剂等影响使得血管组织和骨骼组织在图像显示的强度值(例如灰度值)有部分重叠，或者由于成像结构非常接近，有限的检测器分辨率，血管钙化点，以及介入装置(例如植入的血管支架)使得骨骼和脉管系统结构的识别和分割变得困难起来，大大增加血管提取的难度。

本实施例中，所述选取一层切片图像为沿Z轴方向上的终止层切片图像，z方向可以为垂直轴的方向(冠状面和矢状面形成的相交线的方向)，例如，从被成像物体下部到上部(脚-头方向)、或者上部到下部的方向(头-脚方向)。本实施例中Z轴方向优选为(脚-头方向)。对于包含下肢血管的图像，例如腹主动脉血管，血管半径沿血管走向时逐渐变小，所述沿Z轴方向上的终止层切片图像中为腹主动脉的起始层，利于后续种子点的确定。

本实施例中，基于所述切片图像中像素点的灰度信息，所述灰度信息可以为以包括纹理结构、灰度、平均灰度、强度、颜色饱和度、对比度、亮度等一种或多种的组合。本实施例中所述灰度信息优选为图像中各像素点的灰度(例如CT值)，可以基于基于进行区域生长，获取若干个连通域，例如优选阈值为200。其中，所述阈值的设定可以是全自动、半自动、或手动进行的。例如，可以根据一种或多种运算自动计算或选取阈值。又如，用户或操作者可以通过输入/输出设备的图形用户界面手动选取阈值。

S20：基于所述连通域的圆度，确定种子点。本实施例中以腹主动脉为例，由于腹主动脉空间形态上为粗状的圆柱形，在切片图像中显示为圆形截面，可以根据圆度确认腹主动脉所在的连通域，继而确定所述种子点。种子点的确定方法结合附图2作具体描述。

S21：选取所述切片图像中任意一个连通域，计算所述连通域中所有像素点的边界距离场值；所述距离场值为连通域上的像素点距离最近的边界点的像素宽度，距离边界越远则距离场值越短。本实施例中，所述距离可以为体素宽度，也可以根据距离场公式计算。例如所述连通域中一个像素点的坐标为(x，y),与它最近的边界上的一个的坐标是(i，j)，则点(x，y)的距离场值近似为

S22～23：选取最大边界距离场值为该连通域的半径r，根据圆面积公式计算S＝πr²该半径对应的圆形面积S，并计算所述圆形面积与所述连通域的实际面积S’的比值，所述比值为所述连通域的圆度c。以此类推，获取所述切片图像中所有连通域的圆度。

S24：选取所述切片图像中圆度中与1差值的绝对值的的最小值所在的连通域，即选取圆度最接近1的连通域，圆度越接近于1，说明该连通域的形状越近似于圆形。对于诸如下肢血管的切片图像中，终止层中主要血管为腹主动脉血管，截面图像为较为粗壮的圆形截面，可以根据圆度确认腹主动脉血管所在连通域。然后确定该连通域中最大边界距离场值所在的像素点为所述种子点。

S30：以所述种子点为起点，基于阈值a对所述三维图像进行区域生长，获取第一图像。本实施例中，采用S20步骤确定的种子点为起点，对整个三维图像以阈值a作区域生长，本实施例中，所述阈值a选取灰度值为200，所述第一图像包含大部分血管和骨骼组织构成的多个三维连通域。

S40：继续根据所述种子点基于第一阈值对所述第一图像进行区域生长，获取血管掩膜。本实施例中，获取所述血管掩膜的方法结合图3作具体描述。

S41：计算所述第一图像中各体素点的边界距离场值；本实施例中，计算所述第一图像中各体素点的三维边界距离场，距离边界越远，则该体素点的距离场值越大。本实施例中，所述距离可以为体素宽度，也可以根据距离场公式计算。例如在所述第一图像中的任一各三维连通域中，取体素点的坐标为(x，y，z),与它最近的边界上的一个的坐标是(i，j，k)，则该体素点(x，y，z)的距离场值近似为

S42：根据步骤S20确定的所述种子点为起点，基于所述第一阈值进行区域生长，本实施例中，所述第一阈值选取距离场值大于3.0。由于在所述第一图像中，血管半径沿走向越来越小，而此部分的血管与骨骼相连部分的体素点的距离场值较小，可以通过本步骤控制主要血管不生长至骨骼组织，同时也无法生长细小血管。

S43～S44：对步骤S42中区域生长结果作膨胀操作。由于S42步骤为避免生长到骨骼，区域生长结果获取的是主要血管(例如腹部主动脉)的部分区域，为获取完整的主要血管，需要进行膨胀操作。本实施例中，采用以一定体素点宽度(例如选取4个体素点的宽度)为核对所述区域生长结果进行膨胀，获取所述血管掩膜，所述血管掩膜包含所述第一图像中主要血管组织，如图4所示。

S50：将所述第一图像与血管掩膜相减获取第一减影图像；所述第一减影图像中包含主要骨骼组织，和少部分的细小血管，而主要血管区域已通过相减去除，即将该体素的灰度值设置为背景区域。

S60：粗分割所述第一减影图像，获取初始掩膜；本步骤用于粗提取骨骼组织，下面结合图5作具体描述：

S61：计算所述第一减影图像上所有体素点的边界距离场值；计算方法如前所述，取体素点的坐标为(x，y，z),与它最近的边界上的一个体素点的坐标是(i，j，k)，则该体素点(x，y，z)的距离场值近似为

S62～S65：选取边界距离场值大于阈值b的体素点为骨骼种子点，由于S50步骤获取的减影图像中去除了血管上距离场值较大的点(例如腹主动脉等主要血管)，剩余为较细的血管和主要的骨骼组织，所以减影图像中边界距离场值较大的体素点主要位于骨骼区域，因此可以选取阈值b的体素点为骨骼种子点，根据所述第一减影图像中体素点的灰度信息或边界距离场值进行区域生长，获取所述初始掩膜，即骨骼组织的粗分割结果。

本实施例中，所述阈值b为大于3.0个体素点的宽度；所述体素点的灰度信息述灰度信息可以为以包括纹理结构、灰度、平均灰度、强度、颜色饱和度、对比度、亮度等一种或多种的组合，本实施例中优选体素点为灰度值大于500进行区域生长，用以去除细小血管。或者选取体素点的边界距离场值为大于1.5个体素点的宽度进行区域生长，可以去除贴近骨骼组织的细小血管。

S70：形态学处理所述初始掩膜，获取骨骼掩膜；本步骤用于去除前述步骤中提取初始掩膜时，可能生长到的血管组织，下面结合图6作具体描述：

S71：以阈值c为宽度腐蚀所述初始掩膜，获取若干个连通域；由于前述获取初始掩膜步骤中，基于灰度信息可能会由于造影剂原因而使得血管部分亮度增强(灰度值变大)而导致误提取血管，或者其它原因导致基于距离场值误提取的血管部分，可以先通过腐蚀操作，断开血管与骨骼组织的连接。本实施例中，所述阈值c优选为5个体素点的距离。

S72：选取最大体积的连通域为骨骼区域；因初始掩膜为骨骼组织的粗分割结果，主要血管(例如腹主动脉、椎动脉等主要血管)已通过所述第一减影图像去除，所以经S71腐蚀后获取若干个连通域中，最大体积的连通域可以确定为骨骼区域。

S72：以阈值d为宽度膨胀所述骨骼区域，以恢复因腐蚀步骤而丢失的骨骼组织，获取第一骨骼掩膜，所述第一骨骼掩膜为去除误生长血管的骨骼组织提取结果，如图7a所示；所述阈值d的选取可以与所述阈值c相同，本实施例本实施例中，所述阈值d优选为5个体素点的距离。

S73：由图7a所示，由于骨骼组织在图像中显示为中空结果，内部存在空洞情况，因此需要进行填实操作，可以包括如下步骤：

根据S20获取的所述种子点为起点基于阈值e对所述三维图像进行区域生长，获取第二图像；本实施例中，所述阈值e小于所述第一图像选取的阈值，用以获取更多的骨骼信息，例如由于扫描问题，或者骨骼内部的灰度信息低于骨骼组织等原因导致所述第一图像中骨骼区域存在部分缺失，例如，在CTA中，图像亮度的重叠，成像结构非常接近，有限的检测器分辨率，血管钙化点钙化，以及介入装置(例如植入的血管支架)可使对骨使得骨骼和脉管系统结构的识别和分割变得困难起来。可以通过降低阈值，即以所述阈值e进行区域生长补充缺失的各个信息。接着，逐层填实所述第二图像中各连通域。

本实施例中，逐层填实可以如下操作：选取一层切片图像，基于第一骨骼掩膜的结果，对于背景区域进行区域生长获取背景区域的二值图(例如背景区域的灰度值为1、骨骼区域所在连通域为0)，然后进行取反操作，即背景区域的灰度值设置为0，骨骼区域所在连通域为1，从而达到填实该层切片图像中骨骼区域的目的，依次类推逐层进行填实操作，获取填实后的第一图像。如图7a～7c为逐层填实操作的各阶段结果示意图。

S74：将填实后的第二图像与所述第一骨骼掩膜进行图像叠加，获取所叠加图像，本实施例中，所述图像叠加基于填实后的第二图像与所述第一骨骼掩膜中对应连通域的部分进行叠加，用以填实第一骨骼掩膜中可能存在的空洞部分或者缺失的骨骼边缘。由于叠加图像中提取的连通域可能存在几个体素点距离的空洞，为了进一步完善提取结果，还可以对所述叠加图像进行闭操作。

S75～S76：提取所述叠加图像最大体积的连通域，如前所述该连通域即为去除误提取血管的骨骼区域所在的连通域。对该连通域经闭操作获取所述骨骼掩膜。本实施例中，所述闭操作包括：首先以一定尺寸膨胀该连通域，例如选取尺寸为10个体素点宽度，接着以一定尺寸腐蚀该连通域，以恢复该连通域，达到填实小空洞的目的。需要说明的是，S74步骤获取的叠加图像是为了填充骨骼区域中较为明显(即尺寸较大)的空洞，而本步骤中的闭操作是为了填充细小空洞，例如几个体素点宽度，使获取的所述骨骼掩膜内部完全填实。

S80：将所述第一图像与所述骨骼掩膜相减获取第二减影图像，所述第二减影图像包括若干个连通域；以排除骨骼掩膜中骨骼组织所在的所有体素点，所述若干个连通域包括主要血管断开部分和一些散骨，血管断开部分可以包括一条血管在图像分割过程中所形成的断开的(或者不相连接的)部分。在一些实施例中，一条血管可以具有一个或多个血管断开部分。一个血管断开部分可以具有关于这条血管的一个或多个数据。两个或多个血管断开部分之间可以具有不同数量的血管。

S90：基于血管生长条件连接所述若干个连通域，获取血管提取结果。血管的提取可以根据骨骼掩膜、和/或所述第二减影图像中的血管断开部分(即若干个连通域)。在一些实施例中，血管提取可以根据两个或多个血管断开部分，在所述两个或多个血管断开部分之间生成血管，并采用一种或多种方法，连接所述两个或多个血管断开部分，生成一条或多条血管。下面结合图8作详细说明。

S91：计算所述第二减影图像包括若干个连通域的边界距离场值。

S92～S93：基于最短路径算法，并且排除所述骨骼掩膜中的各个连通域所在区域生长血管中心线，即中心线的走向绕开所述各个掩膜中标记为骨骼组织的区域，获取所述血管提取结果。本实施中，所述最短路径算法可以包含Dijkstra算法、Bellman-Ford算法、A*搜索算法、Floyd-Warshall算法、Johnson算法、Viterbi算法等。以Dijkstra算法为例，该算法中的价值函数可以从补充区域的距离变换以及数据补充位点间的距离推导。如图9所述为本实施例中血管提取结果示意图。

综上所述，本发明提供一种血管提取方法，通过基于所述连通域的圆度确定种子点，根据所述种子点基于阈值进行三维区域生长，获取包含血管和骨骼的第一图像；将所述第一图像与所述骨骼掩膜相减，获取包含有主血管和少部分碎骨的连通域的第二减影图像；基于血管生长条件，排除所述骨骼区域，通过多区域区域生长的方式连接所述血管连通域，有效获取与骨骼粘连的血管。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。