基于广义势场的肝脏管道骨架可视化方法研究

叶 楠，陈国栋

(福州大学物理与信息工程学院，福建福州 350002)

随着计算机图形技术及虚拟仿真技术的不断发展，计算机辅助技术在肝脏疾病诊断和治疗上的应用也越来越广泛［1］。如何能有效表现和分析出血管管道结构已成为一个重要的研究热点。肝脏管道的可视化研究主要包括血管的分割、血管骨架的提取、管道的三维重建以及可视化的显示。其中，骨架作为几何形态的一种很重要拓扑描述方式，能简洁地表示物体的基本形状结构特征，减少物体描述的信息量。在一个血管管道可视化系统中，骨架的表示方式起到了关键性的作用。

目前，在血管管道的应用上，三维的骨架提取算法主要有4种类型:基于细化的方法、基于距离场的方法、基于Voronoi图的方法和基于广义势能场的方法［2］。这些方法存在各自的优缺点［3］:细化法能较好地保留骨架的拓扑结构、居中性和连通性，不足在于不能得到精确骨架点;距离场法能够满足骨架的几何特性，但在物体细小狭窄的地方难以保证骨架的连通性和完整性，且中心性容易受边界的影响;Voronoi图法可以得到较为满意的拓扑结构和几何特性，缺点是在时效上并不乐观，不适合数据量较大的物体;广义势场法［4-6］是从场的方向出发，生成的骨架具有良好的连通性和光滑性，只是在时间复杂度上较高。实际使用时，人们往往根据所需的应用领域选择相应的方法。肝脏血管管道的三维数据集是由一系列的二维切片构成的，数据量较大;管道的拓扑结构呈树形结构，在血管末梢上，血管半径较窄，为了使骨架能够精确表达出肝脏血管管道的拓扑结构、管道的半径等信息，对肝脏血管管道骨架的提取提出了以下要求［7-8］:

1)骨架能够保存原始肝脏血管管道模型的拓扑结构;

2)骨架具有良好的居中特性，骨架线位于原始肝脏血管管道中心线附近;

3)骨架应具有连通性、平滑性和单像素宽;

4)骨架算法应适用在三维空间中;

5)骨架算法对三维物体的位置和噪声不敏感，具有一定的鲁棒性。

针对以上各算法比较结果和骨架要求，本文采用了基于广义势场法的肝脏管道骨架可视化仿真方法。对原始数据预处理分割出血管后，三维重建获得肝脏血管管道的形态结构模型，利用广义势场提取出管道的中心路径关键点以及相应的血管半径，继而跟踪关键点得出管道骨架，然后结合OpenGL实现可视化。本方法优点在于能适用在三维空间中，充分保留了血管管道的拓扑结构，并保证了骨架的平滑度。实验结果证明了该方法的有效性和可行性。

1 预处理工作

1.1 材料准备

美国数字化可视人体数据集来源于美国国立医学图书馆提出的“可视人体计划”工程中采集到的人体截面图像，本文采用该数据集中的肝脏数据作为数据源［9］。

1.2 血管分割

肝脏血管图像的分割是肝脏管道三维重建过程中的第一步，也是至关重要的一步。在诸多分割算法中，区域增长算法是图像分割中简洁而有效的方法，由于血管图像的连通性，以及血泡的封闭性，使得区域增长算法适用于血管分割。其中，连接门限阈值法(Connected Threshold)就是区域增长法中的一种。

连接门限阈值法是采用了注水迭代来实现访问种子像素的邻域。图形学中注水是满水法填充，是用来填充区域的［10］。用户可以提供一个特定的灰度区间［最小阈值，最大阈值］，当种子点的邻近像素的灰度值落入上述区间内时，该算法会将此像素包含入生长区域中，通过反复迭代直至没有新的区域增加为止。

由于采用的数据集具有彩色成分图像，因此，为了增强血管强度，先对数据集进行灰度处理，并进行图像滤波消除噪声，可减少对分割结果的影响。

在本项研究中，按照上述方法，以Visual Studio 2010为开发环境结合C++和ITK开发包对数据进行了血管分割处理，下面为编号为n=1 509的图像分割结果。图1a为原始图像，图1b为灰度图像且去噪后的图像，红色的点即为种子点，图1c为采用［10，33］阈值区间得到的分割结果。

图1 连接阈值法分割肝脏血管图像

1.3 肝脏管道的血管重建

数字化虚拟肝脏管道骨架可视化研究的关键是三维重建，建立具有立体空间、各方位旋转的肝脏管道模型［11］。本项研究在上述血管图像分割的基础上，使用VTK(Visualization ToolKit)开发工具包对图像序列进行了面绘制的三维重建。主要是通过VTK提供的vtkBMPReader读取肝脏血管分割后的一系列二维图形，产生三维信息，并采用vtkContourFilter滤波器抽取等值面为255的血管，而且为了提高渲染的速度和效果，在等值面上利用vtkPolyDataNormals产生一系列的法线，增强平滑感，最后利用窗口进行渲染显示。

本项实验一共采用了VCH-F1的147幅600×630肝脏数据集图像，灰度等级为255。三维重建后的肝脏管道整体结构图如图2所示。

图2 三维重建后的肝脏管道整体结构图

2 基于广义势场的骨架化

文中提出的广义势场法是以静电斥力场为基础，对三维空间中肝脏血管管道骨架的提取。通过生成斥力场，对体素进行归类，检测关键点，直至生成单像素宽的骨架。下面对算法给出更近一步的描述。

2.1 斥力场源的构造

算法基于离散体素的表示，斥力场源是通过在体数据表面边界点上放置同种点电荷来生成的。对于体素归为三类:外部点、内部点和边界点。体素值为0的即为外部点，体素值有非零值即为内部点，若某内部点体素的6邻域中至少有1个值为0，则该体素为边界点。由于肝脏血管管道数据集的生成有人工掺杂、噪声影响以及分割和重建算法自身的不足，因此，体数据中或多或少地存在一些孔洞和内腔。为了提高计算结果的精确性，利用形态学的孔洞填补方法进行处理，可保证体数据中不存在孔洞。

2.2 三维势场值

根据上述的约定，扫描体数据，把边界点坐标都存放于数组B中，现假定内部点中存在一点A，如果有一边界点与该A点的3个方向上的绝对距离都小于等于阈值(这里取值为100)，则认为该边界点对点A有产生场强，若大于该阈值，则可忽略该边界点的作用。点A处的势场值可定义为对该点A有产生场强的边界点C与该点到边界点的欧氏距离R的m次幂成反比，公式表达式如下

式中:FAC为点A处的斥力;CA为从边界点C到点A的方向向量;R为C与A间的距离;m为力的阶数(若m=2即为牛顿斥力)，m的值越大，边界点对内部点的作用就越大，反之，m的值越小，产生的势场就越平滑。由于一体素共有26邻域邻接，所以在点A处的势场值由其26邻接点的平均值决定，最后的势场值可归一化为

式中:数组B存放边界点的坐标。每个边界点上的点电荷对每个内部点产生斥力的总和即称为斥力场。

2.3 关键点与“核骨架”

本部分的输入即为上述归一化后的势力场，输出即为所需的关键点列表。所谓的关键点是斥力场为0处的点，它是矢量场的表示拓扑结果的方式之一，也常作为矢量场可视化的表现内容。计算关键点处势场的Jacobian矩阵特征值和特征向量，根据特征值实部和虚部正负符号不同，关键点可分为三类:当所有特征值实部为负时，该关键点为吸引点;当所有特征值实部为正时，该关键点为排斥点;当特征值的实部和虚部有正有负时，称该关键点为鞍点。核骨架的生成将从鞍点开始，具体步骤如下:

1)遍历所有的鞍点，每个骨架段以鞍点为起始点生成;

2)从鞍点出发，沿着正特征值对应的特征向量为方向，按一定步长前进;

3)若计算出的当前位置和起始位置相同，则算法结束;

4)若下一点是某骨架段A中的点a，则将该点加入到当前段，并且作为段结束的端点，将A段在点a处分开，即A段分成了两段;

5)若下一点是关键点，则直接加入当前段，并且作为段结束的端点;

6)若下一点不是骨架点或关键点，将该点归进当前骨架段，沿着该点的合力方向(即该点的势场值)，继续前进;若前进过程中步数过多，超出规定的步数值，则算法结束。

算法流程图如图3所示。

图3 算法流程图

经过上述的计算，最终输出核骨架的骨架点。

3 骨架可视化及结果分析

利用Visual C++和OpenGL在计算机上实现了上述算法，计算机的配置为Core 2 Duo E4700 2.6 GHz的CPU，2 Gbyte的RAM，NVIDIA GeForce GT 430显卡。对骨架化后的肝脏管道可视化结果如图4所示。

算法骨架提取所需的时间统计信息如表1所示。

从表1中可以看出，在计算三维势场值时需要占用较长的时间，这是因为在对体数据进行读取后，需要遍历各个边界点，边界点的数量多少将决定了算法时间的长短。

4 结论

本文首先对肝脏血管管道做了一系列的前期工作，进行了血管分割和三维重建，继而提出了一种有效的三维骨架化方法。该方法假定在物体边界点上放置同种点电荷作为斥力场源，然后通过计算势场值检测关键点，得出单像素宽的血管骨架点。实验结果表明了该方法的可行性和有效性。

图4 对骨架化后的肝脏管道可视化结果

表1 肝脏血管骨架化统计信息

:

［1］翟朝亮，陈国栋，王娜，等.基于体纹理的肝脏可视化仿真方法研究［J］.电视技术，2012，36(17):169-172.

［2］CORNEA N D，SILVER D，MIN P.Curve-skeleton properties，applications，and algorithms［J］.IEEE Transactions on Visualization and Computer Graphics，2007，13(3):530-548.

［3］CHEN Y，DRECHSLER K，ZHAO W，et al.A thinning-based liver vessel skeletonization method［C］//Proc.Conference on IEEE 2011 International Internet Computing ＆ Information Services(ICICIS).［S.l.］:IEEE Press，2011:152-155.

［4］CHUANG J H，TSAI C H，KO M C.Skeletonisation of three-dimensional object using generalized potential field［J］.IEEE Trans.Pattern Analysis and Machine Intelligence，2000，22(11):1241-1251.

［5］CORNEA N D，SILVER D，YUAN X，et al.Computing hierarchical curve-skeletons of 3D objects［J］.The Visual Computer，2005，21(11):945-955.

［6］马锐，伍铁如.基于广义势场的三维形体多层次线骨架构建［J］.计算机应用，2011，31(1):16-19.

［7］CHEN Y，DRECHSLER K，LAURA C O，et al.A graph description and analysis framework of liver vascular system［J］.Journal of Computational Information Systems，2011，7(6):1998-2005.

［8］王昌.高精度肺部气道树的分割及骨架中心线的提取［D］.合肥:中国科学技术大学，2010.

［9］美国数字化可视人体数据集［EB/OL］.［2013-03-10］.http://vhnet.nlm.nih.gov/.

［10］彭微.连接门限阈值法在肝脏CT图像分割上的应用［J］.咸宁学院学报，2011，31(6):72-73.

［11］方驰华，周五一，黄立伟，等.虚拟中国人女性一号肝脏图像三维重建和虚拟手术的切割［J］.中华外科杂志，2005，43(11):748-752.