1 双树复小波降噪
1.1 双树复小波变换
1.2 双树复小波降噪
1.3 双树复小波降噪的python实现
2.1 经验小波变换
2.2 经验小波降噪
2.3 经验小波降噪的python实现
3.1 测试用例
3.2 降噪结果
3.3 降噪指标
1.1 双树复小波变换
双树复小波变换(Dual-Tree Complex Wavelet Transforms, DT-CWT)由实部树和虚部树的两个并行的实小波变换构成[1]。复小波表示为:
图1 双树复小波变换流程
可以看出,双树复小波变换的分解和重构实现了实部树和虚部树信息的互补,能够保持信号的完全重构性。
1.2 双树复小波降噪
双树复小波变换将原始信号分解成小波系数和尺度系数,也称为细节分量和近似分量,其中噪声成分主要集中在细节分量中(高频部分)。与小波降噪相同,随着分解层数的增加,噪声的能量逐渐衰减,双树复小波分解得到的小波系数变小,确定合理的阈值可以对噪声进行剔除。阈值确定方法详见信号处理基础之噪声与降噪(五) | 小波降噪及python代码实现。
双树复小波降噪的基本流程如下:
(1)信号分解:使用双树复小波变换对信号进行多级分解。DT-CWT 生成两棵树的小波系数,每棵树提供了信号的一半信息,这两棵树合在一起可以提供更完整的频率和方向信息;
(2) 阈值处理:对小波系数进行阈值处理以去除噪声。常用的阈值方法包括软阈值和硬阈值处理。软阈值可以更平滑地处理信号,而硬阈值可能会导致信号失真;
(3)系数重构:使用处理后的小波系数通过逆双树复小波变换重构信号。
1.3 双树复小波降噪的python实现
双树复小波降噪的python实现如下:
import pywt
import numpy as np
import matplotlib.pyplot as plt
def denoised_with_dtwt(signal, wavelet, mode='soft', level=1):
coeffs = pywt.wavedec(signal, wavelet, mode='per')
sigma = np.median(np.abs(coeffs[-1])) / 0.6745
threshold = sigma * np.sqrt(2 * np.log(len(signal)))
new_coeffs = [pywt.threshold(c, value=threshold, mode=mode) for c in coeffs]
new_coeffs[0] = coeffs[0] # 保留近似系数不变
signal_denoised = pywt.waverec(new_coeffs, wavelet, mode='per')
return signal_denoised
2.1 经验小波变换
图2 傅里叶轴的分割
(4)原始信号重构:原始信号可表示为
2.2 经验小波降噪
在经验小波变换的基础上,结合峭度法和相关系数法等方法(详见信号处理基础之噪声与降噪(三) | EMD降噪与VMD降噪及python代码实现)可对信号中的噪声进行去除,基本步骤如下:
2.3 经验小波降噪的python实现
经验小波降噪的python实现如下:
import numpy as np
import matplotlib.pyplot as plt
import ewtpy
def denoised_with_ewt(signal):
ewt, mfb, boundaries = ewtpy.EWT1D(signal)
ewt = ewt.T
ewt_denoised = [component for component in ewt if np.corrcoef(signal, component)[0, 1] >= 0.4]
signal_denoised = np.sum(ewt_denoised, axis=0)
return signal_denoised
import matplotlib.pyplot as plt
import matplotlib
# 创建一个信号
n = 500 # 生成500个点的信号
t = np.linspace(0, 30*np.pi, n, endpoint=False)
s = np.cos(0.1*np.pi*t) + np.sin(0.3*np.pi*t) + np.cos(0.5*np.pi*t) + np.sin(0.7*np.pi*t) # 原始信号
r = 0.5 * np.random.randn(n)
y = s + r # 加噪声
denoised_signal_dtwt = denoised_with_dtwt(y, 'db4', 'soft', level=8)
denoised_signal_ewt = denoised_with_ewt(y)
value_dtwt = denoise_evaluate(s, r, denoised_signal_dtwt)
value_ewt = denoise_evaluate(s, r, denoised_signal_ewt)
print(value_dtwt)
print(value_ewt)
fig = plt.figure(figsize=(16, 12))
plt.subplots_adjust(hspace=0.5)
plt.subplots_adjust(left=0.05, right=0.98, top=0.95, bottom=0.05)
font = {'family': 'Times New Roman', 'size': 16, 'weight': 'normal',}
matplotlib.rc('font', **font)
plt.subplot(5, 1, 1)
plt.plot(t, s, linewidth=1.5, color='b')
plt.title('pure signal', fontname='Times New Roman', fontsize=20)
plt.subplot(5, 1, 2)
plt.plot(t, y, linewidth=1.5, color='b')
plt.title('noised signal', fontname='Times New Roman', fontsize=20)
plt.subplot(5, 1, 3)
plt.plot(t, denoised_signal_dtwt, linewidth=1.5, color='b')
plt.title('DTWT', fontname='Times New Roman', fontsize=20)
plt.subplot(5, 1, 4)
plt.plot(t, denoised_signal_ewt, linewidth=1.5, color='b')
plt.title('EWT', fontname='Times New Roman', fontsize=20)
plt.subplot(5, 1, 5)
plt.plot(t, y-denoised_signal_ewt, linewidth=1.5, color='b')
plt.title('y-EWT', fontname='Times New Roman', fontsize=20)
plt.show()