简介本资源是G-Rilling开发的经验模态分解EMD与希尔伯特-黄变换HHTMATLAB工具箱专为非线性、非平稳信号分析设计面向物理、工程、生物医学等领域的科研人员及信号处理工程师。工具箱完整实现EMD核心算法如cemdc2、emdc、emd等、边界条件处理、IMF可视化emd_visu、hhspectrum、瞬时谱分析及安装/卸载脚本含40个MATLAB函数.m、17个C语言加速模块.c、11个头文件.h及配套示例NSIP2003/SPL2007和工具脚本共72个文件总大小仅98KB轻量高效。目前已有196人学习下载适合需快速部署HHT流程、理解EMD算法细节、复现经典论文结果或开展时频分析研究的中高级用户。1. 不是“装个包就能用”的EMD为什么Rilling版HHT在信号分析中仍被反复引用你下载了package_emd.zip解压后看到emd.m、hht.m、plot_hht.m这几个文件双击运行却报错 “Undefined function ‘emd’”——这不是你的MATLAB版本太低也不是路径没加对而是 Rilling 版 EMD常被误称为 G-Rilling EMD 或 Rilling HHT根本不是标准 MATLAB 工具箱而是一套依赖特定预处理与迭代终止逻辑的手工实现。它不依赖 Signal Processing Toolbox 的emd那是 MathWorks 2018a 后才内置的、基于镜像延拓和Spline插值的现代版本而是严格复现 2003 年 G. Rilling、P. Flandrin 和 P. Gonçalvès 在 IEEE SP Letters 上提出的原始算法以三次样条插值包络、以局部极值为节点、以 IMF 能量比和标准差阈值联合判停。这种实现至今仍在旋转机械故障诊断、脑电微状态识别、地震波瞬时频率追踪等场景中被引用不是因为“更先进”而是因为可复现性高、参数语义明确、每一步迭代都暴露在用户控制之下。适合需要调试 IMF 分解过程、验证 HHT 谱稳定性、或与论文结果做逐行比对的工程师与研究者——如果你只是想快速画个时频谱MathWorks 官方emdhht更省事但若你要写方法学章节、做算法对比实验、或解释为什么某段信号分解出 7 个 IMF 而不是 5 个Rilling 版就是绕不开的基准。2. 从package_emd.zip到可调用函数解压、路径配置与最小验证流程Rilling 版 EMD 的核心是纯 MATLAB 函数无 MEX 编译、无外部依赖但必须满足三个隐性前提MATLAB 版本 ≥ R2007b因使用nargout和结构体字段动态访问、信号长度 100 点避免极值点不足、采样率需已知HHT 频率轴计算必需。下面以 Windows MATLAB R2021b 为例走通从解压到首条 IMF 输出的完整链路。2.1 解压与目录结构确认下载的package_emd.zip解压后应包含以下关键文件注意大小写emd.m主分解函数输入信号x输出 IMF 矩阵imf每行一个 IMF和残余分量residuehht.m希尔伯特谱计算函数输入imf和采样频率fs输出瞬时频率f、时间t、能量谱ampplot_hht.m可视化函数封装hht计算并绘制时频能量图example.m示例脚本含合成信号生成与流程演示README.txt作者说明G. Rilling, 2005强调该实现与原始论文一致提示不要将整个解压目录拖入 MATLAB 路径addpath而应仅添加该目录本身。若同时存在 MathWorks 官方emdMATLAB 会优先调用路径靠前的版本导致which emd返回错误位置。执行restoredefaultpath后再addpath(D:\package_emd)可规避冲突。2.2 最小可运行命令一条信号一个 IMF用合成信号验证是否真正加载成功% 生成测试信号含 20Hz 正弦 50Hz 调幅 噪声 fs 1000; t (0:1/fs:1-1/fs); x sin(2*pi*20*t) 0.5*sin(2*pi*50*t).*cos(2*pi*5*t) 0.1*randn(size(t)); % 调用 Rilling EMD注意必须传入列向量 [imf, residue] emd(x); % 检查输出维度imf 是 N_imf × N_sample 矩阵residue 是列向量 size(imf) % 例如ans [4, 1000] → 得到 4 个 IMF size(residue) % ans [1, 1000]若报错Error using emd: Not enough extrema说明信号极值点过少如全零或直流分量主导。此时需预处理去均值x x - mean(x);Rilling 版未内置去均值必须手动检查长度length(x) 100否则插值不稳定避免饱和若信号含大直流偏移x x - median(x)比mean更鲁棒2.3 参数表emd.m中可直接修改的 5 个关键阈值Rilling 版 EMD 的迭代停止条件由两个嵌套阈值控制全部硬编码在emd.m第 42–49 行。修改前务必备份原文件参数名默认值物理含义修改建议sd0.2标准差准则阈值SD std( (h1-h2)./h1 )降低至0.1可得更精细 IMF但增加迭代次数提高至0.3加速收敛但可能欠分解sdmax0.5SD 准则最大容忍值防止无限循环一般不动仅当信号含强噪声时设为0.6ndmax100单次筛分最大迭代次数处理长信号10⁵点时可增至200避免中途退出tol1e-3能量比准则阈值E₂/E₁ tol对高频成分敏感信号可降至1e-4maxiter1000全局最大筛分次数所有 IMF 总和学术复现实验建议保持1000工程部署可设500注意这些参数不通过函数输入传递必须编辑emd.m源码。例如将sd 0.2;改为sd 0.15;后保存下次调用即生效。这是 Rilling 版与官方emd的根本区别——后者用Name-Value对如SiftRelativeTolerance,0.15传参前者靠代码级定制。3. HHT 谱计算与hht.m的三重校验从 IMF 到时频能量图得到 IMF 后hht.m负责计算每个 IMF 的瞬时频率与幅值并聚合为时频谱。但 Rilling 版hht.m有三个易被忽略的约束直接决定谱图是否可信IMF 必须满足正交性、采样率必须精确、相位展开必须无跳变。跳过任一校验plot_hht画出的“伪彩色图”可能完全失真。3.1 IMF 正交性验证为什么corr(imf(1,:), imf(2,:))不能接近 0 就要重算Rilling EMD 的理论基础是 IMF 的近似正交性各阶分量能量无泄漏。但实际分解中因插值误差与终止阈值IMF 间可能存在相关性。必须在调用hht前验证% 计算 IMF 两两相关系数绝对值矩阵 C abs(corrcoef(imf)); % 提取上三角排除自相关 C_triu triu(C, 1); % 若最大非对角相关系数 0.1说明分解质量差 if max(C_triu(:)) 0.1 warning(IMF 正交性不足HHT 谱可能出现能量混叠); % 应对增大 sd 阈值如 sd0.25重新分解 end3.2hht.m输入要求与关键参数解析hht.m的调用签名是[f, t, amp] hht(imf, fs, option)其中option仅支持plot触发绘图或空默认。但内部有三个隐式参数影响结果fs必须为标量数值若传入fs 1000.0double正常但fs 1000int在某些 MATLAB 版本会触发类型错误强制转为double(fs)imf行数必须 ≥ 2单 IMF 无法计算 Hilbert 变换的相位导数hht.m第 67 行会报错Index exceeds matrix dimensions时间向量t由linspace(0, (N-1)/fs, N)生成N size(imf,2)因此fs的精度直接影响频率轴分辨率% 正确调用假设 imf 为 4×1000 矩阵fs1000Hz [f, t, amp] hht(imf, 1000); % f: 4×1000, t: 1×1000, amp: 4×1000 % 错误调用示例 % [f,t,amp] hht(imf, 1000.0000001); % 频率轴出现非整数 bin后续插值失真 % [f,t,amp] hht(imf(1:2,:), 1000); % 若只取前2阶 IMF但原始信号含高频成分丢失信息3.3plot_hht.m的底层逻辑与可替换绘图方案plot_hht.m默认用imagesc绘制amp矩阵但存在两个缺陷频率轴非线性f是瞬时频率矩阵每行 IMF 的f(i,:)长度与t相同但不同 IMF 的频率范围差异巨大如 IMF1 可能 1–5HzIMF4 达 200–300Hz直接imagesc(t,f,amp)会导致低频区像素挤压能量归一化缺失amp未按 IMF 能量加权高频 IMF 的微弱能量可能被低频 IMF 的强能量淹没替代方案推荐% 手动构建时频谱对每个 IMF 计算其能量占比加权叠加 N_imf size(imf,1); tf_spectrum zeros(length(t), 512); % 频率 bins 数 for i 1:N_imf % 对第 i 个 IMF 的瞬时频率 f(i,:) 和幅值 amp(i,:) 插值到统一频率轴 f_i f(i,:); amp_i amp(i,:); % 构建直方图每个时间点 t(k) 对应频率 f_i(k)能量 amp_i(k)^2 for k 1:length(t) bin_idx round(f_i(k) / (fs/2) * 512); % 映射到 0–512 bins if bin_idx 1 bin_idx 512 tf_spectrum(k, bin_idx) tf_spectrum(k, bin_idx) amp_i(k)^2; end end end % 归一化并绘图 tf_spectrum tf_spectrum / max(tf_spectrum(:)); imagesc(t, linspace(0, fs/2, 512), tf_spectrum); xlabel(Time (s)); ylabel(Frequency (Hz)); colorbar;此方案显式控制频率分辨率、能量权重与可视化范围比plot_hht更符合工程分析需求。4.uninstall_emd不存在清理 Rilling EMD 环境的三种安全方式网络搜索uninstall_emd时多数结果指向删除整个package_emd文件夹——但这只是物理移除MATLAB 的路径缓存path cache和函数句柄function handle可能残留导致后续调用仍命中旧版本。真正的“卸载”需三层清理4.1 清理 MATLAB 路径缓存即使已用rmpath(D:\package_emd)移除路径which emd仍可能返回旧位置因 MATLAB 缓存了函数位置。必须执行% 1. 清除函数位置缓存 clear functions % 2. 强制刷新路径索引 rehash toolboxcache % 3. 验证是否彻底清除 which emd % 应返回 emd not found而非路径提示rehash toolboxcache比rehash path更彻底它重建整个工具箱索引耗时略长但确保无残留。4.2 检查并重置函数重载Overload冲突若曾将emd.m放在某个类的myclass目录下如signal/emd.mMATLAB 会为signal类型对象自动调用该emd。此时which emd可能显示.../signal/emd.m。解决方法进入该signal目录删除emd.m或执行clear classes重置所有类定义4.3 替换为 MathWorks 官方 EMD 的无缝迁移步骤若决定弃用 Rilling 版迁移到官方emd需注意三点差异功能点Rilling 版MathWorks 官方版迁移操作输入信号必须列向量支持行/列向量自动转置删除x x(:)强制转换IMF 输出imf矩阵每行一阶imf结构体含imf,residual,numIMF字段用imf.imf取矩阵imf.residual取残余终止准则sd/tol双阈值SiftRelativeTolerance单阈值将sd0.2映射为SiftRelativeTolerance,0.2HHT 计算hht.m独立函数hht函数直接接受imf结构体hht(imf)可直接调用无需预处理迁移后验证% 官方版最小流程 [imf_struct, ~] emd(x, SiftRelativeTolerance, 0.2); hht_result hht(imf_struct); % 对比 Rilling 版输出维度官方版 imf_struct.imf 与 Rilling 版 imf 应接近5. 用package_emd做轴承故障诊断从振动信号到冲击特征频率的实操技巧Rilling 版 EMD 的不可替代性在于其对瞬态冲击响应的分解保真度——这正是滚动轴承早期故障诊断的核心。当轴承外圈出现微米级裂纹振动信号中会叠加周期性冲击其包络谱峰值对应故障特征频率如 BPFO。而 Rilling EMD 能将冲击成分精准分离至特定 IMF避免官方版因默认平滑插值导致冲击边缘模糊。5.1 冲击信号预处理为什么必须用detrend而非highpass轴承振动信号常含强趋势项如传感器漂移但highpass滤波会引入相位失真破坏冲击时刻定位。Rilling 版要求线性/多项式去趋势% 正确用 detrend 消除趋势保留冲击瞬态 x_detrend detrend(x, linear); % 或 quadratic % 错误highpass 滤波即使 fc10Hz会平滑冲击上升沿 % x_hp highpass(x, 10, fs);5.2 IMF 选择策略如何定位含冲击的 IMF并非所有 IMF 都含故障信息。经验法则是冲击成分通常集中在 IMF2–IMF4假设信号采样率 20kHz故障频率 100–500Hz。验证方法% 计算每个 IMF 的峭度Kurtosis冲击成分峭度显著高于其他 IMF kurtosis_vec zeros(size(imf,1),1); for i 1:size(imf,1) kurtosis_vec(i) kurtosis(imf(i,:)); end [~, idx_impulse] max(kurtosis_vec); % 峭度最大者最可能含冲击 % 提取该 IMF 的包络谱用 hilbert abs imf_impulse imf(idx_impulse,:); analytic hilbert(imf_impulse); envelope abs(analytic); % FFT 包络谱找峰值频率即故障特征频率 f_env (0:length(envelope)-1)*fs/length(envelope); env_fft fft(envelope); plot(f_env(1:end/2), abs(env_fft(1:end/2))); xlabel(Envelope Frequency (Hz));5.3 故障频率验证表BPFO/BPFI 计算与 IMF 匹配轴承故障特征频率公式以深沟球轴承为例外圈故障频率 BPFO f_r * Z/2 * (1 - d/D * cos(α))内圈故障频率 BPFI f_r * Z/2 * (1 d/D * cos(α))其中f_r为转频HzZ为滚动体数d为滚动体直径D为节圆直径α为接触角。轴承参数示例值计算 BPFOIMF 中对应频率范围f_r 30 Hz1800 RPMZ 12,d/D 0.2,α 0°BPFO ≈ 162 Hz查f(idx_impulse,:)的均值频率是否在150–170 Hzf_r 30 HzZ 12,d/D 0.2,α 15°BPFI ≈ 198 Hz若idx_impulse的f均值为195 Hz则内圈故障概率高关键技巧不要只看单个 IMF 的瞬时频率均值而要看其频率分布直方图。健康轴承的 IMF 频率呈窄带分布故障轴承的对应 IMF 会出现双峰基频谐波或宽峰冲击调制这是比均值更可靠的判据。本文还有配套的精品资源点击获取
简介本资源是G-Rilling开发的经验模态分解EMD与希尔伯特-黄变换HHTMATLAB工具箱专为非线性、非平稳信号分析设计面向物理、工程、生物医学等领域的科研人员及信号处理工程师。工具箱完整实现EMD核心算法如cemdc2、emdc、emd等、边界条件处理、IMF可视化emd_visu、hhspectrum、瞬时谱分析及安装/卸载脚本含40个MATLAB函数.m、17个C语言加速模块.c、11个头文件.h及配套示例NSIP2003/SPL2007和工具脚本共72个文件总大小仅98KB轻量高效。目前已有196人学习下载适合需快速部署HHT流程、理解EMD算法细节、复现经典论文结果或开展时频分析研究的中高级用户。1. 不是“装个包就能用”的EMD为什么Rilling版HHT在信号分析中仍被反复引用你下载了package_emd.zip解压后看到emd.m、hht.m、plot_hht.m这几个文件双击运行却报错 “Undefined function ‘emd’”——这不是你的MATLAB版本太低也不是路径没加对而是 Rilling 版 EMD常被误称为 G-Rilling EMD 或 Rilling HHT根本不是标准 MATLAB 工具箱而是一套依赖特定预处理与迭代终止逻辑的手工实现。它不依赖 Signal Processing Toolbox 的emd那是 MathWorks 2018a 后才内置的、基于镜像延拓和Spline插值的现代版本而是严格复现 2003 年 G. Rilling、P. Flandrin 和 P. Gonçalvès 在 IEEE SP Letters 上提出的原始算法以三次样条插值包络、以局部极值为节点、以 IMF 能量比和标准差阈值联合判停。这种实现至今仍在旋转机械故障诊断、脑电微状态识别、地震波瞬时频率追踪等场景中被引用不是因为“更先进”而是因为可复现性高、参数语义明确、每一步迭代都暴露在用户控制之下。适合需要调试 IMF 分解过程、验证 HHT 谱稳定性、或与论文结果做逐行比对的工程师与研究者——如果你只是想快速画个时频谱MathWorks 官方emdhht更省事但若你要写方法学章节、做算法对比实验、或解释为什么某段信号分解出 7 个 IMF 而不是 5 个Rilling 版就是绕不开的基准。2. 从package_emd.zip到可调用函数解压、路径配置与最小验证流程Rilling 版 EMD 的核心是纯 MATLAB 函数无 MEX 编译、无外部依赖但必须满足三个隐性前提MATLAB 版本 ≥ R2007b因使用nargout和结构体字段动态访问、信号长度 100 点避免极值点不足、采样率需已知HHT 频率轴计算必需。下面以 Windows MATLAB R2021b 为例走通从解压到首条 IMF 输出的完整链路。2.1 解压与目录结构确认下载的package_emd.zip解压后应包含以下关键文件注意大小写emd.m主分解函数输入信号x输出 IMF 矩阵imf每行一个 IMF和残余分量residuehht.m希尔伯特谱计算函数输入imf和采样频率fs输出瞬时频率f、时间t、能量谱ampplot_hht.m可视化函数封装hht计算并绘制时频能量图example.m示例脚本含合成信号生成与流程演示README.txt作者说明G. Rilling, 2005强调该实现与原始论文一致提示不要将整个解压目录拖入 MATLAB 路径addpath而应仅添加该目录本身。若同时存在 MathWorks 官方emdMATLAB 会优先调用路径靠前的版本导致which emd返回错误位置。执行restoredefaultpath后再addpath(D:\package_emd)可规避冲突。2.2 最小可运行命令一条信号一个 IMF用合成信号验证是否真正加载成功% 生成测试信号含 20Hz 正弦 50Hz 调幅 噪声 fs 1000; t (0:1/fs:1-1/fs); x sin(2*pi*20*t) 0.5*sin(2*pi*50*t).*cos(2*pi*5*t) 0.1*randn(size(t)); % 调用 Rilling EMD注意必须传入列向量 [imf, residue] emd(x); % 检查输出维度imf 是 N_imf × N_sample 矩阵residue 是列向量 size(imf) % 例如ans [4, 1000] → 得到 4 个 IMF size(residue) % ans [1, 1000]若报错Error using emd: Not enough extrema说明信号极值点过少如全零或直流分量主导。此时需预处理去均值x x - mean(x);Rilling 版未内置去均值必须手动检查长度length(x) 100否则插值不稳定避免饱和若信号含大直流偏移x x - median(x)比mean更鲁棒2.3 参数表emd.m中可直接修改的 5 个关键阈值Rilling 版 EMD 的迭代停止条件由两个嵌套阈值控制全部硬编码在emd.m第 42–49 行。修改前务必备份原文件参数名默认值物理含义修改建议sd0.2标准差准则阈值SD std( (h1-h2)./h1 )降低至0.1可得更精细 IMF但增加迭代次数提高至0.3加速收敛但可能欠分解sdmax0.5SD 准则最大容忍值防止无限循环一般不动仅当信号含强噪声时设为0.6ndmax100单次筛分最大迭代次数处理长信号10⁵点时可增至200避免中途退出tol1e-3能量比准则阈值E₂/E₁ tol对高频成分敏感信号可降至1e-4maxiter1000全局最大筛分次数所有 IMF 总和学术复现实验建议保持1000工程部署可设500注意这些参数不通过函数输入传递必须编辑emd.m源码。例如将sd 0.2;改为sd 0.15;后保存下次调用即生效。这是 Rilling 版与官方emd的根本区别——后者用Name-Value对如SiftRelativeTolerance,0.15传参前者靠代码级定制。3. HHT 谱计算与hht.m的三重校验从 IMF 到时频能量图得到 IMF 后hht.m负责计算每个 IMF 的瞬时频率与幅值并聚合为时频谱。但 Rilling 版hht.m有三个易被忽略的约束直接决定谱图是否可信IMF 必须满足正交性、采样率必须精确、相位展开必须无跳变。跳过任一校验plot_hht画出的“伪彩色图”可能完全失真。3.1 IMF 正交性验证为什么corr(imf(1,:), imf(2,:))不能接近 0 就要重算Rilling EMD 的理论基础是 IMF 的近似正交性各阶分量能量无泄漏。但实际分解中因插值误差与终止阈值IMF 间可能存在相关性。必须在调用hht前验证% 计算 IMF 两两相关系数绝对值矩阵 C abs(corrcoef(imf)); % 提取上三角排除自相关 C_triu triu(C, 1); % 若最大非对角相关系数 0.1说明分解质量差 if max(C_triu(:)) 0.1 warning(IMF 正交性不足HHT 谱可能出现能量混叠); % 应对增大 sd 阈值如 sd0.25重新分解 end3.2hht.m输入要求与关键参数解析hht.m的调用签名是[f, t, amp] hht(imf, fs, option)其中option仅支持plot触发绘图或空默认。但内部有三个隐式参数影响结果fs必须为标量数值若传入fs 1000.0double正常但fs 1000int在某些 MATLAB 版本会触发类型错误强制转为double(fs)imf行数必须 ≥ 2单 IMF 无法计算 Hilbert 变换的相位导数hht.m第 67 行会报错Index exceeds matrix dimensions时间向量t由linspace(0, (N-1)/fs, N)生成N size(imf,2)因此fs的精度直接影响频率轴分辨率% 正确调用假设 imf 为 4×1000 矩阵fs1000Hz [f, t, amp] hht(imf, 1000); % f: 4×1000, t: 1×1000, amp: 4×1000 % 错误调用示例 % [f,t,amp] hht(imf, 1000.0000001); % 频率轴出现非整数 bin后续插值失真 % [f,t,amp] hht(imf(1:2,:), 1000); % 若只取前2阶 IMF但原始信号含高频成分丢失信息3.3plot_hht.m的底层逻辑与可替换绘图方案plot_hht.m默认用imagesc绘制amp矩阵但存在两个缺陷频率轴非线性f是瞬时频率矩阵每行 IMF 的f(i,:)长度与t相同但不同 IMF 的频率范围差异巨大如 IMF1 可能 1–5HzIMF4 达 200–300Hz直接imagesc(t,f,amp)会导致低频区像素挤压能量归一化缺失amp未按 IMF 能量加权高频 IMF 的微弱能量可能被低频 IMF 的强能量淹没替代方案推荐% 手动构建时频谱对每个 IMF 计算其能量占比加权叠加 N_imf size(imf,1); tf_spectrum zeros(length(t), 512); % 频率 bins 数 for i 1:N_imf % 对第 i 个 IMF 的瞬时频率 f(i,:) 和幅值 amp(i,:) 插值到统一频率轴 f_i f(i,:); amp_i amp(i,:); % 构建直方图每个时间点 t(k) 对应频率 f_i(k)能量 amp_i(k)^2 for k 1:length(t) bin_idx round(f_i(k) / (fs/2) * 512); % 映射到 0–512 bins if bin_idx 1 bin_idx 512 tf_spectrum(k, bin_idx) tf_spectrum(k, bin_idx) amp_i(k)^2; end end end % 归一化并绘图 tf_spectrum tf_spectrum / max(tf_spectrum(:)); imagesc(t, linspace(0, fs/2, 512), tf_spectrum); xlabel(Time (s)); ylabel(Frequency (Hz)); colorbar;此方案显式控制频率分辨率、能量权重与可视化范围比plot_hht更符合工程分析需求。4.uninstall_emd不存在清理 Rilling EMD 环境的三种安全方式网络搜索uninstall_emd时多数结果指向删除整个package_emd文件夹——但这只是物理移除MATLAB 的路径缓存path cache和函数句柄function handle可能残留导致后续调用仍命中旧版本。真正的“卸载”需三层清理4.1 清理 MATLAB 路径缓存即使已用rmpath(D:\package_emd)移除路径which emd仍可能返回旧位置因 MATLAB 缓存了函数位置。必须执行% 1. 清除函数位置缓存 clear functions % 2. 强制刷新路径索引 rehash toolboxcache % 3. 验证是否彻底清除 which emd % 应返回 emd not found而非路径提示rehash toolboxcache比rehash path更彻底它重建整个工具箱索引耗时略长但确保无残留。4.2 检查并重置函数重载Overload冲突若曾将emd.m放在某个类的myclass目录下如signal/emd.mMATLAB 会为signal类型对象自动调用该emd。此时which emd可能显示.../signal/emd.m。解决方法进入该signal目录删除emd.m或执行clear classes重置所有类定义4.3 替换为 MathWorks 官方 EMD 的无缝迁移步骤若决定弃用 Rilling 版迁移到官方emd需注意三点差异功能点Rilling 版MathWorks 官方版迁移操作输入信号必须列向量支持行/列向量自动转置删除x x(:)强制转换IMF 输出imf矩阵每行一阶imf结构体含imf,residual,numIMF字段用imf.imf取矩阵imf.residual取残余终止准则sd/tol双阈值SiftRelativeTolerance单阈值将sd0.2映射为SiftRelativeTolerance,0.2HHT 计算hht.m独立函数hht函数直接接受imf结构体hht(imf)可直接调用无需预处理迁移后验证% 官方版最小流程 [imf_struct, ~] emd(x, SiftRelativeTolerance, 0.2); hht_result hht(imf_struct); % 对比 Rilling 版输出维度官方版 imf_struct.imf 与 Rilling 版 imf 应接近5. 用package_emd做轴承故障诊断从振动信号到冲击特征频率的实操技巧Rilling 版 EMD 的不可替代性在于其对瞬态冲击响应的分解保真度——这正是滚动轴承早期故障诊断的核心。当轴承外圈出现微米级裂纹振动信号中会叠加周期性冲击其包络谱峰值对应故障特征频率如 BPFO。而 Rilling EMD 能将冲击成分精准分离至特定 IMF避免官方版因默认平滑插值导致冲击边缘模糊。5.1 冲击信号预处理为什么必须用detrend而非highpass轴承振动信号常含强趋势项如传感器漂移但highpass滤波会引入相位失真破坏冲击时刻定位。Rilling 版要求线性/多项式去趋势% 正确用 detrend 消除趋势保留冲击瞬态 x_detrend detrend(x, linear); % 或 quadratic % 错误highpass 滤波即使 fc10Hz会平滑冲击上升沿 % x_hp highpass(x, 10, fs);5.2 IMF 选择策略如何定位含冲击的 IMF并非所有 IMF 都含故障信息。经验法则是冲击成分通常集中在 IMF2–IMF4假设信号采样率 20kHz故障频率 100–500Hz。验证方法% 计算每个 IMF 的峭度Kurtosis冲击成分峭度显著高于其他 IMF kurtosis_vec zeros(size(imf,1),1); for i 1:size(imf,1) kurtosis_vec(i) kurtosis(imf(i,:)); end [~, idx_impulse] max(kurtosis_vec); % 峭度最大者最可能含冲击 % 提取该 IMF 的包络谱用 hilbert abs imf_impulse imf(idx_impulse,:); analytic hilbert(imf_impulse); envelope abs(analytic); % FFT 包络谱找峰值频率即故障特征频率 f_env (0:length(envelope)-1)*fs/length(envelope); env_fft fft(envelope); plot(f_env(1:end/2), abs(env_fft(1:end/2))); xlabel(Envelope Frequency (Hz));5.3 故障频率验证表BPFO/BPFI 计算与 IMF 匹配轴承故障特征频率公式以深沟球轴承为例外圈故障频率 BPFO f_r * Z/2 * (1 - d/D * cos(α))内圈故障频率 BPFI f_r * Z/2 * (1 d/D * cos(α))其中f_r为转频HzZ为滚动体数d为滚动体直径D为节圆直径α为接触角。轴承参数示例值计算 BPFOIMF 中对应频率范围f_r 30 Hz1800 RPMZ 12,d/D 0.2,α 0°BPFO ≈ 162 Hz查f(idx_impulse,:)的均值频率是否在150–170 Hzf_r 30 HzZ 12,d/D 0.2,α 15°BPFI ≈ 198 Hz若idx_impulse的f均值为195 Hz则内圈故障概率高关键技巧不要只看单个 IMF 的瞬时频率均值而要看其频率分布直方图。健康轴承的 IMF 频率呈窄带分布故障轴承的对应 IMF 会出现双峰基频谐波或宽峰冲击调制这是比均值更可靠的判据。本文还有配套的精品资源点击获取
