一维信号小波去噪原理解析及python实现方式

pip install PyWavelets

def wavelet_denoising(signal):
    # signal: list
    wavelet_basis = 'db8'  # 此处选择db8小波基
    w = pywt.Wavelet(wavelet_basis) 
    maxlevel = pywt.dwt_max_level(len(signal), w.dec_len) # 根据数据长度计算分解层数
    threshold = 0.1  # threshold for filtering
    coeffs = pywt.wavedec(signal, wavelet_basis, level=maxlevel)  # 将信号进行小波分解
    # 使用阈值滤除噪声
    for i in range(1, len(coeffs)):
        coeffs[i] = pywt.threshold(coeffs[i], threshold*max(coeffs[i]))  # 将噪声滤波
        # pywt.threshold 方法有三种阈值模式：硬阈值、软阈值、软硬结合
        # 默认为soft 软阈值模式，其调用格式为
        # pywt.threshold(signal, threshold_value, mode='soft', substitute=0)
        # 处理后值为 signal/np.abs(signal) * np.maximum(np.abs(signal) - threshold_value, 0)
        # 具体可参考pywavelets官方文档的说明
    # 重建信号
    signalrec= pywt.waverec(coeffs, wavelet_basis)
    return signalrec

a = np.loadtxt('C:\\Users\\Administrator\\Desktop\\数据文件\\20201008009.txt', dtype=float)# 20201008009为一个一维列表
index = []
data = []
coffs = []
for i in range(len(a) - 1):
    X = float(i)
    Y = float(a[i])
    index.append(X)
    data.append(Y)
# create wavelet object and define parameters
w = pywt.Wavelet('db8')  # 选用Daubechies8小波
maxlev = pywt.dwt_max_level(len(data), w.dec_len)
print("maximum level is" + str(maxlev))
threshold = 0.1  # Threshold for filtering
# Decompose into wavelet components,to the level selected:
coffs = pywt.wavedec(data, 'db8', level=maxlev)  # 将信号进行小波分解
for i in range(1, len(coffs)):
    coffs[i] = pywt.threshold(coffs[i], threshold * max(coffs[i]))
datarec = pywt.waverec(coffs, 'db8')  # 将信号进行小波重构
mintime = 0
maxtime = mintime + len(data)
print(mintime, maxtime)
plt.rcParams['font.sans-serif']=['SimHei'] #使得标题输出为中文
plt.rcParams['axes.unicode_minus'] = False #使得标题输出为中文
plt.figure(1)
#plt.subplot(3, 1, 1)
plt.plot(index[mintime:maxtime], data[mintime:maxtime])
plt.title("原始信号")
plt.figure(2)
#plt.subplot(3, 1, 2)
plt.plot(index[mintime:maxtime], datarec[mintime:maxtime])
plt.title("小波降噪信号")
plt.figure(3)
#plt.subplot(3, 1, 3)
plt.plot(index[mintime:maxtime], data[mintime:maxtime] - datarec[mintime:maxtime])
plt.xlabel('time (s)')
plt.ylabel('error (uV)')
plt.tight_layout()
plt.show()