当初在一台离心泵的轴承座上采集振动数据时域波形像一坨毛线包络谱里只有转频和几根看不出名堂的谱线。拆检之后才发现外圈已经有一道明显的剥落。问题出在哪不是传感器装错了而是处理思路不对——传统频谱分析对瞬态冲击天然不友好轴承故障产生的那种“哒哒哒”的周期冲击能量被摊在了几千赫兹的宽带范围里单根谱线幅值很小很容易被噪声淹没。后来我把思路从“频域滤波”换成了“稀疏表示”用学习字典做原子匹配、用加权稀疏系数做成分筛选重新处理同一段数据故障特征频率在包络谱里直接就立起来了。这篇文章把我踩过的坑、每个环节为什么这样设计、参数怎么调、代码怎么组织完整拆开讲一遍。适合正在做振动信号处理、设备故障诊断方向或者刚接触稀疏表示但不知道怎么落地到轴承场景的读者。1. 为什么故障冲击在数据里“看不见”信号模型与稀疏表示的基本逻辑1.1 轴承故障信号的构成冲击、噪声、转频三股势力的角力要搞清楚稀疏表示为什么能解决轴承诊断问题先要看清振动信号里到底有什么。传感器贴在轴承座上采集到的信号基本可以拆成三个部分周期性故障冲击成分 s(t)滚动体滚过缺陷时产生的瞬态冲击周期由转速和轴承几何参数决定是我们要找的东西确定性周期成分 c(t)转频及其谐波、齿轮啮合频率、皮带通过频率等它们严格同步于旋转周期噪声与随机干扰 n(t)环境振动、流体噪声、随机脉冲、传感器电噪声等。轴承故障冲击在时域上只占据很短的持续时间占空比很低大部分时间段其实是“安静”的。这种信号在频域里很难看——冲击包含大量频率成分能量弥散到几千赫兹的带宽中单根FFT谱线的幅值往往低于噪声底。而传统做法是先用带通滤波器选出共振带再做包络解调这本身体现了一个思想故障信息不是表现在某几条特征谱线上而是集中在一个瞬态波形里。稀疏表示的底层假设恰好和这个物理事实吻合。它认为信号可以被一组“原子”的少量叠加逼近绝大多数系数接近零。故障冲击是瞬态的、局域的、间歇的正是一种高度稀疏的信号。我们只要让算法从一堆候选波形里挑出最像冲击的那几个原子就可以把故障成分从背景噪声中单独剥离出来。1.2 稀疏表示为什么能“挑”出故障冲击一个抓壮丁的比喻可以把稀疏表示理解成抓壮丁。字典 D 就是一份候选名单里面有正弦波、脉冲、衰减振荡等各种波形的原子。我们手里的振动信号 x 是最终需要凑齐的人头数。普通的最小二乘求解会让每个人都出一份力平均分摊误差结果是谁都不像、谁都有份而稀疏求解要求尽量用最少的壮丁凑齐目标剩下的误差再大也不用每个人去补。匹配追踪类的算法本质就是贪心选择第一轮挑与残差相关性最强的原子第二轮在剔除它贡献的信号以后继续挑下一个。故障冲击因为是突发的瞬态波形总有一个原子与它的相似度极高第一轮就能被选中。而噪声和转频成分没有这种“一击命中”的特性要么是频谱摊得太开要么是波形太平缓难以用单个原子高效表示。对轴承诊断而言这个过程有点像“定向抓取”——把冲击波形直接从混合信号里抠出来。但要说明一个关键前提稀疏表示只在字典与信号高度匹配的时候才有效。如果字典里的原子都是纯正弦波想精确表示一个带衰减的窄带冲击需要十几个原子一起堆稀疏性就无从谈起。这也是为什么后面要引入学习字典。1.3 为什么光靠稀疏还不够普通稀疏的三大软肋我最初直接用标准OMP配合DCT字典处理数据效果有一部分但离“能用”还差得远。普通稀疏求解在本场景下有三个明显软肋第一对噪声不设防。L1范数的凸松弛虽然在数学上保证了稀疏解但实际噪声里总有些片段的能量比较高会被算法当成“有用成分”选进字典系数里。尤其在低信噪比下稀疏解会逐渐退化成对噪声的过拟合。第二对随机冲击没有辨识能力。现场振动信号里经常有随机碰磨、流体脉动造成的非周期冲击它们在时域上和故障冲击长得几乎一样。普通稀疏法没办法区分“有规律的冲击串”和“随机出现的单次冲击”结果会把两者一起重构出来包络谱上多出一堆假谱线。第三对故障冲击的幅度变化不敏感。实际故障从萌生到发展到后期冲击幅值忽大忽小部分小幅值冲击含有重要的早期损伤信息。统一惩罚所有系数会让算法倾向于保留大幅度成分把早期微弱冲击直接当成噪声丢掉。这三点让我意识到单纯的稀疏表示不是终点需要给每个系数额外“加权重”——让算法知道哪些成分应该优先保留哪些成分应该重点压制。2. 加权稀疏是怎么“加权”的两类权重来源与具体设计规则2.1 权重从哪来迭代重加权与物理先验两条路线加权稀疏的核心是把原来求解 min ||α||₁ 的问题改造成 min ||w⊙α||₁其中 w 是一个与系数一一对应的权重向量。权重大的系数被压得更狠权重小的系数更容易被保留。关键问题在于权重 w 怎么设定。我实际用过两条路线各有适用场景。第一条是迭代重加权Iterative Reweighted L1。先跑一次普通稀疏求解得到初始系数 α⁽⁰⁾然后按公式 wᵢ 1/(|αᵢ| ε₀) 更新权重。上一轮系数越大的位置权重就被调得越小下一轮这些成分被保留的代价更低上一轮基本为零的位置权重约等于 1/ε₀数值很大下一轮基本不会被选进去。这样反复迭代三到五次能量会越来越集中到少数几个真正的冲击分量上。ε₀ 是个防除零的小常数太小会数值不稳太大则权重区分度不足。第二条是物理先验加权。滚动轴承的故障特征频率可以按转速和几何参数算出来比如外圈故障特征频率 BPFO、内圈 BPFI、滚动体 BSG、保持架 FTF。既然知道了故障冲击应该出现在哪些频率位置在字典设计时就把这些频率对应原子的权重降低让算法在相同的惩罚代价下优先选择它们。更进一步的变体是“谐波增强”不仅对基频原子降权对二倍频、三倍频、分数谐波附近一小段频带内的原子也做不同程度的降权。2.2 两类权重怎么结合先物理后数据的两阶段策略两条路线不是二选一我推荐两阶段结合。第一阶段先用物理先验加权跑一次稀疏分解得到初步系数目的不是直接出结果而是得到一个“大致在故障频率上有响应”的基准解第二阶段用迭代重加权把上一轮的系数分布作为权重更新的起点。这样既有物理约束兜底——如果转频计算不准或轴承型号搞错了物理先验会将错就错数据分析权重能拉回来一点又有数据自适应能力——真实故障会有许多非理想的边带和谐波单纯靠先验设定覆盖不全。权重向量还需要定期归一化。如果不归一化w 的尺度直接乘到系数上实际效果等价于改变了正则化强度 λ会让算法要么过度稀疏、要么过度放松。我的做法是每轮迭代后让 w 的均值保持在 1 附近这样 λ 的语义始终一致调参时不会因为权重尺度变化而反复试错。2.3 加权对求解流程的影响一个容易被忽略的收敛问题加权稀疏不是简单地把 OMP 的停止条件改一改就行。标准 OMP 每轮选择与残差相关性最强的原子这个逻辑本身不包含权重。要在 OMP 里融入权重常见做法是把问题改造成加权形式把原信号 x 和字典 D 都乘上权重矩阵的平方根然后用标准 OMP 求解变换后的问题最后再把系数映射回去。具体过程是设 W diag(√w)求解 min ||x̃ - D̃α||₂² λ||α||₁其中 x̃ WxD̃ WD得到系数 α 后真正对应的原子权重已经内化到变换里了。实测下来这个处理方式稳定不会破坏 OMP 的选原子逻辑。但我要强调一个收敛上的细节迭代重加权的第 2 轮起如果字典 D̃ 的条件数变差OMP 选原子时会出现相邻轮次来回切换同一个原子组的情况。此时不要盲目加大迭代次数建议每轮迭代后计算重构信号与原始信号的信噪比连续两轮信噪比提升不足 0.5dB 就停。跑过一组对比固定迭代 20 轮和自适应终止的结果几乎一样但自适应终止省了近一半计算时间。3. 学习字典从“选模板”到“学模板”的关键一步3.1 解析字典 vs 学习字典为什么正弦基在这里吃亏前面说 DCT 字典表示冲击需要堆很多原子根源在于解析字典是“通用模板”不是为轴承冲击设计的。DCT 的原子是不同频率的余弦波本质假设是信号由平稳周期成分组成小波包字典稍好一些能描述瞬态但小波基的形状是固定的小波函数和真实故障冲击的衰减振荡形态不一定吻合。学习字典的思路完全不同直接从训练数据里学出一组原子让原子本身逼近待处理信号的典型形态。对轴承故障来说训练样本是从振动信号里切出来的一段段时域波形其中包含大量正常振动、噪声和故障冲击的混合体。学出来的字典原子往往直接就是一个带指数衰减的振荡波形——这不就是故障冲击的天然模板吗使用学习字典还有另一层好处原子形态与故障机制相关。我经常把训练好的字典原子画出来看外圈故障训练出来的原子通常集中在某个高频共振带附近而内圈故障多了明显的调制包络形状。这等于在诊断前先白拿了一层特征可视化信息。3.2 K-SVD 的基本原理与实现要点K-SVD 是目前最常用的字典学习方法名字拆开看K 指字典原子个数SVD 指用奇异值分解更新字典。整个流程是两阶段交替优化稀疏编码阶段固定当前字典用 OMP 对每个训练样本求解稀疏系数字典更新阶段逐个原子更新。更新第 j 个原子时找出所有用到了这个原子的样本计算在这些样本上的残差对残差矩阵做 SVD 分解用最大的左奇异向量作为新原子、对应奇异值更新系数。这个“逐个更新”的设计非常巧妙。它保证每次更新都针对某个原子的专属残差不会像同时更新所有原子那样互相干扰。实际实现里需要特别留两个心一是如果一个原子长期没有被任何样本选中要重新初始化否则字典里会堆积一堆“死原子”二是初始化方式不要用随机矩阵直接从原始信号里随机抽几个训练样本作为初始原子收敛速度会快很多。3.3 分帧、标准化、原子尺度三个容易翻车的细节学习字典不像普通机器学习那样把特征喂进去就行信号处理里的预处理直接决定字典学到的是“故障特征”还是“噪声特征”。第一个细节是分帧。轴承信号很长不能整段丢进 K-SVD。需要按原子长度 N 切帧相邻帧之间有 50% 到 75% 的重叠。重叠率太低会丢掉冲击的相位关系重叠率太高则训练样本冗余度大计算量成倍上升收益不明显。我常用 50% 重叠N 取 128 到 512 个采样点。帧长如果小于冲击衰减时间的一半原子学出来的是“残缺波形”大于冲击间隔原子会包含两三个冲击反而模糊了单个冲击的形态。第二个细节是标准化。每帧必须先减掉均值再做功率归一化。帧不去均值字典里会学出直流分量原子完全没有物理意义不做归一化大幅值帧会主导字典更新小幅值冲击帧的形态信息被压制。第三个细节是原子尺度要与包络解调兼容。原子长度直接决定了后续重构信号能分辨的频率精度。假设采样率 25.6kHz原子长度 256那么单个原子的有效频带宽度约 100Hz。如果故障特征频率的边带间隔小于这个值重构后边带会被抹平。我一般会让原子长度对应的时间至少覆盖个冲击周期长度保证频域分辨率够用。4. 从数据到结论完整诊断流程与代码结构4.1 全流程的六个环节把上述方法串起来一个完整可用的诊断流程如下数据采集与窥窗检查。确认采样率、采集时长、转速信息剔除明显断线、削顶、直流偏置异常的数据段分帧与标准化。按原子长度切帧去均值、归一化功率生成训练样本集字典学习。用 K-SVD 学出字典原子个数按经验公式配置见下文参数表加权稀疏编码。结合物理先验初值和迭代重加权策略对整段信号做稀疏分解重构故障冲击成分。将加权筛选后的系数与原子的线性组合重建为纯净冲击信号包络谱分析与故障定位。对重构信号做 Hilbert 变换取包络再计算 FFT在包络谱中搜索 BPFO、BPFI、BSF、FTF 及其谐波。第六步不是可有可无的验证而是真正的诊断结论输出。稀疏表示本身只是把故障成分从混合信号里分离出来判断“轴承到底坏没坏、坏在哪个元件”最终还是要靠包络谱里的特征频率与理论计算值对照。4.2 Python 代码骨架与核心参数我用的是 Python 实现核心三块字典学习类、加权稀疏编码函数、包络谱提取函数。骨架如下import numpy as np from scipy import signal def standardize(frames): # 去均值 功率归一化 frames frames - frames.mean(axis1, keepdimsTrue) norm np.linalg.norm(frames, axis1, keepdimsTrue) return frames / np.maximum(norm, 1e-12) def ksvd_learn(frames, m128, k8, iters10): # frames: (N_frames, N_atom) n frames.shape[1] D frames[np.random.choice(frames.shape[0], m, replaceFalse)].T X frames.T # n x num_frames for _ in range(iters): # 稀疏编码阶段 coef omp_batch(X, D, k) # 字典更新阶段 for j in range(m): idx np.nonzero(coef[j, :])[0] if len(idx) 0: D[:, j] X[:, np.random.randint(X.shape[1])] continue R X[:, idx] - D coef[:, idx] D[:, j:j1] * coef[j, j:j1] U, S, Vt np.linalg.svd(R) D[:, j] U[:, 0] coef[j, idx] S[0] * Vt[0, :] return D def weighted_omp(x, D, w, k8, eps01e-6): # 加权OMP利用矩阵变换把权重内化 W np.diag(np.sqrt(w)) Dt W D xt W x alpha omp(xt, Dt, k) # 标准OMP return alpha def envelope_spectrum(sig, fs): analytic signal.hilbert(sig) env np.abs(analytic) spec np.abs(np.fft.rfft(env)) freq np.fft.rfftfreq(len(env), 1/fs) return freq, spec实际工程中我会把ksvd_learn里的 OMP 和weighted_omp里的 OMP 统一封装避免重复实现。上面的结构只是为了把三个模块拆开方便理解。4.3 参数设置参考表与调试经验我整理了一张参数参考表方向和量级可以直接复用。参数含义推荐值说明fs采样率12.8kHz ~ 51.2kHz需要覆盖轴承共振带至少含 3~5 倍故障特征频率N原子长度128 ~ 512 点覆盖 3 个以上冲击周期太短丢形态太长丢精度overlap分帧重叠率50%保证冲击跨帧时不丢失m字典原子数64 ~ 256过小表达力不够过大会学过噪声k稀疏度3 ~ 8每帧最多用 k 个原子iterK-SVD 迭代次数8 ~ 15观察重构误差收敛就停ε₀重加权稳定常数1e-4 ~ 1e-6太小数值不稳太大权重不均匀λ稀疏惩罚系数0.05~0.2 倍的 max(Dᵀx)用相对值统一调参不用绝对值调试中最常见的现象是包络谱里故障特征频率不出来排查顺序我建议是先看分帧参数有没有覆盖冲击长度再看稀疏度是否设得太小导致冲击被拆散最后检查加权在迭代过程中权重是否发生了漂移。我自己踩过最深的坑是把 K-SVD 的迭代次数开到了 40结果字典慢慢把噪声也学成了原子重构信号包络谱里到处都是假峰降回 12 次迭代后才恢复正常。5. 实测效果三组方法对比能看出什么5.1 对比设计普通 OMP、加权稀疏、学习字典加权为了评估各环节的实际贡献我在同一组现场数据上跑了三套方案。数据来自一套齿轮箱测试台架上的滚动轴承采样率 25.6kHz采集时长约 10 秒转速约 1500rpm根据轴承几何参数估算外圈故障特征频率 BPFO 约为 103.7Hz。三套方案分别是方案 A固定 DCT 字典 标准 OMP这是许多论文的基线做法方案 B固定 DCT 字典 加权稀疏物理先验 迭代重加权方案 CK-SVD 学习字典 加权稀疏也就是本文完整流程。三套方案都输出重构信号的包络谱我对比的是三个指标BPFO 谱峰幅值、包络谱信噪比谱峰幅值与周边平均噪声之比、可识别的谐波数量。5.2 包络谱结果解读结果差距非常明显。方案 A 的包络谱在 103.7Hz 处确实有谱峰但旁边还有一组幅值更高的干扰峰如果不带先验知识很容易把干扰峰当故障特征误诊率不低。方案 B 的 103.7Hz 谱峰成了全谱最高峰二倍频 207.4Hz 也开始显现包络谱信噪比相比方案 A 提升了约 8dB。方案 C 的包络谱最干净103.7、207.4、311.1 三根谱线清晰可见整体噪声底比方案 A 低了一个量级以上。方案 C 优势的来源值得多说一句。K-SVD 学到的字典原子是一组衰减振荡波形它们的中心频率恰好落在轴承外圈共振带附近加权稀疏在此基础上选择与故障冲击形态最匹配的原子重构信号几乎等于“把故障冲击直接滤了出来”。相比之下方案 B 用 DCT 字典虽然有物理先验加权帮衬但正弦原子描述瞬态冲击的效率太低同样的稀疏度下残差里的故障能量还是不少总有一部分冲击信息留在噪声里出不来。我后来把方案 C 的字典原子可视化发现两件有意思的事一是大部分原子不是单一频率的正弦波而是带有明显衰减包络的窄带振荡和理论上的冲击响应模型完全吻合二是不同原子之间的中心频率差大约几十赫兹相当于 K-SVD 自动做了一个共振带内的频率细分。这种细分是人工带通滤波很难预设的因为滤波器的频带划分是刚性的而信号里的冲击波形是连续变化的。5.3 两个很容易误判的场景第一个场景是内圈故障。内圈故障信号有强烈的调制效应——故障点随内圈旋转载荷区周期性变化冲击幅值被低频调制。如果直接套用上述流程偶尔会发现加权迭代把调制包络的低频成分当成主成分导致重构信号偏向包络而非冲击。我的处理是把帧长加大到覆盖多个调制周期让调制包络作为一个整体被稀疏表示出来再用包络谱定位时就不会丢失调制带来的边带信息。第二个场景是严重变转速。转速波动会让故障特征频率不再是常数而是随时间漂移的曲线。此时把物理先验加权直接用在整段信号上会出问题——先验设定的频率位置已经过时。我的方式是先将信号按转速脉动切成多个准平稳段每段单独计算特征频率和权重。虽然计算量变大了但包络谱里的谱峰不再模糊成一团。还有一个隐藏很深的问题随机冲击。如果数据里有非周期性的随机碰磨加权稀疏偶尔会把单次冲击保留下来。它们在包络谱上不形成整倍数谐波序列就像一根孤立的谱线。识别这种假峰的办法是检查谐波关系故障特征频率一定会伴随二倍频、三倍频而随机冲击孤立出现、不构成序列。我在实际判断时有一条硬性规则——只有基频和高次谐波同时存在且间隔严格等于基频才判为轴承故障。我个人在实际调这套流程时觉得最花时间的不是 K-SVD 的迭代而是加权策略里 ε₀ 和初始权重这两处的小参数。多组实验之后形成的经验是先用解析字典 物理先验跑一遍确认包络谱里存在候选特征频率再切换成学习字典 迭代重加权放大特征由粗到细分两步走。另外一个小技巧是把每轮训练好的字典原子保存成图片观察原子是否逐渐向噪声形态退化这个视觉信号比数值指标更早提醒你参数出了问题。整套方法在常规轴承故障诊断场景里稳定可用如果再往下走把字典学习做成跨工况迁移是本项目最自然的扩展方向。
