小波变换的前因后果（三）

小波去噪是建立在DWT的基础上的，需要进行小波分解、再重构。接上一篇。

小波分析即用Mallat塔式算法对信号进行降阶分解。该算法在每尺度下将信号分解成近似分量与细节分量。近似分量表示信号的高尺度，即低频信息；细节分量表示信号的低尺度，即高频信息。对含有噪声的信号，噪声分量的主要能量集中在小波分解的细节分量中。

二、小波去噪

1、概念

通常情况下，我们在从设备上采集到的信号都是具有一定的噪声的，大多数情况下，可认为这种噪声为高斯白噪声。被噪声污染的信号=干净的信号+噪声。

为什么要使用阈值：由于信号在空间上(或者时间域)是有一定连续性的，因此在小波域，有效信号所产生的小波系数其模值往往较大；而高斯白噪声在空间上(或者时间域)是没有连续性的，因此噪声经过小波变换，在小波阈仍然表现为很强的随机性，通常仍认为是高斯白噪的。那么就得到这样一个结论：在小波域，有效信号对应的系数很大，而噪声对应的系数很小。刚刚已经说了，噪声在小波域对应的系数仍满足高斯白噪分布。如果在小波域，噪声的小波系数对应的方差为sigma，那么根据高斯分布的特性，绝大部分(99.99%)噪声系数都位于[-3*sigma，3*sigma]区间内(切比雪夫不等式， 3sigma准则)。因此，只要将区间[-3*sigma，3*sigma]内的系数置零(这就是常用的硬阈值函数的作用)，就能最大程度抑制噪声的，同时只是稍微损伤有效信号。将经过阈值处理后的小波系数重构，就可以得到去噪后的信号。常用的软阈值函数，是为了解决硬阈值函数“一刀切”导致的影响(模小于3*sigma的小波系数全部切除，大于3*sigma全部保留，势必会在小波域产生突变，导致去噪后结果产生局部的抖动，类似于傅立叶变换中频域的阶跃会在时域产生拖尾)。软阈值函数将模小于3*sigma的小波系数全部置零，而将模大于3*sigma的做一个比较特殊的处理，大于3*sigma的小波系数统一减去3*sigma，小于-3*sigma的小波系数统一加3*sigma。经过软阈值函数的作用，小波系数在小波域就比较光滑了，因此用软阈值去噪得到的图象看起来很平滑，类似于冬天通过窗户看外面一样，像有层雾罩在图像上似的。

比较硬阈值函数去噪和软阈值函数去噪：硬阈值函数去噪所得到的峰值信噪比(PSNR)较高，但是有局部抖动的现象；软阈值函数去噪所得到的PSNR不如硬阈值函数去噪，但是结果看起来很平滑，原因就是软阈值函数对小波系数进行了较大的 “社会主义改造”，小波系数改变很大。因此各种各样的阈值函数就出现了，其目的我认为就是要使大的系数保留，小的系数被剔出，而且在小波域系数过渡要平滑。

如何估计小波域噪声方差sigma的估计，这个很简单：把信号做小波变换，在每一个子带利用robust estimator估计就可以(可能高频带和低频带的方差不同)。 robust estimator就是将子带内的小波系数模按大小排列，然后取最中间那个，然后把最中间这个除以0.6745就得到噪声在某个子带内的方差sigma。利用这个sigma，然后选种阈值函数，就可以去去噪了，在matlab有实现api可使用。

小波阈值去噪过程如下图

在小波分析中经常用到近似和细节，近似表示信号的高尺度，即低频信息；细节表示信号的低尺度，即高频信息。对含有噪声的信号，噪声分量的主要能量集中在小波解的细节分量中。

2、原理

小波阈值去噪的实质为抑制信号中无用部分、增强有用部分的过程。小波阈值去噪过程为：(1)分解过程，即选定一种小波对信号进行n层小波分解；(2)阈值处理过程，即对分解的各层系数进行阈值处理，获得估计小波系数；(3)重构过程，据去噪后的小波系数进行小波重构，获得去噪后的信号。

小波阈值去噪过程

小波分解重构过程

小波分解：X->ca3,cd3,cd2,cd1；小波重构：ca3,cd3,cd2,cd1->X。其中ca为低频信息、近似分量，cd为高频、细节分量。

3、影响降噪效果的因素

（1）小波基的选择

在对信号进行小波分解时需要选择合适的小波基，由于没有任何一种小波基可以对不同类型的信号达到最优的分解效果，因此，如何选择小波基成为小波分解的一个重点。针对现实中的信号，小波基的选择一般要考虑以下几个因素：支撑长度、对称性、消失矩、正则性、相似性。针对一维信号，例如语音信号，通常选择dB小波和sym小波。

（2）分解层数的选择

在对信号进行小波分解时，分解的层数取得越大，则噪声和信号表现的不同特性越明显，越有利于二者的分离，但是分解的层数越大，经过重构的信号失真也会越大，在一定程度上会对信号去噪的效果产生较差的影响。因此，如何选择分解层数以解决信噪分离效果和重构信号失真之间的矛盾呢？

小波分解的频段范围与采样频率有关。若进行N层分解，则各个频段范围为：

假设原始信号X的采样频率为1000Hz，则信号的最大频率为500，对该信号做3层小波分解，则各个频段范围如下图所示。

通常小波分解的频段范围与采样频率有关。若N层分解，则各个频段大小为Fs/2/2^N 。例如：一个原始信号，经历的时间长度为2秒，采样了2000个点，那么做除法，可得出采样频率为1000hz，由采样定理（做除法）得该信号的最大频率为500hz，那么对该信号做3层的DWT，一阶细节的频段为250－500hz，一阶逼近的频段为小于250hz，二阶细节的频段为125－250hz，逼近的频段为小于125hz，三阶细节的频段约为62.5－125hz，逼近的频段为小于62.5hz。对于更多阶的分解也是以此类推的。

（3）具体的分解算法（Mallat塔式算法）

多分辨率分析为正交小波基的构造提供了一种简单方法，同时它还是正交小波变换的快速算法(即 Mallat 算法)基本理论基础。Mallat 算法是由 S.Mallat于1989年提出的，该算法在小波分析中的作用相当于快速傅立叶变换在傅里叶分析中的作用。Mallat 算法由小波滤波器 H、G、和 h、g 对测量的信号进行分解和重构。它的分解算法可以表述为：

分解算法可以用图解的形式表示为：

式中，h，g为时域中的小波重构滤波器，j为分解的层数，若分解的深度(分解的最高层)为J，则j=J−1,J−2,…,1,0，Aj为信号f(t)在第j层的低频部分(近似部分)的小波系数；Dj为信号f(t)在第 j 层的高频部分(细节部分)的小波系数。

这里提供了一维信号下的小波分解代码

"""
将数据序列进行小波分解，每一层分解的结果是上次分解得到的低频信号再分解成低频和高频两个部分。
如此进过N层分解后源信号X被分解为：X = D1 + D2 + … + DN + AN
其中D1,D2,…,DN分别为第一层、第二层到等N层分解得到的高频信号，AN为第N层分解得到的低频信号。
"""
import numpy as np
import matplotlib.pyplot as pltimport pywt
import pywt.dataecg = pywt.data.ecg()
print(len(ecg))
data1 = np.concatenate((np.arange(1, 400),np.arange(398, 600),np.arange(601, 1024)))
x = np.linspace(0.082, 2.128, num=1024)[::-1]
data2 = np.sin(40 * np.log(x)) * np.sign((np.log(x)))mode = pywt.Modes.smoothdef plot_signal_decomp(data, w, title):"""Decompose and plot a signal S.S = An + Dn + Dn-1 + ... + D1"""w = pywt.Wavelet(w)#选取小波函数a = dataca = []#近似分量cd = []#细节分量for i in range(5):(a, d) = pywt.dwt(a, w, mode) #进行5阶离散小波变换ca.append(a)cd.append(d)rec_a = []rec_d = []for i, coeff in enumerate(ca):coeff_list = [coeff, None] + [None] * irec_a.append(pywt.waverec(coeff_list, w)) #重构for i, coeff in enumerate(cd):coeff_list = [None, coeff] + [None] * iif i == 3:print(len(coeff))print(len(coeff_list))rec_d.append(pywt.waverec(coeff_list, w))fig = plt.figure()ax_main = fig.add_subplot(len(rec_a) + 1, 1, 1)ax_main.set_title(title)ax_main.plot(data)ax_main.set_xlim(0, len(data) - 1)for i, y in enumerate(rec_a):ax = fig.add_subplot(len(rec_a) + 1, 2, 3 + i * 2)ax.plot(y, 'r')ax.set_xlim(0, len(y) - 1)ax.set_ylabel("A%d" % (i + 1))for i, y in enumerate(rec_d):ax = fig.add_subplot(len(rec_d) + 1, 2, 4 + i * 2)ax.plot(y, 'g')ax.set_xlim(0, len(y) - 1)ax.set_ylabel("D%d" % (i + 1))# plot_signal_decomp(data1, 'coif5', "DWT: Signal irregularity")
# plot_signal_decomp(data2, 'sym5',
#                   "DWT: Frequency and phase change - Symmlets5")
plot_signal_decomp(ecg, 'sym5', "DWT: Ecg sample - Symmlets5")plt.show()
"""
将数据序列进行小波分解，每一层分解的结果是上次分解得到的低频信号再分解成低频和高频两个部分。
如此进过N层分解后源信号X被分解为：X = D1 + D2 + … + DN + AN
其中D1,D2,…,DN分别为第一层、第二层到等N层分解得到的高频信号，AN为第N层分解得到的低频信号。
"""

二维信号的小波分解代码如下

# 采用小波阈值方法对预先加了标准差为15的高斯噪声图像(256*256)进行小波阈值去噪，
# 其中各层的滤波阈值均为同一值，采用的小波基为sym4小波，分解层数为3，阈值函数为软阈值函数"""
原理：
小波阈值去噪的是指为抑制信号中无用部分、增强有用部分的过程，小波阈值去噪过程为：
1、分解过程，即选定一种小波对信号进行n层小波分解
2、阈值处理过程，即对分解的隔层系数进行阈值处理，获得估计小波系数
3、重构过程，据去噪后的小波系数进行小波重构，获得去噪后的信号影响效果的因素：分解层数、阈值、小波基的选择、阈值函数的选择
"""import numpy as np
import pywt
from cv2 import COLOR_BGR2GRAY, COLOR_RGB2GRAY, cv2, denoise_TVL1
from PIL import Image# img = cv2.imread("lenage15.bmp",COLOR_BGR2GRAY)
img = cv2.imread("lenage15.bmp",0)
print(img.shape)
w = 'sym4' #定义小波基的类型
l = 3 # 变换层次为3
coeffs = pywt.wavedec2(data=img, wavelet=w, level=l) # 对图像进行小波分解
"""
coeffs shape：
4*38*38
"""
threshold = 0.04list_coeffs = []
for i in range(1, len(coeffs)):list_coeffs_ = list(coeffs[i]) # 38*38list_coeffs.append(list_coeffs_)
# print(list_coeffs) # list_coeffs.shape:4*38*38for r1 in range(len(list_coeffs)): # 4for r2 in range(len(list_coeffs[r1])): # 38# 对噪声滤波（软阈值）list_coeffs[r1][r2] = pywt.threshold(list_coeffs[r1][r2], threshold*np.max(list_coeffs[r1][r2]))rec_coeffs = [] # 重构系数
rec_coeffs.append(coeffs[0]) #将原图像的低尺度系数保留进来for j in range(len(list_coeffs)):rec_coeffs_ = tuple(list_coeffs[j])rec_coeffs.append(rec_coeffs_)denoised_img = pywt.waverec2(rec_coeffs, 'sym4')
denoised_img = Image.fromarray(np.uint8(denoised_img))
denoised_img.save("result.bmp")
# 总结：这里能够看到小波去噪对去除高斯噪声的效果还是可以的，但是效果不是很明显。
# 大家可以从阈值的选择以及调整滤波的阈值，甚至将各层的阈值设置为不同的值再去试试。
# 再看看有关的论文关于提升小波阈值效果的实现方法。# 本文只是举个简单的例子来说明小波阈值的处理过程

小波变换的前因后果（三）相关推荐

关于小波变换的一些理解
原地址 http://blog.csdn.net/chenyusiyuan/article/details/1864195 1.小波与傅立叶变换任何学科都是由一门基本学科积累发展起来的,要做到学 ...
基于小波的图像边缘检测,小波变换边缘检测原理
1.什么是"小波神经网络"?能干什么用呀小波神经网络(Wavelet Neural Network, WNN)是在小波分析研究获得突破的基础上提出的一种人工神经网络.它是基于小波 ...
雷达原理---时频分析--3.小波变换-3.1基础知识
文章目录一.短时傅里叶变换的缺陷二.小波变换的优点三.小波变换和傅里叶变换的比较四.小波变换的基础知识(Wavelet Transform,WT) 1. 连续小波变换(Continuous W ...
mne库脑电时频信号分析函数解读：小波变换及方法比较
2023/1/8-2023/1/9 脑机接口学习内容一览: 这一篇博客里,主要研究mne库中的函数mne.time_frequency.tfr_morlet如何完成时频信号分析,提供基本的函数功能.参 ...
从傅里叶(Fourier)变换到伽柏(Gabor)变换再到小波(Wavelet)变换
从傅里叶(Fourier)变换到伽柏(Gabor)变换再到小波(Wavelet)变换本文是边学习边总结和摘抄各参考文献内容而成的,是一篇综述性入门文档,重点在于梳理傅里叶变换到伽柏变换再到小波变换的 ...
读书｜《静心冥想的练习》：体验超越一切理解的平静与快乐
冥想是瑜伽实现入定的一项技法和途径,把心.意.灵完全专注在原始之初之中,最终目的在于把人引导到解脱的境界.瑜伽者通过冥想来制服心灵,并超脱物质欲念:感受到和原始动因直接沟通.通过简单练习冥想,即可帮助 ...
荧光和明场图像融合 matlab,一种用于明场显微成像的多层图像融合算法
一种用于明场显微成像的多层图像融合算法 [技术领域] [0001] 本发明涉及图像处理技术领域,特别涉及一种用于明场显微成像的多层图像融合算法. [背景技术] [0002] 当前对细胞形态表型研宄的 ...
西安工业大学计算机学院科协,卢文科
卢文科教授卢文科,1962年5月出生,西安交通大学博士,东华大学控制科学与工程学科教授.博士生导师,主要从事声表面波.静磁波.小波变换.传感器理论及技术的研究,尤其在小波变换处理器及小波反变换处理器 ...
何恺明暗通道去雾(阅读笔记)
何恺明暗通道去雾文章阅读刘海山,2021年6月24日文献引文信息: He, K. M., Sun, J. & Tang, X. O. Single Image Haze Removal U ...

小波变换的前因后果（三）

小波去噪是建立在DWT的基础上的，需要进行小波分解、再重构。接上一篇。

二、小波去噪

1、概念

2、原理

3、影响降噪效果的因素

小波变换的前因后果（三）相关推荐

最新文章

热门文章