SPIHT算法在DICOM图像压缩中的应用研究

龚正　沈建新

摘 要：DICOM图像文件的储存和传输一直是近代医学图像处理中一个比较热门的研究话题，其中较好的解决方案就是对DICOM图像进行压缩处理。文中通过对DICOM文件的数据结构进行分析，将DICOM文件拆分成图像数据部分和文本信息部分，然后采用9/7提升小波变换对拆分出来的DICOM图像数据进行图像变换，最后对变换后的小波系数进行SPIHT编码完成DICOM图像压缩。主客观对图像的评价表明SPIHT算法在高压缩比时仍能保证图像较好的质量。该方法能够对DICOM图像的存储和传输带来便利也将产生一定的影响。

关键词：医学数字图像通信协议标准；图像压缩；9/7提升小波变换；多级树集合分裂排序

中图分类号：TP391 文献标识码：A

1 引 言

影像储存和通信系统（PACS，全称Picture Archiving and Communication System）广泛应用于我国医院影像管理事业中[1]。PACS主要用来解决CT（计算机断层扫描）、MRI（核磁共振）、US（超声成像）等成像设备产生的图像数据传输和储存的问题。1985年ACR和NEMA联合建立了医学数字成像的DICOM（Digital Imaging and Communications in Medicine）标准[2]。DICOM文件的一个特点就是数据量巨大，其中图像数据占90%以上并且图像的数据量要比普通图像要大得多[3]。文中对DICOM文件的数据结构进行了分析并结合现代图像压缩技术，在确保一定图像质量的情况下，完成大压缩比的DICOM图像压缩工作。

2 DICOM图像压缩原理

DICOM标准的数据结构部分中定义了DICOM如图1所示的文件结构。DICOM图像文件是由文件序言、首缀组成的头文件部分和由多个数据元素组成的数据集部分构成[4]，每个数据元素由标识符（Tag）、数据类型（VR）、数据长度（VL）和数据组成。DICOM图像文件包含了许多图像的附加信息，如设备信息、图像所属患者信息、图像的窗位窗宽大小等，这些信息和图像数据一起存放在各个数据元素中。DICOM数据字典对每个数据元素都分配了一个唯一的Tag，其中图像数据的Tag为（7FE0，0010）[5]。对DICOM图像压缩的本质是将DICOM图像文件中的文本数据和图像数据分离出来并对DICOM文件中的图像数据进行压缩。在对DICOM图像文件解码的过程中，通过对图像数据Tag的定位可以将DICOM文件中的图像数据分离出来并转换成适合图像压缩的BMP文件。

图1 DICOM图像文件结构

此时，DICOM文件的压缩就变成了针对从DICOM文件中分离的出来的BMP文件进行压缩处理，经典图像压缩技术的一般步骤如图2中所示。首先是对原始图像进行变换，然后对变换系数进行量化编码得到压缩图片文件，这时就可以方便快捷的对压缩文件进行存储或者传输。解码器和编码器是一个互逆的过程，通过解码、去量化、反变换可得到压缩重构的图像文件。JPEG标准中使用的是DCT（离散余弦变换）变换，但由于DCT不能很好的处理它自身引起的方块效应和在保证图像质量的前提下不能达到更高的压缩比等原因而逐渐被舍弃[6]。因离散小波变换具有优良的时频特性（能够将低频和高频的信息进行分离）被JPEG2000标准所采纳[7]。

3 图像变换与编码方法

3.1 9/7提升小波变换

小波理论的思想来自伸缩与平移变换，将满足一定条件的小波基本函数Ψ（x）构造成形如（1）公式的构造函数。用尺度因子a将Ψ（x）做伸缩处理，小波函数图像会随a的变小而缩短，随着a的变大而伸长；用b控制Ψ（x）水平移动。

小波变换（Wavelet Transform，简写WT）与FT（Fourier Transform）在某种程度上类似，其定义如公式（2）中所示。小波变换可以看成把一个信号（信号也可以看作一种函数）表示为由多个或者无穷个经过伸缩平移变换而得到的小波信号的叠加。信号作小波分解实际上相当于用具有选择性的滤波器对信号进行滤波处理，它的要点就是对于平稳变化的信息（对应于信号低频部分）可以在大区域上查看；对于快速变化的信息（对应于信号的高频部分）可以在小区域上放大查看。

对于图像信号，可以抽象成二维离散函数F（x，y），二维离散函数的WT就可以通过一维离散行WT和一维离散列WT来实现二维离散一级WT，如图3所示。图像通过一维行WT就可以获得行低频（L）和行高频（H）的小波系数，再经过一次一维列WT就可以将图像的最低频系数集中在左上角LL1处，最高频系数集中在右下角HH1处，此时也就完成了图像信号的一级WT。二级WT就是在原来变换的基础上单独对LL1处的信息进行一次一级WT，多级WT则以此类推就可得到。WT都必须用到了小波函数，而对于二维离散函数的WT则是由小波基来完成。但由于这些小波构造多依赖于FT，构造过程繁杂且算法复杂占用较大内存。因此，1995年，Sweldens W.博士创造性地提出了在时域中采用lifting scheme方法构造小波函数的方法即第二代小波（second generation wavelet）构造方法[8]。这种构造小波方式与第一代小波不同，它不需要通过FT或者基于FT尺度伸缩，而是快速便捷的利用Split，Predict和Update等方法来完成对数字离散信号的小波变换。JEPG2000标准中采用了9/7提升小波和5/3小波提升小波来完成对图像的变换[9]。在确保图片质量的前提下，9/7提升小波变换获得时频特性更好的小波系数，所以这里选用双正交9/7提升小波。

9/7提升小波正变换也是分别通过小波行变换和小波列变换来完成小波变换的，小波行变换和小波列变换的公式相同都是分别通过四步提升步骤（预测1，更新1，预测2，更新2）和一个尺度变换来完成的。具体步骤如下：endprint

①分裂（split）：For[-1≤n≤（N/2）+1] do： { s_0（2n±1）=x（2n±1）； d_0（2n）=x（2n）}；

②预测1（predict1）：For[-1≤n≤（N/2）+1] do： d_1（2n）=d_0（2n）+[α（s_0（2n+1）+s_0（2n-1））]；

③更新1（update1）： For[0≤n≤（N/2）+1] do： s_1（2n-1）=s_0（2n-1）-[β（d_1（2n-2）+d_1（2n））]；

④预测2（predict2）： For[0≤n≤（N/2）] do： d_2（2n）=d_1（2n）-[γ（s_1（2n-1）+ s_1（2n+1））]；

⑤更新2（update2）： For[1≤n≤（N/2）] do： s_2（2n-1）=s_1（2n-1）-[δ（d_2（2n-2）+d_2（2n））]；

⑥尺度变换： For[1≤n≤（N/2）] do：

L（2n-1）= ε s_2（2n-1）；H（2n）=（1/ε） d_2（2n）

表中的相关系数α=-1.58613432， β=-0.05298011854，γ=0.8829110762，δ=0.443506822，ε=1.149604398，s_0和d_0表示原始像素值中的奇数部分与偶数部分。通过一级提升WT后得到了如图1中LL1、LH1、HL1、HH1频率高低不同的子带。图4所示，（a）为经典Lena图片，作为9/7提升小波变换的实验图片；（b）为Lena图片经过一级9/7提升小波变换而得到的分解图。

3.2 SPIHT编码算法

图像经WT后得到的的小波系数在各级子带中存在树状构架，如图5中所示。在这些小波树中，相同空间位置和相同方向上的小波系数有很强的相似性，HH3中的小方格里的的小波系数可以看成是HH2在相应位置小方格里的小波系数的父系。零树编码思想就是将量化后为0的小波系数的所有子孙小波系数都归为0，进而简化小波编码的过程。但是这种算法存在很大缺陷，因为并不

是父系小波系数为0的小波系数其子孙的小波系数都为0，仅仅通过零数编码可能会导致图像的失真比较严重。一种改进的方法是EWZ（嵌入零树编码）算法，它提出将每层小波变换设置一个阈值T（i），如果第i层的小波系数的绝对值大于该阈值，则认为这个小波系数的重要的，相反则认为该小波系数不重要。如果父系的小波系数为非重点小波系数，而且它的全部子孙小波系数都相应的被判为非重点小波系数，故此命名该父系小波系数是关于阈值T的零树根，进而将该父系小波系数的全部子孙小波系数不不列入小波编码的范围。

由于EWZ编码需要设置每层小波的阈值并且要搜索扫描每层小波系数绝对值小于阈值的所有后代，这样就导致了计算量很大，编码速度慢的后果。因此EWZ编码还需要进一步的改进，在目前所有改进的算法中，多级树集合分裂排序（Set Partitioning in Hierarchical Tree，简写SPIHT）算法是较为成功的一种[10]。

SPIHT算法的零树构架中则认定LLN没有子带。SPIHT算法中划分了两种零树，第一种零树有根系数构成，并且该根系数在给定位平面中它的全部子孙都为0；第二种零树与第一种零树相似，不同之处在于它排除了根系数的四个子带。SPIHT算法将这两种零树类型分别标记为类型A和B不重要集合，SPIHT算法还规定了三个有序表用来存放重要信息，分别是LSP（显著系数表）、LIP（不显著系数表）、LIS（不显著子集表）。三个表中的元素都采用（x，y）坐标来标记，C（x，y）表示（x，y）坐标处的小波系数值，H表示所有树根的集合。

具体算法的流程如下：

①初始化：设阈值T=2n，其中

n=log2（max|C（x，y）|）、LSP为空、LIP={（x，y）|（x，y）∈H}、LIS={（x，y）A|（x，y）∈H且都具有孩子}；

②依次对LIP表和LIS表中元素进行排序扫描并编码。对LIP编码只要判断其中的每个元素是否为显著系数，若是则输出1及其符号码、将该系数从LIP表中剪切到LSP表的尾部，否则输出0；对LIS编码的流程比对LIP编码复杂，在判断表中元素属于类型A还是类型B之后还要对其子孙系数进行编码，具体见图6的流程结构图。

③当LIS集合中的最后一个编码完成后，对上一级的得到的LSP表进行精细扫描：对表中每项对应的小波系数进行二进制转换、将该二进制数的第N个重要的位输送到精细位流中。

④完成后将T除以2，重新进入第二步并进行下一级的编码，直到T为某一设定值为止。

4 DICOM图像压缩实现与分析

文中选用了一个DICOM文件Ankle.dcm作为实验文件进行说明。如图7所示，通过VC++语言，根据DICOM文件的数据元素中的标识符、数据类型、数据长度、数据值，来对DICOM文件进行读取，并将DICOM图像文件中的图像数据和文本数据分离出来。图7的左边为DICOM文件Ankle.dcm中分离出来的部分文本信息，图7的右边为DICOM文件Ankle.dcm中的分离出来的格式为.bmp的图像文件。

将分离出来的BMP图像进行双正交9/7提升小波变换，并将小波分解后的系数取整量化、SPIHT编码。表1是在同为三级小波变换，压缩比不同情况下，压缩后图像的PSNR（峰值信噪比）的值。图8为压缩比率为30：1时候的Ankle压缩图像。

图像客观评价的标准有很多种，比较常见的有MSE（Mean square error）、NMSE（Normalized mean squared error）、SNR（Signal Noise Ratio）、PSNR（Peak Signal Noise Ratio）。它们之间有一定的转换关系，这里采用PSNR这个参数作为对图像客观评价的一个标准。如表2中所示，随着压缩比率的增加图像的PSNR值也逐渐降低，事实上，人眼很难观察出PSNR超过30dB的解码图像与原图像的差异[11]，由图8亦可观察出压缩比在30：1时仍能保证图像的清晰度。endprint

5 结 论

文中通过采用了双正交9/7提升小波变换得到较好时频特性的小波系数再结合更为先进的SPHIT算法对小波系数进行了编码进而完成了对从DICOM文件中分离出来的图像信息的压缩，并最终实现了DICOM图像的压缩工作。实验表明这种方法实现了在较高图像质量的条件下可以获得较大压缩比的图像文件，这也在一定程度上提高PACS系统中对图像数据储存和传输的效率。

参考文献

[1] 曹奕雯， 陶蔷. PACS关键技术及国内应用[J]. 医疗卫生装备， 2011，32（1）：77-80.

[2] 梁存升， 冯骥. DICOM标准分析及其应用[J]. 中国医学装备， 2006，3（2）：18-20.

[3] 赵越， 韩滢， 王之琼. DICOM技术在PACS系统中的应用综述[J]. 中国数字医学， 2009，4（2）：57-60.

[4] BIDGOOD W D J， HORII S C，PRIOR F W， et al. Understanding and using DICOM， the data interchange standard for biomedical imaging.[J]. Journal of the American Medical Informatics Association ： JAMIA， 1997，4（3）：199-212.

[5] 李自胜， 肖晓萍， 龚伟. 基于.NET框架的DICOM协议设计与实现[J]. 计算机应用与软件， 2009，26（10）：19-23.

[6] 王海松， 赵杰. 数字图像压缩技术的现状及前景分析[J]. 科技信息， 2010（3）：53-54.

[7] 胡高军， 吴芝路， 高真真， 等. 自适应小波变换及其在JPEG 2000中的应用[J]. 光电子·；激光， 2006，17（9）：1114-1118.

[8] HUI F，LANYING G，JINSHENG X. The Lifting Scheme Based on the Second Generation Wavelets[J]. Wuhan University Journal of Natural Sciences， 2006，11（3）：503-506.

[9] CHEN T，LIN S，LIN Y， et al. JPEG2000 Still Image Coding Quality[J]. Journal of Digital Imaging， 2013，26（5）：866-874.

[10]SAID A，PEARLMAN W A. A new， fast， and efficient image codec based on set partitioning in hierarchical trees[J]. IEEE Transactions on Circuits and Systems for Video Technology， 1996，6（3）：243-250.

[11]SUNG M， KIM H， YOO S K. Clinical Evaluation of Compression Ratios Using JPEG2000 on Computed Radio graphy Chest Images[J]. Journal of Digital Imaging， 2002，12（2）：78-83.endprint