宋 静， 何启泠， 徐 庆， 王海江

(成都信息工程大学电子工程学院，四川成都610225)

0 引言

医学影像是现代医学诊断和辅助治疗最重要的技术手段之一，使医生能够对人体解剖结构及功能的认识更加直观与全面，影像诊断已经广泛应用于医疗诊断、治疗、术前计划和术后监测的各个环节。随着互联网医疗的快速发展，医学影像诊断相关的应用越来越朝着轻量级的方向发展，使医学图像可以在多种应用中展示与浏览而不是依赖于专用的图像浏览软件。通过各种医学图像处理方法帮助医生提高诊断的准确性的同时，也在思考如何更加方便医生进行诊断，如更方便与快捷地远程影像诊断。图像分割是图像处理中的关键性技术，对图像进行进一步分析和处理的前一过程基本都依赖于图像分割。图像分割是指根据图像的不同特性(图像的灰度、颜色、纹理等)分割成不同的区域，这一处理有助于提取目标图像，如某一组织或是病变位置。

在互联网医疗时代，大量的医学图像通过传输完成诊断，交流探讨。在这一过程中对于图像传输的质量和效率的控制直接影响互联网与医疗的结合度以及用户体验。大量的医学图像对诊断有用的区域与图像的实际大小总是存在一定的差值，就会有大量的无效数据加大存储量与降低传输效率。另一方面在三维重建的过程中，图像中的非有效成像区域对于三维可视效果也有极大的影响，改善这两点正是研究的意义所在。ITK主要提供了医学影像图像分割和配准方面的功能，但是在可视化方面还要借助VTK(visualization toolkit)观察结果以及进行交互显示。

图像分割在组织提取或区域提取方面已有很多研究，文中主要从基于图像分割的实际应用方面出发，讨论如何降低数据的传输量，提高图像传输速率，这一点对于开发PACS系统和研究远程医疗具有一定的意义。此外大量的医学影像存在设备干扰影像，在单幅影像中，可能这些干扰影像还不影响医生进行诊断，但是在对这些医学影像进行三维重建或者一些特定处理时，将严重影响成像效果。文中也讨论图像分割在这方面的应用。

1 ITK简介

ITK是一款医学图像处理软件包，是一个开源，跨平台的影像分析扩展软件工具。ITK采用先进的多模态数据分割配准算法，是一个用于处理医学图像的开源软件库，其中有丰富的图像分割与配准的算法程序，并且支持多种开发平台和开发语言，方便开发者进行图像处理方面的研究与应用。

ITK作为一个开源，跨平台的图像分析框架，里面大量的前沿算法实现可以广泛用于图像配准和分割［1］。ITK的主要特征有:(1)ITK提供了通用的方式表示图像(任意维度)和面片(非结构化的meshes)，提供了大量用于进行分割和配准的算法，主要面向医学应用。(2)ITK只提供最小化的文件接口而不提供可视化的用户接口，留由其他库提供。(3)大量使用泛型编程技术。(4)内存模型使用智能指针维护对象的引用计数，使用对象工厂实例化对象。(5)支持多线程并行处理。(6)使用命令/观察者模式进行事件处理。(7)基于数据流的架构进行组织，数据被表示成数据对象，数据对象由处理对象处理，数据对象和处理对象连在一起形成Pipeline。(8)ITK主要采用管道模块结构设计，在对于图像处理后需要将图像文件输出，需要借助VTK显示图像。对于此类需求可通过ITK，VTK联合编程实现。

许多图像处理方法可以在MATLAB中实现，但ITK比MATLAB而言优点主要体现在几个方面:(1)ITK提供了大量的完备的与可灵活调用的图像处理算法，如图像分割，图像配准与滤波等。(2)有效地将ITK的图像读取与处理功能与VTK在图像可视化方面的突出能力结合起来可以实现强大的图像处理功能。(3)ITK作为一个开源的平台，会随着时间不断的丰富并且应用会有更好的延伸性。文中的主要研究内容是基于图像分割的实际应用，重点在应用性的拓展上而不是图像处理算法的实现上，因此ITK平台相对于MATLAB而言是更佳的选择。

2 医学图像分割算法应用

用于诊断的各种扫描图像中，包含大量的信息，如人体器官，病变组织及某些物体。为分析它们的形状、形态和功能，或者对它们进行可视化处理，都需要先将其从所在的情景中提取出来。这一提取过程即为医学图像分割。医学图像分割算法包括阈值分割、区域增长、基于分水岭的分割、FastMarching算法、LevelSet(水平集)等多种分割方法［2］。文中使用的方法主要是区域生长与阈值分割。由于篇幅有限，主要介绍使用到的两种方法以及其在应用过程中所作的调整。为提取目标图像区域，采用先对原图进行二值化处理，再对图像进行区域生长，得到最大的连通区域，除去原图像中的设备干扰成像。

阈值分割就是确定一个或多个阈值，根据图像中像素的灰度值与确定阈值的大小关系进行一个映射，从而达到分割效果［3］。设原图像f(x，y)，分割后的图像为g(x，y)，新得到的图像为二值图像。原图和分割后的图像与阈值关系为

在这一应用中，根据二值图像可以大体确定图像所在的区域，研究常见图像的提取情况，根据二值图像提取目标图像如图1～图6所示(图像均为提取数据在网页显示效果截图)。

图1 原图

图2 二值图像

图3 根据二值图像提取图

图4 原图

图5 二值图像

图6 根据二值图像提取

根据结果发现在提取无干扰成像的图像时，仅根据二值图像来提取就可以达到预想的目的，但是当图像中有干扰成像时提取的情况还不够理想，根据这一情形想到再提取图像前还需要对二值图像作进一步的处理，达到提取的目标数据在保证提取到完整图像的同时要尽可能的减少数据量。通过分析大量医学图像可见医学影像图像的目标总是在图像的中心而不是边缘或角落。把目标当作一个完整的连通区域，干扰图像就可被屏蔽，通过实验发现ITK中的区域生长中的连通阈值提取的方法正好可以解决这一情况。

区域生长的基本思想是提取有相似性质的像素从而构成区域，对每个区域都要选取一个种子点，以和种子像素相似或相同的性质为准则在种子点周围开始生长，合并满足条件的点，然后将这些点再作为种子开始生长，直到没有满足条件的像素点为止停止生长。在这一过程中有3个关键点:分别是种子点的选取;生长准则的确定;区域生长的停止条件［4-5］。

种子点的选择:应用中，由于要求在应用过程中是后台自动处理生成，这一点决定的种子点是不能通过交互来得到的，因此文中选择总是从图像的中心开始生长。这一方法也会产生一个新的问题，即如果图像的中心在二值化的过程中是一个不可生长的像素如何处理，针对这一情况，做了两个帮助生长的处理。首先由根据图像二值化得到的二值图像进行统计，可得到有效像素点数与整个图像中像素点数的比例，根据这一比例在图像中心构造一个有效区域，确保生长顺利进行。这一方法是基于医学影像的目标区域总是在图像的中心。第二个帮助生长的处理是因为在生长过程中发现有一些希望得到的小区域会在生长过程中丢失，于是就采取先对二值图像作膨胀处理，确保在生长过程中不会丢失掉一些边缘信息。表1为中心构造区域大小与图像有效像素与所有像素的比值关系。

表1 有效区域比例与中心构造区域面积关系表

在实验过程中，根据以上介绍的方法可以准确提取出目标图像，并且在无干扰成像的情况下依然可以准确提取出目标图像数据，如图7～图9所示。

图7 原图

图8 连通区域图

图9 根据连通区域图像提取图

3 基于ITK中图像封装类的实现过程

用ITK进行医学图像的处理会比较方便，因为大量的图像处理算法已经在ITK中封装好，在编写程序的过程中只需调用相应的库文件就可加以应用。目的主要是为了讨论如何提取出有效的数据信息，从而可以在网页绘制过程中更快，在传输过程中更快。图10为基于ITK实现这一图像处理过程的流程图。

图10 图像处理过程流程图

3．1 图像二值化处理

用二值化处理过程的目的主要是为后续的处理过程，如膨胀处理，目标连通区域生成［6］。选择Otsu方法来生成二值图像，作为一种常用的灰度图像二值化的方法，用一个基于直方图统计得到的灰度阈值来对图像进行分割，从而把图像分为目标与背景两部分，可以帮助确定目标区域在图像中的所在位置。在ITK中已经把这一方法封装，直接调用方法即可，其阈值是基于统计自动得到，从而可自动对图像进行分割。

用ITK中的itk::OtsuThresholdImageFilter类对输入图像进行二值化处理。其中把SetOutsideValue()设置为0，SetInsideValue()设置为255，分别代表处理后阈值范围外/内的像素值;SetNumber of Histogram Bins()设置为256表示用于计算分割阈值的直方图的灰度级数。

const unsigned int Dimension=2;//定义数据的维数

typedef signed short InputPixelType;//定义读取对象的类型

typedef itk::Image＜InputPixelType，Dimension ＞InputImageType;//定义输入图像的类型

typedef itk::ImageFileReader＜InputImageType＞ReaderType;

ReaderType::Pointer reader=ReaderType::New();//建立一个读取对象

reader-＞SetFileName(InFile);//读取图像文件

typedef itk::OtsuThresholdImageFilter＜InputImageType，InputImageType ＞FilterType;

FilterType::Pointer filter=FilterType::New();//定义FilterType指针

filter-＞SetInput(reader-＞GetOutput());//获得读取的CT图像

const InputPixelType outsideValue=255;//设置分割的低阈值为255

const InputPixelType insideValue=0;//设置分割的高阈值为0

filter-＞SetOutsideValue(outsideValue);

filter-＞SetInsideValue(insideValue);

filter-＞SetNumberOfHistogramBins(256)

filter-＞Update();

3．2 二值化图像膨胀处理

图像膨胀就是用固定结构元素在二值图像上移动，扫描图像中的每一个元素，即用结构元素中的像素与图像像素做“与”，如果相与结果全部为0，则该像素置0，反之置1［7-8］。在这里加入膨胀处理的原因是由二值图像直接得到目标连通区域时可能会丢失一些边缘，通过膨胀处理就可以使目标增大，这样在生成连通区域图时就尽可能的包含图像本身的信息，尽量避免在由二值化到连通区域生成的过程中错误的丢失掉原本属于目标区域的图像。用ITK中的 itk:BinaryDilateImageFilter类对上一过程输出的二值图像进行膨胀。这一过程主要是为的保存目标图像的一些边缘信息。将结构元素的半径设置为5，调用方法为structuringElement.SetRadius()。通过 GrayscaleDilateImageFilter类定义了一个膨胀处理对象grayscaleDilate，因此膨胀处理后的输出为grayscaleDilate-＞GetOutput()。

原图不经过膨胀处理和经过膨胀处理提取的图像效果如图11～13所示。

图11 原图

图12 未经过膨胀处理的提取图

图13 经过膨胀处理的提取图

可看出没经过膨胀处理得到的结果可能会丢失掉一些边缘信息，所以在处理过程中加入了膨胀处理这一过程，提高图像提取时的完整性。

3．3 得到目标图像连通区域

ITK中的ConnectedThreshold分割算法可以得到目标图像的连通区域［9-11］。图像中存在的噪声有可能会降低滤波的效果因此先对图像进行一个平滑处理，接下来选取种子点如上面介绍选取图像的中心点作为种子点。还要通过SetLower()，SetUpper()两个方法来设置图像分割的高阈值和低阈值，这里设置时要特别注意取值范围要包含二值图像的有效值的灰度。SetReplaceValue()设置新的有效值。SetSeed()加载种子点索引。其中size代表原图像的行列数。

typedef float InternalPixelType;//定义用于区域生长的对象类型

typedef itk::Image＜ InternalPixelType，Dimension ＞InternalImageType;//定义用于区域生长的图像类型

typedef itk::ConnectedThresholdImageFilter＜InternalImageType，InternalImageType ＞ ConnectedFilterType;

ConnectedFilterType::PointerconnectedThreshold=ConnectedFilterType::New();

//定义ConnectedFilterType类型指针

connectedThreshold-＞SetInput(caster-＞GetOutput());

//将膨胀处理后输出转换数据类型后作为提取连通区域的输入即caster-＞GetOutput()

connectedThreshold-＞SetLower(125);

connectedThreshold-＞SetUpper(500);

connectedThreshold-＞SetReplaceValue(255.0);

InternalImageType::IndexType index;//定义种子点

index［0］=size［0］ /2;//种子点设置为图像的中心区域

index［1］=size［1］ /2;

connectedThreshold-＞SetSeed(index);

connectedThreshold-＞Update();

3．4 根据连通区域提取目标区域图像数据

通过以上几步处理就可得到目标区域的二值图。扫描这二值图像可以得到图像起止点矩形区域的坐标(x1，y1，x2，y2)，根据这一矩形区域坐标索引从原图中提取数据，提取出来的数据可直接在页面绘制，这一数据为图像的原始数据所以不会影响后续的图像操作如调窗、反显、伪彩等，后续如果还要进行图像分割处理也不会有任何影响。只是通过这一操作节约前端处理的操作时间和传输时间。其中原始像素数据可用迭代法提取，速度会快于遍历提取，主要代码如下:

选取不同医院设备提取的CT，MR(dicom图像)进行实验，处理前后数据量见表2。

表2 图像数据提取前后大小对比

这一处理数据量的变化不是固定不变的，与目标区域和原图的比例紧密相关，从总体看对数据量的减少能起到很大作用。这一处理还可以用在三维模型生成前的图像处理，主要是去除在成像过程中设备成像的干扰。这里对ITK生成三维模型［12-13］就不赘述，实验效果如图14、15所示(由VTK显示的面绘效果)。

图14 未处理过的图像模型

图15 处理过后的图像模型

4 结束语

ITK强大的图像处理功能可以方便开发者快捷的编制程序，提高开发效率。图像压缩可以解决图像传输效率问题，但是在医学图像的处理过程中对图像的无损性要求相对更高，以上介绍的提取数据的方法对于医学影像的目标区域而言始终保存原始数据，因为是从原始图像中提取区域数据，对于单个图像而言可能有更简化的方法提取，但是对于要处理大量医学图像的软件应用层次而言文中介绍的方法具有一定的适用性。无需交互特性使得应用范围可以得到扩展。

［1］ 张芳，李强，王阳萍.ITK及其在医学图像分割中的应用［J］.微计算机信息，2008，24(18):304-306.

［2］ 徐彩云.图像分割算法的研究综述［J］.电脑知识与技术，2014，(11).

［3］ 王磊.多维Otsu方法在图像分割中的研究［D］.济南:山东师范大学，2009.

［4］ 李敏.数字人脑切片图像自动分割算法的初步研究［D］.重庆:重庆大学，2010.

［5］ 邓金莲，李成贵.结合区域生长和 SUSAN的图像分割方法［J］.电光与控制，2012，19(8):38-40.

［6］ 常江.车辆图像局部识别［J］.信息网络安全，2013，(11):94-96.

［7］ 黄颖，杨光琼.结合小波系数的normalized cut分割算法［J］.计算机应用，2011，31(1):182-183.

［8］ 乔诚，李连东，王秀丽.ITK及其在分水岭医学图像分割中的应用［J］.太原科技，2008，(2).

［9］ 牟春洁，张国英.基于区域边界生长的图像分割方法［J］.北京石油化工学院学报，2009，17(4):8-12.

［10］ 吕晓琪，任晓颖，贾东征.基于 ITK，VTK和MFC的DICOM图像读写及显示［J］.中国组织工程研究，2011，15(13):2416-2420.

［11］ 伍云智.基于VTK和ITK算法库的研究与应用［D］.南昌:南昌大学，2010.

［12］ 袁杲，杨玲，朱小波，等.利用ITK和VTK集成实现三维医学图像的分割［J］.计算机应用与软件，2009，26(2):248-250.

［13］ 陈保远.基于ITK和VTK医学图像配准方法的研究与应用［D］.南昌:南昌大学，2012.