MATLAB经验模态分解 (EMD) 信号处理算法实现
项目简介
本项目旨在MATLAB环境中构建一套完整且高效的经验模态分解(EMD)算法工具。作为一种自适应的时频数据分析方法,本实现专门用于解决非线性和非平稳信号的处理难题。
项目的核心在于完全自主实现了EMD算法的筛选(Sifting)过程,不依赖MATLAB工具箱中的现成EMD函数。代码涵盖了从信号极值点提取、三次样条包络拟合、均值去除到IMF分量迭代筛出的全过程。此外,针对EMD算法固有的“端点效应”问题,本项目集成了镜像延拓技术以提高边界处的分解精度。
该工具不仅包含核心算法,还内置了合成测试信号生成器和详细的可视化模块,能够直观展示信号在时域和时频域的分解效果。
主要功能特性
- 完整的EMD核心算法:自主实现了基于三次样条插值的筛选流程,支持多阶固有模态函数(IMF)的自动提取。
- 端点效应处理:采用镜像延拓(Mirror Extension)算法处理信号边界,有效抑制了样条插值在两端产生的发散现象。
- 自适应停止准则:采用了基于柯西收敛准则的标准差(SD)阈值判断,结合最大迭代次数限制,确保筛选过程的收敛性和效率。
- 非平稳信号模拟:内置测试信号生成逻辑,能够合成包含低频、高频、间歇性振荡、线性趋势及噪声的复杂混合信号。
- 多维可视化分析:
*
时域分析:自动绘制原始信号、各阶IMF分量及残余项的波形图。
*
频域特征:基于Hilbert变换计算并展示各主要IMF分量的瞬时频率分布。
系统要求
- MATLAB R2016a 或更高版本(代码主要依赖基础数学运算和
spline、hilbert函数,无需额外工具箱支持)。
使用方法
- 直接运行主程序入口函数,程序将自动执行以下流程:
* 清除工作区环境。
* 生成包含不同频率成分和趋势项的合成测试信号。
* 调用EMD核心算法对信号进行分解。
* 输出分解得到的IMF数量信息。
* 弹出两个图形窗口展示分解结果。
实现细节与逻辑说明
本项目完全基于所提供的源代码实现,具体逻辑流程如下:
1. 信号生成与预处理
程序首先构建一个采样频率为1000Hz,时长为1秒的合成信号,该信号由以下部分叠加而成,用于全面测试算法性能:
- 低频分量:10Hz的正弦波。
- 高频分量:50Hz的余弦波。
- 间歇分量:仅在0.3s至0.7s存在的100Hz调幅信号,用于模拟非平稳突变。
- 趋势项:线性增长的斜坡信号。
- 噪声:随机高斯白噪声。
2. EMD 筛选核心逻辑
算法通过双重循环结构实现分解:
* 不断尝试从当前残差信号中分离IMF,直到残差变为单调函数或极值点不足(小于2个)。
* 设有最大IMF数量限制(默认为10),防止过度分解。
* 从当前数据中通过极值点计算上下包络。
* 计算均值包络并从数据中减去。
* 使用标准差(SD)作为停止准则。当连续两次筛选结果的差异(通过SD衡量)小于预设阈值(0.2)或达到最大筛选次数(100次)时,停止筛选,认定当前分量为通过验证的IMF。
3. 可视化输出
程序执行完毕后生成两组图表:
- 分解结果分析图:动态计算行数,自上而下依次绘制原始信号、所有提取出的IMF分量以及最终的残余项(Residue)。不同IMF分量采用循环配色以便区分。
- 时频特征图:选取前4个主要IMF分量,利用Hilbert变换提取解析信号,计算瞬时频率,并以散点图形式展示频率随时间的变化,直观反映各模态的物理意义(如成功分离出的100Hz、50Hz和10Hz成分)。
关键算法与函数分析
以下是对代码内部关键功能模块的技术分析:
极值搜索 (Extrema Detection)
通过计算信号的一阶差分并分析符号变化来精确识别局部极大值和极小值。该算法能够准确捕捉“左增右减”的极大值点和“左减右增”的极小值点,是后续包络拟合的基础。
镜像延拓 (Mirror Extension)
这是解决EMD端点效应的关键步骤。在进行样条插值前,算法会检查信号两端:
- 如果极值点不仅存在于端点,通过对称原则将靠近边界的内部极值点映射到边界外部。
- 这种延拓人为地在信号两端之外增加了“虚拟”极值点,迫使三次样条插值曲线在通过实际端点时保持合理的趋势,避免了边界处的剧烈摆动或飞逸。
三次样条插值 (Cubic Spline Interpolation)
利用扩展后的极大值和极小值坐标序列,分别拟合出贯穿整个时间轴的上包络线和下包络线。算法中计算均值包络为
(上包络 + 下包络) / 2,这是从原始信号中剥离局部趋势的核心手段。
柯西收敛准则 (Cauchy Convergence Test)
在筛选过程中,计算当前迭代结果与上一次结果之间的归一化平方差(Standard Deviation, SD)。
[ SD = sum frac{(h_{prev} - h)^2}{h_{prev}^2} ]
代码设定阈值为0.2,这在保证IMF正交性和计算效率之间取得了平衡。
瞬时频率估计 (Instantaneous Frequency)
代码利用Hilbert变换构造解析信号
z = hilbert(imf),通过计算相位的微分来获取瞬时频率。为了处理相位卷绕问题,使用了
unwrap 函数,并对结果进行了归一化处理以转换为Hz单位。