CN112734663A - 基于自适应阈值的轮廓波变换医学ct图像降噪方法 - Google Patents
基于自适应阈值的轮廓波变换医学ct图像降噪方法 Download PDFInfo
- Publication number
- CN112734663A CN112734663A CN202011620169.9A CN202011620169A CN112734663A CN 112734663 A CN112734663 A CN 112734663A CN 202011620169 A CN202011620169 A CN 202011620169A CN 112734663 A CN112734663 A CN 112734663A
- Authority
- CN
- China
- Prior art keywords
- image
- noise
- original
- filter
- decomposition
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Withdrawn
Links
- 238000000034 method Methods 0.000 title claims abstract description 40
- 230000009467 reduction Effects 0.000 title claims abstract description 26
- 230000009466 transformation Effects 0.000 title claims abstract description 26
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 35
- 238000012545 processing Methods 0.000 claims abstract description 18
- 230000003044 adaptive effect Effects 0.000 claims abstract description 11
- 230000000694 effects Effects 0.000 claims abstract description 9
- 238000002591 computed tomography Methods 0.000 claims description 64
- 238000000354 decomposition reaction Methods 0.000 claims description 52
- 230000008602 contraction Effects 0.000 claims description 11
- 239000000654 additive Substances 0.000 claims description 9
- 230000000996 additive effect Effects 0.000 claims description 9
- 238000004458 analytical method Methods 0.000 claims description 7
- 238000004364 calculation method Methods 0.000 claims description 3
- 238000009795 derivation Methods 0.000 claims description 3
- 230000014759 maintenance of location Effects 0.000 claims description 3
- 238000005070 sampling Methods 0.000 claims description 3
- 238000010008 shearing Methods 0.000 claims description 3
- 239000000126 substance Substances 0.000 claims description 3
- 230000001629 suppression Effects 0.000 claims description 3
- 238000004088 simulation Methods 0.000 abstract 1
- 230000008569 process Effects 0.000 description 8
- 230000000875 corresponding effect Effects 0.000 description 6
- 238000010586 diagram Methods 0.000 description 6
- 238000003745 diagnosis Methods 0.000 description 4
- 238000002474 experimental method Methods 0.000 description 4
- 206010058467 Lung neoplasm malignant Diseases 0.000 description 3
- 238000013170 computed tomography imaging Methods 0.000 description 3
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 3
- 201000005202 lung cancer Diseases 0.000 description 3
- 208000020816 lung neoplasm Diseases 0.000 description 3
- 238000006243 chemical reaction Methods 0.000 description 2
- 230000002596 correlated effect Effects 0.000 description 2
- 201000010099 disease Diseases 0.000 description 2
- 210000004072 lung Anatomy 0.000 description 2
- 239000002096 quantum dot Substances 0.000 description 2
- 206010028980 Neoplasm Diseases 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 201000011510 cancer Diseases 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000000691 measurement method Methods 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 230000000644 propagated effect Effects 0.000 description 1
- 238000013441 quality evaluation Methods 0.000 description 1
- 238000009877 rendering Methods 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 210000003625 skull Anatomy 0.000 description 1
- 210000004872 soft tissue Anatomy 0.000 description 1
- 238000003325 tomography Methods 0.000 description 1
- 238000011426 transformation method Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
- G06T7/0012—Biomedical image inspection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10081—Computed x-ray tomography [CT]
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Medical Informatics (AREA)
- Radiology & Medical Imaging (AREA)
- Quality & Reliability (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Health & Medical Sciences (AREA)
- General Health & Medical Sciences (AREA)
- Image Processing (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
本发明公开了一种基于自适应阈值的轮廓波变换医学CT图像降噪方法，包括以下步骤：步骤1)获取医学CT图像模型；步骤2)对CT图像使用轮廓波滤波器组分解，得到多分辨率多方向的图像子带高频分量和低频分量；步骤3)根据自适应阈值算法计算轮廓波系数高频分量的阈值；步骤4)使用阈值收缩算法处理高频分量；步骤5)使用三边滤波器处理分解后的低频分量；步骤6)将处理后的低频分量和高频分量进行轮廓波逆变换。本发明使用仿真实验和过去的降噪方法进行比较，在医学CT图像上表现了更好的效果，能更有效地帮助医生通过阅片进行诊断。
Description
技术领域
本发明涉及医学图像去噪领域，尤其涉及医学CT图像，更具体的说，尤其涉及一种适用于医学CT图像的基于自适应阈值的轮廓波变换医学CT图像降噪方法。
背景技术
肺癌是世界上发病几率最高的癌症之一，按发病人数排序，肺癌位居我国恶性肿瘤发病首位，位于男性发病首位和女性发病第二位。在所有的肺癌检测方法中，计算机断层扫描(ComputedTomography，CT)被称为最有效的方法。断层扫描本来是通过环切获得物体断面的技术，在医学中不仅仅被当作断面获得的检查方法，通过多平面图像重建法、多曲面图像重建法、最大值投影法、面渲染等三维图像处理技术来要进行影像诊断。
CT对人体软组织的分辨灵敏度较高，容易受到噪声的影响，降低了低对比度结构的可视性，影响医生的诊断。因此，了解噪声的来源及其在医学CT中的统计分布特征是很有必要的。
CT图像中的噪声会受到多种因素的影响，但主要有以下几个组成部分：
(1)随机噪声：该噪声是通过检测投影中的少量X射线而产生的，它看起来像图像密度的失真，然而图像密度是不知道的，因此产生的噪声称为随机噪声；
(2)量子噪声：X射线的能量，又称统计噪声，X射线是以单个的量子能块的形式传播，因此，这些无限量子的X射线被探测器看到。由于统计学误差，以及测量方法的不同，检测到的X射线量子塔可能会有不同的测量结果。医学CT中的量子噪声是由于检测到的少量X射线量子点可能出现的波动造成的。减少量子噪声影响的唯一方法是增加检测到的X射线量子点的数量，但这样也会增加CT的成本；
(3)电子噪声：当数字电路与模拟电路发生转换时，光电转换也会引入噪声，该噪声即为电子噪声；
(4)舍入噪声：通过模态变换获得X射线衰减减薄，然后送入指定的计算机进行CT图像重建，在计算机处理过程中，由于计算机的字节存储空间有限，通常会对数字信号进行四舍五入，计算出的CT图像会产生一定的四舍五入噪声，该噪声即为舍入噪声。
通过对CT设备的分析，可以看出医学CT中的主要噪声是量子噪声和电子噪声，各种因素相互约束且相互影响。上述多种噪声的相互作用会对医生和放射科医师的诊断造成一定的影响，考虑到为帮助医生更好，更清晰的阅片，从技术角度出发研究，在一定程度下降低各种噪声对医用CT成像的干扰是很有必要的。
发明内容
本发明的目的在于解决现有的医学CT成像中各种噪声的影响对医用CT成像的干扰严重，对医生和放射科医师的诊断造成一定的影响的问题，提出了一种基于自适应阈值的轮廓波变换医学CT图像降噪方法，用于多尺度有向分解图像，并在阈值收缩层面进行处理从而实现图像降噪的方法。
本发明通过以下技术方案来实现上述目的：一种基于自适应阈值的轮廓波变换医学CT图像降噪方法，具体包括以下步骤：
步骤一：获取CT图像模型
首先利用计算机断层扫描获取原始CT图像，原始CT图像接收到的信号氛围灰度信号和噪音信号，其中假设噪声信号为n(x，y)，原始信号为f(x，y)，输出信号为g(x，y)，根据原始CT图像与噪声之间的相互关系，将噪声分成加性噪声和乘性噪声两种形式，其中：
加性噪声与原始CT图像不相关，加性噪声的模型表示为
g(x，y)＝f(x，y)+n(x，y) (1)乘性噪声与原始CT图像相关，乘性噪声的模型表示为
g(x，y)＝f(x，y)*n(x，y) (2)
其中，(x，y)表示原始CT图像中的一串像素点，x，y分别表示原始CT图像中像素点的横坐标和纵坐标；
步骤二：对原始CT图像使用轮廓波滤波器组分解
对步骤一中含有噪声的原始CT图像使用轮廓波滤波器组进行多尺度多方向的分解，轮廓波滤波器组指的是拉普拉斯金字塔滤波器和二维方向滤波器组合而成的双重迭代的滤波器组结构；
由于步骤一中的原始CT图像为医学CT图像通用模型，对公式(1)进行轮廓波变换进行变换，可以得到新的加性噪声模型：
其中，
表示噪声图像的轮廓波系数，
表示无噪图像的轮廓波系数，
表示噪声的轮廓波系数，上标j表示轮廓波变换的分解层数，下标(l，k)，(l，k)∈Z2表示轮廓波域的坐标，J指最大分解层数；
步骤三：计算自适应于阈值算法的阈值
对步骤二种获得的轮廓波系数基于轮廓波变换的自适应降噪算法进行处理，
为了在平衡噪声抑制和特征保留过程中达到所期望的效果，综合经典阈值函数，并添加一个第jth层的自适应参数tj，提出一种新的阈值算法，该阈值算法的表达式如下：
其中，j(＝1，2，...，J)表示轮廓波变换的分解层数，J为最大分解层数，
其中，tj＝2J-j+1，k1＝k2＝0.5(5)
在上公式(4)只有噪声标准差和无噪图像的标准差是未知的；噪声标准差是由轮廓波内高频的轮廓波系数，即式中轮廓波系数
的各个高频子带
绝对值的中值计算得到
因为轮廓波线性变换，因此由式可以得到j层轮廓波方差关系式
由于
和
被建模为零均值的模型，因此轮廓波j层的标准差可有轮廓波系数的标准差计算得到
这里的M＝m×n即是对应轮廓波域内轮廓波系数的总个数，m和n分别是对应轮廓波域内系数的行数和列数；
由公式(8)可得
因此公式(8)和公式(9)便可计算得到噪声标准差σn和无噪声图像的标准差σs，j，此时式中的阈值T也能计算得到；
步骤四：使用阈值收缩算法处理高频分量
在步骤三中得到阈值T后，对比阈值高的轮廓波系数按照收缩算法进行收缩，得到降噪后的轮廓波系数；轮廓波收缩函数利用了贝叶斯估计理论，使用的是对应于均匀型贝叶斯风险函数的最大后验概率，改进的小波收缩算法的推导和步骤如下：
步骤1.1：运用贝叶斯最大后验估计获得小波域中信号的估计；设置
其中，PR(r)服从广义拉普拉斯分布的特殊情况；因此
步骤1.2：根据第二章对噪声模型的分析可以假定噪声的轮廓波系数
服从均值为0的高斯分布：
带入等式，可得
(3)求出ln(P(R|S)(r|s))关于r的一阶导数；为了获取最大后验概率，将随后的代数表达式设为0，由此可得，
其中
是r的估计；
(4)最后，合成的轮廓波收缩公式如下：
步骤五：使用三边滤波器处理分解后的低频分量
对步骤二种产生的低频分量使用三边滤波器处理，具体表达式如下：
其中w′S(x，ξ)是新的权重函数，有空域滤波器wS(x，ξ)、值域滤波器wR(x，ξ)和脉冲权重wI(ξ)组成，可表示为:
w′(x，ξ)＝wS(x，ξ)wR(x，ξ)wI(ξ) (17)
脉冲权重wI(ξ)的表达式为
其中引入fm(x)来估计像素点x是边缘点还是噪点；定义d(x，ξ)，表示点x和ξ之间的像素亮度的绝对差：
d(x，ξ)＝|f(x)-f(ξ)| (19)
由此可得，
其中，gi(x)表示出了当ξ∈Ωx(1)时d(x，ξ)的值外，第ith个最小的值。而m∈[2，7],这里取m＝5；
步骤六：将处理后的低频分量和高频分量进行轮廓波逆变换
步骤四处理过的高频分量和步骤五处理过的低频分量还需要使用轮廓波逆变换换原得到降噪后的CT图像。
进一步的，所述步骤二中对步骤一中含有噪声的原始CT图像使用轮廓波滤波器组进行多尺度多方向的分解包括如下步骤：首先对原始CT图像利用拉普拉斯金字塔滤波器进行多尺度分解，每个级别的拉普拉斯金字塔滤波器分解都会生成原始CT图像的降采样版本以及原始图像与预测之间的差异图像，从而产生带通图像；拉普拉斯金字塔滤波器具在每一级都只生成一个带通图像，同时，每一级的带通图像不具有加扰频率；在原始CT图像进行多尺度分解后，对拉普拉斯金字塔滤波器进行多尺度分解后所有金字塔层在每一级生成的带通图像使用二维方向滤波器组进行多方向分解，二维方向滤波器组由一个通道梅花形滤波器组成。通道梅花形滤波器包含扇形滤波器，一个剪切算子组成；扇形滤波器把图像方向大理成垂直和水平两个方向。后者相当于对图像严格进行重新排序。该分解通过l级二叉树分解有效实现的，该分解会导致2l个具有楔形频率划分的子带。
二维方向滤波器组只捕获图像的方向信息，所以对低频内容处理不善。实际上，低频分量会出没在几个高频子带中。因此，仅二维方向滤波器组不能提供图像的稀疏表示，使用拉普拉斯金字塔滤波器和二维方向滤波器组的组合则可以实现多尺度和多方向分解。组合结果是一个双重迭代的滤波器组结构，被称为轮廓波滤波器组。
在步骤二中仅对含噪图像做分解处理，分解后获得的图像细节高频分量和含噪低频分量会分别在后面的步骤中进行处理。
基于轮廓波变换的自适应降噪算法，主要针对分解后的轮廓波系数进行处理。该步骤是对根据给定的算法规则处理上一步图像信号分解后的轮廓波系数，在该过程中选取的是经典的阈值分析方法。众所周知，阈值函数T深受系数的数量的影响。M越大，阈值T越大，而较大的T值会过滤掉一些有用的信息，导致阈值函数在医学CT图像的噪声去除中失效。
本发明的有益效果在于：
1、本发明把分解得到的高频频分量通过自适应阈值收缩函数滤波，使高频分量的细节更好地保留
2、本发明利用轮廓线段实现了灵活的多分辨率，局部和方向性的图像扩展。
3、本发明克服了轮廓波变换在高维奇异性的无力，应用可以给医生进行更精确的诊断带来帮助。
4、本发明根据多种经典阈值函数，总结优点，提出了新的阈值函数。
5、本发明配合新提出的阈值函数，提出新的阈值收缩算法，有效处理了轮廓波高频分量的系数。
6、本发明使用脉冲权重三边滤波器降低轮廓波低频分量的系数。
7、本发明步骤简便，结构完整，有着完善的理论和实验支持。
附图说明
图1是本发明一种基于自适应阈值的轮廓波变换医学CT图像降噪方法的整体流程示意图。
图2是本发明轮廓波滤波器组原理图。
图3为各种算法应用在Barbara图的(噪声方差＝0.01)的噪声图和实验结果图，图3a为原图，图3b为噪声图，图3c是脊波变换结果图，图3d是曲波变换结果图，图3e时本发明算法结果图。
图4为各种算法应用在颅部医学CT图像(噪声方差＝0.01)的噪声图和实验结果图，图4a为噪声图，图4b为脊波变换图，图4c为曲波变换图，图4d为本发明算法图。
具体实施方式
下面结合附图对本发明作进一步说明：
如图1～2所示，一种基于自适应阈值的轮廓波变换医学CT图像降噪方法，具体包括以下步骤：
步骤一：获取CT图像模型
首先利用计算机断层扫描获取原始CT图像，原始CT图像接收到的信号氛围灰度信号和噪音信号，其中假设噪声信号为n(x，y)，原始信号为f(x，y)，输出信号为g(x，y)，根据原始CT图像与噪声之间的相互关系，将噪声分成加性噪声和乘性噪声两种形式，其中：
加性噪声与原始CT图像不相关，加性噪声的模型表示为
g(x，y)＝f(x，y)+n(x，y) (1)乘性噪声与原始CT图像相关，乘性噪声的模型表示为
g(x，y)＝f(x，y)*n(x，y) (2)
其中，(x，y)表示原始CT图像中的一串像素点，x，y分别表示原始CT图像中像素点的横坐标和纵坐标；
步骤二：对原始CT图像使用轮廓波滤波器组分解
对步骤一中含有噪声的原始CT图像使用轮廓波滤波器组进行多尺度多方向的分解，轮廓波滤波器组指的是拉普拉斯金字塔滤波器和二维方向滤波器组合而成的双重迭代的滤波器组结构；
由于步骤一中的原始CT图像为医学CT图像通用模型，对公式(1)进行轮廓波变换进行变换，可以得到新的加性噪声模型：
其中，
表示噪声图像的轮廓波系数，
表示无噪图像的轮廓波系数，
表示噪声的轮廓波系数，上标j表示轮廓波变换的分解层数，下标(l，k)，(l，k)∈Z2表示轮廓波域的坐标，J指最大分解层数；
步骤三：计算自适应于阈值算法的阈值
对步骤二种获得的轮廓波系数基于轮廓波变换的自适应降噪算法进行处理，
为了在平衡噪声抑制和特征保留过程中达到所期望的效果，综合经典阈值函数，并添加一个第jth层的自适应参数tj，提出一种新的阈值算法，该阈值算法的表达式如下：
其中，j(＝1，2，...，J)表示轮廓波变换的分解层数，J为最大分解层数，
其中，tj＝2J-j+1，k1＝k2＝0.5 (5)
在上公式(4)只有噪声标准差和无噪图像的标准差是未知的；噪声标准差是由轮廓波内高频的轮廓波系数，即式中轮廓波系数
的各个高频子带
绝对值的中值计算得到
因为轮廓波线性变换，因此由式可以得到j层轮廓波方差关系式
由于
和
被建模为零均值的模型，因此轮廓波j层的标准差可有轮廓波系数的标准差计算得到
这里的M＝m×n即是对应轮廓波域内轮廓波系数的总个数，m和n分别是对应轮廓波域内系数的行数和列数；
由公式(8)可得
因此公式(8)和公式(9)便可计算得到噪声标准差σn和无噪声图像的标准差σs，j，此时式中的阈值T也能计算得到；
步骤四：使用阈值收缩算法处理高频分量
在步骤三中得到阈值T后，对比阈值高的轮廓波系数按照收缩算法进行收缩，得到降噪后的轮廓波系数；轮廓波收缩函数利用了贝叶斯估计理论，使用的是对应于均匀型贝叶斯风险函数的最大后验概率，改进的小波收缩算法的推导和步骤如下：
步骤1.1：运用贝叶斯最大后验估计获得小波域中信号的估计；设置
其中，PR(r)服从广义拉普拉斯分布的特殊情况；因此
步骤1.2：根据第二章对噪声模型的分析可以假定噪声的轮廓波系数
服从均值为0的高斯分布：
带入等式，可得
(3)求出ln(P(R|S)(r|s))关于r的一阶导数；为了获取最大后验概率，将随后的代数表达式设为0，由此可得，
其中
是r的估计；
(4)最后，合成的轮廓波收缩公式如下：
步骤五：使用三边滤波器处理分解后的低频分量
对步骤二种产生的低频分量使用三边滤波器处理，具体表达式如下：
其中w′(x，ξ)是新的权重函数，有空域滤波器wS(x，ξ)、值域滤波器wR(x，ξ)和脉冲权重WI(ξ)组成，可表示为:
w′(x，ξ)＝wS(x，ξ)wR(x，ξ)wI(ξ) (17)
脉冲权重wI(ξ)的表达式为
其中引入fm(x)来估计像素点x是边缘点还是噪点；定义d(x，ξ)，表示点x和ξ之间的像素亮度的绝对差：
d(x，ξ)＝|f(x)-f(ξ)| (19)
由此可得，
其中，gi(x)表示出了当ξ∈Ωx(1)时d(x，ξ)的值外，第ith个最小的值。而m∈[2，7],这里取m＝5；
步骤六：将处理后的低频分量和高频分量进行轮廓波逆变换
步骤四处理过的高频分量和步骤五处理过的低频分量还需要使用轮廓波逆变换换原得到降噪后的CT图像。
进一步的，所述步骤二中对步骤一中含有噪声的原始CT图像使用轮廓波滤波器组进行多尺度多方向的分解包括如下步骤：首先对原始CT图像利用拉普拉斯金字塔滤波器进行多尺度分解，每个级别的拉普拉斯金字塔滤波器分解都会生成原始CT图像的降采样版本以及原始图像与预测之间的差异图像，从而产生带通图像；拉普拉斯金字塔滤波器具在每一级都只生成一个带通图像，同时，每一级的带通图像不具有加扰频率；在原始CT图像进行多尺度分解后，对拉普拉斯金字塔滤波器进行多尺度分解后所有金字塔层在每一级生成的带通图像使用二维方向滤波器组进行多方向分解，二维方向滤波器组由一个通道梅花形滤波器组成。通道梅花形滤波器包含扇形滤波器，一个剪切算子组成；扇形滤波器把图像方向大理成垂直和水平两个方向。后者相当于对图像严格进行重新排序。该分解通过l级二叉树分解有效实现的，该分解会导致2l个具有楔形频率划分的子带。
二维方向滤波器组只捕获图像的方向信息，所以对低频内容处理不善。实际上，低频分量会出没在几个高频子带中。因此，仅二维方向滤波器组不能提供图像的稀疏表示，使用拉普拉斯金字塔滤波器和二维方向滤波器组的组合则可以实现多尺度和多方向分解。组合结果是一个双重迭代的滤波器组结构，被称为轮廓波滤波器组。
在步骤二中仅对含噪图像做分解处理，分解后获得的图像细节高频分量和含噪低频分量会分别在后面的步骤中进行处理。
基于轮廓波变换的自适应降噪算法，主要针对分解后的轮廓波系数进行处理。该步骤是对根据给定的算法规则处理上一步图像信号分解后的轮廓波系数，在该过程中选取的是经典的阈值分析方法。众所周知，阈值函数T深受系数的数量的影响。M越大，阈值T越大，而较大的T值会过滤掉一些有用的信息，导致阈值函数在医学CT图像的噪声去除中失效。
对本申请进行实验分析，本发明以头部医学CT图像作为对象，以峰值信噪比PSNR和均方误差MSE作为实验分析的评价指标，其中PSNR越大，表示方法的降噪性能越好，而MSE越大则表明方法的降噪性能欠佳。实验以CT噪声经典图像Barbara图(像素大小512×512)作为自然图像数据输输入到实验中,实验分析整体流程如下图所示，实验通过脊波变换，曲波变换和本发明方法在含噪Barbara图的效果。实验结果图和含噪原图在图3，各种算法应用在头部CT图的效果及原图在图4。
在图像质量评价体系中，可以表示原图像与去噪后图像之间的均方误差，表达式如下
M和N分别别表示二维信号x的长和宽。
近年来，PSNR被证明与其他质量指标相比，在估计图像质量，特别是人类感知到的视频质量时，PSNR的表现较差。具体表达式如下：
下表表一各种算法的在Barbara图和肺部CT图像上的峰值信噪比和均方误差实验数据可以看出，当噪声等级减小时，即σ增大，噪声影响较明显，对降噪算法的性能要求更高。本发明算法在同一噪声等级时，效果优于脊波变换和曲波变换，而当噪声影响增大时，即噪声等级减小，仍能保持较好的降噪效果。
表1不同降噪算法在Barbara图(噪声误差sigma＝0.04)上的PSNR和MSE值
表2 Barbara图不同算法在不同噪声方差sigma下的PSNR/dB值
表3肺部CT图像不同变换方法在不同噪声防方差下的PSNR/dB值
上述实施例只是本发明的较佳实施例，并不是对本发明技术方案的限制，只要是不经过创造性劳动即可在上述实施例的基础上实现的技术方案，均应视为落入本发明专利的权利保护范围内。
Claims (2)
1.一种基于自适应阈值的轮廓波变换医学CT图像降噪方法，其特征在于：具体包括以下步骤：
步骤一：获取CT图像模型
首先利用计算机断层扫描获取原始CT图像，原始CT图像接收到的信号氛围灰度信号和噪音信号，其中假设噪声信号为n(x，y)，原始信号为f(x，y)，输出信号为g(x，y)，根据原始CT图像与噪声之间的相互关系，将噪声分成加性噪声和乘性噪声两种形式，其中：
加性噪声与原始CT图像不相关，加性噪声的模型表示为
g(x，y)＝f(x，y)+n(x，y) (1)
乘性噪声与原始CT图像相关，乘性噪声的模型表示为
g(x，y)＝f(x，y)*n((x，y) (2)
其中，(x，y)表示原始CT图像中的一串像素点，x，y分别表示原始CT图像中像素点的横坐标和纵坐标；
步骤二：对原始CT图像使用轮廓波滤波器组分解
对步骤一中含有噪声的原始CT图像使用轮廓波滤波器组进行多尺度多方向的分解，轮廓波滤波器组指的是拉普拉斯金字塔滤波器和二维方向滤波器组合而成的双重迭代的滤波器组结构；
由于步骤一中的原始CT图像为医学CT图像通用模型，对公式(1)进行轮廓波变换进行变换，可以得到新的加性噪声模型：
其中，
表示噪声图像的轮廓波系数，
表示无噪图像的轮廓波系数，
表示噪声的轮廓波系数，上标j表示轮廓波变换的分解层数，下标(l，k)，(l，k)∈Z2表示轮廓波域的坐标，J指最大分解层数；
步骤三：计算自适应于阈值算法的阈值
对步骤二种获得的轮廓波系数基于轮廓波变换的自适应降噪算法进行处理，
为了在平衡噪声抑制和特征保留过程中达到所期望的效果，综合经典阈值函数，并添加一个第jth层的自适应参数tj，提出一种新的阈值算法，该阈值算法的表达式如下：
其中，j(＝1，2，...，J)表示轮廓波变换的分解层数，J为最大分解层数，
其中，tj＝2J-j+1，k1＝k2＝0.5 (5)
在上公式(4)只有噪声标准差和无噪图像的标准差是未知的；噪声标准差是由轮廓波内高频的轮廓波系数，即式中轮廓波系数
的各个高频子带
绝对值的中值计算得到
因为轮廓波线性变换，因此由式可以得到j层轮廓波方差关系式
由于
和
被建模为零均值的模型，因此轮廓波j层的标准差可有轮廓波系数的标准差计算得到
这里的M＝m×n即是对应轮廓波域内轮廓波系数的总个数，m和n分别是对应轮廓波域内系数的行数和列数；
由公式(8)可得
因此公式(8)和公式(9)便可计算得到噪声标准差σn和无噪声图像的标准差σs，j，此时式中的阈值T也能计算得到；
步骤四：使用阈值收缩算法处理高频分量
在步骤三中得到阈值T后，对比阈值高的轮廓波系数按照收缩算法进行收缩，得到降噪后的轮廓波系数；轮廓波收缩函数利用了贝叶斯估计理论，使用的是对应于均匀型贝叶斯风险函数的最大后验概率，改进的小波收缩算法的推导和步骤如下：
步骤1.1：运用贝叶斯最大后验估计获得小波域中信号的估计；设置
其中，PR(r)服从广义拉普拉斯分布的特殊情况；因此
步骤1.2：根据第二章对噪声模型的分析可以假定噪声的轮廓波系数
服从均值为0的高斯分布：
带入等式，可得
(3)求出ln(P(R|S)(r|s))关于r的一阶导数；为了获取最大后验概率，将随后的代数表达式设为0，由此可得，
其中
是r的估计；
(4)最后，合成的轮廓波收缩公式如下：
步骤五：使用三边滤波器处理分解后的低频分量
对步骤二种产生的低频分量使用三边滤波器处理，具体表达式如下：
其中w′(x，ξ)是新的权重函数，有空域滤波器wS(x，ξ)、值域滤波器wR(x，ξ)和脉冲权重wI(ξ)组成，可表示为：
w′(x，ξ)＝wS(x，ξ)wR(x，ξ)wI(ξ) (17)
脉冲权重wI(ξ)的表达式为
其中引入fm(x)来估计像素点x是边缘点还是噪点；定义d(x，ξ)，表示点x和ξ之间的像素亮度的绝对差：
d(x，ξ)＝|f(x)-f(ξ)| (19)
由此可得，
其中，gi(x)表示出了当ξ∈Ωx(1)时d(x，ξ)的值外，第ith个最小的值。而m∈[2，7]，这里取m＝5；
步骤六：将处理后的低频分量和高频分量进行轮廓波逆变换
步骤四处理过的高频分量和步骤五处理过的低频分量还需要使用轮廓波逆变换换原得到降噪后的CT图像。
2.根据权利要求1所述的基于自适应阈值的轮廓波变换医学CT图像降噪方法，其特征在于：所述步骤二中对步骤一中含有噪声的原始CT图像使用轮廓波滤波器组进行多尺度多方向的分解包括如下步骤：首先对原始CT图像利用拉普拉斯金字塔滤波器进行多尺度分解，每个级别的拉普拉斯金字塔滤波器分解都会生成原始CT图像的降采样版本以及原始图像与预测之间的差异图像，从而产生带通图像；拉普拉斯金字塔滤波器具在每一级都只生成一个带通图像，同时，每一级的带通图像不具有加扰频率；在原始CT图像进行多尺度分解后，对拉普拉斯金字塔滤波器进行多尺度分解后所有金字塔层在每一级生成的带通图像使用二维方向滤波器组进行多方向分解，二维方向滤波器组由一个通道梅花形滤波器组成。通道梅花形滤波器包含扇形滤波器，一个剪切算子组成；扇形滤波器把图像方向大理成垂直和水平两个方向。后者相当于对图像严格进行重新排序。该分解通过l级二叉树分解有效实现的，该分解会导致2l个具有楔形频率划分的子带。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011620169.9A CN112734663A (zh) | 2020-12-30 | 2020-12-30 | 基于自适应阈值的轮廓波变换医学ct图像降噪方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011620169.9A CN112734663A (zh) | 2020-12-30 | 2020-12-30 | 基于自适应阈值的轮廓波变换医学ct图像降噪方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN112734663A true CN112734663A (zh) | 2021-04-30 |
Family
ID=75608500
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011620169.9A Withdrawn CN112734663A (zh) | 2020-12-30 | 2020-12-30 | 基于自适应阈值的轮廓波变换医学ct图像降噪方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112734663A (zh) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114240941A (zh) * | 2022-02-25 | 2022-03-25 | 浙江华诺康科技有限公司 | 内窥镜图像降噪方法、设备、电子装置和存储介质 |
| CN114255182A (zh) * | 2021-12-13 | 2022-03-29 | 郑州轻工业大学 | 基于空间自适应全变差的cs迭代阈值图像去噪重建方法 |
-
2020
- 2020-12-30 CN CN202011620169.9A patent/CN112734663A/zh not_active Withdrawn
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114255182A (zh) * | 2021-12-13 | 2022-03-29 | 郑州轻工业大学 | 基于空间自适应全变差的cs迭代阈值图像去噪重建方法 |
| CN114240941A (zh) * | 2022-02-25 | 2022-03-25 | 浙江华诺康科技有限公司 | 内窥镜图像降噪方法、设备、电子装置和存储介质 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Diwakar et al. | A review on CT image noise and its denoising | |
| Naimi et al. | Medical image denoising using dual tree complex thresholding wavelet transform and Wiener filter | |
| US20060110061A1 (en) | Image processing apparatus, image processing method, storage medium, and program | |
| Sudha et al. | Speckle noise reduction in ultrasound images using context-based adaptive wavelet thresholding | |
| CN109598680B (zh) | 基于快速非局部均值和tv-l1模型的剪切波变换医学ct图像去噪方法 | |
| CN109003232B (zh) | 基于频域尺度平滑Shearlet的医学MRI图像去噪方法 | |
| Thakur et al. | Image quality based comparative evaluation of wavelet filters in ultrasound speckle reduction | |
| Qian et al. | Tree-structured nonlinear filter and wavelet transform for microcalcification segmentation in mammography | |
| CN112734663A (zh) | 基于自适应阈值的轮廓波变换医学ct图像降噪方法 | |
| Diwakar et al. | CT Image noise reduction based on adaptive wiener filtering with wavelet packet thresholding | |
| Kumar et al. | A new locally adaptive patch variation based CT image denoising | |
| Lee et al. | Poisson-gaussian noise reduction for x-ray images based on local linear minimum mean square error shrinkage in nonsubsampled contourlet transform domain | |
| Kavitha et al. | COVID-19 and MRI image denoising using wavelet transform and basic filtering | |
| Gabralla et al. | Denoising CT Images using wavelet transform | |
| AlAsadi | Contourlet transform based method for medical image denoising | |
| CN108846813A (zh) | 基于mfdf分解框架与nsst的医学ct图像去噪方法 | |
| Diwakar et al. | Inter-and intra-scale dependencies-based CT image denoising in curvelet domain | |
| Lazrag et al. | Despeckling of intravascular ultrasound images using curvelet transform | |
| CN109035156B (zh) | 基于dnst的医学ct图像去噪方法 | |
| Diwakar et al. | Blind noise estimation-based CT image denoising in tetrolet domain | |
| Hosseinian et al. | Assessment of restoration methods of X-ray images with emphasis on medical photogrammetric usage | |
| CN108416737A (zh) | 基于dnst的医学ct图像去噪方法 | |
| Diwakar et al. | Internet of medical things: A CT image denoising in tetrolet domain | |
| US20230133074A1 (en) | Method and Apparatus for Noise Reduction | |
| Zhang et al. | A novel denoising method for medical ct images based on moving decomposition framework |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| WW01 | Invention patent application withdrawn after publication | ||
| WW01 | Invention patent application withdrawn after publication |
Application publication date: 20210430 |