
简介面向旋转机械振动分析与故障诊断需求这份资源提供了一套基于MATLAB的阶次分析完整代码适合从事设备状态监测的工程师、科研人员以及信号处理初学者快速入门。资源压缩包共包含3个文件包括主程序脚本、示例数据文件和图形说明文档整体大小仅13.72MB结构清晰、便于下载。已有355人学习利用包内数据可直接复现时域信号采集、预处理、角度域阶次转换及频谱分析的全过程省去自行搭建实验环境的步骤。代码针对变转速运行场景设计读者通过运行示例可理解阶次跟踪如何消除转频波动带来的干扰并将滚动轴承或齿轮故障特征频率与无关成分有效区分。对照图形说明还能直观掌握各处理环节的数据流向和算法要点为后续开展更复杂的振动分析试验提供坚实参考。 做旋转机械振动分析的朋友大概率被同一个问题折磨过设备明明在升速或降速振动频谱却糊成一片原本该是尖峰的特征频率全变成了“小土丘”。网上搜“阶次分析”论文讲原理的不少能直接跑的完整代码却很难找。这篇文章就是把你需要的整条代码链路补完——从转速脉冲提取、角域重采样到阶次谱绘制给出可以在MATLAB里直接运行的完整实现。不管是你手里有一组升速过程的数据还是想把这个方法写进自己的故障诊断程序里照着这篇文章都能落得下来。1. 为什么变转速工况下FFT会“失灵”——阶次分析要解决什么问题1.1 频率模糊的根源信号在时间上根本不平稳很多新手遇到变转速数据第一反应是“采样率不够”“分辨率不够”折腾半天还是没改善。其实问题出在更底层的地方FFT的前提假设是信号在时间上是平稳的而变转速下的旋转机械信号天然不满足这个假设。举个例子一台齿轮箱从1000 RPM升到3000 RPM啮合频率假设24齿会从400 Hz跑到1200 Hz。你把这段升速过程做FFT本质上是在用一系列固定频率的正弦波去拟合一个频率持续变化的信号结果自然是每个频率都被“拉宽”成了频带能量分散到一片峰值幅度大幅下降故障特征淹没在噪声里。这个现象在行业内叫“频率模糊”或者“频谱涂抹”不是数据处理技巧能救的。1.2 阶次以转频为标尺的“归一化频率”阶次Order的定义非常直观转轴每转一圈某个振动事件发生的次数。轴每转一圈齿轮啮合24次就是24阶不平衡引起的振动每转一次是1阶。阶次和Hz之间只差一个实时转速f Order × (RPM / 60)换句话说阶次是把频率用当前转速做了归一化。变转速工况下固定阶次成分的频率虽然随时间漂移但它在每一个转动的角度周期里出现的次数是固定的。这就是阶次分析能绕开非平稳问题的核心逻辑。1.3 阶次分析的总体思路三步走阶次分析的具体实现路径并不神秘核心就三步从键相脉冲或转速计获取瞬时转速曲线根据转速曲线构造“等角度间隔”的时间序列把原始振动信号重采样到角度域再做FFT得到阶次谱。代码层面最麻烦的其实是第2步因为涉及时间域和角度域之间的坐标映射。下面会一步步把每一行的含义讲清楚。2. 核心原理从时间域到角度域的坐标变换变转速信号如何“变平稳”2.1 计算阶次跟踪的数学逻辑先看振动信号在数学上长什么样。假设信号由若干阶次分量组成x(t) Σ A_k · sin(2π · k · N(t) φ_k)其中N(t)是累计转数k就是阶次。令θ 2π·N(t)它是转轴转过的总角度。于是x(θ) Σ A_k · sin(k·θ φ_k)在角度域θ里每个阶次分量都是固定频率k的正弦波。这样一来原本时间上非平稳的信号在角度域重新变成了平稳信号可以放心地做FFT。这个变换过程就是“计算阶次跟踪”Computed Order Tracking也是目前工程上最常用的一种阶次分析方法。2.2 一个容易理解的类比你可以想象转盘上画了一个黑点用普通相机以固定时间间隔拍照转速越快黑点在相邻两帧之间跑得越远画面看起来忽快忽慢。如果换一个“每转10度就拍一张”的同步相机无论转速怎么变黑点在每一帧的位置都基本一致。同步相机的快门逻辑就是阶次分析里的“等角度采样”。给振动信号做等角度重采样本质上就是给测量系统装了一个随转速变化的“虚拟相机”。2.3 每转采样点数连接时间域和角度域的桥梁等角度重采样需要先定义一个关键参数每转采样点数SamplesPerRev。它直接决定了角度域的“采样率”也决定了后续阶次谱能分析到多少阶。这里要记住两个公式角度采样间隔Δθ 2π / SamplesPerRev最高可分析阶次Order_max SamplesPerRev / 2如果已知信号的采样率fs和最高转速RPM_max那么SamplesPerRev fs / (RPM_max / 60)。这是个很实用的估算公式比如采样率5120 Hz最高转速3000 RPM每秒最高50转那每转大约能采102.4个点最高分析阶次约51阶。如果想看到更高的阶次要么提高采样率要么降低分析的最高转速范围。3. 完整代码实现从键相脉冲到阶次谱的一站式方案3.1 代码结构总览整个阶次分析流程我拆成了四个文件来组织方便你直接套用到自己的数据上文件功能OrderAnalysis_Demo.m主脚本模拟一组升速数据并跑通全流程extractRPM.m从键相脉冲时刻计算瞬时转速并拟合转速曲线angularResample.m核心时间域到角度域的重采样orderSpectrum.m对角度域信号做FFT并绘制阶次谱实际工程中你只需要把主脚本里的模拟信号换成从数据采集仪导出的真实信号即可。为了先验证算法正确性我用一组已知阶次成分的模拟数据做自检这个习惯建议保留。3.2 模拟信号的生成首先保证“标准答案”在手第一段代码生成模拟数据转频从10 Hz线性升到30 Hz即600 RPM到1800 RPM振动信号包含1阶和2.5阶分量再加一点噪声。这组数据的好处是“标准答案”已知跑完阶次谱后能在1阶和2.5阶的位置看到清晰的峰值。%% 参数设置 fs 5120; % 采样率 5120 Hz T 10; % 时长 10 秒 t (0:fs*T-1) / fs; % 时间序列 % 线性升速转频从 10 Hz 到 30 Hz rps 10 (30-10) * t / T; % 累计转数对转频积分 revs 10*t (30-10) * t.^2 / (2*T); % 累计角度弧度 theta_t 2 * pi * revs; % 键相脉冲时刻每当累计角度跨过 2*pi 的整数倍 target_theta 2*pi * (1:floor(max(revs))); pulse_t interp1(theta_t, t, target_theta, linear); % 振动信号1阶 2.5阶 噪声 vib 1.0 * sin(2*pi*1.0*revs) ... 0.6 * sin(2*pi*2.5*revs) ... 0.05 * randn(size(t));这里有个细节值得说一下键相脉冲的时刻是通过interp1(theta_t, t, target_theta)反插出来的。因为角度是单调递增的这比写成循环去查找过零点的方式快得多代码也更简洁。3.3 键相脉冲处理与转速曲线提取决定后续精度的地基从键相脉冲提取转速不难难的是把转速曲线弄得足够平滑。原始脉冲间隔的倒数就是瞬时转速但直接微分噪声很大必须做拟合或平滑。% 每个脉冲间隔对应转轴转过一圈 dT diff(pulse_t); rps_inst 1 ./ dT; % 转速点的时间戳取脉冲间隔中点 rps_t_mid (pulse_t(1:end-1) pulse_t(2:end)) / 2; % 多项式拟合阶数取 5 通常够用 p polyfit(rps_t_mid, rps_inst, 5); rps_smooth polyval(p, t); % 再计算对应的累计角度后面重采样要用 theta_fit 2 * pi * cumtrapz(t, rps_smooth);用polyfit做全局拟合的好处是能把脉冲抖动平均掉得到一条光滑的转速曲线。需要注意的是多项式阶数过高会过拟合脉冲噪声阶数过低则拟合不了剧烈变速5到7阶是一个经验区间。如果你手里的数据变速太剧烈改成spline插值或者分段拟合会更好。3.4 角域重采样整篇文章最核心的几行代码这里是阶次分析的“心脏”。思路是目标角度均匀划分然后用插值反查出每个目标角度对应的时刻再在那个时刻对振动信号取值。% 每转采样点数决定了阶次分辨率和最高分析阶次 SamplesPerRev 512; dtheta 2 * pi / SamplesPerRev; theta_axis 0 : dtheta : max(theta_fit); % 反插已知 theta_fit 对应 t求 theta_axis 对应的时刻 t_resample interp1(theta_fit, t, theta_axis, pchip); % 对振动信号做插值得到角度域信号 vib_order interp1(t, vib, t_resample, pchip); % 末尾可能出现 NaNtheta_axis 超出 theta_fit 范围直接截掉 valid_idx ~isnan(t_resample); vib_order vib_order(valid_idx);这里两个interp1都要解释一下。第一个interp1(theta_fit, t, theta_axis, pchip)是把“角度-时间”关系反转过来查询属于反插第二个interp1(t, vib, t_resample, pchip)是在原始时间波形上取值。插值方法我推荐pchip而不是spline两者都很光滑但spline在转速剧烈变化时会产生过冲导致重采样后的波形在局部出现虚假的“毛刺振荡”而pchip不会超过原始数据的取值范围对工程信号更安全。3.5 阶次谱计算角度域FFT与坐标轴换算角度域信号已经满足平稳性剩下的就是标准FFT流程。唯一需要留神的是横坐标的单位换算。L length(vib_order); w hann(L, periodic); % 加窗抑制频谱泄漏 Y fft(vib_order .* w); A abs(Y(1:floor(L/2))) / L * 2; % 阶次分辨率 SamplesPerRev / L order_res SamplesPerRev / L; order_axis (0:floor(L/2)-1) * order_res; % 绘制阶次谱 figure; plot(order_axis, A); xlim([0 10]); xlabel(Order); ylabel(Amplitude); title(Order Spectrum); grid on;跑完这段阶次谱上应该在1阶和2.5阶的位置出现两个明显的峰值1阶幅值约为1.02.5阶幅值约为0.6。如果能看到这两个峰恭喜整条链路已经通了现在可以换成你自己的真实数据。3.6 完整代码文件组织方式上面几段代码是分步拆解。实际使用我把它们封装成了三个函数文件加一个主脚本主脚本的结构大致如下%% OrderAnalysis_Demo.m 主脚本 % 1. 加载/生成数据模拟或实测 % 2. 提取转速曲线 % [rps_smooth, theta_fit] extractRPM(t, pulse_t); % 3. 角域重采样 % [vib_order, order_axis] angularResample(t, vib, theta_fit, SamplesPerRev); % 4. 计算阶次谱并绘图 % orderSpectrum(vib_order, SamplesPerRev);这样拆的好处是修改互不影响你想换转速提取算法只动extractRPM想改插值方法只动angularResample。实测项目里临时改需求会很频繁这种结构能省下不少返工时间。4. 阶次分辨率与最高分析阶次两个绕不开的参数与工程取值4.1 阶次分辨率取决于总转数而不是采样率很多刚从传统FFT转过来的人会习惯性地想“提高采样率就能提高分辨率”。这话在阶次分析里只对了一半。阶次分辨率的公式是ΔOrder SamplesPerRev / N 1 / 总转数也就是说阶次分辨率只和参与分析的总转数有关。想要分辨相隔0.1阶的两个相邻阶次分量需要至少10转的数据想分辨0.01阶就需要100转。转速上升越快、数据越短阶次分辨率就越差。这跟传统FFT里“频率分辨率1/时长”的规律是一脉相承的。4.2 最高分析阶次欠采样的坑最高分析阶次受限于每转采样点数的一半也就是角度域的Nyquist频率。设SamplesPerRev 512那么最高只能可靠分析256阶。但这里有个容易踩的坑SamplesPerRev是重采样时你指定的固定值如果原始信号的采样率不够高、转速又很快那么在高转速段会出现实际采样点数不足的情况导致隐性的欠采样混叠。工程上事先要估算SamplesPerRev必须小于 fs / max(RPS)。比如采样率5120 Hz、最高转速30转/秒SamplesPerRev绝不能超过170。安全的做法是留20%以上裕量取128或更低。4.3 两个工程上的经验取值参数推荐取值说明每转采样点数128~1024取决于要分析的阶次范围和原始采样率插值方法pchip平滑且无过冲优先选它转速拟合阶数5~7变速平滑时用5剧烈变速时用7窗函数Hannperiodic抑制泄漏兼顾主瓣宽度这些数值都不是拍脑袋定的我之前分别用阶次谱的“谱峰尖锐程度”和“幅值误差”做过交叉验证落在这些区间时结果最稳。5. 踩坑实录四个实际问题和排查思路5.1 键相脉冲抖动导致转速曲线毛刺现场键相脉冲并不像模拟数据那么干净。光电传感器装在联轴器上对齿盘加工误差和安装偏心都非常敏感电涡流传感器对轴面划痕也很敏感。曾经遇到过一个现场脉冲间隔在相邻两个点之间相差接近一倍转速曲线直接出现周期性的尖刺重采样后的阶次谱在低阶背景上隆起一片。排查思路是这样的先画出rps_inst对时间的关系看到周期性毛刺后基本确定是机械误差而不是随机干扰。解决方法是先剔除异常脉冲间隔超出中值30%的点直接去掉再做多项式拟合或者spline平滑。注意一定要先剔野值再拟合否则拟合会被一两个极端值拉歪。5.2 interp1插值出现NaN阶次谱开头或结尾缺一大截这个问题几乎每个第一次写阶次分析代码的人都会遇到。interp1默认不进行外插凡是超出已知数据范围的查询点都会返回NaN。在构造theta_axis时如果取到max(theta_fit)的边缘末尾就可能出现几个NaN点FFT遇到NaN会整体输出NaN阶次谱就废了。最简单的处理就是我在3.4里写的构造完t_resample之后用isnan把无效点截掉。另一种做法是theta_axis 0 : dtheta : max(theta_fit) - dtheta;从源头避免超界。这里提个醒截掉末尾点会损失最后小半转的数据但对总转数较大的分析没有实质影响。5.3 阶次谱幅值标定错误加窗之后忘了修正很多人拿到阶次谱发现幅值大得离谱根因是FFT幅值标定不对。单边谱幅值必须乘2/N这是标准操作但加了Hann窗之后窗函数的相干增益是0.5单根谱线幅值会额外衰减一半。换句话说加Hann窗后要恢复真实幅值修正因子不是2/N而是4/N。实际操作中如果只是做故障诊断、关注相对幅值变化统一加窗后即使不修正也不影响判断但如果要做声功率级计算或者需要跟理论值对比这个4/N的修正必须正确处理。5.4 急变速情况下重采样信号失真升速太猛比如几秒内从500 RPM拉到5000 RPM转速拟合稍微偏差一点角度路径就会积出明显的相位误差重采样后的波形在局部出现拉伸或压缩。这时候阶次谱上的谱峰虽然还在但形状左右不对称幅值也偏低。我的经验是急变速数据不要指望一次全局多项式拟合能行改成基于脉冲间隔的局部线性转速估计配合滑动窗口平滑往往更稳。如果条件允许直接用增量式编码器替代每转一个脉冲的键相传感器角度基准精确度能提升一个量级。最后分享一点实操体会阶次分析这方法本身不新难就难在数据处理链路长任何一个环节出问题都会让结果面目全非。我现在的调试顺序固定不变先画转速曲线确认脉冲没丢再看角域波形有没有毛刺最后才看阶次谱。每一次重采样后的中间结果都会存成图扫一眼虽然繁琐但能省下大量排查时间。你把这套代码跑通后建议把模拟信号里的阶次改一改比如加一个0.5阶的松旷特征再去观察阶次谱的变化会比你直接上真实数据更容易建立直觉。本文还有配套的精品资源点击获取