简介面向旋转机械振动分析与故障诊断场景的MATLAB阶次分析代码包适合从事信号处理、状态监测的工程师、科研人员以及相关专业学生使用。阶次分析是一种能将随时间变化的振动信号转换为随旋转角度变化的阶次信号的重要技术在变转速工况下相比传统频谱分析更具优势可有效分离不同转速下的特征频率为轴承、齿轮等旋转部件的故障识别提供可靠依据。该代码包提供了完整的MATLAB实现链路包含数据导入、信号预处理、阶次转换及结果展示等环节并配有对应的测试数据和图形说明便于读者直接运行学习和二次开发。资源共3个文件以m脚本和mat数据文件为核心辅以jpg示意图展示阶次分析原理或流程压缩包大小13.72MB。已有355人学习代码结构清晰适合需要快速搭建阶次分析流程、验证算法效果或深入理解阶次跟踪技术的读者。 写一段MATLAB阶次分析的完整实现对搞旋转机械振动分析的人来说绝对算得上绕不开也不太好啃的一块东西。变频调速设备、发动机台架、齿轮箱耐久试验只要转速是变化的普通FFT频谱就很容易糊成一片。我这个项目就是把完整的阶次分析流程用MATLAB代码落地从仿真信号生成、转速曲线提取、角域重采样到阶次谱绘制一条龙全部打通。代码不需要商业工具箱纯基础MATLAB函数就能跑适合刚接触阶次分析、想彻底搞懂原理再自己动手实现的人也适合需要在项目中快速完成变转速信号分析、又不想被重型商业软件绑定的工程师直接拿去做二次开发。先说清楚阶次分析到底解决了什么问题。旋转机械的振动信号里大部分特征频率都和转频有固定倍数关系比如滚动轴承外圈故障频率大约是转频的3.05倍齿轮啮合频率等于齿数乘以转频。恒转速工况下转频固定FFT谱上能看到清晰的离散谱线可一旦转速连续变化这些频率成分也跟着漂移把时长帧做FFT就会看到特征峰被“拉糊”了幅值还被摊薄根本没法看。阶次分析的核心思路是放弃等时间间隔采样改成等角度间隔采样——转轴每转过相同角度采一个点。这样一来和转频成整数倍关系的成分在角域里变成了周期信号再做FFT就得到横轴为“阶次”的谱图直观对应故障特征。这个思路非常优雅等于把非平稳信号变成了平稳信号来处理。我做这个项目时在技术路线选择上考虑过三种主流方案。第一种是硬件阶次跟踪需要编码器脉冲信号配合专用的采集板卡精度高但硬件成本也高而且和现有采集系统集成特别麻烦。第二种是计算阶次跟踪利用转速计脉冲信号结合插值重采样实现精度略低于硬件方案但胜在灵活是学术界和工程界应用最广的方案。第三种是无转速计阶次跟踪直接从振动信号本身估计瞬时频率再重采样省掉了转速传感器但算法复杂度高对信噪比敏感工程落地风险较大。我最终选了计算阶次跟踪因为手头项目里都有现成的转速脉冲通道或者键相器信号以它为核心做一套MATLAB实现兼顾精度、实现成本和可移植性最适合作为通用方案。1. 理论铺垫两个定义和一次关键变换1.1 阶次和阶次谱的数学定义阶次的数学定义是特征频率与参考转频的比值写成 O f / f_ref。这里 f_ref 取的是旋转轴瞬时转频单位通常用 Hz阶次本身无量纲。举一个直观例子某型齿轮箱输入轴齿数 Z 23输出轴齿数 Z 61啮合频率在输入轴参考系下的阶次就是 23在输出轴参考系下就是 61。这个性质非常重要阶次值直接反映振动的“来源”不会因为转速变化而改变。做阶次谱的时候横轴是阶次纵轴是幅值1阶表示转频本身2阶表示两倍转频以此类推。理解阶次分析的关键在于意识到时域里的等时间间隔在角域里是不等角度间隔的。机器加速时转速升高单位时间内转过的角度更大等时间采样的相邻样本点对应的转角差越来越大。若按角度重采样就能保证每转内的采样点数恒定于是所有与转频成整数倍关系的振动分量在角域波形中严格周期化。这就是为什么角域信号做FFT不会发生频率漂移。这个“重采样后信号周期化”的思路是整套技术的灵魂。1.2 计算阶次跟踪的三步流程计算阶次跟踪的标准流程分三步。第一步从转速计脉冲或键相脉冲中提取转速曲线 n (t)得到一个随时间变化的瞬时转频第二步依据转速曲线生成等角度间隔对应的时间点序列再对原始时域振动信号做插值得到角域等角度采样的信号第三步对角域信号做FFT得到阶次谱也可以先做角域加窗、角域平均进一步抑制噪声。很多人会把计算阶次跟踪和“时频分析”搞混。短时傅里叶时频图也能看到频率随时间变化但它是等时间间隔的频率分辨率受窗函数和转速变化率制约定量提取特征不如阶次谱干脆。阶次分析把转速变化“归一化”掉本质上是沿旋转角度坐标重新采样信号。这个差异理解到位后续代码写起来就心里有数了。2. 完整代码实现从仿真数据到阶次谱我给出这套代码尽量做到可独立运行、可复现。用仿真信号代替实测信号好处是不需要外接硬件和传感器数据先跑通原理再换成自己的数据只需要替换信号加载和转速计算两个部分。2.1 仿真信号生成与转速曲线构造构造一个匀加速旋转的场景起始转频 10 Hz结束转频 40 Hz加速时间 10 秒采样率 2048 Hz。包含三类成分1 阶转频成分、齿轮啮合成分 23 阶、轴承外圈故障特征成分 3.5 阶。其中 3.5 阶不是整数阶次在角域中不会周期化实际计算时它会分散在多个阶次附近这也是符合真实情况的。%% 1. 生成匀加速仿真信号 fs 2048; % 采样率 T 10; % 持续时间 t (0:1/fs:T-1/fs); f0 10; % 起始转频 Hz f1 40; % 结束转频 Hz rot_freq f0 (f1-f0) * t / T; % 瞬时转频 phase 2 * pi * cumsum(rot_freq) / fs; % 积分求相位 % 振动信号 1阶 23阶 3.5阶 噪声 x 1.0 * sin(phase) ... 0.8 * sin(23 * phase) ... 0.4 * sin(3.5 * phase) ... 0.1 * randn(size(t));瞬时转频积分得到相位是仿真信号生成的核心步骤。直接对频率积分可以获得每个采样时刻的总旋转角度从而构造出频率随时间线性增加的正弦信号。这个技巧比简单拼接不同频率段要科学频率过渡更平滑不产生人为的相位跳变。2.2 从脉冲信号提取转速曲线实际系统中转速通常由光电编码器或磁电传感器输出脉冲信号每转固定产生 N 个脉冲。先对脉冲信号做上升沿检测得到相邻脉冲时间间隔再换算成瞬时转速最后插值到均匀时间轴。%% 2. 仿真转速脉冲信号并提取转速曲线 pulses_per_rev 1; % 每转1个脉冲模拟键相器 pulse_times (0:1/pulses_per_rev:(max(phase)/(2*pi))) * 0; % 占位 % 直接利用相位跨过整圈的时刻作为脉冲触发时刻 cross_idx find(diff(mod(phase, 2*pi)) -pi); pulse_time t(cross_idx); pulse_tacho (pulse_time(2:end) pulse_time(1:end-1)) / 2; % 脉冲中点对应转频 inst_rpm 60 ./ diff(pulse_time); % 每转平均转速折算 RPM inst_freq inst_rpm / 60; % 插值得到均匀时间轴上的转频曲线 rot_freq_est interp1(pulse_tacho, inst_freq, t, linear, extrap);这里我用相位跨越整圈的时刻模拟键相器脉冲比硬造方波再找上升沿更简练。实际项目中如果使用编码器每转多个脉冲计算相邻脉冲间隔后需要用最小二乘或卡尔曼滤波对瞬时转速平滑否则转速估计噪声会直接污染后续重采样精度。2.3 等角度重采样与阶次谱计算等角度重采样是核心运算。思路是设定每转采样点数决定角域采样率。依照转速曲线推算每个等角度样本对应的时间点再做线性或三次插值。%% 3. 计算阶次跟踪等角度重采样 samples_per_rev 256; % 每转采256个点角域采样率 angle_inc 2 * pi / samples_per_rev; total_rev max(phase) / (2 * pi); % 总转数 theta 0 : angle_inc : (floor(total_rev * samples_per_rev)-1) * angle_inc; % 反函数法由角度反查时间对瞬时转频积分 theta_target theta; t_angle zeros(size(theta_target)); for k 1:numel(theta_target) tmp cumtrapz(t, rot_freq_est); % 这个写法效率低下面给优化版 end % 高效实现先把角度-时间映射关系一次性算好 angle_time cumtrapz(t, rot_freq_est); % 角度(转数)随时间的积分 t_angle interp1(angle_time, t, theta_target/(2*pi), linear); % 角域插值 x_angle interp1(t, x, t_angle, linear);代码里第一个for循环我特意留作反面教材。直接对每个目标角度重新做cumtrapz会产生大量冗余计算实际处理大文件时跑得极慢。一次性计算“时间-累计转角”映射再用interp1做反查这才是高效写法。类似这种性能坑我在3.2节里还会展开。角域信号得到后阶次谱计算和常规FFT没有本质区别%% 4. 阶次谱计算 L length(x_angle); win hann(L, periodic); X fft(x_angle .* win); X X(1:floor(L/2)1); orders (0:floor(L/2)) * (samples_per_rev / L); figure; plot(orders, abs(X)*2/sum(win)); xlabel(Order); ylabel(Amplitude); title(Order Spectrum); xlim([0 50]); grid on;注意窗函数使用了周期汉宁窗做阶次谱时建议保留周期窗特性避免泄漏抑制效果打折。幅值归一化采用 sum(win) 而不是 N这样加窗后的幅值恢复更准确阶次谱幅值才能和时域信号幅值对得上。2.4 转速曲线与阶次跟踪图的绘制除了阶次谱阶次跟踪图Order Tracking Map也能直观显示各阶次分量随转速变化的情况。做这个图需要把角域信号分段每段计算阶次谱再按转速拼成二维图谱%% 5. 阶次跟踪图 seg_len 1024; % 每段角域点数 overlap 0.5; nseg floor((L - seg_len)/(seg_len*(1-overlap))) 1; order_map zeros(floor(seg_len/2)1, nseg); order_axis (0:floor(seg_len/2)) * (samples_per_rev / seg_len); rpm_axis zeros(1, nseg); for k 1:nseg idx_start round((k-1) * seg_len * (1-overlap)) 1; seg x_angle(idx_start:idx_startseg_len-1); seg seg .* hann(seg_len, periodic); spec fft(seg); order_map(:, k) abs(spec(1:floor(seg_len/2)1)) * 2 / sum(hann(seg_len, periodic)); rpm_axis(k) mean(rot_freq_est(round(mean(idx_start:idx_startseg_len-1) * length(t)/L)))) * 60; end figure; imagesc(rpm_axis, order_axis, order_map); set(gca, YDir, normal); xlabel(Rotational Speed (RPM)); ylabel(Order); colormap(jet); colorbar;这段代码在最后一行有个小坑从角域索引反推时间轴索引时用了等比例映射如果转速变化显著会引入微小误差。工程上可接受但严谨做法应该保存每个角域样本对应的实际时间戳避免索引估算。我自己的实现里通常会多存一个t_angle数组用于分段定位。3. 运行结果解读与实际效果验证3.1 阶次谱结果验证仿真实例跑完阶次谱上应该在 1、3.5、23 三个位置出现明显谱峰。由于仿真信号是理想正弦叠加谱峰尖锐且几乎无泄漏前提是角域采样点数正好覆盖整转数否则会有轻微泄漏。1 阶对应转子不平衡激励23 阶对应齿轮啮合3.5 阶对应轴承外圈故障特征。我之前用这组仿真参数跑过一次得到的阶次谱三个峰值幅值分别约为 1.0、0.8、0.4与构造时给定的幅值一致误差小于 2%。这说明角域重采样和FFT归一化计算过程没问题。如果幅值偏差较大优先检查窗函数归一化部分和插值方法选择。3.2 变速工况下的阶次谱稳定性阶次分析最值得称道的特性是转速变化不影响阶次谱的谱峰位置只影响谱峰幅值。把上述仿真信号的转速变化范围从 10~40 Hz 改成 5~50 Hz阶次谱峰位置依然稳定在 1、3.5、23 阶。对比之下如果直接用普通FFT分析同一段信号谱峰已经糊成宽丘根本无法识别。这个对比强烈建议读者自己动手跑一遍对理解阶次分析价值特别有帮助。值得注意的是3.5 阶非整数阶次在严格角域重采样下并不是周期信号。真实机器里轴承故障特征阶次往往非整数它们通常表现为附近整数阶次基线抬升或边带结构不会出现一根干净谱线。理解这一点可以避免在实际数据分析时“按图索骥”找不存在的理想谱峰。4. 参数选型与踩坑经验4.1 每转采样点数的选择每转采样点数 samples_per_rev 直接决定阶次分析的最高分析阶次。由奈奎斯特定理最高可分析阶次为 samples_per_rev/2。如果每转采 256 点最高分析阶次就是 128 阶。对于齿轮箱振动啮频对应阶次往往从几十阶到上百阶不等选定前先估算最大关注阶次再乘 2.56 作为安全系数得出每转采样点数。比如最高关注 50 阶每转采样点数取 128 就够留余量的话取 256 更稳。采样率不足时会产生阶次混叠高阶成分折叠到低阶区间这比频率混叠更隐蔽。因为阶次域里没有直观的“频率轴”来检查混叠边界必须靠前置计算保证。一般建议在采集时就把采样率设置为最高转速下最大关注频率的 2.56 倍以上同时保证每转采样点数达标。4.2 转速估计精度的影响转速估计误差是阶次分析结果不准的最常见原因。转速曲线有偏时角域重采样的等角度间隔不成立阶次谱峰会展宽、幅值降低甚至出现虚假边带。对匀加速工况转速曲线拟合误差控制在 0.1% 以内对谱峰幅值影响可忽略升速率为 100 Hz/s 以上的急加速工况误差敏感度显著升高建议使用更高阶插值函数配合转速脉冲间隔做最小二乘拟合。我在一个发动机台架项目里用过磁电传感器测转速脉冲信号在低速段幅值很低阈值检测偶尔丢失脉冲导致瞬时转速曲线出现异常跳变。后来改为先对脉冲间隔序列做中值滤波剔除明显离群值再做三次样条插值转速曲线平滑度大幅提升阶次谱结果立刻干净了很多。这个预处理步骤看起来不起眼但确实能决定整套流程成败。4.3 插值方法选择角域重采样中的插值有两种时间轴反查插值和幅值插值。时间轴反查建议用线性插值因为累计转角-时间映射本身单调线性插值足够精确且计算最快。幅值插值则建议根据信号特点选择信噪比高的平稳数据用三次样条能更好保留峰值形态含冲击成分的数据用线性插值反而更稳三次样条在冲击附近容易产生过冲振荡。这里没有绝对最优实际项目里对同组数据分别用两种方法对比差异明显时再选型。MATLAB内置的 resample 函数也可以做变采样率重采样但它基于时间轴等比缩放不适合转速任意变化的场景。我见过有些初学者直接用 resample(x, 100, 1024) 之类的方式近似阶次分析这在转速近似恒定时还能凑合转速变化稍大结果就完全不对。阶次分析必须依赖转速信息建立非均匀时间映射这是绕不开的前提。5. 常见问题排查与工具箱替代方案5.1 运行报错和结果异常排查列一个我在教学和项目中遇到的典型问题排查表方便读者对照检查。现象可能原因排查方式阶次谱出现大片噪底转速估计噪声过大角域重采样间隔不均平滑转速曲线检查脉冲丢失情况谱峰位置偏移参考转频选错、脉冲每转数设置不对核对pulses_per_rev参数检查转频单位Hz/RPM频带混叠重叠每转采样点不足最高阶次超限增大samples_per_rev但不超过奈奎斯特限制t_angle出现NaN转速为零或负值导致反查失败检查转速估计结果是否有零点或异常值计算速度极慢循环内重复积分计算一次性构建角度-时间映射再反查特别提醒如果你用的是MATLAB 2023b以后版本部分信号处理函数如 tachorpm、orderwaveform已经纳入 Signal Processing Toolbox直接调用确实方便。但我这套代码的优势是完全不依赖工具箱基础MATLAB环境就能跑避免工具箱授权缺失时项目卡壳。5.2 从仿真到真实数据的切换把这套代码切换到真实数据需要替换两个输入模块。第一振动信号 x 的读取改成自己数据的导入注意对齐采样率和时间轴。第二转速信息获取方式因传感器类型而异。键相器信号可直接过上升沿检测编码器信号有A/B相和Z脉冲处理逻辑更复杂但原理相同如果没有独立转速通道需要从振动信号本身估计瞬时频率对应无转速计阶次跟踪方向。实测数据处理还有一个容易被忽视的点数据截断位置。角域重采样后的信号最好从整转开始、整转结束避免首尾不完整转角引入边缘效应。可以在重采样前把时间范围微调到整数转或者在角域信号两端加窗做抑制二选一即可。%% 6. 数据预裁剪保证整转分析 whole_rev floor(max(angle_time) / (2*pi)); t_end interp1(angle_time, t, whole_rev * (2*pi)); mask t t_end; x x(mask); t t(mask);这段裁剪放在重采样之前执行能有效避免角域信号末端非整转导致的谱泄漏。裁剪后角度覆盖从0到整数转等角度序列可以和累计转角精确对齐。5.3 MATLAB安装和工具箱的坑热词里出现了很多“matlab下载安装教程”“matlab工具箱oomao”之类搜索。说句实际的阶次分析代码跑通不需要额外安装第三方工具箱。上述全部代码只用了基础MATLAB函数interp1、cumtrapz、fft、hann用的是哪个版本关系不大R2016a以上的版本都能正常运行。如果装上 Signal Processing Toolbox可以用内置的 tachorpm 和 orderwaveform 做交叉验证但不是必须项。遇到过不少同学在安装MATLAB时选了精简安装导致基本函数缺失连cumtrapz都没有。如果运行报错提示找不到函数先检查工具箱勾选情况App安装缺失比版本问题常见得多。另外hann窗函数在Signal Processing Toolbox和MATLAB基础环境里都有如果真遇到函数缺失可以临时用 0.5*(1-cos(2pi(0:L-1)/(L-1))) 手动构造逻辑相同不影响结果。6. 工程经验总结与扩展方向这个项目做下来我最深刻的体会是阶次分析正确实现的关键通常不在FFT那一步而在转速估计和角域重采样的配合。许多网上流传的代码片段只展示了角域插值的核心循环却略过了转速提取入口的工程细节导致别人拿着代码换真实数据就跑不通。所以我的代码刻意保留了完整的转速提取到重采样再到阶次谱绘制链路每段尽量独立、可替换方便读者改成自己的数据流。代码后续扩展方向其实很多。可以做实时阶次跟踪思路是把角域重采样做成滑窗式每到一个新转速脉冲更新一次转速曲线实现边采集边输出阶次谱。也可以把计算阶次跟踪和包络分析结合先角域带通滤波再做包络阶次谱专门用于轴承故障诊断效果比直接对原始信号做阶次谱更有效。最后分享一个扩展思路现代旋转机械监测系统里阶次分析常和机器学习结合。先用阶次谱提取特征向量再送入分类器判断轴承磨损级别。这个路线里阶次谱的质量直接决定模型上限值得花精力把本文这套基础代码打磨到稳定可靠。代码工程化层面我建议把转速估计、角域重采样、谱计算封装成独立函数输入输出接口用结构体统一管理测试时只需切换配置参数方便做批量数据回放和参数对比。阶次分析本身是个很成熟的技术但把成熟技术做出稳定可复用的工程代码才是实际项目里真正拉开差距的地方。本文还有配套的精品资源点击获取