前阵子帮一个产线朋友处理减速机振动数据，他第一时间把FFT频谱发过来，问“为什么频谱上找不到外圈故障的边带”？其实这是很多刚接触滚动轴承故障诊断的人都会卡住的地方——不是FFT算错了，而是选错了分析频段。滚动轴承早期故障产生的冲击能量非常小，它的特征被调制在高频共振带上，直接对原始信号做FFT，那点冲击成分早被淹没在宽频噪声里了。要解开这个结，我们就得用谱峭度先把整个频带“扫描”一遍，锁定故障冲击最集中的共振频带，再对该频带做包络谱分析。这篇就是基于MATLAB快速谱峭度+包络谱整套流程的完整记录，里面包含我实际调试时踩过的坑和总结出的参数选择经验，适合正在做故障诊断课程设计、毕业论文，以及刚入门设备状态监测的工程师参考。
1. 为什么要用谱峭度和包络谱：轴承故障信号的“隐身术”与“现形记”
1.1 直接做FFT为什么看不到故障频率
很多教材里有个典型说法：轴承外圈故障振动信号中会出现BPFO及其谐波成分。于是大家拿到数据后第一时间做FFT，期望在频谱里看到一根根清晰的谱线，结果往往是整个低频段乱糟糟一片，连转频都找不仔细，更别提BPFO了。这不是FFT代码写错了，而是实际轴承振动信号的结构远比理想模型复杂。
轴承故障的物理过程是：滚动体滚过缺陷时产生一个极短的冲击脉冲，这个脉冲的频带非常宽，会激励起轴承座、端盖、传感器安装部位的固有共振。我们实际采集到的信号，可以近似看作一串窄带高频衰减振荡，振荡频率即系统共振频率，重复频率即故障特征频率。用调幅信号的语言来说，就是“高频载波被低频脉冲调制”。直接对这样的信号做FFT，我们能看到的频谱主体集中在共振频率附近，而故障特征频率由于调制效应会以边频带形式出现，但这些边频带的幅值相比共振峰小得多，很容易被噪声掩盖，在频谱上看就是“什么都没有”。
1.2 包络分析的本质：把高频共振信号搬回低频
包络谱分析（Envelope Spectrum Analysis）的思路很直接：既然故障信息藏在调制包络里，我们就把包络提取出来，再对包络做FFT。这样调制边带被解调回故障特征频率及其谐波，在低频段形成清晰的谱线。实现上最常用的工具是Hilbert变换，得到信号的解析信号后取模，得到的包络信号再FFT，就是包络谱。对于包络分析，关键前置步骤是带通滤波——只有把分析频带对准包含故障冲击能量的共振频带，解调出来的信号才足够干净。
1.3 谱峭度扮演的角色：自动寻找最佳共振频带
既然包络分析需要一个“合适的频带”，那问题就变成了：这个频带在哪？传统做法是用经验反复尝试，比如在2kHz到10kHz之间几个频带挨个试，看哪个包络谱效果最好。试过的人都知道，这非常依赖经验，而且不同转速、不同工况下共振频带会改变，换一台设备就得重新试。
谱峭度（Spectral Kurtosis）就是来解决这个问题的。峭度是四阶统计量，用来衡量信号偏离高斯分布的程度。正常轴承运行时的振动信号接近高斯分布，峭度接近0；一旦出现冲击性故障，信号中会出现大量“尖峰”，峭度显著升高。谱峭度则是把峭度概念推广到每个频率点，计算每个频带上的峭度值，然后画成一张“频率-峭度”图。峭度最大的频带，就是包含最多冲击能量、最适合做包络分析的频带。
快速谱峭度（Fast Kurtogram）由Antoni在2007年提出，通过1/3-二叉树结构高效划分频带，计算每个子带的谱峭度值，形成一张二维色图。这张图可以直观地告诉我们：故障冲击能量集中在哪个频率范围、哪个频带宽度最合适。有了它，包络谱分析就从一个“试错游戏”变成了一个有明确依据的工程流程。
2. 快速谱峭度的核心逻辑：从峭度的数学定义到1/3二叉树频带划分
2.1 峭度的直觉理解与数学形式
先聊透峭度。对一段时域信号x，其峭度定义为：
K = E[(x - μ)⁴] / σ⁴
其中μ是均值，σ是标准差。正态分布的峭度刚好是3，所以工程上常用超额峭度K-3，把高斯信号归零。峭度大说明信号中存在大量远离均值的“尖峰”。你可以把它想象成一个“尖峰探测器”：一个班的学生体重，如果大家都差不多，分布接近正态；如果突然有极少数人胖得离谱，分布就会出现明显的“重尾”，峭度就会明显上升。
轴承故障信号里，滚动体周期性滚过缺陷产生的冲击脉冲，就是这些“极少数胖得离谱的学生”。峭度越大，说明冲击成分在信号中占有的能量比例越高。这也是为什么很多状态监测系统直接把时域峭度作为早期报警指标——简单、快、无需任何先验频率信息。
2.2 谱峭度：在每个频带上分别计算峭度
时域峭度只能反映整个信号的整体冲击性，无法告诉我们“冲击能量集中在哪个频带”。谱峭度的改进正是针对这一点：先把信号做短时傅里叶变换（STFT），在每个频率点f上得到一段复值时间序列X(t, f)，然后计算这个序列的四阶累积量：
K(f) = E{|X(t,f)|⁴} / (E{|X(t,f)|²})² - 2
这里的核心思想是：如果某个频带内信号只是平稳高斯噪声，那么X(t,f)的幅度分布近似瑞利分布，K(f)约等于0；如果该频带内叠加了周期性冲击引起的共振衰减振荡，X(t,f)的幅度会出现明显的闪变，K(f)就会显著增大。所以谱峭度本质上是在逐个频带上做“冲击性检验”，把所有可能包含故障信息的频带都筛出来。
2.3 快速算法的关键：1/3二叉树结构与工程意义
直接在每个频率点都做一次四阶矩估计，计算量非常大，不适合工程应用。Antoni的快速谱峭度算法巧妙地用带通滤波器组代替逐点计算，采用1/3二叉树结构把频带逐层细分：第0层是整个频带，第1层分为两个子带，第2层每个子带再分两个，依此类推；同时在每层之间插入1/3倍频程的频带划分，保证频率分辨率与带宽的合理平衡。每个子带的中心频率和带宽确定后，计算滤波输出信号的复包络峭度，最终生成Kurtogram色图。
这个结构有两个工程上的直接好处。第一是计算速度极快，实测一段10秒、采样率25.6kHz的信号，在普通笔记本上完成6层分解加绘图的耗时通常不超过两秒，完全满足离线分析和在线快速诊断的需要。第二是频带信息直观，色图的横轴是频率，纵轴是分解层数（对应不同带宽），颜色代表峭度值。算法自动推荐峭度最大的子带，但工程师也能通过色图看到次优频带、多频带，为后续人工判断提供依据。
2.4 快速谱峭度工具箱的基本用法
在MATLAB中做快速谱峭度，可以直接使用John Antoni公开的Kurtogram工具箱。我自己的用法是把工具箱文件全部加入当前工程目录，然后用下面这段代码调用核心函数：
% x: 振动加速度信号 % fs: 采样频率 % level: 最大分解层数，经验取4~6 [K, level_f, freq_f] = fast_kurtogram(x, fs, level); % 找到峭度最大的子带 [Kmax, idx] = max(K(:));其中level_f和freq_f分别记录每个子带的分解层数（对应带宽）和中心频率。拿到Kmax对应的索引后，就可以换算最优带通滤波器的上下限频率。需要注意的是，Kurtogram工具箱输出的是特征频率向量与峭度矩阵，使用时务必检查矩阵维度与频带边界对应关系，别在索引换算时弄错。
3. MATLAB实现全流程：从振动数据到最优共振频带的完整代码
3.1 工程问题描述与参数准备
我们以一套典型的滚动轴承故障诊断为场景：电机转速约为1800 rpm（转频fr=30 Hz），加速度传感器安装在轴承座上，采样频率fs=25600 Hz，采样时长10秒。被测轴承为深沟球轴承6205，节圆直径D=39.04 mm，滚动体直径d=7.94 mm，滚动体个数Z=9，接触角α=0°。
先用轴承特征频率公式算出各故障频率，作为后续判别的基准：
%% 轴承参数与特征频率计算 fs = 25600; % 采样率 Hz fr = 30; % 转频 Hz D = 39.04e-3; % 节圆直径 m d = 7.94e-3; % 滚动体直径 m Z = 9; % 滚动体数量 alpha = 0; % 接触角 rad BPFO = Z * fr / 2 * (1 - d/D * cos(alpha)); % 外圈故障特征频率 BPFI = Z * fr / 2 * (1 + d/D * cos(alpha)); % 内圈故障特征频率 BSF = D / (2*d) * fr * (1 - (d/D * cos(alpha))^2); % 滚动体故障特征频率 FTF = fr / 2 * (1 - d/D * cos(alpha)); % 保持架故障特征频率代入数值后得到：BPFO约107.5 Hz，BPFI约162.5 Hz，BSF约70.2 Hz，FTF约11.9 Hz。这些频率在后面对照包络谱谱线时非常关键，尤其是BPFO和BPFI，它们之间相差大约55 Hz，与转频的调制关系也能帮助我们区分内圈还是外圈故障。
3.2 加载信号并通过快速谱峭度定位最优频带
读入振动数据后，我习惯先绘制时域波形的粗览图，确认数据中没有明显的仪器饱和、断点或异常冲击，再进入谱峭度计算。这样能避免把传感器敲击、噪声脉冲误当成轴承故障信号。
%% 加载信号（示例代码，按实际文件格式调整） load('bearing_fault.mat'); % 内含 x 和 fs t = (0:length(x)-1) / fs; %% 快速谱峭度计算 level = 6; [K, level_f, freq_f] = fast_kurtogram(x, fs, level); %% 找到最优频带 [Kmax, idx_max] = max(K(:)); [lv_opt, f_opt] = ind2sub(size(K), idx_max); % 根据level_f和freq_f矩阵换算最优频带的中心频率和带宽 bw_opt = fs / 2^(level_f(lv_opt, f_opt) + 1); % 简化估算，具体需按工具箱输出换算 fc_opt = freq_f(lv_opt, f_opt); % 最优频带中心频率 f_low = fc_opt - bw_opt/2; f_high = fc_opt + bw_opt/2; fprintf('最优频带: %.1f Hz ~ %.1f Hz, 峭度: %.2f\n', f_low, f_high, Kmax);这里有一个容易犯的坑：level_f和freq_f在不同版本工具箱中的结构略有差异，有的记录的是子带边界频率，有的记录中心频率，直接套用索引可能得到错误频带。我自己的做法是先打印几个频带的freq_f数值，与色图横轴对比确认，再决定如何换算带宽。如果嫌麻烦，也可以用kurtogram函数直接输出Kurtogram图，图上标注的最优频带参数是经过换算的，可以直接读取。
3.3 用带通滤波器提取共振信号并用Hilbert求包络谱
得到最优频带后，下一步就是对原始信号做带通滤波、提取包络、计算包络谱。滤波器我推荐用Butterworth加filtfilt零相位滤波，避免普通滤波引起的相位失真干扰冲击位置判断。滤波阶数4到6阶足够，太高会把共振频带削出振铃，太低则阻带衰减不足。
%% 带通滤波：截取最优共振频带 ord = 4; [b, a] = butter(ord, [f_low, f_high]/(fs/2), 'bandpass'); x_f = filtfilt(b, a, x); %% 希尔伯特变换提取包络 env = abs(hilbert(x_f)); %% 包络谱计算 N = length(env); env = env - mean(env); % 去除直流分量 EnvSpec = abs(fft(env)); f_axis = (0:N/2-1) * fs / N; EnvSpec = EnvSpec(1:N/2); %% 绘制包络谱与特征频率标注 figure; plot(f_axis, EnvSpec); xlim([0, 500]); hold on; for k = 1:4 xline(k*BPFO, '--r', sprintf('%dBPFO', k)); xline(k*BPFI, '--g', sprintf('%dBPFI', k)); end xlabel('频率 (Hz)'); ylabel('幅值'); title('最优频带包络谱'); legend('包络谱', 'BPFO谐波', 'BPFI谐波'); grid on;运行后如果故障特征明显，包络谱中会在BPFO及其2倍、3倍频处出现清晰谱峰，对应的就是我们最关心的故障指示。这里强调一个操作细节：包络谱的纵轴是线性幅值还是对数幅值，会影响对低频小幅值谱线的判断。我一般先用线性坐标看主峰，再用对数坐标观察微弱的高次谐波，两者配合能更全面评估故障严重程度。
3.4 没有工具箱时的替代方案：分频带峭度扫描
有的情况下不方便下载第三方工具箱，此时可以自写一个简化的分频带峭度扫描程序，作为快速谱峭度的轻量替代。思路不复杂：把整个频率范围均匀划分成若干候选频带，用带通滤波器逐个提取后计算时域峭度，峭度最高的频带就是近似最优频带。这个方法的频率分辨率不如1/3二叉树精细，但对于早期故障已经足够。
%% 简化版分频带峭度扫描 nband = 48; % 候选频带数量，可调 bw = fs / nband; % 每个频带的宽度 kurt_band = zeros(1, nband); for k = 1:nband fL = (k-1) * bw; fH = k * bw; [bb, aa] = butter(6, [fL, fH]/(fs/2), 'bandpass'); y = filtfilt(bb, aa, x); kurt_band(k) = kurtosis(y); end [~, best] = max(kurt_band); fL_best = (best-1) * bw; fH_best = best * bw; fprintf('简化扫描最优频带: %.0f ~ %.0f Hz\n', fL_best, fH_best);但这个简化方法有一个明显局限：频带宽度是固定的，如果真实共振频带比划分带宽窄或宽，峭度值会被平均稀释。所以我的建议是，先用粗划分（比如32个频带）看大致位置，再用更细的划分在候选位置附近加密扫描，结合包络谱效果综合判断。
4. 谱峭度+包络谱联合诊断：从谱线到故障类型的判别方法
4.1 特征频率与转频的调制关系
拿到包络谱后，如何判断是外圈、内圈还是滚动体故障？最直接的依据是看谱峰是否出现在对应的理论特征频率及其谐波处。但在实际谱图中，还需要关注以下几点：
外圈故障的特征是BPFO及其谐波处有谱峰，且谱峰处通常伴随转频fr的边带。由于外圈固定安装于轴承座，载荷分布相对稳定，所以BPFO频段的谱线通常尖锐、能量集中，边带相对较弱。
内圈故障则不同。内圈随转子一起旋转，故障点周期性进入和退出承载区，导致冲击幅度周期性变化，因此BPFI处不仅有谱峰，其两侧还会出现明显的转频fr及其倍频调制边带。换句话说，“BPFI谱峰带fr边带”是内圈故障的典型指纹。
滚动体故障的包络谱特征相对复杂。由于滚动体自转频率较高，且故障点有时与内圈接触、有时与外圈接触，冲击幅度不稳定，BSF处谱峰往往较矮，且伴随保持架特征频率FTF或转频的边带。实际诊断中BSF附近噪声偏大，需要结合更高阶谐波综合判断。
保持架故障的特征频率FTF通常低于10 Hz，且幅值很小，容易被频谱低频段的高能量漏过掩盖，因此包络谱分析对保持架故障的灵敏度有限，这里不做重点。
4.2 内圈故障和外圈故障的判别实例
我拿一组典型的内圈故障数据举例。谱峭度选出的最优频带在4.5kHz到7kHz附近，该频带包络谱在162.5 Hz处出现绝对主峰，同时在107.5 Hz处几乎看不到谱线。乍一看像是单纯的内圈故障（BPFI），但仔细观察162.5 Hz主峰两侧，能看到间隔约30 Hz的一对边带，也就是132.5 Hz和192.5 Hz。转频是30 Hz，所以这个模式与“内圈故障、周期性进入承载区调制”完全吻合，判断为内圈故障无疑。
换一组外圈故障数据，包络谱主峰出现在107.5 Hz处，而且2倍频215 Hz、3倍频322.5 Hz处都有可辨识的谱峰，主峰两侧的边带减弱到几乎不可见。即使不看时域波形和包络波形，单凭这个谱图也能做出外圈故障的判断。
我的体会是，特征频率绝对值重要，但谱峰之间的“相对关系”——谐波序列是否完整、边带间距是否等于转频、主峰与边带的幅值比——才是区分故障类型的关键。只看单一谱线位置很容易误判。
4.3 多频带对比：为什么只看最优频带不够
工程中我常遇到一种情况：Kurtogram色图上不只一个高峭度区域，最优频带包络谱很干净，但次优频带包络谱能揭示更多信息。举一个具体例子：某次分析中，算法选出的最优频带在6kHz附近，包络谱显示出明显的BPFO峰；但另一个峭度次高的频带在2.5kHz附近，包络谱中除了BPFO，还出现了明显的啮合频率边带。后查得知该设备齿轮箱与轴承座之间存在结构耦合，齿轮故障的冲击也调制到了这一频带。
这意味着谱峭度选出的最优频带，反映的是“冲击能量最集中的频带”，不一定是“唯一包含诊断信息的频带”。因此我的做法是：在自动推荐的基础上，对Kurtogram色图中峭度前3到5名的频带都做一次包络谱，把多频带包络谱并排显示对比。如果所有高峭度频带都指向同一故障频率，结论可靠性会大幅提高；如果不同频带指向不同故障，则应警惕多重故障或来自其他部件的干扰。
5. 实测中的参数陷阱与调参经验：那些会让结果“翻车”的细节
5.1 数据长度不足导致频率分辨率不够
包络谱的频率分辨率取决于包络信号时长，Δf = fs/(N_band)，其中N_band是包络信号长度。假设信号时长只有0.25秒，那么包络谱频率分辨率约4 Hz。对转频30Hz、BPFO约107.5Hz的轴承来说，4 Hz的分辨率勉强够用，但如果故障特征频率之间相差很小（比如BSF与某个边带仅差2 Hz），就完全无法区分。
我给初学者的建议是：至少采集20到50个故障周期。外圈故障周期为1/BPFO约0.0093秒，50个周期约0.47秒，所以1秒以上的数据对BPFO级别的故障足够；但如果要观察转频边带，尤其是保持架故障（FTF约12Hz），则需要至少5到10秒数据才能清晰分辨。压箱底的经验是：故障诊断宁可多采几秒数据，也不要为了节省存储空间缩短采样时长，后期发现分辨率不够再补采，成本更高。
5.2 转速波动对谱线模糊化的影响
如果设备在分析期间存在明显转速波动，故障冲击的重复周期不再严格恒定，包络谱中的谱峰会被“抹开”，幅值降低、带宽变宽，甚至难以识别。这对谱峭度算法本身是致命打击——它假设故障频率稳定，一旦转速漂移，同一频带上计算的峭度值会被多段不同频率的冲击平均掉，导致最优频带选择失准。
遇到转速波动工况，我会优先使用阶比跟踪或重采样技术，先通过键相脉冲把信号从时间域映射到角度域，再做谱峭度和包络分析。如果没有键相信号，至少应检查一段数据的瞬时转速（比如用转速计或基于转频的STFT估计），如果波动超过1%，就要谨慎解释包络谱结果。
5.3 最优频带接近Nyquist频率时的“假阳性”
快速谱峭度有时会把最优频带选在接近fs/2的高频区。这时候要警惕，因为它可能不是轴承共振，而是信号采集链路的伪迹或结构噪声。例如传感器的固有谐振峰值较高，或电缆屏蔽不佳引入的电磁干扰，都会产生高频冲击性成分，拉高高频段峭度。
排查方法很简单：先看时域波形，如果“冲击”出现的时间间隔完全随机、与轴转周期毫无关联，大概率不是轴承故障；再看原始频谱中该频带附近是否有工频及其谐波，如果有，应考虑电源干扰。我处理过一批现场数据，谱峭度把最优频带定位到12kHz，包络谱却没有任何轴承特征频率，反而出现了100Hz电源纹波的谐波，最终确认是传感器线缆屏蔽层破损。这类假阳性在实验室里不常见，但在现场非常常见。
5.4 带通滤波参数不当导致的信号失真
很多人在选择带通滤波器时随手用bandpass函数，默认阶数不够高，结果带外噪声泄漏严重；也有人把阶数拉到20以上，结果滤波后的信号出现明显的振铃现象，包络谱中产生大量虚假伪峰。
我推荐的原则是：Butterworth滤波器阶数4~8阶，使用filtfilt做零相位滤波。中心频带较窄时适当提高阶数到8，较宽时4到6阶足够。滤波完成后必须对比滤波前后的时域波形，确保冲击序列依旧保持原有的周期性。如果发现滤波后冲击之间出现了额外的“尾巴”或振荡，多半是阶数过高或频带边界选择不当，需要调整。
5.5 峭度值下降不等于故障消除：晚期故障的识别问题
一个鲜为人知的坑是：随着故障从早期到中晚期发展，时域峭度和谱峭度反而可能下降。原因是故障严重时冲击变成连续的大面积接触，不再是单一“尖峰”，同时多个冲击叠加后的分布趋于平滑，峭度回落到接近正常水平。也就是说，如果你只在峭度值最高的时段采集数据，可能恰好错过最紧急的报警窗口。
这就是为什么我不建议单独依赖谱峭度做趋势监测，而应该同时跟踪包络谱特征频率处的幅值、时域峭度、有效值（RMS）等多维指标。实际项目中，我习惯记录每个时刻的最优频带峭度值和对应特征频率幅值，把它们存入历史数据库，观察趋势变化。如果峭度值下降但BPFO幅值仍在上升，说明故障在恶化但冲击性减弱，此时更应保持警惕。
6. 程序扩展思路：从单段信号诊断到批量状态监测
6.1 批处理脚本：一次处理几十个数据文件
实际项目里很少只分析一段信号，最常见的需求是把一批文件全部跑一遍，输出一个汇总表。我的做法是写一个批处理函数，外层用dir遍历文件夹，内层调用前面实现的单段分析函数，最后把故障频率幅值、最优频带、峭度值写入Excel。
files = dir('data/*.mat'); results = table(); for i = 1:length(files) load(fullfile(files(i).folder, files(i).name)); % 调用谱峭度+包络谱分析函数 [bpfo_amp, bfi_amp, kurt_max, f_low, f_high] = bearing_diag(x, fs, fr, D, d, Z); results = [results; table({files(i).name}, bpfo_amp, bfi_amp, kurt_max, f_low, f_high)]; end writetable(results, 'diagnosis_results.xlsx');批处理的关键在于容错。个别文件可能因传感器松动、数据损坏等原因产生异常波形，如果不加保护，整个循环会中断，前面的结果全部丢失。所以我的循环里必加try-catch，异常文件单独记录，继续处理后续文件。
6.2 与包络谱结合构建健康趋势指标
单次诊断只能回答“现在是否故障”，回答不了“故障如何发展”。要构建趋势，需要把每次分析得到的特征频率幅值、峭度、RMS等指标串联成时间序列。以BPFO幅值为例，健康状态时几乎为零，早期故障出现后开始升高，中期可能保持稳定或波动，严重时达到峰值甚至下降。以RMS为例，故障初期变化不大，中后期快速上升。
实际监测项目里，我会把多个指标放进同一个趋势图中，用不同颜色区分，并叠加报警阈值。阈值设定方法有两种：一是基于健康样本统计的3σ准则，二是基于行业标准经验值。前者更客观，但前提是故障前有多组健康数据；后者更省事，但需要根据具体设备调整。
6.3 与其他诊断方法的融合
谱峭度+包络谱是非常有效的方法，但并非万能。现场诊断时，我还会配合时域统计指标（RMS、峭度、峰值因子）、频段能量比、倒频谱、深度学习的频谱图像识别等方法交叉验证。谱峭度的最大价值在于“降维”——把宽频段的搜索问题变成少数几个候选频带，大幅节省人工试错时间。而包络谱的价值在于“精判”——在选定的频带内，给出物理意义明确的特征频率谱线。两者结合，构成了我目前最顺手的一套分析流程。
已跑通的这套MATLAB程序，我已经在实验室故障注入数据和现场采集数据上反复验证过。我自己优化后的流程是：先用Kurtogram快速定位候选频带，再用多频带包络谱交叉验证，最后结合时域特征核对。整个过程从加载数据到输出诊断结论，通常不到5分钟。最后分享一个小技巧：把fast_kurtogram的色图和包络谱图保存在同一个PDF中，诊断报告就能非常直观地向非专业人员展示“为什么在这个频带分析”“故障频率在哪里”——这比单纯贴一张频谱图有说服力得多。
