python的先进制造技术工业场景模拟第三十一篇:读取数控铣削振动采集数据，按粗加工，半精，精加工分组，对比振动幅值差异。

周四上午数控车间精加工区域。这批铝合金壳体粗铣→半精铣→精铣三道工序现场工艺员老马把振动传感器数据线从机床主轴拆下来粗加工时振动就大半精好了一些精加工时我以为很稳了结果表面出了振纹。我想知道每道工序的振动到底差多少、是不是精加工那刀振动其实没控制好。我点开采集软件导出的 CSV。这表里有什么老马问。每条记录一行有时间戳、主轴转速、进给、切削深度、振动加速度三轴数据我指着屏幕问题是全部混在一起没有按工序分组。粗加工切深 3mm半精 1mm精加工 0.3mm振动特征肯定不一样。你得先分组再对比。我就想干一件事老马说把粗铣、半精铣、精铣三组的振动幅值分别算出来对比均值和峰值画个箱线图再看看三组之间有没有显著差异。如果精加工振动其实比半精大那振纹就解释通了——不是刀具问题是工艺参数没匹配好。比如粗加工组振动均值 4.2 m/s²峰值 12.8半精均值 2.1峰值 6.3精加工均值 1.8峰值 5.7我接话如果精加工峰值 5.7 接近半精的 6.3说明精加工并没有真正静下来。按理说切深从 1mm 降到 0.3mm振动应该再降一半才对。没降到位要么是主轴径向跳动大要么是夹具刚性不够。对老马点头还有我想看振动幅值的分布形态——是正态还是偏态有没有异常尖峰如果有个别采样点振动突然飙到 15 m/s²那可能是刀具崩刃或者撞刀前兆。最好还能用统计检验确认三组振动差异是不是真的显著别光看均值差零点几就下结论。用 pandas 按工序分组numpy 算 RMS 和峰值matplotlib 画箱线图分布叠加scipy 做单因素 ANOVA 和事后 Tukey 检验scikit-learn 做异常振动检测networkx 建工序-刀具-振动等级的关联图我开工程数据自包含合成一批三工序振动数据下载就能跑。敲了行原型groups df.groupby(process)[vibration_rms].agg([mean,std,max])f_stat, p_val scipy.stats.f_oneway(*[g[vibration_rms] for _,g in df.groupby(process)])完整版用 OOP 封好我说一个类管数据加载一个类按工序分组一个类算振动指标一个类做统计检验一个类做异常检测一个类建关联图一个类出图。输出分组统计表、ANOVA 结果、Tukey 对比、异常清单存 results/。老马凑近看那以后看报告粗铣 RMS 4.23 m/s²半精 2.08精铣 1.91。ANOVA p0.001三组差异显著。但 Tukey 说半精和精铣之间 p0.31差异不显著——精加工振动没比半精好多少异常检测标出 3 个尖峰在精加工段查了是刀柄锥面有切屑导致径向跳动。换刀柄清洁后重测精铣 RMS 降到 0.87。对我接话振动不是一个数是分组对比统计检验异常定位。数字孪生里要建加工质量预测模型振动分组特征就是核心输入特征之一。一、实际应用场景真实痛点场景设定数控铣削车间主轴加装三轴加速度传感器每道工序粗铣/半精铣/精铣连续采集振动加速度数据。工艺员需要量化评估每道工序的振动水平判断精加工是否真正受控定位异常振动源优化切削参数。现场原话叙事化不是我们不想控振动老马说是以前没有数据。操作工凭耳朵听——这个声音不对然后凭经验降转速或者减切深。可耳朵不能写进工艺卡。上次客户投诉表面振纹我查了切削参数、刀具磨损、夹具刚性都看起来正常。后来装了振动传感器才发现精加工段有几个尖峰RMS 其实跟半精差不多。没人注意到因为以前根本没量化过。还有参数优化问题老马补充粗加工我想开快点但怕振动太大伤主轴。半精想找振动-效率平衡点。精加工必须稳。这三道工序的振动容忍阈值应该不一样但现在全靠一把尺子量。核心矛盾全工序混合振动流水 与 按工序分组对比 统计显著性检验 异常尖峰定位 工序-刀具关联分析 之间的断层。二、痛点分析映射到滨州职业学院《先进制造技术》课程模型《先进制造技术》模块 本篇痛点对应数控加工与CAD/CAM技术切削参数优化、加工稳定性 工序级振动分组对比与参数匹配先进制造技术基础信号分析、加工质量评估 振动幅值统计与异常检测智能制造与数字孪生加工过程状态监控 振动特征作为质量预测输入FMS与先进生产管理设备健康管理、预防性维护 异常振动预警与根因追溯先进制造新模式数据驱动工艺迭代 统计检验指导参数优化决策一句话总结我们需要一个数控铣削振动按工序分组对比分析程序用pandas 分组聚合numpy 算 RMS/峰值scipy 做 ANOVATukeyscikit-learn 异常检测networkx 建工序-刀具关联图实现从混合振动流到工序分组量化 显著性检验 异常定位。三、核心逻辑讲解大白话3.1 问题本质把振动想成开车过路面把铣削加工想成开车走三段路* 粗铣 走烂泥路切深大、进给快车颠得厉害* 半精铣 走碎石路切深中等颠得轻一些* 精铣 走高速路切深小应该很平* 振动幅值 车身颠簸程度加速度越大越颠* RMS 平均颠簸感不是某个瞬间是整段路的综合感受* 峰值 突然一个大坑崩刃、撞刀前兆* ANOVA 三段路的颠簸程度到底是不是真的不一样还是只是随机波动* Tukey 具体哪两段路差异显著烂泥路和高速路肯定不同但碎石路和高速路呢3.2 业务逻辑 → 代码映射导入振动采集数据│▼ VibeLoader (pandas)读取 CSVtimestamp, process, tool_id, spindle_speed,feed_rate, cut_depth, acc_x, acc_y, acc_z校验合成振动合成幅值│▼ ProcessGrouper (pandas)按工序分组粗加工(cut_depth≥2) / 半精(0.5≤cut_depth2) / 精加工(0.5)或按 process 列直接分组│▼ VibeCalculator (numpy)振动指标RMS √(Σ(acc²)/N)峰值 max(|acc|)峰峰值 max - min峭度 Σ((acc-mean)/std)⁴/N│▼ VibeStatistics (scipy)统计检验单因素 ANOVA: 三组 RMS 是否有显著差异Tukey HSD: 两两对比Shapiro-Wilk: 每组正态性检验│▼ VibeAnomalyDetector (scikit-learn)异常检测IsolationForest 找振动尖峰或 3σ 规则标出超限点│▼ VibeGraph (networkx)工序-刀具-振动等级关联图节点工序/刀具/振动等级边关联, 边权该组合的平均RMS│▼ VibeVisualizer (matplotlib)可视化1. 三组箱线图(RMS对比)2. 三组 RMS 分布叠加直方图3. 时域波形图(标注异常点)4. 工序-刀具关联图5. ANOVATukey 结果文本表│▼ SyntheticVibeGenerator (numpy)合成数据三工序×多刀具RMS 粗半精精含异常尖峰(崩刃模拟)3.3 为什么不能只看均值视角 问题整体均值 2.8 m/s² 掩盖精加工其实 4.0 的问题看峰值 一个尖峰拉高整体误判分组RMS 每道工序的真实振动水平ANOVA 差异是真的还是随机波动3.4 分析前后对比维度 凭耳朵听 本程序振动水平 声音大 RMS 定量分组工序对比 感觉 ANOVA 显著性异常尖峰 事后发现 实时检测定位参数优化 拍脑袋 统计驱动维护预警 坏了再说 振动趋势预警四、OOP 代码实现4.1 项目结构cnc_vibe_analysis/├── cnc_vibe_analysis/│ ├── __init__.py│ ├── vibe_loader.py # 振动数据加载│ ├── process_grouper.py # 工序分组│ ├── vibe_calculator.py # 振动指标计算│ ├── vibe_statistics.py # 统计检验(scipy)│ ├── vibe_anomaly.py # 异常检测(sklearn)│ ├── vibe_graph.py # 工序-刀具关联图(networkx)│ ├── visualizer.py # 可视化│ └── synthetic_data.py # 合成数据├── tests/│ ├── __init__.py│ └── test_vibe_analysis.py├── results/│ ├── vibe_boxplot.png│ ├── vibe_dist_overlay.png│ ├── vibe_timeseries.png│ ├── vibe_network.png│ ├── process_stats.csv│ ├── anova_result.txt│ ├── anomaly_points.csv│ └── vibe_report.txt└── run_vibe_analysis.py4.2 核心源码detailssummary/summary数控铣削振动数据加载器import pandas as pdfrom pathlib import Pathfrom typing import Optionalclass VibeLoader:加载振动采集CSVdef __init__(self, filepath: str vibe_data.csv,encoding: str utf-8):self.filepath Path(filepath)self.encoding encodingself._raw: Optional[pd.DataFrame] Nonedef load(self) - pd.DataFrame:if not self.filepath.exists():raise FileNotFoundError(f文件不存在: {self.filepath})self._raw pd.read_csv(self.filepath, encodingself.encoding)rename {}for tgt, al in {timestamp: [timestamp, 时间, time],process: [process, 工序, proc],tool_id: [tool_id, 刀具, tool],spindle_speed: [spindle_speed, 主轴转速, rpm],feed_rate: [feed_rate, 进给, feed],cut_depth: [cut_depth, 切深, depth],acc_x: [acc_x, 加速度X, ax],acc_y: [acc_y, 加速度Y, ay],acc_z: [acc_z, 加速度Z, az],}.items():if tgt not in self._raw.columns:for a in al:if a in self._raw.columns:rename[a] tgtbreakself._raw self._raw.rename(columnsrename)# 数值列num_cols [spindle_speed, feed_rate, cut_depth,acc_x, acc_y, acc_z]for c in num_cols:if c in self._raw.columns:self._raw[c] pd.to_numeric(self._raw[c], errorscoerce)# 合成幅值if all(c in self._raw.columns for c in [acc_x, acc_y, acc_z]):self._raw[amplitude] (self._raw[acc_x]**2 self._raw[acc_y]**2 self._raw[acc_z]**2).pow(0.5).round(4)# 工序列处理if process not in self._raw.columns and cut_depth in self._raw.columns:self._raw[process] self._raw[cut_depth].apply(lambda d: 精加工 if d 0.5else (半精加工 if d 2.0 else 粗加工))self._raw[process] self._raw.get(process, 未知).fillna(未知).str.strip()return self._raw.reset_index(dropTrue)/detailsdetailssummary/summary工序分组import numpy as npimport pandas as pdfrom typing import Optional, Dictclass ProcessGrouper:按工序分组振动数据def __init__(self):passdef group(self, df: pd.DataFrame) - Dict[str, pd.DataFrame]:返回 {工序名: 子DataFrame}groups {}for proc, g in df.groupby(process):groups[proc] g.reset_index(dropTrue)return groupsdef group_order(self) - list:标准工序顺序return [粗加工, 半精加工, 精加工]def ordered_groups(self, df: pd.DataFrame) - list:按标准顺序返回 [(工序名, DataFrame)]grps self.group(df)result []for p in self.group_order():if p in grps:result.append((p, grps[p]))# 加其他非标准工序for p, g in grps.items():if p not in self.group_order():result.append((p, g))return result/detailsdetailssummary/summary振动指标计算 (numpy)import numpy as npimport pandas as pdfrom typing import Optionalclass VibeCalculator:计算RMS/峰值/峭度等振动指标def __init__(self):passdef rms(self, signal: np.ndarray) - float:均方根值if len(signal) 0:return 0.0return round(float(np.sqrt(np.mean(signal**2))), 4)def peak(self, signal: np.ndarray) - float:峰值(绝对值最大)if len(signal) 0:return 0.0return round(float(np.max(np.abs(signal))), 4)def peak_to_peak(self, signal: np.ndarray) - float:峰峰值if len(signal) 0:return 0.0return round(float(np.max(signal) - np.min(signal)), 4)def kurtosis(self, signal: np.ndarray) - float:峭度(冲击特征)if len(signal) 4:return 0.0mu np.mean(signal)sd np.std(signal, ddof1)if sd 0:return 0.0return round(float(np.mean(((signal - mu) / sd) ** 4)), 4)def process_stats(self, df: pd.DataFrame) - pd.DataFrame:按工序计算所有指标rows []for proc, g in df.groupby(process):amp g.get(amplitude, pd.Series([0]))if amplitude not in g.columns:# 用三轴合成if all(c in g.columns for c in [acc_x,acc_y,acc_z]):amp (g[acc_x]**2 g[acc_y]**2 g[acc_z]**2).pow(0.5)else:continuerows.append({process: proc,count: len(g),rms: self.rms(amp.values),peak: self.peak(amp.values),peak_to_peak: self.peak_to_peak(amp.values),kurtosis: self.kurtosis(amp.values),mean: round(float(np.mean(amp.values)), 4),std: round(float(np.std(amp.values, ddof1)), 4),max_tool: g[tool_id].mode().iloc[0] if tool_id in g.columns and len(g[tool_id].mode())0 else N/A,})return pd.DataFrame(rows).sort_values(rms, ascendingFalse).reset_index(dropTrue)/detailsdetailssummary/summary振动统计检验 (scipy)import numpy as npimport pandas as pdfrom scipy import statsfrom typing import Dict, Optional, Listclass VibeStatistics:ANOVA Tukey 正态性检验def __init__(self, alpha: float 0.05):self.alpha alphadef normality_test(self, df: pd.DataFrame,value_col: str amplitude) - pd.DataFrame:每组正态性检验rows []for proc, g in df.groupby(process):vals g[value_col].dropna().valuesif len(vals) 3:rows.append({process: proc, shapiro_stat: np.nan,shapiro_p: np.nan, is_normal: False})continueif len(vals) 5000:vals vals[:5000]stat, p stats.shapiro(vals)rows.append({process: proc,shapiro_stat: round(float(stat), 4),shapiro_p: round(float(p), 4),is_normal: bool(p self.alpha),})return pd.DataFrame(rows)def anova(self, df: pd.DataFrame,value_col: str amplitude) - Dict:单因素ANOVAgroups []for _, g in df.groupby(process):vals g[value_col].dropna().valuesif len(vals) 1:groups.append(vals)if len(groups) 2:return {f_statistic: np.nan, p_value: np.nan,significant: False}f_stat, p_val stats.f_oneway(*groups)return {f_statistic: round(float(f_stat), 4),p_value: round(float(p_val), 4),significant: bool(p_val self.alpha),}def tukey_hsd(self, df: pd.DataFrame,value_col: str amplitude) - pd.DataFrame:Tukey HSD 两两对比(简化版: 用 pairwise t-test Bonferroni)procs list(df[process].unique())rows []from itertools import combinationsn_comparisons len(list(combinations(procs, 2)))for a, b in combinations(procs, 2):a_vals df[df[process] a][value_col].dropna().valuesb_vals df[df[process] b][value_col].dropna().valuesif len(a_vals) 2 or len(b_vals) 2:continuestat, p stats.ttest_ind(a_vals, b_vals, equal_varFalse)# Bonferroni correctionp_adj min(1.0, p * n_comparisons)rows.append({group_a: a, group_b: b,mean_a: round(float(np.mean(a_vals)), 4),mean_b: round(float(np.mean(b_vals)), 4),t_statistic: round(float(stat), 4),p_value_raw: round(float(p), 6),p_value_bonf: round(float(p_adj), 6),significant: bool(p_adj self.alpha),})return pd.DataFrame(rows).sort_values(p_value_bonf).reset_index(dropTrue)/detailsdetailssummary/summary振动异常检测 (scikit-learn 3σ)import numpy as npimport pandas as pdfrom sklearn.ensemble import IsolationForestfrom typing import Optionalclass VibeAnomalyDetector:检测振动尖峰和异常模式def __init__(self, method: str both,contamination: float 0.05,sigma_threshold: float 3.0,random_state: int 42):self.method methodself.contamination contaminationself.sigma_threshold sigma_thresholdself.random_state random_statedef detect_3sigma(self, df: pd.DataFrame,value_col: str amplitude) - pd.DataFrame:3σ 规则异常检测out df.copy()out[z_score] np.nanout[is_anomaly_3sigma] Falsefor proc, g in df.groupby(process):idx g.indexvals g[value_col].valuesif len(vals) 3:continuemu np.mean(vals)sd np.std(vals, ddof1)if sd 0:continuez np.abs((vals - mu) / sd)out.loc[idx, z_score] zout.loc[idx, is_anomaly_3sigma] z self.sigma_thresholdreturn outdef detect_isolation_forest(self, df: pd.DataFrame,value_col: str amplitude) - pd.DataFrame:IsolationForest 异常检测out df.copy()iso IsolationForest(contaminationself.contamination,random_stateself.random_state)X out[[value_col]].fillna(0)pred iso.fit_predict(X)out[is_anomaly_if] pred -1out[if_score] iso.decision_function(X)return outdef detect(self, df: pd.DataFrame,value_col: str amplitude) - pd.DataFrame:综合检测out self.detect_3sigma(df, value_col)out self.detect_isolation_forest(out, value_col)out[is_anomaly] out[is_anomaly_3sigma] | out[is_anomaly_if]return outdef summary(self, df: pd.DataFrame) - pd.DataFrame:异常摘要sub df[df[is_anomaly]]if sub.empty:return pd.DataFrame()rows []for proc, g in sub.groupby(process):rows.append({process: proc,anomaly_count: len(g),anomaly_ratio: round(len(g) / max(1, len(df[df[process]proc])), 4),max_amplitude: round(float(g[amplitude].max()), 4),mean_amplitude: round(float(g[amplitude].mean()), 4),})return pd.DataFrame(rows).sort_values(anomaly_count, ascendingFalse).reset_index(dropTrue)/detailsdetailssummary/summary工序-刀具-振动等级关联图 (networkx)import networkx as nximport pandas as pdimport numpy as npfrom typing import Optionalclass VibeGraph:建工序-刀具-振动等级三层关系网def __init__(self):self.G nx.Graph()def build(self, df: pd.DataFrame) - nx.Graph:self.G.clear()for proc, pg in df.groupby(process):pnode fP:{proc}self.G.add_node(pnode, ntypeprocess)for tool, tg in pg.groupby(tool_id):tnode fT:{tool}self.G.add_node(tnode, ntypetool)self.G.add_edge(pnode, tnode,weightround(float(tg[amplitude].mean()), 2))# 振动等级rms np.sqrt(np.mean(tg[amplitude].values**2))if rms 1.5:level 低elif rms 3.0:level 中else:level 高lnode fL:{level}self.G.add_node(lnode, ntypelevel)self.G.add_edge(tnode, lnode, weightround(rms, 2))return self.Gdef high_vibe_tools(self) - pd.DataFrame:高振动刀具清单rows []for n, d in self.G.nodes(dataTrue):if d.get(ntype) tool:ws [self.G[n][nb][weight] for nb in self.G.neighbors(n)]levels [nb for nb in self.G.neighbors(n)if self.G.nodes[nb].get(ntype) level]rows.append({tool: n.split(:)[-1],avg_vibe: round(sum(ws)/len(ws), 2) if ws else 0,level: levels[0].split(:)[-1] if levels else 未知,})df pd.DataFrame(rows)if df.empty:return dfreturn df.sort_values(avg_vibe, ascendingFalse).reset_index(dropTrue)/detailsdetailssummary/summary可视化 (matplotlib)import numpy as npimport pandas as pdimport matplotlib.pyplot as pltfrom pathlib import Pathplt.rcParams[font.sans-serif] [SimHei, DejaVu Sans]plt.rcParams[axes.unicode_minus] Falseclass VibeVisualizer:def __init__(self, results_dir: str results):self.results_dir Path(results_dir)self.results_dir.mkdir(exist_okTrue)def boxplot(self, df, value_colamplitude):三组箱线图procs df[process].unique()fig, ax plt.subplots(figsize(10, 6))data [df[df[process]p][value_col].values for p in procs]bp ax.boxplot(data, labelsprocs, patch_artistTrue,showmeansTrue, meanlineTrue)colors [#E74C3C, #F39C12, #27AE60, #3498DB, #9B59B6]for i, patch in enumerate(bp[boxes]):patch.set_facecolor(colors[i % len(colors)])patch.set_alpha(0.7)ax.set_ylabel(振动幅值 (m/s²))ax.set_title(各工序振动幅值箱线图, fontsize14, fontweightbold)ax.grid(axisy, alpha0.3)plt.tight_layout()plt.savefig(self.results_dir / vibe_boxplot.png,dpi150, bbox_inchestight)plt.close()def dist_overlay(self, df, value_colamplitude):分布叠加直方图procs df[process].unique()fig, ax plt.subplots(figsize(11, 6))colors [#E74C3C, #F39C12, #27AE60, #3498DB, #9B59B6]for i, p in enumerate(procs):vals df[df[process]p][value_col].dropna()if len(vals) 0:ax.hist(vals, bins30, densityTrue, alpha0.5,colorcolors[i % len(colors)], labelp, edgecolorwhite)ax.set_xlabel(振动幅值 (m/s²))ax.set_ylabel(密度)ax.set_title(各工序振动幅值分布叠加, fontsize14, fontweightbold)ax.legend()ax.grid(alpha0.3)plt.tight_layout()plt.savefig(self.results_dir / vibe_dist_overlay.png,dpi150, bbox_inchestight)plt.close()def timeseries(self, df, value_colamplitude,anomaly_colis_anomaly, max_points2000):时域波形(标注异常)利用AI解决实际问题如果你觉得这个工具好用欢迎关注长安牧笛