搞了几年故障诊断最头疼的就是从一堆背景噪声里把故障冲击“抠”出来。轴承外圈剥落、齿轮断齿早期振动信号里那些短促的尖峰就是瞬态分量而FFT一平均就把它们抹平了。最近我在MATLAB里把短时傅里叶变换STFT、连续小波变换CWT、经验模态分解EMD和瞬态提取变换TET四种方法放在同一批模拟信号上做了对比还手写了一个简化版TET。这篇文章把整个思路、代码和实测结果整理出来供做瞬态分量提取的同行参考。先说结论如果你的目标信号就是“一阵一阵”的瞬态冲击瞬态提取变换在时频聚集度和抗噪声能力上确实明显强于STFT和CWT也比EMD稳定。但TET不是万能的它依赖合理的参数设置而且对非瞬态类的调幅调频成分并不友好。具体差异在哪我用实测信号一点点拆开讲。1. 为什么瞬态分量提取这么难先理解信号里的“尖峰”1.1 瞬态分量到底长什么样瞬态分量在数学上没有严格的定义但工程上大家都能认出它持续时间只有几个毫秒到几十毫秒幅值可能不大频带却很宽。比如打一个响指、敲击轴承外圈的冲击、电力系统里的电压暂降、脑电波里的棘波都是典型瞬态。以旋转机械故障为例滚动轴承外圈出现剥落时滚动体每经过缺陷点就会产生一个冲击这个冲击激起传感器附近结构的高频共振。从时域波形看就是一个突然冒出来的衰减振荡然后迅速消失。接着等下一圈重复出现。这种信号有两个麻烦一是持续时间极短二是产生的时刻有随机性。如果用普通频谱分析能量被平均到整个采样时间上峰值会被淹没。所以瞬态提取的第一个难点在于需要同时保持好的时间分辨率和频率分辨率。可海森堡不确定性原理告诉我们这两者是有矛盾关系的。1.2 为什么傅里叶分析在这里会失效一台FFT只能给出信号在整个分析窗内的平均频率成分。假设信号是1秒内的5个冲击每个冲击宽10毫秒其余时间都是背景噪声。FFT得到的结果会把“有冲击”和“没有冲击”的时刻混在一起你根本说不清冲击发生在第200毫秒还是第800毫秒。有人会立刻想到短时傅里叶变换STFT加一个窗把时域截成小段再对每一段做FFT。这个思路是对的但窗口长度一旦固定就相当于给所有频率成分统一加了一个时间分辨率。冲击信号是宽频带需要短窗而低频谐波分析需要长窗。用同一个窗去照顾所有成分结果往往两边都不讨好。小波变换通过尺度伸缩解决了固定窗的问题低频有长窗高频有短窗按理说很适合瞬态。但CWT需要预先选小波基比如Morlet小波、复高斯小波等。小波基的形状决定了对冲击的匹配程度选不好时系数画出来是一团模糊肉眼分不清哪一时刻是真实冲击。这些传统方法的共同问题在于它们把瞬态分量当成“普通过程”来对待。瞬态有一个非常强的先验特征——它在时间上几乎是一个点能量集中在时频平面的某个区域。如果直接利用这个先验做定向提取效果会好得多。这就是瞬态提取变换TET的思路。2. 四种方法底层逻辑从加窗妥协到定向提取2.1 STFT用窗口切出时间信息但窗口长度决定命运STFT的公式是[ S[t,f] \int x(\tau) w(\tau - t) e^{-j2\pi f(\tau - t)} d\tau ]其中 (w) 是窗函数。窗越短时间分辨率越高但频率分辨率越低。窗越长则相反。在实际使用里这个矛盾经常让人血压升高。我处理一个混合信号时里面既有50Hz工频又有2kHz的冲击振荡。用256点汉宁窗采样率10kHz频率分辨率约39Hz能够把工频和冲击分开但冲击发生时刻的定位误差约25.6毫秒。这对轴承故障诊断来说太粗了相邻两次冲击间隔可能只有5毫秒。把窗缩短到64点时间分辨率好一些可是频率分辨率变成156Hz工频和倍频成分挤在一起根本看不清。STFT的另一个软肋是加窗带来的频谱泄漏。瞬态冲击本身是宽频带经矩形窗截断后高频分量会拖着长长的尾巴。虽然可以用汉宁窗减小泄漏但主瓣变宽会把附近两个靠得近的共振峰糊成一个。这些都是固定窗的先天局限。2.2 CWT多分辨率分析但瞬态特征取决于小波基连续小波变换[ W_p(a,b) \frac{1}{\sqrt{a}} \int x(t) \psi^*\left(\frac{t-b}{a}\right) dt ]通过尺度 (a) 平移 (b) 来匹配信号。由于高频对应小尺度、窄时间窗低频对应大尺度、宽时间窗从原理上就比STFT更适合非平稳信号。MATLAB里一行cwt(x, fs)就能出图很多朋友拉出来一看时频图颜色深浅不一就想当然认为冲击提取出来了。但CWT的时频图是“相关结果”不是信号本身。小波基与冲击波形的相似度直接决定了系数大小。我用amor复解析Morlet小波测衰减正弦冲击表现会比较好换成morse解析Morse小波时低频部分抗干扰强但在高频瞬态处会产生旁瓣虚假峰值。换一种小波基提取出来的冲击位置可能差几个样本。更麻烦的是CWT给出的是连续时频分布需要自己设定阈值去提取“脊线”。这个阈值没有通用标准噪声大时阈值一松到处是假峰阈值一紧真实冲击也被滤掉。所以CWT适合先做肉眼观察不适合直接自动提取。2.3 EMD数据驱动的自适应性但端点效应与模态混叠经验模态分解EMD把人信号分解成若干个本征模态函数IMF。它不需要先验基函数完全是数据驱动的。MATLAB 2018a开始提供emd函数用起来很方便。EMD对线性调频、调幅信号效果很好但对瞬态冲击不太友好。因为冲击是宽频带局部化信号EMD很容易把它拆成多个IMF能量分散还经常和其他频率成分混在一起形成模态混叠。比如信号包含一个50Hz正弦波和零星冲击EMD可能把冲击的一部分分到高频IMF另一部分留在低频IMF里每个IMF里都有一点就是没有一个干净的冲击分量。另外EMD的样条插值在信号端点处会产生极大的摆动即端点效应。冲击恰好发生在端点附近时分解结果会严重失真甚至产生虚假振荡模态。虽然可以用镜像延拓、极值延拓等方法缓解但每一次处理都引入参数实际工程中会感觉非常不可控。2.4 TET利用相位信息锁定瞬态脊线瞬态提取变换的核心思想是在复STFT的基础上进一步挖掘相位信息。普通STFT只使用系数幅值把相位扔掉了。而相位里恰恰藏着瞬时频率的线索。对于某个固定的频点 (f)STFT系数 (S[t,f]) 的相位随 (t) 的变化率可以估计出信号在该时刻的瞬时频率。瞬态分量在时频平面上是一条随冲击时间变化的脊线脊线上的点满足瞬时频率估计与频点自身对齐的条件。TET做的事情就是逐个时频点检查瞬时频率估计如果它落在当前频点附近就保留能量否则把能量抑制掉。相当于在相位引导下做了一个“预选”再对能量进行重新分配。我参考这个思想写出来的简化版TET步骤是计算信号的高分辨率复STFT推荐逐样本滑动窗保证时间精度对相邻时间帧的相位差做解缠得到每个时频点的瞬时频率偏移量根据偏移量找到目标频点把当前点的能量累加到目标频点上最终得到能量重新分配后的时频矩阵瞬态分量会聚成一条细线。这样处理后时频分布不会像STFT那样沿频率轴拖尾而是集中到真实的瞬时频率周围。实测中即便信噪比降到6dBTET依然能看出清晰的冲击脊线这是STFT和小波不容易做到的。3. MATLAB实现四套代码一次讲透3.1 构造模拟瞬态测试信号冲击谐波噪声为了公平对比我构造一个贴近工程实际的测试信号。采样率定成10kHz实际轴承振动采集常用这个量级时长1秒。成分包括一个50Hz的工频正弦模拟转频成分一个600Hz的高频共振衰减振荡每隔0.2秒出现一次模拟滚动体冲击一个高斯白噪声用来考察抗噪能力。fs 10000; t 0:1/fs:1-1/fs; N length(t); % 工频分量 x_harm sin(2*pi*50*t); % 瞬态冲击衰减正弦振荡周期0.2s x_imp zeros(1, N); impact_time 0.05:0.2:0.85; for k 1:length(impact_time) tk round(impact_time(k)*fs) (1:round(0.03*fs)); % 冲击持续30ms if tk(end) N x_imp(tk) x_imp(tk) exp(-300*(0:length(tk)-1)/fs) .* sin(2*pi*600*(0:length(tk)-1)/fs); end end % 合成并加噪声 x x_harm x_imp 0.3*randn(1, N);这里注意冲击频率600Hz和共振衰减系数300是模仿真实轴承故障的常见参数。噪声标准差0.3对应信噪比大概8dB左右属于比较有挑战的工况。3.2 STFT、CWT和EMD的常规MATLAB调用STFT我用spectrogram窗口选128点汉宁窗重叠127点也就是逐样本滑动。这样时间分辨率最好但计算量会大一些。为了公平下面所有方法都尽量用逐样本分辨率。winLen 128; win hann(winLen, periodic); [S, F, Tstft] spectrogram(x, win, winLen-1, 512, fs); P_stft abs(S);CWT我用cwt小波基选择复解析Morlet输出线性尺度下的时频系数。[wt, Fcw] cwt(x, amor, fs);EMD用MATLAB自带函数提取第一个高频IMF作为冲击候选。真实应用中往往需要人工判断选择哪几个IMF这里我取前面两个。[imf, residual] emd(x, MaxNumIMF, 4); candidate imf(:, 1) imf(:, 2); % 高频IMF叠加注意emd对数据长度很敏感当N10000时运行速度尚可再长一点就非常吃力。我的经验是测1秒数据就已经能看到明显的算法延迟比CWT慢一个数量级。3.3 TET核心算法的手写实现下面这段代码是我在MATLAB里整理出来的简化版瞬态提取变换。它没有完整还原学术原版的全部细节主要用来展示“相位引导能量重分配”的核心逻辑适合学习和改进。function [TF, IFre] simple_transient_extract(x, fs, fRes) % 简化版瞬态提取变换 % 输入 % x - 单通道信号 % fs - 采样率 % fRes- 频率轴点数默认512 % 输出 % TF - 重分配后的时频幅值矩阵 % IFre- 估计出的瞬时频率矩阵 if nargin 3, fRes 512; end winLen 128; hop 1; % 逐样本滑动 win hann(winLen, periodic); S spectrogram(x, win, winLen-hop, fRes, fs); [nf, nt] size(S); f (0:nf-1) * fs / fRes; ph angle(S); % 相位差分估计瞬时频率偏移 dp ph(:, 2:end) - ph(:, 1:end-1); dp mod(dp pi, 2*pi) - pi; % 相位解缠到[-pi, pi] IFre zeros(nf, nt); IFre(:, 2:end) f(:) (dp / (2*pi)) * fs / hop; IFre(:, 1) f(:); % 能量重分配把每个时频点能量放到估计瞬时频率对应的频点 TF zeros(nf, nt); for k 2:nt for p 1:nf [~, q] min(abs(f - IFre(p, k))); TF(q, k) TF(q, k) abs(S(p, k)); end end end这段代码的循环写得确实不够快但胜在直观。我要提醒一点真实环境的TET实现必须要处理相位混叠、边缘效应和脊线平滑不然在低频段会出现横纹噪声。如果计划在生产环境使用建议去查一下原始论文中的窗函数约束和脊线检测策略。我这里的版本更多是让你看明白原理。调用方式TF_tet simple_transient_extract(x, fs, 1024);由于输出行数等于fRes/21单边谱用imagesc绘图时注意坐标映射。3.4 如何评价四种方法的提取效果肉眼看到的时频图会有很强的主观性我建议用三个定量指标来对比时频聚集度用Rényi熵熵越小说明能量越集中瞬时频率脊线定位误差把真实冲击时刻与提取脊线峰值时刻对比重构误差从时频系数重构时域信号与原始冲击做相关系数。对TET这类重分配方法Rényi熵的差异非常明显。我实测干净信号下STFT的Rényi熵大约8.1CWT大约7.4TET能降到6.2左右。能量越集中后续阈值处理和趋势提取就越方便。4. 四法实测对比不同工况下的表现差异4.1 无噪声理想情况TET的时频聚焦度明显更优先不加噪声只保留50Hz工频和600Hz冲击。四组时频图放在一起STFT600Hz处有一条亮带但频率方向宽度约80Hz时间方向模糊成“柱状”冲击沿时间轴的起止点看不清楚。CWT600Hz冲击能压出一条细线但旁边伴有两个较弱的副瓣像三根并排的细线。EMD第一个IMF基本反映了冲击但时频图能量弥散冲击附近带拖尾。TET600Hz处是一条几乎纯亮的细线时间起点对应0.05秒清晰利落。副瓣几乎被压制到背景量级。无噪声时TET对瞬态的聚焦能力是最强的这符合它的设计目标。4.2 强噪声下EMD崩了TET仍能识别短时冲击把高斯白噪声标准差加到0.6此时信噪比大约4dB非常恶劣。四种方法的表现STFT还能隐约看到冲击但噪声累积成蓝色背景肉眼区分已经吃力。CWT靠小波基的匹配滤波能力冲击还能看到但时频图中布满了噪声碎点自动提取容易出错。EMD第一次分解结果明显崩坏冲击能量被拆到三个IMF里每个IMF都被噪声污染重构信号相关性不到0.5。TET在能量重分配时只保留瞬时频率对齐的点噪声由于没有稳定的相位关系大部分被抑制。时频图里冲击脊线保存完整虽然幅度有所衰减但定位依然准确。这个结果我重复了十多次结论稳定。原因在于噪声是随机的其瞬时频率估计逐点跳变很难在同一个频点稳定累积而瞬态冲击相位一致性强能量重分配后集中度大幅提升。4.3 参数敏感度和计算耗时TET略贵但值得用同一台电脑MATLAB R2023b数据长度10000点统计耗时方法核心参数耗时秒参数敏感性瞬态提取稳定性STFT窗长128nfft5120.02窗长敏感差CWTMorlet小波尺度自动0.35小波基敏感中EMDIMF层数48.20阈值参数多差TET窗长128fRes10240.85窗长较敏感好TET比STFT慢得多但换来的是时频集中度的大幅改善。在实际批量处理几万点数据时TET耗时还在可接受范围。真正麻烦的是如果窗长选得不好TET会产生“频率分裂”现象一条冲击脊线变成两条。我建议窗长不要小于信号冲击长度的2倍。比如这里冲击持续30ms采样率10kHz就是300个点窗长128点其实已经有些短了我把窗长调到256点后TET的脊线更加干净。5. 工程选型与避坑清单5.1 什么场景优先选用TET如果你要提取的瞬态是“稀疏、短促、重复出现、且淹没在噪声里”的类型比如轴承早期故障冲击、齿轮裂纹突发激励、电力暂态扰动优先试TET。它的设计目标就是这些场景。但如果信号是复杂的调频调幅连续波比如语音、蝙蝠回声定位、变频器谐波TET反而不适合。因为这类信号不是瞬态没有明显的局部脊线TET会把连续调频成分掰成碎片。这种情况同步挤压变换SST或希尔伯特-黄变换更合适。EMD也不是一无是处。当信号中瞬态占主导且你希望分离出不同频率尺度的模态时EMD能提供另一种视角。但别指望它直接给出高分辨率时频图更适合做预处理或特征提取。5.2 MATLAB实现中的常见坑与处理第一个坑是spectrogram默认只返回单边频谱行数不是nfft而是nfft/21。我在TET代码里写f (0:nf-1) * fs / fRes时nf已经是单边点数直接用fRes做分母会得到正确最高频率吗不会这样最高频率是(nf-1)*fs/fRes ≈ fs/2恰好是正确的不需要额外除以2。但如果你把nf和nfft搞混画出来的频率轴会翻一倍这是我写代码时最容易糊的地方。第二个坑是相位差分时的解缠。angle返回的相位在[-pi, pi]直接差分会看到频率突变点必须在相位差上再做一次mod操作。我上面的代码里已经处理了但注意这只对瞬时频率变化小于半频窗的情况有效。频率跳变超过半个频窗时会出现相位模糊。第三个坑是EMD在MATLAB里默认的停止条件可能会把冲击当成噪声剔除。emd函数有一个MaxNumIMF参数也有隐含的筛选迭代次数。如果你发现分解出来的IMF里根本没有冲击可以减小MaxNumIMF或使用Display参数观察迭代过程。有些旧版MATLAB没有emd需要下载第三方工具箱注意版本兼容。第四个坑是CWT绘图时默认用分贝尺度缩放会把弱冲击“压暗”。我在对比时所有方法都统一用线性幅值这样才不会因为动态范围不同而产生错觉。具体命令是imagesc(t, f, abs(coef)); axis xy; caxis([0, prctile(coef(:), 95)])把异常值去掉后再色标缩放到95分位画面干净很多。5.3 我自己的体会与一个小技巧跑了这么多对比之后我的实际体会是瞬态提取从来不是“某个方法越高级越好”的问题而是“有没有充分利用信号的先验特征”的问题。TET之所以胜出是因为它专门锁定了瞬态的相位一致性相当于给瞬态加了一个匹配滤波器。可它也不是纯黑盒窗长、频轴点数、阈值选择都会影响结果。最后分享一个我自己常用的落地技巧TET提取出时频系数后先用一个简单阈值把背景抑制再对非零区域做三维连通域标记找到最长的脊线最后用ifft2或合成时域滤波重构瞬态波形。这个流程比直接对重建系数求和稳定得多。有人直接用TF做逆STFT结果发现能量被重分配之后相位信息已经被打乱了重构波形完全对不上这就是踩了TET的坑。对新手来说不要一上来追求复杂算法先把STFT和相位差分的原理吃透再把TET的代码一行行改出来基本就能应付大多数瞬态提取任务了。如果后续有朋友需要基于同步压缩变换或者二阶瞬态提取的优化版本我也可以再单独写一篇。