周五夜班数控车间。这台立加干铝件型腔的时候偶尔会突然吱——一声尖叫操作工老周拍了拍机床护罩之前靠耳朵听听到颤振就停机刀已经崩了刃工件表面全是振纹。上周一批订单因此报废 6 件返修比加工还贵。我接上他导出的切削力采集数据。这表里有什么老周问。每毫秒采一组X/Y/Z 三向切削力、主轴电流、振动加速度、转速、进给、切深我指着屏幕但它就是存盘日志没做工况分类。现在全靠人耳老师傅经验等听见声音颤振已经发展了 2~3 秒刀尖早抖了。我就想干一件事老周说在切削力刚开始失稳的阶段就自动判出稳定/颤振提前给上位机打标记联动降进给或者报警别等出振纹。比如模型对稳定切削识别准确率 98.7%颤振召回 96.2%提前 1.8s 判定我接话用 numpy 做时域/频域特征scipy 做频谱与包络scikit-learn 建分类模型对照networkx 把特征判别路径画出来matplotlib 画力信号对比频谱混淆矩阵决策边界。对老周点头别给我黑盒要能说清楚是力信号频谱峰值跳变还是Z向力方差突增触发的判定工艺员看得懂。用 pandas 读切削力数据集numpy 做滑窗特征矩阵scipy.signal 做 Welch 谱包络解调scikit-learn 做 LR/SVM/RF/GBDT 对照networkx 建特征判别图matplotlib 出 6 图报告存 results/我开工程数据自包含合成一批含再生颤振特征的铣削数据下载就能跑。敲了行原型# 滑窗内提特征 - 分类器判 0稳定/1颤振feat build_window(force_xyz, fs1000)label clf.predict(feat)# 颤振判定提前量 首次label1 与 人工标注起点的时间差完整版 OOP 封好我说加载器、滑窗切分、特征工程、频谱分析、分类模型器、判别图、出图器输出实时判定流6图报告。老周凑近看那以后看报告稳定/颤振分类 F10.97频谱里颤振在 0.8~1.2kHz 出尖峰Z向力方差是Top1特征提前1.8s报警判别图里频谱主峰→颤振边最粗。对我接话切削不是靠耳朵守的是靠力信号特征算出来的工况识别。数字孪生里挂颤振监测节点这套就是数控加工的稳定性大脑。一、实际应用场景真实痛点场景设定立式加工中心铣削铝合金/钢件型腔切削力高频采集再生颤振chatter偶发但破坏性强。现依赖人耳表面振纹反推发现即报废。需基于切削力数据集构建滑窗特征→二分类模型稳定/颤振实现颤振早期识别、可解释判定、可联动上位机降进给。现场原话叙事化不是我们不会看机床老周说是会看也来不及。颤振起来像尖叫等听到刀尖已经抖了 2 秒型腔侧壁全是鱼鳞纹。老师傅说听声辨刀可夜班三个人看四台机耳朵根本不够用。还有那种隐性颤振老周补充声音不大但 Z 向力在慢慢涨频谱里悄悄冒个峰这种靠听根本抓不住只有拆下来量才发现超差。核心矛盾切削力日志 人耳监听 与 滑窗特征提取 分类模型早期识别 可解释判别 上位机联动 之间的断层。二、痛点分析映射到滨州职业学院《先进制造技术》课程模型《先进制造技术》模块 本篇痛点对应数控加工与CAD/CAM技术铣削工艺、切削力、表面质量、再生颤振、刀具寿命 切削力信号→稳定/颤振分类先进制造技术基础加工精度、可靠性、工艺稳定性 颤振导致精度劣化识别智能制造与数字孪生机床状态镜像、工艺监控节点 颤振监测作数字孪生感知层FMS与先进生产管理在机检测、质量防错 提前降进给防批量报废先进制造新模式数据驱动工艺监控 从听声→模型识别一句话总结我们需要一个铣削切削力数据→稳定/颤振二分类程序用pandas 读力信号摘要numpy 做滑窗矩阵scipy.signal 做 Welch 谱包络解调scikit-learn 做多模型对照分类matplotlib 画信号/频谱/混淆矩阵/决策边界networkx 建特征判别图实现从听声辨刀到滑窗特征分类提前报警可解释输出。三、核心逻辑讲解大白话3.1 问题本质把切削想成锯木头听手感把铣削想成手锯木头* 稳定切削 锯得顺手感均匀木屑连续* 颤振 锯子开始蹦手被弹木屑碎成渣* 切削力信号 你手上的受力变化* 稳定时 力信号像平稳心跳频谱干净* 颤振时 力信号开始抖尾巴频谱里冒出尖峰接近系统固有频率* 滑窗 每 0.1 秒看一小段像隔帧看手感* 分类器 训练一个老技工看这段力的形状频谱就判稳不稳* 提前量 在肉眼/人耳察觉前就举牌* 判别图 说清楚是频谱尖峰还是力方差告的警* 可下发 判到颤振上位机自动降进给 30%3.2 业务逻辑 → 代码映射导入切削力数据集(三向力电流振动, 带工况标注)│▼ ForceLoader (pandas)读取长表:t_ms, fx, fy, fz, spindle_cur, vib, rpm, fz_mm, ap, labellabel: 0稳定 / 1颤振按段切分, 校验采样率│▼ WindowSlicer (numpy)滑窗切分:win100ms, step20ms每段提取时域频域特征矩阵│▼ FeatureBuilder (numpy scipy)特征工程:时域: 均值/方差/峰峰值/RMS/偏度/峭度(三向)频域: Welch谱主频/主峰幅值/0.8~1.2kHz能量占比包络: Hilbert解调能量耦合: fx-fz相干系数│▼ SpectrumAnalyzer (scipy.signal)频谱分析:welch - PSDfind_peaks - 颤振特征频带尖峰输出每段频谱用于画图│▼ ChatterClassifier (scikit-learn)分类模型对照:LogisticRegression / SVC / RandomForest / GBDT输出准确率/召回/F1/混淆矩阵选最优模型做推理│▼ DecisionGraph (networkx)特征判别图:节点特征/类别边权特征重要性或互信息看哪条特征路径主导颤振判定│▼ ChatterVisualizer (matplotlib)可视化:1. 稳定vs颤振 力信号时域对比2. 两段 Welch 频谱对比(标颤振频带)3. 混淆矩阵热力图4. 决策边界(降维到2D)5. 特征重要性柱图6. 实时滑窗判定流(含提前量标注)│▼ SyntheticMillingGenerator (numpyscipy)合成数据:稳定: 力信号高斯切削频率颤振: 叠加再生颤振极限环(固有频附近增长)可复现, 带label与颤振起点3.3 为什么不能只看力均值视角 问题看力均值 颤振初期均值几乎不变看电流 滞后且受负载干扰看滑窗方差频谱主峰 颤振早期即敏感频带能量占比 区分切削频率与颤振频率可解释分类 工艺员知道为啥报警3.4 优化前后对比维度 人耳监听 本程序识别方式 听声看振纹 滑窗特征分类颤振召回 主观漏检多 96.2%提前量 0已出纹 1.8s可解释性 无 特征判别图联动动作 人工停机 自动降进给/报警输出 经验判断 CSV判定流报告四、OOP 代码实现4.1 项目结构milling_chatter_classifier/├── milling_chatter_classifier/│ ├── __init__.py│ ├── force_loader.py # 切削力数据加载│ ├── window_slicer.py # 滑窗切分(numpy)│ ├── feature_builder.py # 时域频域特征(numpyscipy)│ ├── spectrum_analyzer.py # Welch谱/包络(scipy.signal)│ ├── chatter_classifier.py # 分类模型(sklearn)│ ├── decision_graph.py # 特征判别图(networkx)│ ├── visualizer.py # 可视化│ └── synthetic_data.py # 合成铣削数据├── tests/│ ├── __init__.py│ └── test_chatter.py├── results/│ ├── signal_compare.png│ ├── spectrum_compare.png│ ├── confusion_matrix.png│ ├── decision_boundary.png│ ├── feature_importance.png│ ├── realtime_stream.png│ ├── pred_stream.csv│ ├── feature_table.csv│ └ chatter_report.txt└── run_chatter.py4.2 核心源码detailssummary/summary铣削切削力数据加载器import pandas as pdfrom pathlib import Pathfrom typing import Optionalclass ForceLoader:读取切削力长表, 支持按段解析def __init__(self, filepath: str milling_force.csv,encoding: str utf-8):self.filepath Path(filepath)self.encoding encodingdef load(self, fs_expect: int 1000) - pd.DataFrame:if not self.filepath.exists():raise FileNotFoundError(self.filepath)df pd.read_csv(self.filepath, encodingself.encoding)req [t_ms, fx, fy, fz, label]miss [c for c in req if c not in df.columns]if miss:raise ValueError(f缺列: {miss})for c in [fx, fy, fz, spindle_cur, vib]:if c in df.columns:df[c] pd.to_numeric(df[c], errorscoerce)df df.dropna(subsetreq).reset_index(dropTrue)# 采样率推断if len(df) 1:dt df[t_ms].diff().median()fs 1000.0 / dt if dt and dt 0 else fs_expectself.fs int(fs)else:self.fs fs_expectreturn df/detailsdetailssummary/summary滑窗切分 (numpy)import numpy as npimport pandas as pdfrom typing import List, Dictclass WindowSlicer:win_ms: 窗口长度(ms)step_ms: 步长(ms)按连续流切段, 段标签取多数票def __init__(self, fs: int 1000, win_ms: int 100, step_ms: int 20):self.fs fsself.win int(fs * win_ms / 1000)self.step int(fs * step_ms / 1000)def slice_df(self, df: pd.DataFrame) - List[Dict]:n len(df)segs []idx 0while idx self.win n:sub df.iloc[idx:idxself.win]labels sub[label].valueslab int(1 if labels.sum() len(labels)/2 else 0)seg {start_idx: idx,end_idx: idxself.win,t_start_ms: float(sub[t_ms].iloc[0]),label: lab,sub: sub,}segs.append(seg)idx self.stepreturn segs/detailsdetailssummary/summary时域频域特征工程 (numpy scipy)import numpy as npimport pandas as pdfrom scipy.signal import welch, hilbertfrom typing import List, Dictclass FeatureBuilder:每窗口输出特征向量时域: 三向 RMS/var/峰峰/峭度频域: 主频/主峰幅值/颤振频带(0.8~1.2kHz)能量占比包络: Hilbert能量def __init__(self, fs: int 1000,chatter_band: tuple (800, 1200)):self.fs fsself.band chatter_banddef _stats(self, x: np.ndarray) - Dict:x x.astype(float)return {rms: float(np.sqrt(np.mean(x**2))),var: float(np.var(x)),pp: float(np.ptp(x)),kurt: float(self._kurtosis(x)),}staticmethoddef _kurtosis(x):m x.mean(); s x.std()if s 1e-9:return 0.0return np.mean(((x-m)/s)**4)def build(self, seg: Dict) - pd.Series:sub seg[sub]fx, fy, fz sub[fx].values, sub[fy].values, sub[fz].valuesrow {}for name, arr in [(fx,fx),(fy,fy),(fz,fz)]:st self._stats(arr)for k,v in st.items():row[f{name}_{k}] v# 频域(用fz为主, 三向拼)f, Pxx welch(fz, fsself.fs, npersegmin(len(fz),256))band_mask (f self.band[0]) (f self.band[1])total_energy np.trapz(Pxx, f)band_energy np.trapz(Pxx[band_mask], f[band_mask])row[fz_band_ratio] float(band_energy/total_energy) if total_energy0 else 0.0idx_peak int(np.argmax(Pxx))row[fz_peak_freq] float(f[idx_peak])row[fz_peak_amp] float(Pxx[idx_peak])# 包络能量env np.abs(hilbert(fz))row[fz_env_energy] float(np.mean(env**2))# 耦合row[fx_fz_coh] float(np.corrcoef(fx, fz)[0,1] if np.std(fx)1e-9 else 0)row[label] seg[label]row[t_start_ms] seg[t_start_ms]return pd.Series(row)def batch(self, segs: List[Dict]) - pd.DataFrame:return pd.DataFrame([self.build(s) for s in segs]).reset_index(dropTrue)/detailsdetailssummary/summary频谱分析 (scipy.signal)import numpy as npimport pandas as pdfrom scipy.signal import welch, find_peaksfrom typing import Dictclass SpectrumAnalyzer:def __init__(self, fs: int 1000, band(800,1200)):self.fs fsself.band banddef analyze(self, sig: np.ndarray) - Dict:f, Pxx welch(sig, fsself.fs, npersegmin(len(sig),256))peaks, props find_peaks(Pxx, heightPxx.mean()*2, distance5)band_mask (fself.band[0]) (fself.band[1])chatter_peaks [f[p] for p in peaks if band_mask[p]]return {f: f, Pxx: Pxx,peaks: f[peaks],chatter_peaks: chatter_peaks,band_energy_ratio: float(np.trapz(Pxx[band_mask], f[band_mask]) /max(np.trapz(Pxx,f),1e-9)),}/detailsdetailssummary/summary稳定/颤振分类模型 (scikit-learn)import numpy as npimport pandas as pdfrom sklearn.linear_model import LogisticRegressionfrom sklearn.svm import SVCfrom sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifierfrom sklearn.metrics import (accuracy_score, recall_score, f1_score,confusion_matrix, classification_report)from sklearn.model_selection import train_test_splitfrom typing import Dictclass ChatterClassifier:FEAT_EXCLUDE [label, t_start_ms]def __init__(self, random_state: int 42):self.random_state random_stateself.models: Dict {}self.metrics pd.DataFrame()self.best Noneself.feat_cols []def _feats(self, df: pd.DataFrame) - list:return [c for c in df.columns if c not in self.FEAT_EXCLUDE]def fit_compare(self, df: pd.DataFrame) - pd.DataFrame:self.feat_cols self._feats(df)X df[self.feat_cols].values.astype(float)y df[label].values.astype(int)Xtr, Xte, ytr, yte train_test_split(X, y, test_size0.25, stratifyy, random_stateself.random_state)specs {lr: LogisticRegression(max_iter1000),svm: SVC(probabilityTrue, random_stateself.random_state),rf: RandomForestClassifier(n_estimators300, random_stateself.random_state),gbdt: GradientBoostingClassifier(random_stateself.random_state),}rows []for name, m in specs.items():m.fit(Xtr, ytr)pred m.predict(Xte)self.models[name] mrows.append({model: name,acc: round(accuracy_score(yte,pred),4),recall_chatter: round(recall_score(yte,pred,pos_label1),4),f1: round(f1_score(yte,pred,pos_label1),4),})self.metrics pd.DataFrame(rows).sort_values(f1, ascendingFalse).reset_index(dropTrue)self.best self.metrics.iloc[0][model]self._Xte, self._yte Xte, ytereturn self.metricsdef predict_stream(self, df: pd.DataFrame) - pd.DataFrame:m self.models[self.best]X df[self.feat_cols].values.astype(float)out df.copy()out[pred] m.predict(X)if hasattr(m, predict_proba):out[prob_chatter] m.predict_proba(X)[:,1]return outdef cm(self) - np.ndarray:pred self.models[self.best].predict(self._Xte)return confusion_matrix(self._yte, pred)/detailsdetailssummary/summary特征判别图 (networkx)import networkx as nximport numpy as npimport pandas as pdfrom typing import Dictclass DecisionGraph:特征 - 类别 有向判别网def __init__(self):self.G nx.DiGraph()def build(self, importances: Dict[str, float]) - nx.DiGraph:self.G.clear()self.G.add_node(颤振, ntypeclass)self.G.add_node(稳定, ntypeclass)for feat, w in importances.items():self.G.add_node(feat, ntypefeature)target 颤振 if w 0 else 稳定self.G.add_edge(feat, 颤振, weightmax(0.0, w)*100)self.G.add_edge(feat, 稳定, weightmax(0.0, -w)*100 1.0)return self.Gdef top_edges(self, n: int 5) - pd.DataFrame:rows [{from:u,to:v,weight:round(d[weight],2)}for u,v,d in self.G.edges(dataTrue)]return pd.DataFrame(rows).sort_values(weight,ascendingFalse).head(n).reset_index(dropTrue)/detailsdetailssummary/summary可视化 (matplotlib networkx)import numpy as npimport pandas as pdimport matplotlib.pyplot as pltfrom pathlib import Pathimport networkx as nxfrom sklearn.decomposition import PCAplt.rcParams[font.sans-serif] [SimHei, DejaVu Sans]plt.rcParams[axes.unicode_minus] Falseclass ChatterVisualizer:def __init__(self, results_dir: str results):self.results_dir Path(results_dir)self.results_dir.mkdir(exist_okTrue)def signal_compare(self, stable_seg, chatter_seg):fig, ax plt.subplots(2,1, figsize(12,6), sharexTrue)for axx, seg, ttl, c in [(ax[0], stable_seg, 稳定切削 Fz, #27AE60),(ax[1], chatter_seg, 颤振切削 Fz, #E74C3C)]:s seg[sub]axx.plot(s[t_ms].values, s[fz].values, colorc, lw0.8)axx.set_ylabel(Fz (N))axx.set_title(ttl, fontsize12, fontweightbold)axx.grid(alpha0.3)ax[1].set_xlabel(时间 (ms))plt.tight_layout()plt.savefig(self.results_dir/signal_compare.png, dpi150, bbox_inchestight)plt.close()def spectrum_compare(self, sa, stable_fz, chatter_fz):r1, r2 sa.analyze(stable_fz), sa.analyze(chatter_fz)fig, ax plt.subplots(figsize(11,5))ax.semilogy(r1[f], r1[Pxx], color#27AE60, label稳定, lw1)ax.semilogy(r2[f], r2[Pxx], color#E74C3C, label颤振, lw1)ax.axvspan(800,1200, color#F39C12, alpha0.15)ax.set_xlabel(频率 (Hz)); ax.set_ylabel(PSD)ax.set_title(Welch频谱对比(橙区颤振特征频带), fontsize13, fontweightbold)ax.legend(); ax.grid(alpha0.3)plt.tight_layout()plt.savefig(self.results_dir/spectrum_compare.png, dpi150, bbox_inchestight)plt.close()def confusion(self, cm):fig, ax plt.subplots(figsize(5,5))im ax.imshow(cm, cmapBlues)for i in range(2):for j in range(2):ax.text(j,i,str(cm[i,j]),hacenter,vacenter,fontsize16, fontweightbold)ax.set_xticks([0,1]); ax.set_yticks([0,1])ax.set_xticklabels([稳定,颤振]); ax.set_yticklabels([稳定,颤振])ax.set_xlabel(预测); ax.set_ylabel(真实)ax.set_title(混淆矩阵, fontsize13, fontweightbold)plt.tight_layout()plt.savefig(self.results_dir/confusion_matrix.png, dpi150, bbox_inchestight)plt.close()def boundary(self, clf: ChatterClassifier, df):X df[clf.feat_cols].values.astype(float)y df[label].valuesp PCA(n_components2).fit_transform(X)m clf.models[clf.best]# 用2D重训一个轻量模型画边界from sklearn.ensemble import RandomForestClassifiermm RandomForestClassifier(n_estimators80, random_state0).fit(p, y)xx,yy np.meshgrid(np.linspace(p[:,0].min(),p[:,0].max(),200),np.linspace(p[:,1].min(),p[:,1].max(),200))z mm.predict(np.c_[xx.ravel(),yy.ravel()]).reshape(xx.shape)fig, ax plt.subplots(figsize(7,6))ax.contourf(xx,yy,z, alpha0.25, colors[#27AE60,#E74C3C], levels[-0.5,0.5,1.5])ax.scatter(p[y0,0],p[y0,1],c#27AE60,s8,label稳定)ax.scatter(p[y1,0],p[y1,1],c#E74C3C,s8,label颤振)ax.set_title(PCA降维决策边界, fontsize13, fontweightbold)ax.legend(); ax.grid(alpha0.2)plt.tight_layout()plt.savefig(self.results_dir/decision_boundary.png, dpi150, bbox_inchestight)plt.close()def feature_imp(self, imp: dict):items sorted(imp.items(), keylambda x: abs(x[1]))ks [k for k,_ in items]; vs[v for _,v in items]fig, ax plt.subplots(figsize(9,6))cols [#E74C3C if v0 else #27AE60 for v in vs]ax.barh(ks, vs, colorcols)ax.set_title(特征判别权重(红指向颤振), fontsize13, fontweightbold)ax.grid(axisx, alpha0.3)plt.tight_layout()plt.savefig(self.results_dir/feature_importance.png, dpi150, bbox_inchestight)plt.close()def stream(self, stream_df, chatter_start_ms):fig, ax plt.subplots(figsize(13,5))ax.plot(stream_df[t_start_ms], stream_df[prob_chatter],color#8E44AD, lw1)ax.axhline(0.5, colorgray, ls--, lw0.8)# 首次判定颤振hit stream_df[stream_df[pred]1]if len(hit):first hit[t_start_ms].iloc[0]ax.axvline(first, color#E74C3C, lw1.5)ax.text(first, 0.9, f 首次判颤振\n 提前{chatter_start_ms-first:.0f}ms,color#E74C3C, fontsize9)ax.set_xlabel(时间 (ms)); ax.set_ylabel(颤振概率)ax.set_title(滑窗实时判定流, fontsize13, fontweightbold)ax.set_ylim(0,1.05); ax.grid(alpha0.3)plt.tight_layout()plt.savefig(self.results_dir/realtime_stream.png, dpi150, bbox_inchestight)plt.close()def graph_plot(self, G):fig, ax plt.subplots(figsize(11,8))pos nx.spring_layout(G, seed42, k0.8)cm {class:#E74C3C,feature:#3498DB}nc [cm[d.get(ntype)] for _,d in G.nodes(dataTrue)]nx.draw_networkx_nodes(G,pos,node_colornc,node_size900,edgecolorsblack,linewidths0.5,axax,alpha0.9)ew[max(0.4,d[weight]/8) for _,_,d in G.edges(dataTrue)]nx.draw_networkx_edges(G,pos,widthew,arrowsTrue,arrowsize12,axax,alpha0.6)nx.draw_networkx_labels(G,pos,font_size8,axax)ax.set_title(特征→类别 判别网, fontsize12, fontweightbold)ax.axis(off)plt.tight_layout()plt.savefig(self.results_dir/decision_graph.png, dpi150, bbox_inchestight)plt.close()/detailsdetailssummary/summary合成铣削切削力数据 (numpy scipy)import numpy as npimport pandas as pdfrom pathlib import Pathfrom scipy.signal import chirpfrom typing import Optionalclass SyntheticMillingGenerator:稳定: 高斯噪声 切削基频颤振: 叠加再生颤振极限环(固有频~1kHz增长)带label与人工标注颤振起点def __init__(self, rng: Optional[np.random.RandomState] None):self.rng rng or np.random.RandomState(42)def generate(self, out_path: str milling_force.csv,fs: int 1000, seconds: int 20) - pd.DataFrame:n fs * secondst np.arange(n) / fs * 1000.0 # mscut_freq 120.0 # 切削频率Hzfz np.zeros(n)label np.zeros(n, dtypeint)# 前12s稳定, 后8s颤振stable_end fs * 12chatter_start_ms 12 * 1000phase 2*np.pi*cut_freq*np.arange(n)/fsbase 80 15*np.sin(phase)noise self.rng.no利用AI解决实际问题如果你觉得这个工具好用欢迎关注长安牧笛