CN112995637B

CN112995637B - 一种基于三维离散小波变换的多剖面医学图像压缩方法

Info

Publication number: CN112995637B
Application number: CN202110258764.0A
Authority: CN
Inventors: 朱江; 刘芷伊; 李建奇
Original assignee: Xiangtan University
Current assignee: Xiangtan University
Priority date: 2021-03-10
Filing date: 2021-03-10
Publication date: 2023-02-28
Anticipated expiration: 2041-03-10
Also published as: CN112995637A

Abstract

本发明提出了一种基于三维离散小波变换的多剖面医学图像压缩方法。首先，对多张医学剖面图像进行多层三维离散小波变换得到若干子带；然后，子带按照小波树结构规律以及从低到高的频率顺序进行空间上的拼接操作，得到小波系数矩阵；然后，对小波系数矩阵系数矩阵进行三维MORTON扫描，得到编码顺序；最后，按照编码顺序对小波系数矩阵进行三维嵌入式零树结合游程编码，获得压缩编码数据流。本发明的优点在于高效、快速地压缩同一个器官的多张剖面图像，在保持图像高信噪比的同时达到客观的压缩比，有效解决二维图像帧间冗余问题，提高了多剖面图像的数据传输效率。

Description

一种基于三维离散小波变换的多剖面医学图像压缩方法

技术领域

本发明涉及一种基于三维离散小波变换的多剖面医学图像压缩方法，尤其涉及三维离散小波变换和三维嵌入式零树编码。

背景技术

在医学领域，CT、MRI等医学图像被广泛用于辅助医疗诊断。近年来，随着互联网技术、计算机技术的蓬勃发展，使远程医疗诊断成为了可能。为了对异地患者第一时间实施诊断，医生需要及时获得患者的医学图像。然而，与自然场景图像相比，医学图像对图像质量有更高的要求，相应地，其数据量也更大。此外，人体器官通常需要若干张剖面图才能完整地反映器官的病理情况，导致待传输的数据量急剧增加。因此，如何在不影响图像质量的前提下，对患者的多剖面医学图像进行压缩，达到减少传输数据量，提高图像传输速度的目的，成为了医学图像压缩处理领域的难题。

在图像压缩领域，基于离散余弦变换(Discrete Cosine Transform,DCT)的JPEG(Joint Photographic Experts Group)标准依然是如今主流的压缩算法，其对应的JPEG图像格式具有文件小、利于网络传输等优点而得到普及应用。但是JPEG图像质量一般，不适用于高分辨率要求的医学图像压缩。随着基于离散小波变换(Discrete Wavelet Transform，DWT)的JPEG2000标准提出，基于DWT的图像压缩算法得到了广泛的关注，并且在多种图像上展现出了优越的压缩性能。然而，主流压缩算法仅能对单张图像压缩，难以满足同时压缩多张多剖面医学图像的要求。

发明内容

针对上述问题，考虑到基于DCT变换的压缩算法不利于恢复高质量医学图像，并且要在能得到高质量恢复图像的前提下还能达到高压缩比，发明了一种基于三维离散小波变换的多剖面医学图像压缩方法，旨在弥补DCT变换压缩算法带来的精度损失的同时去除三维医学图像帧间冗余来提高压缩性能。具体内容如下：

S100.读取2^l张多剖面医学图像图源，其中l∈N₊；将所有图像转换为灰度图像，大小设定为2^m×2ⁿ，其中m,n∈N₊，将变换后的图像储存到同一矩阵中，得到大小为2^m×2ⁿ×2^l的三维图像数据矩阵Image；对Image进行p层三维离散小波变换，其中p应满足p≤min(m,n,l)，最终得到(7×p)+1个一维子带：LLL_p、LLH_i、HLL_i、HLH_i、LHL_i、LHH_i、HHL_i、HHH_i，其中(i＝1,2,...,p)；

S200.对S100中获得的所有子带，按照由低频到高频，低频子带集中在左上角，高频子带集中在右下角的规律进行拼接，得到大小为2^m×2ⁿ×2^l的小波系数矩阵MAT；层间频率T遵循的规律为T_{第i层子带}＜T_{第i-1层子带}，子带频率t遵循的规律为t_LLL＜t_LLH＜t_LHL＜t_LHH＜t_HLL＜t_HLH＜t_HHL＜t_HHH；

S300.对小波系数矩阵MAT进行三维嵌入式零树结合游程编码，获得最终压缩编码数据流。

本发明具有如下优点：

1.本发明能够一次性处理多张人体器官的图像，高效、快速地压缩同一个器官的多张剖面图像，有效解决二维图像帧间冗余问题，实现了多剖面图像的高效数据传输；

2.本发明能够消除图像方块失真效应通过合理地量化图像数据矩阵，解决传统JPEG标准下量化表不够灵活问题。

附图说明

图1是本发明的流程图；

图2是本发明处理的医学图像图源(CT图)；

图3是三维离散小波变换子带变化示意图，具体如下：图3(a)、图3(b)、图3(c)分别为一层三维离散小波变换子带示意图、三层三维离散小波变换子带示意图、小波树结构示意图；

图4是三维MORON扫描示意图，具体如下：图4(a)为子带中扫描顺序，图4(b)为子带间扫描顺序；

图5是嵌入式零树编码主扫描流程图；

图6是三层三维离散小波变换后图像视觉变化示意图，具体如下：图6(a)、图6(b)、图6(c)分别是第一层、第二层、第三层三维离散小波变换图像变化；

图7是图像分块结果；

图8是第(39，16)个像素子块原图；

图9是像素子块通过本发明后并恢复的视觉效果。

具体实施方式

实施本发明，以8张人脑CT扫描图为图源，MATLAB为实现环境为例进行说明。8张脑部的CT图像如图2所示，CT图像是由一系列灰度像素按照一定尺寸排列而成，保存格式为bmp，从文件中直接读取会以RBG形式储存，开始压缩之前，应将其转换为灰度图像，转换后的每张CT图像数据矩阵为512*512，以uint8类型储存。

结合图1，本发明具体实施步骤如下：

S100.读取8张多剖面医学图像图源，将所有图像转换为灰度图像，大小设定为512*512，将变换后的图像储存到同一矩阵中，得到大小为512*512*8的三维图像数据矩阵Image，对Image进行3层三维离散小波变换，得到22个一维子带LLL₃、LLH_i、HLL_i、HLH_i、LHL_i、LHH_i、HHL_i、HHH_i，其中(i＝1,2,3)；

S110.对输入大小为512*512*8的图像数据第k个通道第j列的所有行取出一维数据S_i,j,k，其中i∈(1,512)，j∈(1,512)，k∈(1,8)；

S120.选取haar小波作为小波基，对取出的一维数据S_i,j,k进行一维离散小波变换，得到高频低频系数，计算公式如下：

式中，f(t)为输入信号，φ(t)和

为离散小波基，c为低频近似系数，d为高频近似系数，j₀为任意一起始点，j≥j₀；

S130.将c存入低频带，d存入高频带，返回低频子带L和高频子带H；

S140.分别对L、H的第k个通道第i行的所有列取出一维数据S_i,j,k，其中i∈(1,256)，j∈(1,512)，k∈(1,8)，依次执行步骤S120-S130得到次子带LL、LH、HL、HH；

S150.对LL、LH、HL、HH的第i行第j列的所有通道取出一维数据S_i,j,k，其中i∈(1,256)，j∈(1,256)，k∈(1,8)，依次执行步骤S120-S130到最终子带LLL₁、LLH₁、HLL₁、HLH₁、LHL₁、LHH₁、HHL₁、HHH₁，然后执行S160；

S160.取出最低频子带LLL₁作为输入，重复执行步骤S110-S150，得到子带LLL₂、LLH₂、HLL₂、HLH₂、LHL₂、LHH₂、HHL₂、HHH₂；

S170.S160.取出最低频子带LLL₂作为输入，重复执行步骤S110-S150，得到子带LLL₃、LLH₃、HLL₃、HLH₃、LHL₃、LHH₃、HHL₃、HHH₃，最终共得到22个子带，其中第i层子带大小为2^9-i×2^9-i×2^3-i；以三层三维离散小波变换为例，图像变化情况如图6所示。

S200.对S100中获得的所有子带，按照如图3(b)的规律进行拼接，得到大小为512*512*8的小波系数矩阵MAT；层间频率T遵循的规律为T_{第i层子带}＜T_{第i-1层子带}，层内频率t遵循的规律为t_LLL＜t_LLH＜t_LHL＜t_LHH＜t_HLL＜t_HLH＜t_HHL＜t_HHH。

S300.对小波系数矩阵MAT进行三维嵌入式零树结合游程编码，获得最终压缩编码数据流；

S310.对系数矩阵MAT进行三维MORTON扫描，子带内部执行如图4(a)所示的二维MORON扫描，将MAT转换为一维系数序列，子带间扫描顺序为LLL、HLL、LHL、HHL、LLH、HLH、LHH、HHH，具体扫描顺序如图4(b)；

S320.设定编码次数E＝7，解码次数D＝7；

S330.求出初始阈值

其中max为MAT中最大值，得到阈值列表T＝{T₀,T₁,...,T_E}，其中

S340.主扫描对S310输出的一维系数序列的每一个元素类型进行判断，输出“P”(正重要系数)、“N”(负重要系数)、“T”(零树)、“Z”(孤立零)四种符号类型，并将判断结果存入列表CodeList；其中，P表示元素绝对值大于当前阈值，且大于0的系数；N表示元素绝对值大于当前阈值，但小于0的系数；T表示元素小于当前阈值，且其子代均为非重要系数的系数；Z表示元素小于当前阈值，但其子代里有重要系数的系数；每完成一次主扫描编码都需要更新扫描列表、重要系数列表、已扫描列表；

S350.构造数目为N＝[2T(0)-T(d)]/T(d)的量化器，其中，第n个量化器中储存值的为A＝T(d)*(n+0.25)和B＝T(d)*(n+0.75)，分别用“0”和“1”表示，式中，T(0)为初始阈值，T(d)为当前阈值；

S360.进行辅扫描：首先确定重要系数所选择的量化器，量化器的编号记为rl，定义如下

rl＝floor(imprtvalue(j)/T(d)) (4)

其中，imprtvalue(j)表示第j个重要系数，

用重要系数减去量化器起始值rl*T(d)，再与当前阈值T(d)的二分之一相比较，大于

则选择量化值B，输出量化符号“1”，否则选择A，输出量化符号“0”，将判断的量化结果存入列表QuantiFlagList；

S370.主扫描输出系数类型列表CodeList和辅扫描输出量化符号列表QuantiFlagList，CodeList里为S330的判断结果，数据类型为“P”、“N”、“T”、“Z”，QuantiFlagList里为S360的输出结果，数据类型为“0”、“1”；

S380.把CodeList、QuantiFlagList转换为一维数组，分别进行游程编码，输出最终编码数据流；

S390.以如图8所示的第(39，16)个像素子块遍历压缩算法，并通过逆过程解码，其恢复效果对比如图9所示。

Claims

1.一种基于三维离散小波变换的多剖面医学图像压缩方法，至少包含以下几个步骤：

S100.读取2^l张多剖面医学图像图源，其中l∈N₊；将所有图像转换为灰度图像，大小设定为2^m×2ⁿ，其中m，n∈N₊；将变换后的图像储存到同一矩阵中，得到大小为2^m×2ⁿ×2^l的三维图像数据矩阵Image，对Image进行p层三维离散小波变换，其中p≤min(m，n，l)，最终得到(7×p)+1个一维子带：LLL_p、LLH_x、HLL_x、HLH_x、LHL_x、LHH_x、HHL_x、HHH_x，其中(x＝1，2，...，p)，三维离散小波变换的具体实现为三个一维小波变换组合；

其中，S100包括S110-S130：

S110.对输入大小为2^m×2ⁿ×2^l的图像数据第k个通道第j列的所有行取出一维数据S_i，j，k，进行一维离散小波变换，其中i∈(1，2^m)，j∈(1，2ⁿ)，k∈(1，2^l)，得到子带L、H；

S120.分别对L、H的第k个通道第i行的所有列取出一维数据S_i，j，k，其中i∈(1，2^m-1)，j∈(1，2ⁿ)，k∈(1，2^l)，执行一维离散小波变换，得到次子带LL、LH、HL、HH；

S130.对LL、LH、HL、HH的第i行第j列的所有通道取出一维数据S_i，j，k，其中i∈(1，2^m-1)，j∈(1，2^n-1)，k∈(1，2^l)，执行一维离散小波变换，得到最终子带LLL₁、LLH₁、HLL₁、HLH₁、LHL₁、LHH₁、HHL₁、HHH₁，其中LLL₁，进入下一层三维离散小波变换。

S200.遵循小波树的结构规律对S100中获得的所有子带进行空间上的拼接操作，其中低频子带集中在左上角，高频子带集中在右下角，得到大小为2^m×2ⁿ×2^l的小波系数矩阵MAT；层间频率T遵循的规律为T_{第x层子带}＜T_{第x-1层子带}，子带频率t遵循的规律为t_LLL＜t_LLH＜t_LHL＜t_LHH＜t_HHLL＜t_HLH＜t_HHL＜t_HHH；

其中，S300包括S310-S380：

S310.对系数矩阵MAT进行三维MORTON扫描，子带内部执行二维MORON扫描，将MAT转换为一维系数序列，子带间扫描顺序为LLL、HLL、LHL、HHL、LLH、HLH、LHH、HHH；

S320.设定编码次数E，解码次数D；

S330.求出初始阈值

其中max为MAT中最大值，得到阈值列表T＝{T₀，T₁，...，T_E}，其中

rl＝floor(imprtvalue(y)/T(d)) (4)

其中，imprtvalue(y)表示第y个重要系数，

则选择量化值B，输出量化符号″1″，否则选择A，输出量化符号″0”，将判断的量化结果存入列表QuantiFlagList；

S370.主扫描输出系数类型列表CodeList和辅扫描输出量化符号列表QuantiFlagList；

S380.把CodeList、QuantiFlagList转换为一维数组，分别进行游程编码，输出最终编码数据流。