经验模态分解（EMD）python实现及示例程序

最新推荐文章于 2024-07-22 18:34:50 发布

Captain_Koucha

最新推荐文章于 2024-07-22 18:34:50 发布

阅读量812

点赞数 5

文章标签： python 开发语言信号处理

本文链接：https://blog.csdn.net/f2ancisz/article/details/139331218

版权

什么是经验模态分解（EMD）？

EMD是一种数据分析方法，特别适用于处理非线性和非平稳信号。它能将复杂的信号分解成若干个简单的本征模态函数（Intrinsic Mode Functions，IMFs）和一个残差。

为什么需要EMD？

在现实生活中，许多信号并不是线性和稳定的。例如，天气数据、地震波、金融市场数据等。传统的分析方法（如傅里叶变换）对这些信号的处理效果有限。EMD可以更好地处理这些复杂信号。

EMD的基本原理

EMD的核心思想是将复杂信号分解成若干个更简单的成分。具体步骤如下：
寻找局部极值：对原始信号，找到所有的局部极大值和极小值。
构建包络线：通过这些极值点，使用样条曲线构建上下包络线。
计算均值：上下包络线的均值即为当前信号的局部趋势。
去趋势：将原始信号减去局部趋势，得到一个新的信号。
重复步骤：对新的信号重复上述步骤，直到得到的信号满足IMF的条件：信号在零线两侧的极值点数目和过零点数目相差不超过一个。
通过上述步骤，可以将复杂信号分解为若干个IMF和一个残差。

EMD的优缺点

优点：
适应性强：不需要预设基函数，适合处理非线性、非平稳信号。
局部性好：可以捕捉信号的局部特征，保留细节信息。
缺点：
边界效应：在信号两端可能产生误差，需要特别处理。
模式混淆：不同IMF之间可能出现混叠现象，影响分解效果。

总结

经验模态分解（EMD）是一种强大的信号处理方法，特别适合分析非线性和非平稳信号。它通过迭代去趋势的方式，将复杂信号分解为若干个本征模态函数（IMF）和一个残差，使得我们可以更好地理解和分析信号的特性。

示例程序：

先贴代码：

import matplotlib.pyplot as plt
import numpy as np
import rasterio
from PyEMD import EEMD

file_path = 'your file path '
with rasterio.open(file_path) as src:
    image = src.read(1) 
    transform = src.transform
fig, ax = plt.subplots()
ax.imshow(image, cmap='gray')
plt.title("Click two points to determine a line, then press Enter")
points = []
def onclick(event):
    if len(points) < 2:
        points.append((event.xdata, event.ydata))
        ax.plot(event.xdata, event.ydata, 'ro')
        fig.canvas.draw()
cid = fig.canvas.mpl_connect('button_press_event', onclick)

def onenter(event):
    if event.key == 'enter':
        fig.canvas.mpl_disconnect(cid)
        plt.close(fig)
fig.canvas.mpl_connect('key_press_event', onenter)
plt.show()

if len(points) < 2:
    print("Two points were not selected.")
else:
    (x1, y1), (x2, y2) = points
    num_points = int(np.hypot(x2 - x1, y2 - y1))
    x, y = np.linspace(x1, x2, num_points), np.linspace(y1, y2, num_points)
    x, y = x.astype(int), y.astype(int)
    gray_values = image[y, x]

    # 可视化灰度剖面
    fig, ax = plt.subplots()
    ax.plot(np.arange(len(gray_values)), gray_values, label="value")
    ax.set_title("Gray Level")
    ax.set_xlabel("Point Index")
    ax.set_ylabel("Gray Level")
    ax.legend()
    plt.show()

    # EMD分解
    emd = EEMD()
    IMFs = emd.eemd(gray_values)
    num_IMFs = IMFs.shape[0]

    # 求IMF分量
    fig, axes = plt.subplots(num_IMFs, 1, figsize=(10, 2*num_IMFs))
    for i in range(num_IMFs):
        axes[i].plot(IMFs[i], label=f"IMF {i+1}")
        axes[i].legend()
    plt.tight_layout()
    plt.show()

    # 求方差贡献率
    variances = np.var(IMFs, axis=1)
    total_variance = np.sum(variances)
    variance_contributions = variances / total_variance * 100

    # 可视化方差贡献率
    fig, ax = plt.subplots()
    ax.bar(np.arange(1, num_IMFs+1), variance_contributions)
    ax.set_title("Variance Contribution Rate of Each IMF")
    ax.set_xlabel("IMF")
    ax.set_ylabel("Variance Contribution Rate (%)")
    plt.show()

使用方法：
0. 填自己的图片路径，然后把没装的包装了。

以我的图片为例，选两个目标点构成要截取的线，然后回车确认：
得到该线段上的灰度剖面：
关掉这个窗口（可以自己写别的触发），会直接弹出分解结果，以及模式贡献量：

注：模式分解是一种数据驱动的方法，没有预设基函数，不同信号特性会导致不同的分解结果，即使同一个信号，在不同噪声条件下分解结果也不同。示例程序暂时没有对图像进行严格的降噪和去趋势等处理，用户可以结合自己的方法，用小波过一遍或者滤波一下。或者做一些合并IMF之类的操作，以便得到更清晰的信号特征。