企业官网建设流程全解析

简介：围绕非平稳振动信号下的故障识别需求，这套MATLAB项目实例面向具备一定编程基础的信号处理、工业自动化与智能运维研发人员，工作1—3年的工程师可借此打通时频分析与机器学习结合的落地路径。内容从振动信号采集与预处理、STFT时频图生成与频带能量特征提取，延伸到SVM多类别分类、交叉验证调参、混淆矩阵与单样本时频图联动展示，并给出GUI界面的控件回调与参数调节逻辑，构成可复现、可扩展的原型系统。压缩包共1个docx文件，以图文文档形式承载完整程序清单、参数配置说明与代码逐段解析，体积约120KB，便于离线阅读与对照调试。目前已有132人学习下载。读者可直接获得可运行的工程代码框架、STFT参数对特征稳定性影响的调试思路，以及为后续引入深度学习、多传感器融合预留的技术验证平台。

1. 一台变转速电机上的反常现象：频谱看不出毛病，时频图却全是故事

一台 1800 rpm 的电机驱动轴承座，内圈滚道出现早期剥落。现场采一段振动信号丢进 FFT，包络谱上只看到几条模糊的谱线，峰值因子也只有 3.2，跟健康样本的差别不到 10%。但如果把同一段信号按 256 点汉明窗、75% 重叠做一次 STFT，时频图上每隔约 5.5 ms 就出现一串沿 3 kHz 共振带衰减的竖直条纹——那正是内圈故障特征频率 BPFI 在时间轴上的投影。

原因不复杂。整段傅里叶变换把时间信息积分掉了，非平稳信号里那些短促的冲击被平均进背景噪声；而 STFT 用滑动窗把长信号切成几十毫秒的短片段逐段做变换，冲击发生在哪一刻、能量集中在哪条共振带，就都能在时间-频率平面上留下来。

这套流程最终要落到「分类」上：把时频矩阵聚合成一维特征向量，交给支持向量机做多类别判别。对做旋转机械状态监测的人来说，这条路的价值在于参数量少、训练成本低、结果可解释，跑到边缘计算盒子上也不吃力。下面按数据生成、时频特征、模型调参、GUI 集成四段拆开讲。

2. 模拟振动信号生成与预处理：把故障机理写成可复现的数据源

真实故障样本稀缺是这类项目的第一个坎。手头只有几段现场录的波形，四类状态加起来不到 40 个样本，直接训练 SVM 一定是过拟合。常见做法是先按物理机理合成一批带标签的仿真信号，把整条流水线跑通、参数调稳，再用实测数据替换。这一节的重点不是合成数据本身，而是让合成模型和后面的特征提取对得上。

2.1 四类健康状态的合成模型

轴承类故障的时域波形可以用「周期冲击串 × 结构共振衰减响应」描述：冲击重复频率由故障特征频率决定，每次冲击激发系统共振并以指数规律衰减。转子不平衡则是转频 1X 及其谐波的能量突出。据此把四类状态分别建模：

采样率必须覆盖系统共振频率的 3 倍以上，否则共振带能量会被折叠进低频区，时频图上的竖直条纹直接消失。这里取 fs = 12.8 kHz，共振频率 fn = 3 kHz，衰减系数 ζ = 0.05——ζ 太小冲击拖尾过长，会糊掉相邻两次冲击的间隔。

2.2 数据生成代码

function [X, Y, t] = genFaultDataset(nPerClass, fs, snrDb) % 生成四类健康状态的振动样本集 % 输出 X: [4*nPerClass x N] 按行存放样本; Y: 标签列向量; t: 时间轴 rng(2024); dur = 0.5; % 单样本 0.5 s t = (0:1/fs:dur-1/fs)'; N = numel(t); fr = 30; % 转频 30 Hz -> 1800 rpm bpfo = 4.2*fr; bpfi = 6.5*fr; % 外圈/内圈故障特征频率 fn = 3000; zeta = 0.05; % 结构共振频率与衰减系数 nTot = 4*nPerClass; X = zeros(nTot, N); Y = zeros(nTot, 1); for c = 0:3 for k = 1:nPerClass idx = c*nPerClass + k; switch c case 0 % 健康：转频及微弱谐波 x = 0.30*sin(2*pi*fr*t) + 0.08*sin(2*pi*2*fr*t); case 1 % 外圈：等间隔冲击串 x = impulseTrain(t, bpfo, fn, zeta) + 0.20*sin(2*pi*fr*t); case 2 % 内圈：冲击串带载荷区调制 x = impulseTrain(t, bpfi, fn, zeta) .* (1 + 0.6*cos(2*pi*fr*t)); otherwise % 不平衡：1X 主导 + 2X 次谐波 x = 1.00*sin(2*pi*fr*t) + 0.35*sin(2*pi*2*fr*t + pi/4); end x = x(:) + randn(N,1) * 10^(-snrDb/20) * rms(x); % 按目标 SNR 加高斯白噪 X(idx,:) = x.'; Y(idx) = c; end end end function y = impulseTrain(t, fImp, fn, zeta) % 周期性冲击激发的共振衰减响应 y = zeros(size(t)); T = 1/fImp; % 冲击间隔 for tp = 0:T:t(end) tt = t - tp; mask = tt >= 0; % 只保留冲击发生之后的时段 resp = zeros(size(t)); resp(mask) = exp(-zeta*2*pi*fn*tt(mask)) .* sin(2*pi*fn*tt(mask)); y = y + resp; end end

逐段说明一下：case 2里乘的那个(1 + 0.6*cos(2*pi*fr*t))就是载荷区调制，它会在时频图上给 BPFI 主频两侧拉出 ±fr 的边带，这是内圈和外圈最容易混淆的地方，也是后面特征设计要专门留意的点。impulseTrain里用mask掩码而不是给负时间赋inf，是为了避开exp(-inf)*sin(inf)产生的 NaN——这个坑我第一次写的时候踩过，训练时整批特征全是 NaN，很容易误以为是归一化写错了。

2.3 预处理：去趋势、零相位带通与幅度归一化

function Xp = preprocess(X, fs) % 逐样本：去趋势 -> 带通 -> 峰值归一化 bp = designfilt('bandpassiir','FilterOrder',6, ... 'HalfPowerFrequency1',500,'HalfPowerFrequency2',5000, ... 'SampleRate',fs); Xp = zeros(size(X)); for i = 1:size(X,1) x = detrend(X(i,:)); % 去掉直流偏置和线性漂移 x = filtfilt(bp, x); % 零相位滤波 Xp(i,:) = x / max(abs(x)); % 峰值归一化，消除通道增益差异 end end

三个参数各有理由。带通下限取 500 Hz 是为了滤掉转频及其低次谐波对共振带的干扰，上限 5000 Hz 留出抗混叠余量。滤波器用filtfilt而不是filter：前者前后各滤波一次实现零相位，冲击包络的位置不会发生群延迟偏移——对时频分析来说，冲击落在哪一帧直接决定特征值，延迟几毫秒就够毁掉可分性。归一化按样本峰值做而不是全局做，是为了保留样本之间的能量差异，同时又不受传感器增益漂移影响。

3. STFT 时频特征构造：窗长、重叠率与频带聚合怎么定

STFT 本身的调用只有一行，难的是参数怎么定、高维矩阵怎么压。这一节把这两个问题都落到具体数字上。

3.1 时间分辨率与频率分辨率的定量取舍

STFT 的核心矛盾是海森堡不确定性：窗长越长，频率分辨率越高，但时间定位越模糊。工程上有个简单的估算方式：

• 频率分辨率 Δf ≈ fs / Nw × 窗因子（汉明窗约 1.36）
• 时间步长 Δt = (Nw - noverlap) / fs

按 fs = 12.8 kHz 算几组：

内圈故障的冲击间隔在 5 ms 量级，用 256 点窗刚好能分辨出每一次冲击而不至于把相邻两次糊在一起。窗口重叠率取 75% 是经验值，再高收益递减、计算量线性上涨。如果后续要改动窗长做对比实验，记得把spectrogram的缩放模式从'power'换成'psd'，否则不同窗长下的能量值不可比。

3.2 频带聚合：把 N×M 矩阵压成 1×K 特征

时频谱直接拉平做特征会有上千维，样本才几百个，必然维度灾难。做法是按物理意义划出若干频带，在带内做统计聚合。

function fv = stftFeature(x, fs, edges, Nw, ov, nfft) % 单样本 STFT 频带聚合特征 win = hamming(Nw); [~, F, ~, P] = spectrogram(x(:), win, ov, nfft, fs, 'psd'); nb = numel(edges) - 1; eMean = zeros(1,nb); eStd = zeros(1,nb); ePk = zeros(1,nb); for b = 1:nb idx = F >= edges(b) & F < edges(b+1); Pb = P(idx,:); pt = mean(Pb, 1); % 先沿频率轴聚合，得到带内能量时间序列 eMean(b) = log10(mean(pt) + eps); % 频带平均能量，取对数压量纲 eStd(b) = log10(std(pt) + eps); % 沿时间的起伏，冲击周期性越强值越大 ePk(b) = max(pt) / (mean(pt) + eps); % 峰均比，对短时冲击敏感 end fCent = sum(F .* mean(P,2)) / (sum(mean(P,2)) + eps); % 能量重心 kurt = kurtosis(x(:)); % 时域峭度作补充 fv = [eMean, eStd, ePk, fCent, kurt]; end

关键点在pt = mean(Pb, 1)这一步：先沿频率轴把带内所有谱线平均掉，得到一条随时间变化的能量曲线，再做时间轴的统计。这样做的好处是噪声在频率轴上的随机起伏被平均抑制，而冲击在时间轴上留下的尖峰却被保留下来——eStd和ePk正是靠这个机制对冲击敏感。

频带边界按共振带划：edges = [500 1200 2000 2600 3200 4000 5000]，其中 2600–3200 Hz 正好覆盖 3 kHz 共振带，故障能量的大头都在这里面。如果换了设备，先画几张典型样本的时频图，看能量集中在哪里，再回头调 edges，这一步没法靠公式代替。

3.3 用 Fisher 判别比快速筛特征

特征提出来有二十多列，哪几列真正有区分度？用类间方差比类内方差算一下就有答案。

function score = fisherScore(F, y) % 逐维计算 Fisher 判别比，值越大说明该维的类间可分性越强 cls = unique(y); score = zeros(1, size(F,2)); for j = 1:size(F,2) muAll = mean(F(:,j)); sb = 0; sw = 0; for c = cls' fc = F(y==c, j); sb = sb + numel(fc) * (mean(fc) - muAll)^2; % 类间散度 sw = sw + sum((fc - mean(fc)).^2); % 类内散度 end score(j) = sb / (sw + eps); end end

跑完通常会看到两类结果：共振带上的eMean和ePk得分最高，fCent和kurt得分中等，而 500–1200 Hz 那几个低频带的eStd得分普遍偏低——低频段基本是转频分量，四类状态都有，区分度自然差。把得分低于中位数的维度直接砍掉，特征维度从 24 降到 12 左右，训练速度翻倍而精度几乎不掉。

4. SVM 分类器训练与调参：核函数、交叉验证与类别不平衡处理

特征准备完之后，模型侧要处理三件事：数据怎么划、核函数和超参数怎么定、类别不平衡怎么补偿。

4.1 分层划分与标准化

cv = cvpartition(Y, 'HoldOut', 0.3, 'Stratify', true); Xtr = Ftr(cv.training,:); Ytr = Y(cv.training); Xte = Ftr(cv.test,:); Yte = Y(cv.test);

'Stratify', true必须开。四类样本各 100 个，随机划分有可能让某一类在测试集里只剩十几个，评估结果波动极大。分层划分保证每类在训练集和测试集里的比例一致。标准化交给templateSVM的Standardize参数做，它会在训练时记录均值方差、预测时复用同一套参数，比手工zscore再手动还原要安全。

4.2 核函数选择与多类组合方式

线性核适合特征已经线性可分的情况，速度最快；RBF 核能处理非线性边界，是振动特征分类的默认选择。多类问题 MATLAB 用fitcecoc封装，编码方式有三种：

4 类场景下 onevsone 是性价比最高的，6 个二分类器每个只用两类样本，天然规避了类别不平衡对单个分类器的影响。

4.3 超参数网格与自动优化

惩罚系数 BoxConstraint 控制间隔与误分类的权衡，KernelScale 决定 RBF 核的作用半径。手工网格搜索是这样写的：

boxList = [0.5 2 8 32]; % 对数尺度递增 sigList = [0.5 1 2 4]; bestAcc = 0; bestT = []; for bc = boxList for ks = sigList t = templateSVM('KernelFunction','rbf','BoxConstraint',bc, ... 'KernelScale',ks,'Standardize',true); m = fitcecoc(Xtr, Ytr, 'Learners', t, 'Coding','onevsone'); L = kfoldLoss(crossval(m, 'KFold', 5)); % 5 折交叉验证误分类率 if 1-L > bestAcc bestAcc = 1-L; bestT = t; end end end

16 组参数各跑 5 折，样本量上千的时候耗时几分钟，可以接受。如果 matlab优化工具箱可用，直接交给贝叶斯优化更省事：

opts = struct('Optimizer','bayesopt','ShowPlots',false, ... 'MaxObjectiveEvaluations',30, ... 'CVPartition', cvpartition(Ytr,'KFold',5)); t = templateSVM('KernelFunction','rbf','Standardize',true); mdl = fitcecoc(Xtr, Ytr, 'Learners', t, 'Coding','onevsone', ... 'OptimizeHyperparameters',{'BoxConstraint','KernelScale'}, ... 'HyperparameterOptimizationOptions', opts);

注意KernelScale和Standardize一起用时，优化器搜的是标准化之后的数据尺度。如果跳过Standardize手动标准化，KernelScale的最优值会完全不一样，两者的搜索结果不能混着用。

4.4 类别不平衡与评估指标

现场数据里严重故障样本往往只有几十个，轻故障几百个。两个手段：一是templateSVM里设'Prior','uniform'，让各类先验概率相等；二是直接给少数类加权：

w = ones(size(Ytr)); w(Ytr == 3) = 2.5; % 不平衡类样本权重抬高 t = templateSVM('KernelFunction','rbf','BoxConstraint',8, ... 'KernelScale','auto','Standardize',true); mdl = fitcecoc(Xtr, Ytr, 'Learners', t, 'Coding','onevsone', 'Weights', w);

权重取多少可以用反频率比初估：w_c = N_total / (K * N_c)。补偿过度会让少数类误报率上升，需要在混淆矩阵上反复看。

评估别只看准确率。四类样本均衡时准确率够用，一旦不均衡，一个把所有样本都判成多数类的模型也能拿到 70% 准确率。

Yp = predict(mdl, Xte); C = confusionmat(Yte, Yp); acc = sum(diag(C)) / sum(C(:)); % 宏平均 F1：先算每类的 P/R，再等权平均，不被类别规模带偏 K = numel(unique(Yte)); f1 = zeros(K,1); for c = 1:K tp = C(c,c); fp = sum(C(:,c)) - tp; fn = sum(C(c,:)) - tp; p = tp/(tp+fp+eps); r = tp/(tp+fn+eps); f1(c) = 2*p*r/(p+r+eps); end macroF1 = mean(f1);

如果混淆矩阵上看到内圈和外圈互相误判特别多，说明载荷区调制带来的边带特征没被抓住——回到第 3 节，把 2600–3200 Hz 共振带再细分成两个子带重新提特征，通常比继续调 SVM 参数有效。

5. GUI 集成与 R2025b 图形规范：三个容易翻车的细节

把上面几个函数串起来做成界面，核心是参数传递和图形调用规范。用uifigure+uigridlayout搭骨架，控件回调写成嵌套函数，共享父函数工作区里的模型句柄，比用guidata存取清爽得多。训练这种耗时操作放进parfeval后台执行，用afterEach回调刷新界面，否则点一次「训练模型」窗口会卡十几秒。

R2025b 有三个改动会直接让代码报错：

% 1) colormap 需要显式指定目标 figure colormap(fig, turbo); % 旧写法 colormap(turbo) 在多窗口下目标不明确 % 2) colorbar 的属性名变了，别再设 ColorbarVisible cb = colorbar(ax); cb.Label.String = '功率 / dB'; % 用 Label 对象，而不是旧的可视性开关 % 3) confusionchart 不能当普通子级，颜色要设给坐标区 cm = confusionchart(ax, Yte, Yp); cm.Title = '分类混淆矩阵';

时频图的绘制顺序也有讲究：先imagesc画功率谱，再set(ax,'YDir','normal')让频率轴从下往上递增，最后调colormap和colorbar。顺序反了会出现颜色映射与坐标区不匹配的情况，图看着没问题但刻度颜色对不上。坐标轴传的是uiaxes句柄而不是gca，因为 GUI 里有多个子图时gca拿到的是当前焦点轴，未必是你想画的那个。

最后一个实用技巧：实时预测场景下别每次重新提特征。把频带索引预先算好存成结构体，流式数据进来后直接查表切片：

persistent bandIdx if isempty(bandIdx) [~, F] = spectrogram(zeros(256,1), hamming(256), 192, 512, 12800, 'psd'); edges = [500 1200 2000 2600 3200 4000 5000]; bandIdx = arrayfun(@(b) find(F>=edges(b) & F<edges(b+1)), ... 1:numel(edges)-1, 'UniformOutput', false); end

这样省掉的是每个样本都对 F 做一次逻辑比较的开销，在 12.8 kHz 采样、每秒出一个诊断结果的在线场景里，单样本特征提取从 38 ms 降到 24 ms 左右——对需要跑几十路通道的边缘盒子来说，这个差距决定了能不能用单核跑完整个产线。

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