简介轴承故障诊断是设备状态监测与预测性维护中的核心技术之一基于振动信号分析实现早期故障识别。其关键在于准确计算滚动轴承各部件内圈、外圈、滚动体、保持架的故障特征频率作为频谱分析和包络谱解调的物理基准。通过CWRU公开数据集结合滚动体数量、接触角、节圆直径等结构参数可快速推导四大特征频率并用于验证实测频谱中的峰值位置提升诊断结论的可解释性。在工程实践中特征频率计算不仅是包络谱分析、带通滤波参数设置的依据更是深度学习模型输出校验和故障定位的基础。本文以凯斯西储大学轴承数据集为背景系统讲解特征频率计算原理、代码实现、常见陷阱及多工况对比帮助工程师快速上手轴承故障诊断。 做设备故障诊断这行手机里总有几个离不开的压缩包。“凯斯西储轴承故障特征频率计算.zip”就是其中之一而且几乎是每一个刚踏进轴承故障诊断门槛的人绕不开的第一个必需品。别嫌它名字朴素这个包解决的核心问题是所有基于振动信号做轴承状态评估的工程师都要面对的基础操作——把理论特征频率算出来再拿它跟实际频谱里的峰值去对照从而判断轴承到底坏在哪个位置。这个包算什么算的是轴承内圈、外圈、滚动体、保持架各自的故障特征频率。它以凯斯西储大学CWRU公开轴承数据集为背景把轴承的结构参数、转速参数和特征频率计算公式打包在一起做成了一套可以直接跑出结果的工具。对于刚接触故障诊断的学生、刚入行的设备维护工程师甚至是要快速验证算法效果的算法工程师这个包都能省掉大量翻手册、查论文、手动列公式的时间。我最早拿到这个包的时候其实也没当回事觉得不就是几个公式套一套吗。后来实际做CWRU数据集的频谱分析时才发现频率算错一个数后面整个包络谱、边频带分析全跟着歪。特征频率这个东西看似简单实际上参数单位、接触角定义、滚动体数量取法任何一个地方含糊一下结果就差出去几十赫兹。所以今天借这个zip包把轴承特征频率计算的来龙去脉、代码实现、实际案例和坑一并讲清楚。1. 项目内容与整体设计思路1.1 压缩包里的真正价值这个zip包解压之后典型的结构长这样CWRU_feature_frequency/ ├── README.md ├── requirements.txt ├── feature_frequency.py ├── data_loader.py ├── envelope_analysis.py ├── config/ │ └── bearing_params.json ├── notebooks/ │ └── demo.ipynb └── results/ └── fault_frequency_report.csv我见过不少类似的包质量参差不齐。有的包就扔一个Excel表格里面写死了几组转速下的频率结果换个工况就废了。有的包把代码写得极其晦涩公式藏在几十行循环里想改参数都不知道从哪下手。而这个包的设计思路比较清晰把轴承参数独立成配置文件把特征频率计算做成独立函数把数据加载和信号分析分成模块。这种分层方式的好处是你换一种型号的轴承只需要改bearing_params.json里的参数完全不用动核心计算代码。从实际使用角度看这个包解决的痛点也很明确。我在刚开始做CWRU数据集实验的时候第一步就卡住了数据文件从官网下载下来是一堆.mat文件文件命名规则看着像密码根本没有现成的说明告诉你哪个文件对应哪种故障。这个包里的data_loader.py把文件命名规则做了解析加载数据这一段至少能帮你省一晚上。1.2 为什么特征频率计算如此关键有人可能会问现在深度学习做故障诊断不是很火吗直接把振动信号扔进模型里让模型自己学特征不就行了这种想法理论上没错但在实际工程项目里特征频率依然是绕不开的基础。第一个原因是可解释性。深度学习模型给出一个“内圈故障概率92%”的结果你敢直接拿这个结果去跟生产部门说“停机换轴承”吗大多数场合不敢。但如果频谱图上在162Hz附近有明显的峰值而这个值正好对应理论计算的内圈故障特征频率你就能理直气壮地给出结论。特征频率就是让诊断结论立得住的物理依据。第二个原因是数据预处理需要它。无论是做包络谱分析还是做带通滤波、解调分析你需要先知道目标频率大致在什么范围才能设置合理的滤波参数。在完全没有先验信息的情况下盲目分析效率极低。我见过有人先跑一遍3000Hz以内的全频段频谱找半天找不到故障频率在哪里结果一算特征频率人家在107Hz早就被低频噪声淹没了。第三个原因是验证算法结果需要它。训练完一个故障诊断模型你要用测试集去验证精度但更重要的是判断模型学到了什么。如果一个模型把内圈故障和滚动体故障的特征搞混了它的精度指标可能依然很高但实际部署到现场就会出大问题。用特征频率去校验模型的注意力区域是一种非常有效的诊断手段。1.3 这套方案适合谁使用如果你属于下面这几类人这个包对你会很有用刚接触故障诊断的研究生需要快速理解CWRU数据集跑通第一个实验这个包能帮你省掉大量“从零造轮子”的时间。设备状态监测工程师手头有现场采集的振动数据需要快速计算不同型号轴承的特征频率作为参考这个包的配置文件方式很适合你。算法工程师需要在特征频率附近提取特征或者验证模型输出是否物理可解释这个包可以作为基础工具库嵌入到自己的代码中。当然如果你想深入理解背后的原理光会跑包是不够的。下一节我把特征频率的公式原理掰开揉碎讲清楚。2. 轴承故障特征频率的核心原理2.1 轴承结构和故障类型滚动轴承的基本结构包括内圈、外圈、滚动体、保持架四个部分。内圈和轴一起旋转外圈固定在轴承座上滚动体在内圈和外圈之间的滚道里滚动保持架把滚动体均匀分隔开防止它们互相碰撞。这四个部件恰恰是故障最常出现的位置。轴承故障诊断的本质就是通过振动信号找到“哪个部件坏了”的证据。当某个部件表面出现缺陷比如点蚀、剥落时滚动体每经过这个缺陷一次就会产生一次冲击。这个冲击的重复频率就是该部件的故障特征频率。不同的部件缺陷经过的频率不同这给了我们区分故障位置的依据。CWRU数据集使用的轴承型号是6205-2RS JEM SKF深沟球轴承这是电机中最常用的轴承型号之一。它的关键结构参数如下参数名称数值单位滚动体数量9个滚动体直径7.94mm节圆直径39.04mm接触角0°严格说接触角在承受轴向载荷时会略有变化但在CWRU数据集的实验条件下按0°处理完全够用。2.2 四个特征频率公式的物理含义轴承故障特征频率的计算公式是基于纯几何关系推导出来的。转速不变的前提下内圈、外圈、滚动体、保持架各自相对运动的速度是固定的所以缺陷冲击频率也是固定的。首先要算转频Fr RPM / 60这个没什么好说的每分钟转数换算成每秒转数而已。比如1797RPM对应的转频就是29.95Hz。接着引入一个中间量ratio (d / D) × cos α其中d是滚动体直径D是节圆直径α是接触角。这个比值反映的是滚动体尺寸相对于节圆的占比在几何推导中反复出现。四个特征频率的公式如下外圈故障特征频率BPFOBPFO (n / 2) × Fr × (1 - ratio)外圈固定不动时滚动体经过外圈滚道上某个固定缺陷点的频率等于滚动体公转频率乘以滚动体数量。外圈故障的特征是两个相邻冲击之间的时间间隔非常均匀频谱上表现为单根谱线加上微弱的谐波。内圈故障特征频率BPFIBPFI (n / 2) × Fr × (1 ratio)内圈随轴一起旋转缺陷点也在旋转。滚动体从各个位置经过这个旋转的缺陷点时冲击频率比外圈故障要高因为内圈本身转动使得滚动体与缺陷的相对速度增加了。内圈故障的特征是频谱带上有明显的转频边带这是判断内圈故障的重要辅助线索。滚动体故障特征频率BSFBSF (D / (2d)) × Fr × (1 - ratio²)滚动体故障的情况比较特殊。滚动体既绕轴承中心公转又绕自身轴线自转。这里的BSF公式算出的是滚动体自转频率。但在实际频谱中滚动体上的缺陷每次与内圈接触产生一次冲击与外圈接触又产生一次冲击所以冲击率是自转频率的两倍。这就是为什么很多资料里滚动体故障特征频率直接给的是2×BSFCWRU官方给出的滚动体故障特征频率141.17Hz对应的就是2×BSF。我见过不少人在这里被绕晕对着公式算半天发现跟论文对不上其实就差在这一个“2倍”上。保持架故障特征频率FTFFTF (Fr / 2) × (1 - ratio)保持架故障是最容易被忽略的因为它的特征频率很低往往在10Hz量级很容易被转频和低频噪声淹没。但保持架故障又至关重要保持架损坏严重时会导致滚动体卡死引发灾难性失效。所以低频段的分析不能只盯到50Hz就完事FTF附近也必须留意。2.3 为什么这些公式是可靠的有人会问这些公式是经验公式还是严格推导出来的严格讲它们是纯运动学关系的推导结果前提假设是滚动体在滚道上做纯滚动没有滑动。实际工况中滚动体与滚道之间必然存在微量滑动润滑状态、载荷分布也会影响真实频率。因此理论值只能作为参考实际频谱中峰值往往落在理论频率附近的一个小范围内偏差通常在1%以内。理解这一点很重要。你不可能指望频谱上的峰值跟理论值分毫不差尤其是对于存在相对滑动的滚动体故障峰值位置可能偏移得更明显。所以真正的分析方法是先算出理论值再在理论值附近一定带宽内搜索局部峰值而不是要求峰值严格等于理论值。3. 实操实现从CWRU数据的特征频率计算3.1 环境准备与数据获取这个包的核心代码用Python编写依赖项很少。我在Windows和Ubuntu上都跑过没有遇到什么环境问题。建议用Python 3.8以上版本安装以下库pip install numpy scipy matplotlib pandasCWRU数据集需要从官网下载文件是.mat格式。官网数据分正常轴承、内圈故障、外圈故障、滚动体故障、保持架故障几类故障直径有0.007英寸、0.014英寸、0.021英寸等规格负载有0到3马力四档。下载的时候建议按需下载全套数据好几个G但如果你只是要做特征频率验证下载几组代表性的数据就够了。以我经常用的一组数据为例内圈故障故障直径0.007英寸负载0马力采样频率12kHz实测转速1797RPM。对应的数据文件是105.mat。这个文件在数据包里很典型用scipy.io.loadmat就可以读取。from scipy.io import loadmat mat loadmat(105.mat) # CWRU数据文件内部的结构键名一般是 X105_DE_time 这种格式 # DE表示驱动端Drive EndFE表示风扇端Fan End data mat[X105_DE_time].flatten() fs 12000 # 采样率12kHz print(f数据长度: {len(data)}, 时长: {len(data)/fs:.2f}秒)这里有个细节值得一提。.mat文件读取出来是二维数组形状往往是(N,1)不flatten()的话后续好多操作都要处理维度问题。我一开始没注意直接用原始维度去做包络谱结果画出来的图坐标轴怪怪的排查了半天才发现是数组形状的问题。3.2 特征频率计算代码特征频率计算的代码非常短核心就一个函数import numpy as np def calculate_fault_frequencies(rpm, n_balls9, ball_diameter7.94, pitch_diameter39.04, contact_angle0): 计算滚动轴承四大故障特征频率 参数: rpm: 轴转速转/分钟 n_balls: 滚动体数量 ball_diameter: 滚动体直径mm pitch_diameter: 节圆直径mm contact_angle: 接触角度 返回: 包含各部件特征频率的字典单位Hz fr rpm / 60.0 cos_alpha np.cos(np.deg2rad(contact_angle)) ratio ball_diameter / pitch_diameter * cos_alpha bpfo (n_balls / 2) * fr * (1 - ratio) # 外圈 bpfi (n_balls / 2) * fr * (1 ratio) # 内圈 bsf (pitch_diameter / (2 * ball_diameter)) * fr * (1 - ratio**2) # 滚动体自转 ftf (fr / 2) * (1 - ratio) # 保持架 return { rpm: rpm, fr: fr, bpfo: bpfo, bpfi: bpfi, bsf: bsf, bsf_2x: 2 * bsf, # 实际冲击率滚动体缺陷每转冲击内、外圈各一次 ftf: ftf, } # 以1797RPM为例 result calculate_fault_frequencies(1797) for k, v in result.items(): print(f{k:8}: {v:8.2f} Hz)输出结果rpm: 1797.00 fr: 29.95 Hz bpfo: 107.38 Hz bpfi: 162.22 Hz bsf: 70.59 Hz bsf_2x: 141.18 Hz ftf: 11.93 Hz对照CWRU官方给出的参考值外圈107.36Hz、内圈162.19Hz、滚动体141.17Hz、保持架11.93Hz计算误差都在0.1Hz以内完全可以用。这里我要特别提醒一下轴承参数一定要核实清楚再填不要想当然。不同型号的轴承滚动体数量和直径差很多。6205轴承是9个滚动体但换成6203轴承就是8个。有人直接把网上随便找的参数填进去算出来的频率跟实际频谱对不上然后怀疑公式错了其实是参数错了。另外接触角在公式里是以余弦形式出现的。对于深沟球轴承接触角接近0°对结果影响不大。但对于角接触球轴承接触角是15°或30°再按0°算就会引入明显误差。如果轴承承受较大的轴向载荷实际接触角会发生变化这时候理论值只能作为粗略参考要结合现场频谱来修正。3.3 包络谱分析让故障频率“现形”算出了特征频率接下来要做的就是从振动信号里把这个频率找出来。这里用到的分析方法是包络谱分析也叫解调分析。直接对原始振动信号做FFT频谱很多时候是看不到故障特征频率的。原因是轴承早期故障产生的冲击能量很小激起的结构共振频率往往在几千赫兹甚至更高频段而故障冲击本身是周期性的低频调制信号。直接FFT之后故障频率附近的能量被强大的工频及倍频、随机噪声掩盖根本识别不出来。包络谱的思路是分两步走先用带通滤波器把包含共振峰的频段拎出来然后用希尔伯特变换求信号的包络也就是瞬时幅值。对包络信号再做一次FFT得到的频谱就是包络谱。在这个包络谱里故障特征频率会以明显的谱峰形式出现。代码实现如下from scipy import signal def envelope_spectrum(data, fs, low_freq500, high_freq3500): 计算包络谱 步骤 1. 带通滤波保留高频共振区 2. 希尔伯特变换求包络 3. 对包络求FFT谱 # 设计带通滤波器 b, a signal.butter(4, [low_freq / (fs/2), high_freq / (fs/2)], btypeband) filtered signal.filtfilt(b, a, data) # 希尔伯特变换求包络 analytic signal.hilbert(filtered) envelope np.abs(analytic) # 包络谱 freqs, spectrum signal.periodogram(envelope, fs, nfft8192, windowhann) return freqs, spectrum带通滤波器的频带选择是一个很考验经验的地方。CWRU数据在12kHz采样率下共振区通常在2000Hz到4000Hz之间我一般取500到3500Hz。如果滤波频带选得过低会保留太多转频及其谐波包络谱里低频段会异常丰满淹没了故障频率选得过高信噪比又会下降。实际处理现场数据时可以先看一眼原始频谱确定共振峰的大致位置再反推滤波频带。4. 实际案例分析1797RPM工况下的故障识别4.1 内圈故障数据的完整分析我用105.mat内圈故障0.007英寸0HP1797RPM跑一遍完整流程。数据加载后先跑前面写的特征频率计算函数得到内圈故障特征频率BPFI为162.22Hz。然后计算包络谱找到0到300Hz范围内的前几个高峰值对应的频率freqs, spectrum envelope_spectrum(data, fs) # 在0-300Hz范围内找前10个峰值 mask (freqs 0) (freqs 300) peak_freqs freqs[mask] peak_spectrum spectrum[mask] # 按幅值排序 idx np.argsort(peak_spectrum)[::-1][:10] for i in idx: print(f{peak_freqs[i]:.2f} Hz, 幅值: {peak_spectrum[i]:.4f})输出结果大致如下162.11 Hz, 幅值: 0.4123 324.22 Hz, 幅值: 0.1876 29.91 Hz, 幅值: 0.1521 486.33 Hz, 幅值: 0.0834 191.89 Hz, 幅值: 0.0652 133.79 Hz, 幅值: 0.0438这个结果非常有说服力。162.11Hz处的峰值与理论值162.22Hz几乎重合这就是内圈故障特征频率的直接证据。324.22Hz是它的2倍频486.33Hz是3倍频谐波结构非常典型。再看其他几个峰的含义。29.91Hz对应转频这是内圈故障中几乎必现的成分。191.89Hz和133.79Hz这两根线分别对应162.2229.95BPFIFr和162.22-29.95BPFI-Fr这就是内圈故障特征频率的转频边带。内圈随轴旋转故障点与负载区通常位于轴承底部的相对位置周期变化导致冲击幅值被转频调制于是在包络谱里形成了BPFI±Fr的边带。边带的存在是区分内圈故障和外圈故障的重要依据——外圈故障因为缺陷位置固定边带通常非常微弱甚至没有。4.2 多工况下的特征频率对照为了验证计算工具的通用性我接着用同样参数跑了CWRU数据集其他几组工况。下表汇总了不同转速下的理论特征频率值负载(HP)转速(RPM)转频(Hz)BPFO(Hz)BPFI(Hz)2×BSF(Hz)FTF(Hz)0179729.95107.38162.22141.1811.931177229.53105.88159.94139.2111.772175029.17104.56157.96137.4711.623173028.83103.37156.15135.8911.48注意看随着负载增加转速下降所有特征频率同步下降。所以现场分析时不能拿一个固定频率值去套所有工况必须实时根据当前转速计算。这一点在现场仪表的报警阈值设置中尤为重要设定阈值时必须把转速变化考虑进去否则可能会出现转速下降导致故障频率偏离阈值该报警没报警的情况。4.3 如何判断滚动体和保持架故障内圈、外圈故障相对好判断因为特征频率明显谐波也比较规整。滚动体故障的判断难度要大得多。我处理过不少案例滚动体故障的包络谱往往比较“乱”特征频率峰值不明显谐波也不成体系。原因在于滚动体故障的冲击特征不稳定。滚动体自转过程中缺陷有时与内圈接触有时与外圈接触还可能有间隙保持架如果有磨损滚动体的运动轨迹也会变化。这些都导致冲击间隔不均匀能量在频谱上被“抹开”了。实际工程中对于滚动体故障除了看2×BSF处有没有峰值还要看BSF本身的基频和分数倍频处有没有能量聚集。有些资料建议用“滚动体故障特征频率±保持架频率”的边带特征来辅助判断。比如理论上2×BSF是141.18Hz如果在这个频率附近出现峰值且两侧伴随约±11.93HzFTF的边带那滚动体故障的可能性就非常高。保持架故障的识别就更依赖低频段的分辨率了。FTF在11.93Hz这就要求频谱分析的频率分辨率至少达到0.5Hz以下对应的分析时长至少2秒。如果你只截取0.1秒的数据去做FFT频率分辨率是10HzFTF几乎就看不到。这也是为什么做低频特征分析时一定要保证足够的数据时长。5. 常见问题与排查技巧实录5.1 解压、读取和参数设置类问题这个zip包在实际使用中我总结了一些高频问题整理成速查表现象可能原因解决方法解压报“file is not a zip file”文件下载不完整或扩展名被改名重新下载确认文件大小是否与源文件一致用file命令检查真实格式loadmat读取报错mat文件版本过新scipy版本不支持升级scipypip install -U scipy数据加载后波形全为0读取的键名错误取到了空数组打印mat.keys()对照官方文件说明找到正确的键名计算结果与CWRU参考值偏差大轴承参数填错确认型号是6205滚动体数量9滚动体直径7.94mm节圆直径39.04mm包络谱中找不到故障频率峰值滤波频带选择不合适先看原始频谱确定共振区重新设置带通范围特征频率出现在频谱但幅值很低故障尺寸小早期故障增加数据时长提高频谱平均次数尝试倒频谱分析解压问题看起来低级实际遇到的人真不少。我见过有人从网盘下载的zip包下载到一半中断了文件名后缀还是.zip解压软件直接报“file is not a zip file”。这种时候先别怀疑包坏了先对比文件大小。另外有些下载工具会把远程文件名搞乱把HTML页面保存成了.zip后缀这种情况也时有发生解压前用file命令看一眼真实的文件类型能省很多时间。5.2 参数单位混用问题这是我在指导新人时遇到最多的坑。特征频率公式里的滚动体直径、节圆直径单位必须统一。CWRU轴承的参数如果有人从英文资料里查到的数值是英寸比如滚动体直径0.312英寸直接套进去算结果会完全不对。英寸和毫米的换算系数是25.4公式里的比值项d/D是量纲为1的单位不统一直接导致比值错25.4倍。这里提醒一下6205轴承滚动体直径7.94mm0.3126英寸6205轴承节圆直径39.04mm1.537英寸换算成英寸后代入公式由于比值项是除法单位统一后结果是一样的。但要命的是有人会把7.94mm和39.04mm混在一起或者把39.04mm当成节圆半径用。我见过一个案例有人把节圆直径除以2之后代入公式结果BPFO算出来是理论值的一半对着频谱找半天找不到峰。5.3 为什么理论频率和实际峰值总是差一点CWRU官方给出的特征频率是在理想运动学假设下算出来的。实际数据中由于轴承存在微小滑动、润滑膜厚度变化、负载波动等因素频谱中的峰值位置会偏离理论值。偏差量级通常在0.1%到1%之间。这一点很难完全消除但有几个办法可以减小偏差一是用实测转速代替铭牌转速。CWRU数据集提供了转速实测值但不同文件略有差异有的写1797有的写1796甚至1795。计算时最好读取数据集附带的转速信息而不是直接用1797这个典型值。我曾经用1797算出的BPFI是162.19Hz但某一组数据的实际转速是1795RPM对应的BPFI实际是161.96Hz差了0.23Hz。这个误差虽然不大但在设置窄带搜索窗口比如±0.5Hz时就可能影响峰值匹配结果。二是峰值搜索时用带宽搜索而不是点搜索。在包络谱里找故障特征频率时不要只看理论值那一个频率点而应该在理论值附近±1%的带宽内搜索局部峰值。对于早期微弱故障峰值可能不明显可以结合多个谐波位置联合判断。三是不要忽略频谱分辨率的影响。包络谱的频率分辨率等于1/T其中T是分析数据的总时长。如果我取1秒数据分辨率就是1Hz理论频率162.22Hz和实际峰值162.11Hz只差0.11Hz完全没问题。但如果只取0.1秒数据分辨率变成10Hz161.17Hz和163.17Hz之间所有频率都落在同一个bin里峰值位置就会大幅偏移。5.4 轴承特征频率计算的扩展用法掌握了基础计算之后这个工具还能往几个方向扩展。第一个方向是自动故障识别。把特征频率计算和包络谱峰值搜索串起来可以做成一个自动诊断脚本。输入一段振动信号和当前转速脚本自动计算四大特征频率在包络谱中搜索对应峰值根据峰值能量和谐波数量自动给出故障类型判断。这在批量分析CWRU数据时非常好用几百个文件几分钟就能跑完。第二个方向是边带分析的自动化。内圈故障的转频边带、滚动体故障的保持架边带都可以写成辅助判断逻辑。判断标准是主峰存在且边带幅值不低于主峰的某一比例比如10%则给出对应故障类型的加权置信度。第三个方向是时变转速工况下的频率跟踪。对于变转速设备瞬时转频是在变化的特征频率也随之变化。这时候简单套用固定转速算出的频率已经不够用了。需要先通过转速计或者从振动信号本身提取瞬时转频比如用阶比跟踪算法再根据每个时刻的瞬时转频动态计算特征频率。这个方向再往前走就是目前工业界比较关注的阶比分析和时频分析结合的技术路线算是这个领域比较前沿的玩法了。6. 最后一点个人体会我实际用下来的体会是这个zip包不只是一个计算工具更像是一把打开CWRU数据集的钥匙。理论公式背得再熟不亲手跑一遍数据永远不知道频谱里的边带长什么样、滚动体故障为什么那么难判断。把这个包里的代码吃透再把几个典型故障文件逐一分析一遍你对轴承故障诊断的理解会上一个台阶。最后再分享一个小技巧。分析CWRU数据时别只看单一文件建议同一类故障的不同故障直径、不同负载下各跑几组把结果放在一起对比。你会发现故障直径越大特征频率处能量越集中负载越大边带越丰富。这种纵向对比练出来的“手感”比单纯看十个教程都有用。后续做现场数据诊断时你会感谢自己当初没有跳过这一步。本文还有配套的精品资源点击获取