首页/文章/ 详情

机械故障诊断信号的幅域分析 - 幅值概率密度函数 |在CWRU轴承数据上实战

1年前浏览733
脑袋chegnz

背景及摘要

很多初学者刚接触故障诊断可能觉得很简单,套用深度学习模型进行训练,分类准确率达到99%即可。

在写论文时,这样的确没问题。但是在工作或者在解决实际问题时,几乎是用不上。

最近在读时献江老师等人的《机械故障诊断及典型案例解析》这本书时,发现其内容真的很丰富,故障诊断方法讲解的详细全面,且有实际丰富的案例,更加容易消化理解各方法特点及其应用场景。


信号波形是某种物理量随时间变化的关系。机械振动信号处理的基本方法有幅域分析时域分析频域分析


在幅值上的各种统计处理通常称为幅域分析仅对波形的幅值用简单的方法进行统计分析, 如计算波形的最大值平均值有效值,或研究时域波形幅值的概率分布形式等,


信号在时间域内的变换或统计分析称为时域分析频域分析确定信号的频率结构,即信号中包含哪些频率成分,分析结果以频率为自变量各种物理量的谱线或曲线。


不同的分析方法是从不同的角度观察、分析信号,使信号处理的结果更加利于故障分析与诊断。


本节将介绍幅域分析方法中的幅值概率密度函数,包括其数学公式定义及CWRU轴承数据实战案例。


       
目录        
     

1.  幅值概率密度函数定义

2.  基于CWRU数据代码实战准备和数据介绍

   2.1 导入包

   2.2 定义CWRU数据读取函数

3.  轴承内圈故障幅值概率密度函数案例分析

   3.1 时域图绘制

   3.2 编程思路分析

4.  封装成一个plt_amp_prob_density_fun()函数

   4.1 轴承滚动体轴承故障幅值概率密度函数案例分析

   4.2 轴承外圈轴承故障幅值概率密度函数案例分析

   4.3 正常状态轴承故障幅值概率密度函数案例分析

5.  总结

6.  课后思考

1、幅值概率密度函数定义

随机信号的幅值概率密度函数表示信号的幅值落在某一个指定区间内的概率,幅值概率密度函数提供了随机信号沿幅值域分布的信息,是随机信号的主要统计特性参数之一。在图1中,              值落在                                之间的时间为              ,其总的观测时间为                ,则出现频率可以用              的值确定,当              趋于无穷大时,这一比值就趋于              值落在              到x+Δx之间的概率

              

当              趋于零时,就得到该点的幅值概率密度函数:

            

图1 时域波形及幅值概率密度函数

典型信号的时域波形和幅值概率密度函数如图2所示。根据随机过程理论,随机信号的幅值概率密度函数符合正态分布规律,而确定性信号如简谐信号的幅值概率密度函数则呈盆形曲线,如图2(a)所示。一般故障信号是随机信号简谐信号混合体,所以信号幅值概率密度函数的正态分布曲线上端出现盆型漏斗时[图2(b)],往往预示着系统存在故障征兆


2典型信号的时域波形和幅值概率密度函数

可以采用幅值概率密度函数直接进行设备状态的检测与诊断。图3是新旧机床变速箱的噪声分布规律,可见,新旧两个变速箱的分布规律有着明显的差异。这是因为机床齿轮箱中的零件由于磨损等原因使得配合间隙增大,在机床噪声中会出现大幅值周期性冲击成分,使得噪声信号的方差增加,分散度加大,甚至使曲线顶部变平或出现局部的凹形。

3 新旧机床变速箱噪声分布规律


2、基于CWRU数据代码实战准备和数据介绍

2.1 导入包





import scipy.io as scioimport numpy as npimport matplotlib.pyplot as pltfrom scipy import signal, fftpack, stats
       











from matplotlib import rcParams
config = {    "font.family": 'serif', # 衬线字体    "font.size": 10, # 相当于小四大小    "font.serif": ['SimSun'], # 宋体    "mathtext.fontset": 'stix', # matplotlib渲染数学字体时使用的字体,和Times New Roman差别不大    'axes.unicode_minus': False # 处理负号,即-号}rcParams.update(config)
       

2.2 定义CWRU数据读取函数

使用CWRU轴承数据进行分析,选取了4个mat文件,包括内圈故障、外圈故障、滚动体故障和正常状态。


文件名解释(以“1730_12k_0.007-InnerRace”为例):

1730:转速,单位rpm

12k:采样频率,Hz

0.007:故障大小,单位inch,1inch=25.4mm

InnerRace:代表为内圈故障


更多CWRU数据集介绍见

故障诊断学习|CWRU数据集介绍 第1期

故障诊断学习|CWRU数据集介绍 第2期

故障诊断学习|CWRU数据集介绍 第3期

故障诊断学习|CWRU数据集介绍 第4期


下面定义一个mat文件数据读取函数












def data_acquision(FilePath):    """    fun: 从cwru mat文件读取加速度数据    param file_path: mat文件绝对路径    return accl_data: 加速度数据,array类型    """    data = scio.loadmat(file_path)  # 加载mat数据    data_key_list = list(data.keys())  # mat文件为字典类型,获取字典所有的键并转换为list类型    accl_key = data_key_list[3]  # 获取'X108_DE_time'    accl_data = data[accl_key].flatten()  # 获取'X108_DE_time'所对应的值,即为振动加速度信号,并将二维数组展成一维数组    return accl_data
       

3、内圈故障幅值概率密度函数案例分析


           

3.1 时域图绘制








file_path = r'E:/研究生/pytorch/CSDN代码/fault_diagnosis_signal_processing/第4篇-包络谱/1730_12k_0.007-InnerRace.mat'xt = data_acquision(file_path)  #从mat文件获取振动加速度数据##--------绘制时域图-------##plt.figure(figsize=(12,3))plt.plot(xt)plt.show()print(xt)
       




[ 0.22269856  0.09323776 -0.14651649 ... -0.36125573  0.31138814  0.17055689]
       

面对这10s的时域信号,如何实现获取绝对幅值概率密度函数呢?

3.2 编程思路分析

  • step1:确定幅值绝对值最大值max_value

  • step2:在幅值区间[-max_value, max_value]内等分划分成n个小区间,区间长度为interval_len

  • step3:在第i个幅值小区间[-max_value+i * interval_len, max_value+(i+1) * interval_len]内,统计落入该区间的幅值个数count_num

  • step5:将n个区间统计的个数count_num形成一个列表count_num_list

  • step6:绘制柱状图,即得到幅值概率密度图













def interval_num_count(data, low, high):    '''    fun: 统计一维数据data落入某一个区间[low, high]内的数量    param low: 区间下限    param high: 区间上限    return count_num: 落入某一个区间[low, high]内的数量    '''    count_num = 0    for i in range(len(data)):        if data[i]>low and data[i]<high:            count_num += 1    return count_num
       











n = 10  # 分段数量xt = xt - np.mean(xt)                     # 去直流分量(也叫零均值化处理)max_value = np.abs( xt[np.argmax( np.abs(xt) )] )   # 确定幅值绝对值最大值max_valuecount_num_list = [] for i in range(n):    interval_len = max_value*2/n          # 在幅值区间[-max_value, max_value]内等分划分成n个小区间,区间长度为interval_len    low = -max_value + i*interval_len     # 第i个幅值小区间[-max_value+i * interval_len, max_value+(i+1) * interval_len]的下限    high = -max_value + (i+1)*interval_len    # 第i个幅值小区间[-max_value+i * interval_len, max_value+(i+1) * interval_len]的上限    count_num = interval_num_count(data=xt, low=low, high=high)  # 统计落入该区间的幅值个数count_num    count_num_list.append(count_num)      # 将n个区间统计的个数count_num形成一个列表count_num_listplt.bar(x=range( len(count_num_list) ), height=count_num_list) # 绘制柱状图,即得到幅值概率密度图
       

是不是已成雏形,有正态分布形状了。只是原始区间分成了10份,看不太出来


下面分成100份












n = 100xt = xt - np.mean(xt)                     # 去直流分量(也叫零均值化处理)max_value = np.abs( xt[np.argmax( np.abs(xt) )] )   # 确定幅值绝对值最大值max_valuecount_num_list = [] for i in range(n):    interval_len = max_value*2/n          # 在幅值区间[-max_value, max_value]内等分划分成n个小区间,区间长度为interval_len    low = -max_value + i*interval_len     # 第i个幅值小区间[-max_value+i * interval_len, max_value+(i+1) * interval_len]的下限    high = -max_value + (i+1)*interval_len    # 第i个幅值小区间[-max_value+i * interval_len, max_value+(i+1) * interval_len]的上限    count_num = interval_num_count(data=xt, low=low, high=high)  # 统计落入该区间的幅值个数count_num    count_num_list.append(count_num)      # 将n个区间统计的个数count_num形成一个列表count_num_listplt.bar(x=range( len(count_num_list) ), height=count_num_list) # 绘制柱状图,即得到幅值概率密度图
       

这下有正态分布那味了。不得不说这个形状还挺优美的,看着挺舒服。

4、封装成一个plt_amp_prob_density_fun()函数


           


















def plt_amp_prob_density_fun(data, n):    '''    fun: 绘制幅值概率密度函数    param data: 输入数据,1维array    param n: 分割成段数的数量    return: 绘制幅值概率密度函数    '''    max_value = np.abs( xt[np.argmax( np.abs(xt) )] ) #    count_num = []    for i in range(n):        interval = max_value*2/n              # 区间长度为interval_len        low = -max_value + i*interval         # 区间下限        high = -max_value + (i+1)*interval    # 区间上限        count = interval_num_count(data=xt, low=low, high=high)  # 统计落入该区间的幅值个数        count_num.append(count)    plt.bar(x=range( len(count_num) ), height=count_num)  # 绘制柱状图    plt.show()
       

4.1 滚动体故障轴承幅值概率密度函数案例分析

此时的正态分布形状比较瘦小,峭度K值大于3。K值大于3时,说明信号中冲击成分幅值增大。

4.2 外圈故障轴承幅值概率密度函数案例分析








file_path = r'E:/研究生/pytorch/CSDN代码/fault_diagnosis_signal_processing/第4篇-包络谱/1730_12k_0.007-OuterRace3.mat'xt = data_acquision(file_path)plt.figure(figsize=(12,3))plt.plot(xt)plt.show()n = 100      # 设定分成份数plt_amp_prob_density_fun(data=xt, n=n)
       

在顶部出现了第1节介绍时说的盆型漏斗。

4.3 正常状态轴承幅值概率密度函数案例分析








file_path = r'E:/研究生/pytorch/CSDN代码/fault_diagnosis_signal_processing/第4篇-包络谱/1730_48k_Normal.mat'xt = data_acquision(file_path)plt.figure(figsize=(12,3))plt.plot(xt)plt.show()n = 100      # 设定分成份数plt_amp_prob_density_fun(data=xt, n=n)
       

看着与标准正态分布挺像的。

5、总结


           

与正常状态轴承的幅值概率密度函数相比,其故障轴承状态有两种情况

  • 第一种情况:概率密度函数偏瘦小,此时峭度K大于3

  • 第二种情况:概率密度函数顶部出现盆型漏斗,该现象预示存在故障


       

       

       

       

6、课后思考


       

该方法的不足:

采集的样本点数数量与幅值概率密度函数成正比。因此对采样样本数量有有一定要求。

除了这个不足,你还发现哪些不足呢?


该程序存在两个不足:

  • 第一个是纵坐标没有归一化,因为概率之和应该为1。

  • 第二个是横坐标应该[-max_value, max_value],也就是概率密度函数图应该分布在纵坐标两边。

发动你的小脑袋来完善这个程序吧。


注明:

本文内容摘抄自《机械故障诊断及典型案例解析》,仅供学习交流,有侵权,烦请联系处理!

代码链接:https://github.com/HappyBoy-cmd/fault_diagnosis_signal_processing


参考文献
[1] 《机械故障诊断及典型案例解析,第2版,时献江等,2020年第一次印刷》

编辑:李正平

审核、校对:张泽明,张勇

如需转载,请后台联系作者

说明:部分图片来源网络,若有侵权,烦请联系处理

来源:故障诊断与python学习

振动python理论电机渲染
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2023-06-22
最近编辑:1年前
故障诊断与python学习
硕士 签名征集中
获赞 70粉丝 66文章 143课程 0
点赞
收藏
未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈