CN109360162A - 基于bitonic滤波与分解框架的医学超声图像去噪方法 - Google Patents

Abstract

基于bitonic滤波和分解框架的医学超声图像去噪方法，包括如下步骤：步骤1)将图像通过分解框架分解；步骤2)构造一个rank滤波器；步骤3)对C1,C3分量进行形态学开、闭操作；计算开、闭误差并对其进行平滑；步骤5)对经过开闭运算后的C1，C3分量加权求和得到bitonic滤波后的分量Crec1，Crec3；步骤6)对分量Crec1，Crec3进行逆变换得到最后的去噪图像。本发明通过实验分析与去噪领域中的几种算法进行了对比，有效的应用在医学超声去噪领域；通过框架分解超声图像得到三个分量，再将bitonic滤波应用在其中两个分量上。通过了大量的实验数据对比，提出了基于bitonic滤波和分解框架的医学超声图像去噪方法，能够更好的有帮助医师的分析诊断。

Description

基于bitonic滤波与分解框架的医学超声图像去噪方法

技术领域

本发明涉及一种医学超声图像去噪方法。

背景技术

随着科技的发展，在医学成像领域，超声成像、CT、MRI等成像技术已应用于医学临床诊断中。超声医学，是声学、医学、光学及电子学相结合的学科。凡研究高于可听声频率的声学技术在医学领域中的声医学。涉及的内容广泛，在预防、诊断、治疗疾病中有很高的价值。

超声成像是利用超声束扫描人体，通过对反射信号的接收、处理，以获得体内器官的图象。常用的超声仪器有多种：A型(幅度调制型)是以波幅的高低表示反射信号的强弱，显示的是一种“回声图”。M型(光点扫描型)是以垂直方向代表从浅至深的空间位置，水平方向代表时间，显示为光点在不同时间的运动曲线图。以上两型均为一维显示，应用范围有限。B型(辉度调制型)即超声切面成像仪，简称“B超”。是以亮度不同的光点表示接收信号的强弱，在探头沿水平位置移动时，显示屏上的光点也沿水平方向同步移动，将光点轨迹连成超声声束所扫描的切面图，为二维成像。至于D型是根据超声多普勒原理制成.C型则用近似电视的扫描方式，显示出垂直于声束的横切面声像图。近年来，超声成像技术不断发展，如灰阶显示和彩色显示、实时成像、超声全息摄影、穿透式超声成像、超声计并机断层圾影、三维成像、体腔内超声成像等。

在过去的几十年中，人们大量地研究了在保持图像的主要特征(边缘、纹理、颜色、对比度等)的同时去除图像噪声的问题，本发明使用医学超声图像为研究对象，由于斑点噪声的存在严重影响了超声图像的质量，导致了超声医学图像质量较差。斑点噪声在图像上表现为空间域内相关的形状各异的小斑点，它将掩盖那些灰度差别很小的图像特征。传统滤波方法往往会破坏超声图像原有的图像特征结构。对于临床医生而言，斑点噪声对他们的准确诊断造成了很大的干扰，即使是滤波后的图片，也会因为滤波破坏的图像结构信息影响他们的判断，特别是对于经验不是很丰富的医生造成的影响更大。因此，从临床应用的角度出发，需要研究对超声医学图像保持图像结构信息的去噪方法，为医生做出更准确的诊断提供技术支持，降低人工诊断的风险。现有的滤波方法在实验阶段一般是针对噪声比较严重的情况能得到较高的信噪比为目标，这样的滤波方法固然能得到较好的滤波效果，但是往往要耗费大量的时间，而现有的快速的滤波方法如同中值滤波、高斯滤波得到的滤波效果很一般。真实的医学图像的噪声实际上很少会有像实验假设时那样严重的噪声，所以对于医学超声图像而言，研究医学超声图像的快速的并有一定保护图像结构信息的去噪方法具有非常重要的意义。

发明内容

为了克服传统滤波方法在保留图像结构信息上的不足，以及为了克服现有滤波方法需要消耗大量时间的不足，本发明提供了一种基于图像分解框架与双调滤波(bitonicfilter)医学超声图像的去噪算法，用于解决医学超声图像的去噪。

现有技术中，许多经典的滤波方法已经在图像滤波方面发挥了重大的作用，但是这些方法往往会破坏图像原有的例如图像边缘等结构信息，并且极少有能兼顾时间代价和保证滤波效果的滤波方法。本发明采用了最新的图像框架分解技术与最新的基于图像形态学滤波的双调滤波(bitonic filter)对医学超声图像进行滤波，此方法具有速度快、去噪明显的，保留图像结构信息更强的医学超声图像去噪声方法，最后通过仿真验证了方法的可行性与优化的效果。

本发明的优点是：提出了一种基于bitonic滤波与图像分解框架的医学图像去噪算法，图像分解框架克服了传统去噪方法保留图像结构能力的不足，bitonic滤波解决了传统快速滤波算法的不足，并具备一定的保留图像结构以及图像边缘信息的能力，将其应用于医学超声图像去噪能更好的保护图像边缘信息，更快地完成图像滤波，给医师的诊断提供了方便。

为使本发明的目的、技术方案和优点更加清晰，下面就对本发明的技术方案作进一步描述，基于bitonic滤波与分解框架的医学超声图像去噪方法，步骤如下：

步骤1)将图像通过分解框架分解。

对超声图像I进行分解，分解后能得到三个分量C1,C2,C3，分解的公式如下：

分量C1与C3是将图像不同的信息分开来，分别包含了图像的细节和近似信息，C2始终为0分量，其中算子K的定义公式如下：

上式中λ为平滑参数，I_x为图像I在(x,y)点关于x的偏导，I_y为图像I在(x,y)点关于y的偏导，为图像I在(x,y)点的梯度，经实验得出λ最合适的值为0.001，在I_x与I_y都为零的点，K设为单位矩阵。上式中，图像I位于(x,y)点的梯度值计算公式如下：

步骤2)构造一个rank滤波器。

我们将先构造一个将要用在下面的步骤中的rank滤波器，该滤波器的作用是将图像的局部区域所有的像素值排序后选定需要的位置的像素代替图像局部区域中心的像素值，其计算公式如下：

rank_w,c(x)＝c^thcentile{x_i},i∈w (4)

上式中参数w为局部窗口大小，参数c为排序后选取的序号，局部窗口大小w的选择影响着图片被滤波器平滑的程度，而参数c大小的选择是决定了对图像细节的保留程度，一般最好的效果是c选在10的时候，较低的c可以以一些轻微的非线性失真为代价来保持孤立的精细细节，20或更高的值对修复含有椒盐噪声的图像十分有效，不过也会更多的丢失图像原本的细节部分。本发明中对C1,C3两个分量w的值选为3，对于C1分量，c的值选为10，而对于C3分量，c的值选为20。

步骤3)对C1,C3分量进行形态学开、闭操作。

用步骤2)中构造的rank滤波对两个分量进行形态学开操作和闭操作来得到一个开变量Open和闭变量Close，计算开变量公式如下：

Open_w,c(x)＝rank_w,100-c(rank_w,c(x)) (5)

计算闭变量的公式如下：

Close_w,c(x)＝rank_w,c(rank_w,100-c(x)) (6)

式(5)，(6)中变量x即输入的图像，参数w，c即为步骤2)中介绍的rank滤波器的参数。将分量图像C1,C3作为上式中的变量x代入上式可得开变量OC1、OC3,和闭变量CC1,CC3。

步骤4)计算开、闭误差并对其进行平滑。

计算出经过开、闭操作后的图片与原图相减后取绝对值作为开、闭操作误差值，分别为开误差EOpen和闭误差EClose，并将这两个变量与适当长度的高斯核进行卷积，得到平滑后的开、闭误差，高斯核的长度与前面选定的rank滤波的窗口大小相关，开误差计算式如下：

EOpen(x)＝|G(x-Open_w,c(x))| (7)

闭误差的计算式如下：

EClose(x)＝|G(Close_w,c(x)-x)| (8)

上式中x为输入图像，函数G(.)为高斯平滑函数。将步骤3)中得到的开变量OC1、OC3,和闭变量CC1,CC3带入式(6)、(7)，得到分别对应于C1,C3变量的平滑后的开误差EOC1、EOC3和闭误差ECC1、ECC3。

步骤5)对经过开闭运算后的C1，C3通过分量开、闭误差加权求和公式得到bitonic滤波后的分量Crec1，Crec3。

根据bitonic滤波的最后一步的计算公式：

将以上步骤中计算得到的开闭变量OC1,OC3,CC1,CC3，和开闭误差EOC1,EOC3,ECC1,ECC3代入上式中，可得到分量C1,C3经过bitonic滤波后的分量Crec1,Crec3。

步骤6)对分量Crec1，Crec3进行逆变换得到最后的去噪图像。

经过滤波之后的各分量可经过逆变换合并成最终的去噪图像，逆变换的公式如下：

上式中，K是步骤1)中由图像梯度和偏导构成的算子，将其与滤波后的分量中的像素值组成的矩阵相乘，得到包含最终去噪图像像素的矩阵，公式是将每个图像中的某点像素代入，计算出合成后的该点的像素值，最终得到的去噪图像的像素值即为上式中的I(x,y)。

本发明具有以下优点：

1.本发明使用分解框架和bitonic滤波，将快速无需任何先验条件的bitonic滤波器应用在经分解框架分解之后的图像分量上，能得到更好的医学去噪效果的同时能减少计算的时间消耗。

2.本发明中bitonic滤波过程十分快速，能更好的适应于实时性要求较强的情况下。

3.本发明结构简洁，并且使用了更加新颖的理论。

附图说明

图1a为医学超声原图，图1b和图1c分别为经框架分解后的医学超声图的C1,C3分量图；

图2a,2b分别为经过bitonic滤波后的两个分量Crec1，Crec3的图；

图3为本发明整体步骤流程图；

图4为案例分析整体流程；

图5a～5e为各种算法在经典图Lena(σ_n＝20)上实验结果的比较，其中图5a是原图，图5b是噪声图，图5c是中值滤波算法效果图，图5d是高斯滤波算法效果图，图5e是本发明算法效果图；

图6a～图6d为各种算法在医学超声图像上实验结果的比较，其中图6a是真实医学超声原图，图6b是中值滤波算法效果图，图6c是高斯滤波算法效果图，图6d是本发明算法效果图；

具体实施方式

以下结合附图对本发明做进一步说明。

基于bitonic滤波与分解框架的医学超声图像去噪方法，步骤如下：

步骤1)将图像通过分解框架分解。

对超声图像I进行分解，分解后能得到三个分量C1,C2,C3，分解的公式如下：

分量C1与C3分别包含了图像的细节和近似信息，C2为0分量，其中算子K的定义公式如下：

上式中λ为平滑参数，I_x为图像I在(x,y)点关于x的偏导，I_y为图像I在(x,y)点关于y的偏导，为图像I在(x,y)点的梯度，经实验得出λ最合适的值为0.001，在I_x与I_y都为零的点，K设为单位矩阵。

图像1b，1c显示了医学超声图像经框架分解后的两个分量图。

步骤2)构造一个rank滤波器。

我们将先构造一个将要用在下面的步骤中的rank滤波器，该滤波器的作用是将图像的局部区域所有的像素值排序后选定需要的位置的像素代替图像局部区域中心的像素值，其计算公式如下：

rank_w,c(x)＝c^thcentile{x_i},i∈w (4)

上式中参数w为局部窗口大小，参数c为排序后选取的序号，本发明中w的值选为3，c的值是10。

步骤3)对C1,C3分量进行形态学开、闭操作。

用步骤2)中构造的rank滤波对两个分量进行形态学开操作和闭操作来得到一个开变量Open和闭变量Close，计算开变量公式如下：

Open_w,c(x)＝rank_w,100-c(rank_w,c(x)) (5)

计算闭变量的公式如下：

Close_w,c(x)＝rank_w,c(rank_w,100-c(x)) (6)

式(4)，(5)中变量x即输入的图像，参数w，c即为步骤2)中介绍的rank滤波器的参数。

步骤4)计算开、闭误差并对其进行平滑。

计算出经过开、闭操作后的图片与原图的误差值，分别为开误差EOpen和闭误差EClose，并将其与适合长度的高斯核进行卷积，得到平滑后的开、闭误差。高斯核的长度与前面选定的rank滤波的窗口大小相关，开误差计算式如下：

EOpen(x)＝|G(x-Open_w,c(x))| (7)

闭误差的计算式如下：

EClose(x)＝|G(Close_w,c(x)-x)| (8)

上式中x为输入图像，函数G(.)为高斯平滑函数。

步骤5)对经过开闭运算后的C1，C3分量加权求和得到bitonic滤波后的分量Crec1，Crec3。

根据bitonic滤波的最后一步的计算公式：

可得到分量C1,C3经过bitonic滤波后的分量Crec1,Crec3。图2a,2b为经过bitonic滤波后的两个分量Crec1，Crec3的图。

步骤6)对分量Crec1，Crec3进行逆变换得到最后的去噪图像。

经过滤波之后的各分量可经过逆变换合并成最终的去噪图像，逆变换的公式如下：

最终得到的去噪图像即为上式中的I。

本发明整体步骤流程图如图3所示。

案例分析

本发明通过以具体的医学超声图像以及去噪领域的经典图像lena为对象，通过bitonic滤波与分解框架结合得到快速且效果不错的去噪结果，同时通过与现有技术对比来表现了本发明的优越性。案例分析的流程图如图4所示。

本发明使用峰值信噪比(PSNR)来变现图像重构后的质量，PSNR定义如下：

其中N表示图像中的像素数目，表示弗罗比尼斯范数，255是像素可以在灰度图像中获得的最大值。PSNR数值越大，去噪效果越好。

为了判别算法对图像结构的保护程度，还用结构相似度(SSIM)作为性能评判指标，SSIM的定义如下：

其中，SSIM值越大表明算法对图像结构的保护能力越好。

本发明的实验的硬件参数为CPU：酷睿i5-4210U双核主频1.70GHz 2.40GHz,运行内存：3.67GB。软件使用的是在微软windows7 64位操作系统下运行的MATLAB2014a。本实验采用真实医学超声图像和经典Lena图作为输入数据，可以进行有效的对比实验，案例分析整体流程图如图4。实验通过对比中值滤波，高斯滤波和本文方法三种快速去噪算法的性能，以及与FFST和BM3D这两种去噪效果较好的算法比较消耗的时间来证明本发明的优越性。各种算法应用在图Lena的实验效果图如图5，各种算法应用在医学超声图的实验效果图如图6所示。

表1，2中可看出，在经典图像Lena去噪时，本发明的去噪效果以及对图像结构的保护都优于另外两种传统的快速图像算法。在表3中可以看出，本发明的去噪时间不会随噪声方差的变化而产生较大变化，并且始终保持较快的处理速度。

表1：Lena图在不同去噪算法在不同噪声的PSNR/dB值

算法	σ<sub>n</sub>＝10	σ<sub>n</sub>＝20	σ<sub>n</sub>＝30	σ<sub>n</sub>＝40
					本发明算法	33.4349	30.2967	27.7748	25.7819
高斯滤波	30.0090	29.9064	27.2113	25.0507
					中值滤波	32.4387	28.5124	25.6260	23.3708

表2：Lena图在不同去噪算法在不同噪声的SSIM值

算法	σ<sub>n</sub>＝10	σ<sub>n</sub>＝20	σ<sub>n</sub>＝30	σ<sub>n</sub>＝40
					本发明算法	0.8921	0.7713	0.6486	0.5405
高斯滤波	0.8826	0.7408	0.6037	0.4898
					中值滤波	0.8477	0.6662	0.5185	0.4069

表3：医学超声图在不同去噪算法的运行时间s

算法	本方法	FFST	BM3D
				时间	1.2206	4.3286	3.0051

本说明书实施例所述的内容仅仅是对发明构思的实现形式的列举，本发明的保护范围不应当被视为仅限于实施例所陈述的具体形式，本发明的保护范围也及于本领域技术人员根据本发明构思所能够想到的等同技术手段。

Claims (1)

1.基于bitonic滤波与分解框架的医学超声图像去噪方法，步骤如下：

步骤1)将图像通过分解框架分解；

对超声图像I进行分解，分解后能得到三个分量C1,C2,C3，分解的公式如下：

分量C1与C3是将图像不同的信息分开来，分别包含了图像的细节和近似信息，C2始终为0分量，其中算子K的定义公式如下：

上式中λ为平滑参数，I_x为图像I在(x,y)点关于x的偏导，I_y为图像I在(x,y)点关于y的偏导，▽I为图像I在(x,y)点的梯度，经实验得出λ最合适的值为0.001，在I_x与I_y都为零的点，K设为单位矩阵；上式中，图像I位于(x,y)点的梯度值计算公式如下：

步骤2)构造一个rank滤波器；

将先构造一个将要用在下面的步骤中的rank滤波器，该滤波器的作用是将图像的局部区域所有的像素值排序后选定需要的位置的像素代替图像局部区域中心的像素值，其计算公式如下：

rank_w,c(x)＝c^thcentile{x_i},i∈w (4)

上式中参数w为局部窗口大小，参数c为排序后选取的序号，局部窗口大小w的选择影响着图片被滤波器平滑的程度，而参数c大小的选择是决定了对图像细节的保留程度，一般最好的效果是c选在10的时候，较低的c可以以一些轻微的非线性失真为代价来保持孤立的精细细节，20或更高的值对修复含有椒盐噪声的图像十分有效，不过也会更多的丢失图像原本的细节部分；对C1,C3两个分量w的值选为3，对于C1分量，c的值选为10，而对于C3分量，c的值选为20；

步骤3)对C1,C3分量进行形态学开、闭操作；

用步骤2)中构造的rank滤波对两个分量进行形态学开操作和闭操作来得到一个开变量Open和闭变量Close，计算开变量公式如下：

Open_w,c(x)＝rank_w,100-c(rank_w,c(x)) (5)

计算闭变量的公式如下：

Close_w,c(x)＝rank_w,c(rank_w,100-c(x)) (6)

式(5)，(6)中变量x即输入的图像，参数w，c即为步骤2)中介绍的rank滤波器的参数；将分量图像C1,C3作为上式中的变量x代入上式可得开变量OC1、OC3,和闭变量CC1,CC3；

步骤4)计算开、闭误差并对其进行平滑；

计算出经过开、闭操作后的图片与原图相减后取绝对值作为开、闭操作误差值，分别为开误差EOpen和闭误差EClose，并将这两个变量与适当长度的高斯核进行卷积，得到平滑后的开、闭误差，高斯核的长度与前面选定的rank滤波的窗口大小相关，开误差计算式如下：

EOpen(x)＝|G(x-Open_w,c(x))| (7)

闭误差的计算式如下：

EClose(x)＝|G(Close_w,c(x)-x)| (8)

上式中x为输入图像，函数G(.)为高斯平滑函数；将步骤3)中得到的开变量OC1、OC3,和闭变量CC1,CC3带入式(6)、(7)，得到分别对应于C1,C3变量的平滑后的开误差EOC1、EOC3和闭误差ECC1、ECC3；

步骤5)对经过开闭运算后的C1，C3通过分量开、闭误差加权求和公式得到bitonic滤波后的分量Crec1，Crec3；

根据bitonic滤波的最后一步的计算公式：

将以上步骤中计算得到的开闭变量OC1,OC3,CC1,CC3，和开闭误差EOC1,EOC3,ECC1,ECC3代入上式中，可得到分量C1,C3经过bitonic滤波后的分量Crec1,Crec3；

步骤6)对分量Crec1，Crec3进行逆变换得到最后的去噪图像；

经过滤波之后的各分量可经过逆变换合并成最终的去噪图像，逆变换的公式如下：

上式中，K是步骤1)中由图像梯度和偏导构成的算子，将其与滤波后的分量中的像素值组成的矩阵相乘，得到包含最终去噪图像像素的矩阵，公式是将每个图像中的某点像素代入，计算出合成后的该点的像素值，最终得到的去噪图像的像素值即为上式中的I(x,y)。