EMD降噪python代码

以下是一些使用EMD进行降噪的Python代码示例：

1. 使用PyEMD库实现EMD降噪

from PyEMD import EMD
import numpy as np
import matplotlib.pyplot as plt

# 生成带噪声的信号
t = np.linspace(0, 1, 1000)
s = np.sin(2*np.pi*10*t) + np.sin(2*np.pi*20*t) + np.random.randn(len(t))*0.5

# 使用EMD降噪
emd = EMD()
IMFs = emd(s)
residual = IMFs[-1]
denoised = s - residual

# 绘制结果
plt.subplot(2, 1, 1)
plt.plot(t, s)
plt.title('Original Signal')
plt.subplot(2, 1, 2)
plt.plot(t, denoised)
plt.title('Denoised Signal')
plt.show()

2. 使用EMD算法手动实现降噪

import numpy as np
import matplotlib.pyplot as plt

def emd(signal, max_imf=None):
    # 停止条件
    if max_imf and max_imf <= 0:
        return [signal]
    # 计算均值
    mean = np.mean(signal)
    # 极值点检测
    min_idx = []
    max_idx = []
    for i in range(1, len(signal)-1):
        if signal[i] > signal[i-1] and signal[i] > signal[i+1]:
            max_idx.append(i)
        elif signal[i] < signal[i-1] and signal[i] < signal[i+1]:
            min_idx.append(i)
    if not min_idx or not max_idx:
        return [signal]
    # 构造上下包络线
    upper = np.interp(range(len(signal)), max_idx, signal[max_idx])
    lower = np.interp(range(len(signal)), min_idx, signal[min_idx])
    # 计算平均包络线
    mean_env = (upper + lower) / 2
    # 计算当前IMF
    imf = signal - mean_env
    # 如果当前IMF为常量，则停止分解
    if np.allclose(imf, np.zeros_like(imf)):
        return [signal]
    # 分解下一层IMF
    rest = emd(imf, max_imf-1) if max_imf else emd(imf)
    return rest + [imf]

# 生成带噪声的信号
t = np.linspace(0, 1, 1000)
s = np.sin(2*np.pi*10*t) + np.sin(2*np.pi*20*t) + np.random.randn(len(t))*0.5

# 使用EMD降噪
IMFs = emd(s)
residual = IMFs[-1]
denoised = s - residual

# 绘制结果
plt.subplot(2, 1, 1)
plt.plot(t, s)
plt.title('Original Signal')
plt.subplot(2, 1, 2)
plt.plot(t, denoised)
plt.title('Denoised Signal')
plt.show()