☰
轴承故障振动信号时频分析:STFT与CWT参数标定实战
2026/10/31 11:17:40 网站建设 项目流程

简介:围绕凯斯西储大学轴承故障数据集展开的信号时频分析文档,面向从事机械故障诊断、振动信号处理与设备健康监测的学习者和研究人员,帮助其理解正常、内圈、滚珠、外圈故障信号在时频域中的差异。包内共1个docx文件,约1.05MB,内容集中呈现数据集介绍、短时傅里叶变换与连续小波变换的分析流程及示例代码,便于按章节查阅和复现。已有1584人学习,说明该主题在故障诊断入门与实验参考中具有较高关注度。文档结合驱动端、风扇端和基座振动数据,对比STFT在重叠比例0.5、不同尺度下的时频图,并讨论morl、cmor1-1、cmor1.5-2、cgau8等小波函数及多尺度CWT效果,最终给出参数选择思路。读者可据此掌握时频图判读、尺度权衡和故障特征辨识方法,并参考代码搭建自己的分析流程。

1. 从一张时频图开始:为什么先做信号分析再谈故障分类

拿到一份轴承振动数据集,很多人的第一反应是直接上深度学习,把一段 1024 点信号丢进网络等一个准确率数字。我在做这套轴承故障数据集信号时频分析方案时走的是反过来的路:先用短时傅里叶变换和连续小波变换把原始振动变成一张能看懂的时间-频率图,让内圈、滚珠、外圈三种故障的冲击节奏摆到明面上,再谈分类。这条路慢,但它能回答一个绕不开的问题——你的模型到底在看什么。适合三类人:课程设计要交完整时频分析流程的、毕设开题需要先把数据摸熟的、以及已经在做故障诊断但被准确率卡住的。整套流程围绕一个公开的轴承故障基准数据集展开,从读文件、切窗、参数标定一路走到出图和验收。

2. 数据准备与通道取舍:三个加速度计位置和 1024 点滑窗的取法

2.1 实验台构成与三个测量位置的物理含义

这套基准数据的采集平台由一台 1.5 kW(2 马力)电机、一个扭矩传感器/译码器、一个功率测试计和一套电子控制器组成。电机是振动源,扭矩传感器夹在电机与负载之间负责记录转速和负载,功率测试计吸收负载并测量输出功率。真正决定你能分析什么的是加速度计的安装位置:驱动端、风扇端、基座,三个位置测的是同一台机器在不同传力路径上的响应,不是同一份数据的三个副本。

驱动端的加速度计贴着电机壳体靠近联轴器的一侧,信号里电机转子的旋转基频、传动系统的啮合成分和轴承故障冲击占主导,轴承早期缺陷在这条通道上最敏感。风扇端离轴承座远一些,风扇叶片的通过频率和气流扰动会混进来,风扇自身的不平衡、轴承异常更容易在这里显形。基座测的是整机振动,传递路径最长,高频冲击被结构和安装刚度滤掉一部分,剩下的是电机不平衡、底脚松动这类低频、全局性的问题。

理解这一点,选通道就有依据了:做内圈、滚珠、外圈的故障辨识,驱动端是首选;如果要做多通道融合或者做故障定位,才需要把风扇端和基座一起拉进来。我一般会先看驱动端,把结论跑通,再回头验证另外两路能不能提供互补信息,而不是一上来就三路拼一起。

提示:通道不是越多越好。三路信号采样点相同但物理含义不同,混在一起做特征会给后续分类引入不必要的维度,答辩时也很难解释每一维在干什么。

2.2 MAT 文件结构与通道键名的对应关系

数据以 MAT 格式分发,每个文件对应一种健康状态和一种负载工况。用loadmat读进来拿到的是一个字典,键名由编号、位置、信号类型三段组成,比如X097_DE_time里DE代表驱动端、time代表时域波形。这一步最容易翻车的地方是凭记忆写键名——不同故障件、不同负载下编号会变,写死一个键名跑通了一个文件,换一个就报 KeyError。

稳妥的做法是先打印键名核对,再进流程。下面这段代码把四类样本(正常、0.021 英寸内圈、滚珠、外圈故障)读进来,统一压平成一维并只保留同一段窗口,保证后面所有对比都建立在相同条件上。

import numpy as np from scipy.io import loadmat import matplotlib.pyplot as plt # 中文字体与负号,避免标题变成方块 plt.rcParams["font.family"] = "Microsoft YaHei" plt.rcParams["axes.unicode_minus"] = False # 四类样本:文件名 -> 驱动端时域波形键名 files = { "正常": ("0_0.mat", "X097_DE_time"), "内圈": ("21_1.mat", "X209_DE_time"), "滚珠": ("21_2.mat", "X222_DE_time"), "外圈": ("21_3.mat", "X234_DE_time"), } WIN = 1024 # 单样本长度,四类必须完全一致 signals = {} for label, (path, key) in files.items(): mat = loadmat(path) # 返回 dict,先 print(mat.keys()) 核对键名 if key not in mat: # 键名写错时立刻暴露,而不是带着空数据往下跑 raise KeyError(f"{path} 中不存在 {key},可用键:{list(mat.keys())}") sig = mat[key].reshape(-1) # 列向量压平成一维 signals[label] = sig[:WIN] # 只取前 1024 点,起点统一 # 原始波形四联图,先看冲击是否肉眼可辨 fig, axes = plt.subplots(2, 2, figsize=(18, 9)) for ax, (label, sig) in zip(axes.ravel(), signals.items()): ax.plot(sig, linewidth=0.8) ax.set_title(label) ax.set_xlabel("采样点") ax.set_ylabel("幅值") plt.tight_layout() plt.show()

reshape(-1)不能省,MAT 里存的是列向量,直接切片会得到二维形状,后面送进stft或cwt会报维度错误。WIN=1024是我在这个场景下习惯取的长度:驱动端常见采样率在 12 kHz 量级,1024 点大约覆盖几十毫秒,足够包含若干个完整的故障冲击周期,又不会长到把工况变化混进来。真正的采样率一定要以数据集附带的说明文件为准,这一点后面避坑章节会专门说。

2.3 1024 点滑窗与样本切分

只取前 1024 点做对比图没问题,但要形成可复现的分析结论,得把整段信号切成等长样本。滑窗的两个参数是窗长和步长,步长小于窗长就有重叠。做时频分析的可视化对比,我一般用 50% 重叠,这样相邻时频图之间的过渡连续,冲击不会正好卡在窗口边界被切掉一半。

切分时有三条硬规矩:四个类别的窗长必须相同;起点索引必须对齐,不能正常信号从头取、故障信号从中间取;每类保留的样本数尽量一致,否则后续任何统计对比都带着采样偏差。这三条听起来像常识,但我在评审别人的课设时,十份里有三份栽在第二条上——图的差异其实来自数据段不同,跟故障类型没关系。

样本切完要顺手记一张清单:类别、文件名、键名、窗长、重叠率、样本数。这张表在写报告时直接就是数据说明章节的骨架,答辩被问到"你用了多少样本"也不用现编。数据准备这一步做扎实,后面的 STFT 和 CWT 才是在同一把尺子上比较,否则调参数调到最后,你分不清是参数变好了还是数据换了。

3. 短时傅里叶变换的参数标定:窗口 16、32、64 的对比与 32 的由来

3.1 时间-频率折中为什么必须在数据上试

短时傅里叶变换的思路很直白:用一个中心位于时间 α 的窗函数 g(t−α) 在信号上滑动,在每个窗内做傅里叶变换,得到该时刻的局部频谱。窗越宽,窗内包含的周期数越多,频率分辨率越高,但时间上的定位就越模糊;窗越窄,时间定位越准,频率轴上的谱线糊成一片。这个折中不是靠公式推出来的,它取决于你的信号里故障冲击的间隔和你要抓的特征频率,所以必须在真实数据上试。

试的方法就是固定重叠率,只变窗口长度,把同一段内圈故障信号切成不同尺度的时频图并排看。原文的做法是重叠率取 0.5,比较 16、32、64 三个尺度,观察大尺度频率分辨率更高、小尺度时间分辨率更高,最后在四类信号的辨识度上权衡,选了 32。这个结论我不打算照抄,但方法值得照搬——参数不是选最漂亮的图,是选四类之间差异最稳定的那一档。

3.2 重叠率 0.5 与 nperseg/noverlap 的落地代码

scipy.signal.stft的接口里,nperseg是窗长,noverlap是重叠的采样点数而不是比例。很多人把 0.5 直接填进noverlap,得到的图看着没问题,实际上是 1024 点窗只重叠半个点,等于没重叠。正确写法是先按比例算样本数再传进去。

from scipy.signal import stft import numpy as np import matplotlib.pyplot as plt overlap_ratio = 0.5 # 重叠比例,全流程统一 target = signals["内圈"] # 固定用同一段数据做尺度对比 fig, axes = plt.subplots(2, 2, figsize=(18, 10), dpi=100) for ax, win in zip(axes.ravel(), (16, 32, 64, 128)): noverlap = int(win * overlap_ratio) # 比例换算成样本数 f, t, Z = stft(target, fs=12000, nperseg=win, noverlap=noverlap) ax.pcolormesh(t, f, np.abs(Z), shading="gouraud") ax.set_title(f"窗口 {win} - 内圈") # 标题由变量生成,避免和图不符 ax.set_xlabel("时间 / s") ax.set_ylabel("频率 / Hz") plt.tight_layout() plt.show()

fs填真实采样率,频率轴才是物理频率;如果只做相对比较,填 1.0 也能出图,但那样纵轴没有单位,报告里没法写。noverlap用int()向下取整,窗长为 16 时重叠 8 点,窗长为 1024 时重叠 512 点,比例始终是 0.5。shading="gouraud"让色块之间平滑过渡,看冲击条纹比默认的 flat 更清楚。

注意:标题里的窗口值一定要用 f-string 从变量取,不要手写。原文的图注写"尺度 16-内圈"而代码里window_size=32,这类不一致在答辩现场被指出来很尴尬。

3.3 尺度对比实验的记录表与选型结论

光看图容易凭印象下结论,把观察结果记成表,选型才有依据。下面是我按这个场景整理的判定维度,数值是定性描述,你换成自己数据跑一遍填进去即可。

窗口长度时间分辨率频率分辨率内圈冲击条纹滚珠与正常区分度计算量
16高,冲击起止清晰低,谱线粘连条纹细密但边界糊一般最小
32中中条纹与谱线均可辨较好小
64较低较高条纹变粗,拖尾明显好中
128低,冲击被抹平高细节丢失一般大

从这张表能看出为什么 32 是个不错的落点:时间分辨率还没有差到看不出冲击间隔,频率分辨率又足以把故障特征频率从转频谐波里分出来。窗口再小,滚珠故障那种幅值低、间隔不规则的冲击容易被噪声淹没;窗口再大,冲击的时间位置被平滑成一片,反而看不出周期性。

选完参数还有一步收尾:把 32 这个结论在四类数据上各跑一遍,确认正常、内圈、滚珠、外圈四张图两两之间都有肉眼可辨的差异。如果某两类图像几乎一样,说明单靠 STFT 的时频图不足以支撑分类,这时候要么换分析手段(下一章的 CWT),要么在特征层面补统计量,而不是继续微调窗口长度硬凑。

4. 连续小波变换的小波选型:morl、cmor1-1、cmor1.5-2、cgau8 与尺度序列

4.1 复数小波为何比实小波更适合振动信号

连续小波变换用一族由母小波伸缩平移得到的函数去匹配信号,得到的是时间-尺度平面上的系数。振动信号的特点是冲击衰减快、频带分布宽,用实小波分析时系数会出现正负振荡,取模之后相位信息全丢,幅值图上经常看到正负交替的条纹,读图的人要靠经验判断哪个是真正的冲击包络。

复数小波在这件事上有天然优势:它同时给出实部和虚部,取模得到的是包络幅值,相位随时间的变化也保留下来。Morlet 类小波在时域和频域都有较好的局部化性质,处理振动信号时域波形和频谱形状都接近,所以这类母小波在轴承故障分析里出镜率很高。cmor是复数 Morlet 的变体,名字里的两个数字分别对应中心频率和带宽参数,调这两个数就是在调小波的时频形状。

原文列出了四种候选:morl(Morlet)、cmor1-1、cmor1.5-2、cgau8(复数高斯小波)。其中cgau系列适合近似高斯型信号,cmor1.5-2的带宽参数更大,时域上更紧凑,对瞬态冲击的响应更尖锐。这四种小波在同一段内圈故障信号上、尺度都设 128 的条件下做对比,结论是cgau8和cmor1.5-2对不同故障的辨识度都较高,再综合考虑频率与时间分辨率,选cmor1.5-2继续往下做。

4.2 中心频率与尺度序列 cparam 的构造

小波变换里真正决定分析频段的是尺度序列,而尺度序列不是随便给一个等差数列就行。常用的做法是先取母小波的中心频率fc,再按cparam = 2 * fc * totalscal算出常数项,最后用cparam / np.arange(totalscal, 0, -1)生成从大到小排列的尺度序列。这样生成的尺度在对数意义上分布更均匀,低频段尺度密、高频段尺度疏,符合小波"低频看细节、高频看轮廓"的特性。

import pywt import numpy as np totalscal = 128 # 尺度总数 wavenames = ["cgau8", "morl", "cmor1-1", "cmor1.5-2"] scales_dict = {} for name in wavenames: fc = pywt.central_frequency(name) # 母小波中心频率 cparam = 2 * fc * totalscal # 尺度序列常数项 scales_dict[name] = cparam / np.arange(totalscal, 0, -1) print(name, "fc =", round(fc, 4), "首尾尺度 =", round(scales_dict[name][0], 2), round(scales_dict[name][-1], 2))

fc这一行千万别漏。原文代码里morl那一段把fc2 = pywt.central_frequency(wavename2)注释掉了,后面又拿cparam2 = 2 * fc2 * totalscal去算,跑起来直接 NameError;就算上一轮循环里残留了同名变量,尺度序列也是错的,图能出来但频段完全对不上。这类错误不会报错得很难看,它会安安静静给你一张错的图,这是最麻烦的一种。

totalscal决定尺度序列的密度,不是直接决定分辨率。128 意味着在关心的频段上切 128 份,想看得更细可以加到 256,代价是计算量翻倍、图上的等高线更密反而不易读。

4.3 四种小波与四档尺度的对比实验

小波函数定下来之后,还要再扫一遍尺度。原文固定cmor1.5-2,在滚珠故障数据上比较 32、64、128、256 四档,得到的结论是尺度越大越关注低频特征,尺度越小越关注高频特征,这正好印证了小波变换的双尺度性质。落到实操就是一句话:想抓高频冲击用较小尺度,想看低频调制和转频成分用较大尺度。

import matplotlib.pyplot as plt fs = 12000 # 真实采样率,频率轴才有物理意义 sampling_period = 1.0 / fs name = "cmor1.5-2" target = signals["滚珠"] base = scales_dict[name] # 以 128 个尺度为基准 fig, axes = plt.subplots(2, 2, figsize=(18, 10), dpi=100) for ax, n in zip(axes.ravel(), (32, 64, 128, 256)): # 按目标尺度数重新采样尺度序列,保持 cparam 不变 scales = base[::max(1, len(base) // n)][:n] coef, freqs = pywt.cwt(target, scales, name, sampling_period=sampling_period) amp = np.abs(coef) # 复数系数取模得到包络幅值 ax.contourf(np.arange(len(target)) / fs, freqs, amp, levels=40, cmap="jet") ax.set_title(f"{name} 尺度 {n} - 滚珠") ax.set_xlabel("时间 / s") ax.set_ylabel("频率 / Hz") plt.tight_layout() plt.show()

sampling_period决定freqs的换算,填错的话纵轴频率整体缩放,图还是那张图,数全错。pywt.cwt在不同版本里这个参数有的是位置参数、有的支持关键字,如果报 unexpected keyword argument,就改成第三个位置参数传。contourf的levels控制等高线密度,40 左右在冲击特征和噪声底之间平衡得比较好,太少会把弱冲击抹平,太多满屏都是噪声纹理。

四档尺度扫完,我的取舍逻辑是这样的:32 档能看出冲击的密集程度,但故障之间的差异主要体现在高频段的能量分布;256 档低频结构清楚,代价是冲击的时间位置被拉宽。128 档在两边的折中比较均衡,配合cmor1.5-2的时域紧支撑,四类故障在包络幅值图上分得开,所以我沿用 128 作为默认尺度,只有在专门分析高频共振带时才降到 32 或 64。

5. 避坑与排查:时频分析里最容易翻车的五个点

5.1 四类样本取的不是同一段,图上的差异全是假的

现象:正常和三种故障的时频图看着差异明显,但换一段数据重跑,差异就没了。

原因:正常信号从第 0 点取,故障信号从中间某点取,或者窗长不一致。振动信号在不同时间段幅值本来就有波动,你看到的"故障特征"其实是数据段本身的差异。

解决:把起点索引和窗长写成常量,四类共用;切完样本后统一做一次归一化或统一记录均方根值,任何一类明显偏离就回头查。我习惯在样本清单里加一列起始索引,出问题直接对表,不用重新读数据猜。

5.2 采样率被当成窗长填进去,频率轴失去物理意义

现象:时频图纵轴数字很大或很小,和已知的转频、故障特征频率对不上。

原因:把 1024 这个窗长当成了采样频率传给fs或sampling_period。原文的 CWT 代码里sampling_rate = 1024,在这份脚本里它其实是窗口长度,不是采集卡的采样率。

解决:采样率只从数据集说明文件取,一次写进配置文件,全流程引用同一个变量。填完做一次自检:正常信号的时频图上,转频应该落在预期位置,如果差一个数量级,说明采样率填错了。

5.3 小波函数名写错或中心频率漏算,图能出来但结论错

现象:换一个小波函数,图看起来差不多,但尺度与频率的对应关系变了;或者直接报 NameError。

原因:cmor1.5-2这类带参数的小波名,数字位置、分隔符写错会被当成另一个小波;fc漏算则cparam用的是上一轮的残留值。

解决:把「算中心频率 → 算 cparam → 生成尺度序列 → 做变换」封装成一个函数,每次调用都从零算起,不依赖循环外变量。写完后打印每个小波的首尾尺度做交叉核对。

5.4 参数和图注不同步,答辩时被当场戳穿

现象:图注写"尺度 16",代码里window_size = 32;或者四张子图的标题是手写死的一串固定文字。

原因:改代码时只改了参数没改标题。这在赶进度的课设里非常常见。

解决:所有标题、文件名、图注都用变量拼接,禁止手写。出图后抽一张反查:把参数打印到终端,和图上标题逐字比对,不一致就说明有硬编码残留。

5.5 图好看但分不开类,问题出在可视化而不是分析

现象:时频图配色漂亮,看上去很有"专业感",但四类之间没有稳定差异,调cmap、调对比度都没用。

原因:一直在调可视化参数,没有动分析参数。动态范围不统一时,同一套参数下正常信号的噪声底会被拉亮,故障特征被压暗,看起来就像"分不开"。

解决:固定色标的取值范围,比如统一用分位数归一化到 0 到 1;所有图用同一个cmap和同一个levels;先在灰度图上确认结构差异存在,再上配色。配色是给人看的,分析结论必须在灰度图上就成立。

6. 验收与答辩:怎么证明这套时频方案经得起追问

这套东西做完,交出去之前我自己会过一遍验收清单,四项全过才算收工。清单里的"留证据"那一列很关键——它决定了你答辩时是翻笔记还是翻图。

检查项合格判据留证据的方式
数据一致性四类样本窗长、起点、样本数完全一致样本清单表,含文件名、键名、起始索引
STFT 参数窗口长度有对比实验支撑,不是拍脑袋定16/32/64/128 四联图 + 判定表
CWT 参数小波函数与尺度都有对比记录,cparam可复算四种小波四联图 + 尺度四联图
结论可复现换一台机器、换一个随机种子,图一致固定参数脚本 + 采样率常量

答辩现场最常被问的三个问题,回答方式其实都藏在前面的实验里。第一个是"为什么 STFT 窗口选 32 不选 64"——不要答"试出来的",要答"64 的频率分辨率更高,但内圈冲击在时间上被拉宽,滚珠故障那种弱冲击的间隔读不出来,32 是两类特征都能保住的最小窗口"。第二个是"CWT 为什么不用 morl 而用 cmor1.5-2"——答"cgau8 和 cmor1.5-2 在四类上的辨识度接近,但 cmor1.5-2 的带宽参数更大,时域更紧凑,冲击包络更尖锐,配上 128 尺度在低频调制上更清楚"。第三个最狠:"你的时频图能证明什么"——这时候要能指着图说出具体现象,比如外圈故障的冲击间隔对应哪个特征频率、内圈的调制边带落在哪个频段,而不是说"看起来不一样"。

有个血泪经验:我第一次交这套分析的时候,四张图全是手工调过对比度的,单独看每张都很清楚,并排放在一起就发现正常信号的噪声底比故障信号亮一截,被判了"图的差异可能来自显示设置"。从那以后我每次出图都强制走一遍统一色标范围,先出灰度图确认结构,再上配色,标题一律变量拼接。这套流程看着啰嗦,但它让每一张图都经得起"你把这张图换台机器再跑一遍"的追问。希望帮到你。

本文还有配套的精品资源,点击获取

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询