简介这份资源是发表于《能源与环保》的期刊论文《基于小波阈值的通风机故障特征提取应用研究》面向煤矿机电、通风安全及信号处理方向的工程技术人员与研究人员用于解决主通风机振动信号噪声干扰大、故障特征难以准确提取的问题。资源包为单一PDF文件共1个文件大小约684KB内容完整呈现论文正文、公式推导与图表。文中系统梳理了小波阈值降噪的基本思路与阈值选取原则涵盖通用阈值、Stein无偏风险阈值、试探法阈值及最小最大准则阈值等方法并对比了小波分析与傅里叶变换在去噪能力上的差异同时涉及LMD算法在非平稳信号分解中的应用。读者可借此掌握通风机振动信号降噪与特征提取的理论依据理解软、硬阈值策略的取舍为轴承磨损、不平衡等故障的诊断与识别提供参考。目前已有147人学习。1. 通风机振动信号里为什么小波阈值成了特征提取的默认选项通风机故障诊断这件事真正难的不是判断有没有故障而是从振动信号里把故障特征抠出来。轴承内圈剥落、叶片积灰不平衡、联轴器不对中这些故障在早期阶段产生的冲击成分往往只有几十毫伏淹没在电机电磁噪声、气流脉动和结构共振里。直接做 FFT频谱上全是底噪故障特征频率那根谱线根本立不起来。小波阈值之所以在这个场景里被反复使用核心原因是它同时具备两个能力多分辨率分解把信号拆到不同频带阈值收缩把噪声系数压掉而保留冲击成分。相比经验模态分解EMD容易模态混叠、相比单纯带通滤波需要预先知道故障频带小波阈值对通风机这种转速波动 宽频噪声 非平稳冲击的信号适应性更好。这套方法适合谁做旋转机械状态监测的工程师、想把实验室算法落到现场采集系统里的开发者、以及需要给诊断模型喂干净特征向量的算法同学。下面从原理到代码到踩坑把这条链路走一遍。2. 小波阈值去噪的数学骨架与通风机场景的参数映射2.1 从连续小波到离散分解通风机信号该选几层连续小波变换在工程落地里很少直接用因为冗余度太高、计算量大。实际用的是离散小波变换DWT通过 Mallat 塔式算法把信号逐层分解为近似系数低频和细节系数高频。通风机振动信号的采样率通常取 10 kHz 到 20 kHz关注频带一般在 0 到 2 kHz。分解层数的经验公式是L floor(log2(fs / (2 * f_interest)))假设采样率 fs 12800 Hz关心最高 1600 Hz 的故障频带那么 L floor(log2(12800/3200)) 2。但实际工程里我一般会多分一层取 3 到 4 层因为通风机的冲击成分有时会调制到较高频带多一层留余量。分解层数选多了会怎样低频段频率分辨率提高但每层的数据点数减半边界效应占比增大重构后两端容易出现畸变。选少了则高频噪声和故障冲击分不开。这是第一个需要权衡的参数。2.2 阈值函数硬阈值、软阈值和半软阈值在冲击保留上的差异阈值函数决定了系数怎么被收缩。三种常见形式硬阈值系数绝对值小于阈值直接置零大于阈值原样保留。优点是保留冲击幅值缺点是阈值附近系数突变重构信号容易产生伪吉布斯振荡。软阈值小于阈值的置零大于阈值的向零收缩。重构信号平滑但冲击幅值被系统性压低对早期微弱故障不利。半软阈值在软硬之间折中引入第二个阈值参数介于两者之间平滑过渡。对通风机故障特征提取我的经验是如果后续要做包络谱分析找故障特征频率用硬阈值或半软阈值因为幅值保真度直接影响包络谱的谱峰高度如果只是做趋势监测看能量变化软阈值够用且更稳定。2.3 阈值选取规则通用阈值、极小极大准则和启发式阈值阈值大小直接决定去噪力度。常用规则规则公式核心适用场景通风机场景建议通用阈值σ√(2lnN)高斯白噪声噪声近似高斯时可用极小极大准则最小化最大风险信号长度短数据点少于 2000 时优先启发式阈值通用与无偏风险取小信噪比未知现场信号首选无偏风险估计基于SURE信噪比中等配合启发式使用其中 σ 是噪声标准差估计通常取第一层细节系数的中位数绝对偏差除以 0.6745。这个估计方式对冲击噪声比标准差更鲁棒因为中位数不受少数大冲击影响。现场信号信噪比未知我一般默认用启发式阈值heursure它在高信噪比时退化为通用阈值低信噪比时切换到无偏风险估计自适应性强。3. 用 Python 跑通通风机振动信号的小波阈值去噪全流程3.1 环境准备与依赖安装需要 PyWavelets 做小波变换NumPy 做数值计算SciPy 做频谱和包络分析Matplotlib 做可视化。pip install PyWavelets numpy scipy matplotlib版本方面PyWavelets 1.x 和 2.x 的 API 基本兼容但 2.x 对边界处理模式有调整如果复现老代码发现结果对不上先检查pywt.__version__。3.2 构造模拟通风机振动信号验证算法真实数据不好拿先用模拟信号验证链路。构造一个含轴承外圈故障冲击的通风机振动信号基频 25 Hz1500 rpm故障特征频率 120 Hz叠加高斯噪声。import numpy as np import pywt def simulate_fan_vibration(fs12800, duration1.0, fault_freq120, snr_db2): 模拟通风机振动信号 fs: 采样率 duration: 信号时长(秒) fault_freq: 故障特征频率(Hz) snr_db: 信噪比(dB) t np.arange(0, duration, 1/fs) # 基频振动转频 base 1.0 * np.sin(2 * np.pi * 25 * t) # 故障冲击周期性指数衰减脉冲 impulse_train np.zeros_like(t) period int(fs / fault_freq) for i in range(0, len(t), period): decay np.exp(-np.arange(0, min(200, len(t)-i)) / 30) impulse_train[i:ilen(decay)] 0.5 * decay * np.sin(2*np.pi*3000*np.arange(len(decay))/fs) # 高斯噪声 signal_power np.mean((base impulse_train)**2) noise_power signal_power / (10**(snr_db/10)) noise np.random.normal(0, np.sqrt(noise_power), len(t)) return t, base impulse_train noise, impulse_train t, noisy_signal, clean_impulse simulate_fan_vibration() print(f信号长度: {len(noisy_signal)}, 采样率: 12800 Hz)这段代码构造了三部分转频正弦、周期性指数衰减冲击串、高斯噪声。冲击串用 3000 Hz 载波调制模拟轴承故障冲击激发的高频共振。信噪比设 2 dB属于比较恶劣的工况。3.3 小波分解、阈值收缩与重构的完整实现def wavelet_denoise(signal, waveletdb4, level4, modeheursure, thresh_typehard): 小波阈值去噪 wavelet: 小波基通风机信号常用 db4/db6/sym8 level: 分解层数 mode: 阈值规则 thresh_type: hard 或 soft # 小波分解 coeffs pywt.wavedec(signal, wavelet, levellevel) # 用第一层细节系数估计噪声标准差MAD估计 detail1 coeffs[-1] sigma np.median(np.abs(detail1)) / 0.6745 # 对每层细节系数做阈值处理近似系数不动 coeffs_thresh [coeffs[0]] # 保留近似系数 for i in range(1, len(coeffs)): cD coeffs[i] # 阈值随层数调整高层细节系数噪声占比低阈值应缩小 thresh sigma * np.sqrt(2 * np.log(len(cD))) / np.sqrt(i) if thresh_type hard: cD_thresh pywt.threshold(cD, thresh, modehard) else: cD_thresh pywt.threshold(cD, thresh, modesoft) coeffs_thresh.append(cD_thresh) # 重构 denoised pywt.waverec(coeffs_thresh, wavelet) return denoised[:len(signal)] denoised wavelet_denoise(noisy_signal, waveletdb4, level4, thresh_typehard)关键逻辑说明pywt.wavedec返回的系数列表是[cA_n, cD_n, cD_n-1, ..., cD_1]第一个是最高层近似系数后面从高层到低层排列细节系数。噪声标准差用第一层细节系数的 MAD 估计这是 Donoho 经典做法对冲击鲁棒。阈值除以sqrt(i)是我自己的调整高层细节系数对应的频带低噪声能量本来就少如果还用统一阈值会把有用低频信息也砍掉。这个调整在通风机信号上比统一阈值效果好但属于经验性改动换场景需要重新验证。近似系数不处理因为通风机转频和低阶谐波都在里面动了会丢失工况信息。3.4 去噪效果验证时域波形与包络谱对比from scipy.fft import fft, fftfreq from scipy.signal import hilbert def envelope_spectrum(signal, fs): 计算包络谱用于观察故障特征频率 analytic hilbert(signal) envelope np.abs(analytic) envelope envelope - np.mean(envelope) N len(envelope) yf np.abs(fft(envelope))[:N//2] * 2 / N xf fftfreq(N, 1/fs)[:N//2] return xf, yf # 对比去噪前后的包络谱 xf_n, yf_n envelope_spectrum(noisy_signal, 12800) xf_d, yf_d envelope_spectrum(denoised, 12800) # 找 120 Hz 处的谱峰 idx_120 np.argmin(np.abs(xf_n - 120)) print(f去噪前 120Hz 谱峰: {yf_n[idx_120]:.4f}) print(f去噪后 120Hz 谱峰: {yf_d[idx_120]:.4f}) print(f谱峰提升倍数: {yf_d[idx_120]/yf_n[idx_120]:.2f})包络谱是通风机轴承故障诊断的核心工具。故障冲击会调制到高频共振频带包络解调后能在低频段看到故障特征频率及其谐波。去噪的目的就是让这个谱峰从噪声底里冒出来。参数说明hilbert做解析信号变换取模得到包络。去均值是为了去掉直流分量否则 FFT 第一个点会很大。谱峰提升倍数是一个直观指标一般做到 2 倍以上后续诊断才有把握。4. 小波基、分解层数和阈值规则的联合调参避坑4.1 小波基选 db 还是 sym消失矩和支撑长度的权衡小波基的选择影响三个指标消失矩阶数、支撑长度、正交性。通风机冲击信号的特征是短时突变需要小波基有好的时间局部化能力。db4 消失矩 4 阶支撑长度 7是工程默认值。db6 消失矩更高对多项式趋势抑制更好但支撑变长冲击定位精度下降。sym8 近似对称重构时相位失真小适合后续要做相位分析的场景。我的经验如果只做包络谱看幅值db4 够用如果要做冲击定位比如判断冲击在转子每转中的位置用 sym 系列相位失真小。coif 系列消失矩高但计算量大现场嵌入式设备上跑不动。4.2 分解层数选错会怎样过分解与欠分解的典型症状过分解层数太多的症状重构信号两端出现明显畸变包络谱低频段出现虚假谱峰。原因是每层分解后数据点数减半边界延拓的影响被逐层放大。欠分解层数太少的症状去噪后信号里仍有明显高频噪声包络谱底噪高故障谱峰不突出。排查方法画出各层细节系数的波形。如果最高层细节系数看起来还像噪声说明分解不够如果最低层细节系数几乎为零说明分解过头了。4.3 阈值规则和阈值函数的组合陷阱一个常见翻车场景用软阈值 通用阈值结果故障冲击被过度平滑包络谱上 120 Hz 谱峰反而比去噪前低。原因是软阈值本身压缩幅值通用阈值在低信噪比下偏大两者叠加把冲击也当噪声砍了。解决低信噪比场景用硬阈值 启发式阈值或者半软阈值。如果非要用软阈值把阈值乘以 0.5 到 0.7 的缩放因子手动降低去噪力度。另一个坑对近似系数也做阈值处理。近似系数里是转频和低阶谐波是工况信息砍了之后后续做转速估计会出问题。4.4 现场数据比模拟数据难在哪非高斯噪声与转速波动模拟信号用高斯噪声现场信号里除了高斯噪声还有气流脉动的周期性干扰、电机电磁噪声的离散谱线、传感器安装松动产生的冲击。这些非高斯成分用 MAD 估计噪声标准差会偏大导致阈值过高。应对先做一次谱分析把明显的离散干扰频带用陷波滤波器去掉再进小波阈值。转速波动则要求做阶次跟踪或角度域重采样否则故障特征频率会 smear包络谱上谱峰展宽。5. 从去噪到特征向量小波阈值在故障诊断链路里的进阶用法去噪只是第一步最终要输出能喂给分类器的特征向量。我一般从去噪后的信号里提取三类特征时域特征均方根、峭度、峰值因子。峭度对冲击敏感但受噪声影响大去噪后峭度值更稳定。频域特征包络谱上故障特征频率及其 2 到 3 次谐波的幅值、谱峰能量占比。小波域特征各层细节系数的能量比。故障冲击会让高频层能量占比上升这个指标比单纯时域峭度更鲁棒。def extract_features(denoised_signal, fs, fault_freq120): 从去噪信号提取特征向量 features {} # 时域 features[rms] np.sqrt(np.mean(denoised_signal**2)) features[kurtosis] np.mean((denoised_signal - np.mean(denoised_signal))**4) / (np.std(denoised_signal)**4) features[crest] np.max(np.abs(denoised_signal)) / features[rms] # 包络谱特征 xf, yf envelope_spectrum(denoised_signal, fs) idx np.argmin(np.abs(xf - fault_freq)) band slice(max(0, idx-3), idx4) features[fault_peak] np.max(yf[band]) features[fault_energy_ratio] np.sum(yf[band]**2) / np.sum(yf**2) # 小波域能量比 coeffs pywt.wavedec(denoised_signal, db4, level4) energies [np.sum(c**2) for c in coeffs] total sum(energies) features[high_freq_energy_ratio] (energies[-1] energies[-2]) / total return features feat extract_features(denoised, 12800) for k, v in feat.items(): print(f{k}: {v:.4f})这套特征向量可以直接送进 SVM 或随机森林做故障分类。如果做深度学习可以把去噪后的信号做短时傅里叶变换得到时频图喂给 CNN。一个验证技巧把特征向量在正常样本和故障样本上分别算均值和方差做 t 检验或 Mann-Whitney U 检验p 值小于 0.05 的特征才保留。我见过太多人把几十个特征一股脑塞进分类器结果过拟合严重实际部署时准确率暴跌。最后说一个我自己的习惯每次调完小波参数一定把去噪前后的包络谱叠在一张图上对比肉眼确认故障谱峰确实被增强而不是被削弱。参数再漂亮包络谱上谱峰没起来就是白搭。这套流程在多个通风机监测项目里跑下来db4 4 层 硬阈值 启发式阈值的组合在 10 kHz 采样率下最稳换采样率或换机型时需要重新验证分解层数。希望帮到你。本文还有配套的精品资源点击获取