xyz机械故障诊断大作业

作业名称：FFT滚动轴承故障诊断
院系：机械工程系
学号：
姓名：
指导教师：
西南交通大学峨眉校区
摘要
滚动轴承是旋转机械的主要损失之一，在以往的动检工作中，我们对旋转机械滚动轴承强烈震动原因分析不足，不能满足设备维修工作的需要。

所以要定期对旋转机械进行动态监测，根据所测数据做出诊断分析，及时发现滚动轴承强烈震动情况。

通过FFT方法分析轴承的信号图，对滚动轴承振动的产生原因进行深入分析，不断总结经验，提高故障分析能力，掌握造成滚动轴承强烈振动的原因，及时消除振动，为设备安全提供可行性措施。

关键词：轴承振动FFT
一、快速傅里叶变换(FFT)的原理及公式
非周期性连续时间信号x(t)的傅里叶变换可以表示为
式中计算出来的是信号x(t)的连续频谱。

但是，在实际的控制系统中能够得到的是连续信号x(t)的离散采样值x(nT)。

因此需要利用离散信号x(nT)来计算信号x(t)的频谱。

有限长离散信号x(n)，n=0，1，…，N-1的DFT定义为：
可以看出，DFT需要计算大约N2次乘法和N2次加法。

当N较大时，这个计算量是很大的。

利用WN的对称性和周期性，将N点DFT分解为两个N／2点的 DFT，这样两个N／2点DFT总的计算量只是原来的一半，即(N／2)2+(N／2)2=N2／2，这样可以继续分解下去，将N／2再分解为N／4点 DFT等。

对于N=2m 点的DFT 都可以分解为2点的DFT，这样其计算量可以减少为(N／2)log2N次乘法和
Nlog2N次加法。

图1为FFT与DFT-所需运算量与计算点数的关系曲线。

由图可以明显看出FFT算法的优越性。

将x(n)分解为偶数与奇数的两个序列之和，即
x1(n)和x2(n)的长度都是N／2，x1(n)是偶数序列，x2(n)是奇数序列，则
其中X1(k)和X2(k)分别为x1(n)和x2(n)的N／2点DFT。

由于X1(k)和X2(k)均以N／2为周期，且WN k+N/2=-WN k，所以X(k)又可表示为：
上式的运算可以用图2表示，根据其形状称之为蝶形运算。

依此类推，经过m-1次分解，最后将N点DFT分解为N／2个两点DFT。

图3为8点FFT的分解流程。

FFT算法的原理是通过许多小的更加容易进行的变换去实现大规模的变换，降低了运算要求，提高了与运算速度。

二、滚动轴承运行时的实时数据分析
总共给出了8组通过现场测试得到的滚动轴承运行时的实时数据，其中，2组正常轴承数据（normal1.mat、normal2.mat），2组内圈故障数据（inner-race1.mat、inner-race2.mat），2组外圈故障数据（outer-race1.mat、outer-race2.mat），2组滚动体故障数据（ball1.mat、ball2.mat）。

以上8组信号采样频率均为12000Hz。

从给定的8组数据中，任意选取1组正常数据，2组故障数据（故障类型不同），利用matlab软件编程，对所选取的信号进行分析。

图像
（1）
（2）
（3）
分析信号频谱图，正常数据图1的频率分布在1000HZ左右，故障组图2及图3的频率分布在0-4000HZ，相比故障组1不同的是，在3000HZ左右分布的比较多。

通过频率的分析，正常组因为频率分布较低因此满足轴承的使用，故障组图2和故障组图3是不能满足使用要求的。

分析信号时域图，正常组图1的幅值基本在-0.2到0.2之间且以0对称分布，故障组图2的幅值基本在-0.5到0.5之间且以0对称分布，故障组图3的幅值基本分布在-0.1到0.1之间且以0对称分布。

因为在滚动轴承实际使用中，若振动的幅值太大，有可能会发生事故，因此要求滚动轴承的幅值要尽量小。

结论：利用FFT方法，通过对滚动轴承故障诊断分析，找出滚动轴承在运行当中的各种规律，对比正常工作与非正常工作时的时域图及频谱图，明白了滚动轴承正常工作时频率比较稳定，幅值也较稳定，所以根据滚动轴承的时域图及频谱图便可检测其是否正常工作。

参考文献:
[1] 陆爽,田野.滚动轴承故障特征识别的时频分析方法研究[J].机床与液压, 2005(6): 183-185.
[2] 李军伟,韩捷.小波包-双谱分析和Hilbert-双谱分析的滚动轴承故障诊断方法对比研究[J].中国工程机械学报,2005(3):297-301.
[3] 仇学青,张鑫.滚动轴承故障诊断研究的国内现状与发展方向[J].煤矿机械, 2007(6): 6-8.
[4] 高成,董长虹,郭磊,等.Matlab小波分析与应用[M].北京:国防工业出版社, 2007:6-7.
附：编程程序
图1：
(1)x=X097_BA_time;%信号数组
subplot(2,1,1);
plot(x);%时域波形
xlabel('时间序列');
ylabel('幅值');
title('信号时域图');
fs=12000;%采样频率
N=length(x);
n=0:N-1;
y=fft(x,N);%进行fft变换
m=abs(y(1:N/2))*2/N;%求信号的真实幅值f=n*fs/N; %进行对应的频率转换subplot(2,1,2)
stem(f(1:N/2),m(1:N/2));%绘出频谱图xlabel('频率/Hz');
ylabel('幅值');
title('信号频谱图');
grid on;
图2：
(2)x=X278_BA_time;%信号数组
subplot(2,1,1);
plot(x);%时域波形
xlabel('时间序列');
ylabel('幅值');
title('信号时域图');
fs=12000;%采样频率
N=length(x);
n=0:N-1;
y=fft(x,N);%进行fft变换
m=abs(y(1:N/2))*2/N;%求信号的真实幅值
f=n*fs/N; %进行对应的频率转换
subplot(2,1,2)
stem(f(1:N/2),m(1:N/2));%绘出频谱图
xlabel('频率/Hz');
ylabel('幅值');
title('信号频谱图');
grid on;
图3：
(3)x=X294_BA_time;%信号数组
subplot(2,1,1);
plot(x);%时域波形
xlabel('时间序列');
ylabel('幅值');
title('信号时域图');
fs=12000;%采样频率
N=length(x);
n=0:N-1;
y=fft(x,N);%进行fft变换
m=abs(y(1:N/2))*2/N;%求信号的真实幅值f=n*fs/N; %进行对应的频率转换subplot(2,1,2)
stem(f(1:N/2),m(1:N/2));%绘出频谱图xlabel('频率/Hz');
ylabel('幅值');
title('信号频谱图');
grid on;。