燕山大学
课 程 设 计 说 明 书
题目：机械振动信号分析及故障报警
学院（系）： 电气工程学院
年级专业： 10级仪表3班
学 号： 1001030201
学生姓名：
指导教师：
教师职称：
电气工程学院《课程设计》任务书
课程名称： “单片机原理及应用——数字信号处理”课程设计
院（系）：电气工程学院 基层教学单位： 自动化仪表系
说明：1、此表一式四份，系、指导教师、学生各一份，报送院教务科一份。
2、学生那份任务书要求装订到课程设计报告前面。
目录
第一章 摘要
第二章 总体设计方案
第三章 基本原理
第四章 MATLAB界面设计
第五章 各模块设计及程序
第六章 设计心得及总结
参考文献
第一章摘要
机械振动信号分析是现代机械故障诊断的一个有效方法。在诸多信号分析的手段中，小波分析与傅氏变换相结合的方法得到广泛应用。因为这种方法更适合于提取微弱机械振动的特征信号。
但是与其他分析工具一样，小波分析工具有自己的特点，如果不能正确使用，反而会影响对信号的正确分析。从本质上说，小波分析是用小波函数与被被分析的信号函数做一系列的互相关运算，因此选用小波函数不当会引起分析的误差或误判。
第二章 总体设计方案
对机械振动信号进行采样，把采样的数据进行时域和频域上的分析，包括FFT，功率谱，倒谱分析。提取时域波形指标如均值、峰峰值、峭度、偏度、脉冲因数等。以一种指标为标准，分析振动信号产生的变化。本次课设利用matlab软件，实现对机械振动信号时频域的分析以及故障的判断。因为频域分析特征值的提取较麻烦，这里我们用其中一种参数的计算量为标准来判断是否发生故障。
第三章 基本原理
3.1小波变换
与Fourier变换相比，小波变换是空间（时间）和频率的局部变换，因而能有效地从信号中提取信息。通过伸缩和平移等运算功能可对函数或信号进行多尺度的细化分析，解决了Fourier变换不能解决的许多困难问题。小波变换联系了应用数学、物理学、计算机科学、信号与信息处理、图像处理、地震勘探等多个学科。数学家认为，小波分析是一个新的数学分支，它是泛函分析、Fourier分析、样调分析、数值分析的完美结晶；信号和信息处理专家认为，小波分析是时间—尺度分析和多分辨分析的一种新技术，它在信号分析、语音合成、图像识别、计算机视觉、数据压缩、地震勘探、大气与海洋波分析等方面的研究都取得了有科学意义和应用价值的成果。 信号分析的主要目的是寻找一种简单有效的信号变换方法，使信号所包含的重要信息能显现出来。小波分析属于信号时频分析的一种，在小波分析出现之前，傅立叶变换是信号处理领域应用最广泛、效果最好的一种分析手段。傅立叶变换是时域到频域互相转化的工具，从物理意义上讲，傅立叶变换的实质是把这个波形分解成不同频率的正弦波的叠加和。正是傅立叶变换的这种重要的物理意义，决定了傅立叶变换在信号分析和信号处理中的独特地位。傅立叶变换用在两个方向上都无限伸展的正弦曲线波作为正交基函数，把周期函数展成傅立叶级数，把非周期函数展成傅立叶积分，利用傅立叶变换对函数作频谱分析，反映了整个信号的时间频谱特性，较好地揭示了平稳信号的特征。
小波变换是一种新的变换分析方法，它继承和发展了短时傅立叶变换局部化的思想，同时又克服了窗口大小不随频率变化等缺点，能够提供一个随频率改变的时间一频率窗口，是进行信号时频分析和处理的理想工具。它的主要特点是通过变换能够充分突出问题某些方面的特征，因此，小波变换在许多领域都得到了成功的应用，特别是小波变换的离散数字算法已被广泛用于许多问题的变换研究中。从此，小波变换越来越引起人们的重视，其应用领域来越来越广泛。
3.2 傅里叶变换
有限长序列可以通过离散傅里叶变换(DFT)将其频域也离散化成有限长序列.但其计算量太大,很难实时地处理问题,因此引出了快速傅里叶变换(FFT). 1965年，Cooley和Tukey提出了计算离散傅里叶变换（DFT）的快速算法，将DFT的运算量减少了几个数量级。从此，对快速傅里叶变换（FFT）算法的研究便不断深入，数字信号处理这门新兴学科也随FFT的出现和发展而迅速发展。根据对序列分解与选取方法的不同而产生了FFT的多种算法，基本算法是基２DIT和基２DIF。FFT在离散傅里叶反变换、线性卷积和线性相关等方面也有重要应用。
快速傅氏变换（FFT），是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。
设x(n)为N项的复数序列，由DFT变换，任一X（m）的计算都需要N次复数乘法和N-1次复数加法，而一次复数乘法等于四次实数乘法和两次实数加法，一次复数加法等于两次实数加法，即使把一次复数乘法和一次复数加法定义成一次“运算”（四次实数乘法和四次实数加法），那么求出N项复数序列的X（m）,即N点DFT变换大约就需要N^2次运算。当N=1024点甚至更多的时候，需要N2=1048576次运算，在FFT中，利用WN的周期性和对称性，把一个N项序列（设N=2k,k为正整数），分为两个N/2项的子序列，每个N/2点DFT变换需要（N/2）2次运算，再用N次运算把两个N/2点的DFT变换组合成一个N点的DFT变换。这样变换以后，总的运算次数就变成N+2（N/2）2=N+N2/2。继续上面的例子，N=1024时，总的运算次数就变成了525312次，节省了大约50%的运算量。而如果我们将这种“一分为二”的思想不断进行下去，直到分成两两一组的DFT运算单元，那么N点的DFT变换就只需要Nlog2N次的运算，N在1024点时，运算量仅有10240次，是先前的直接算法的1%，点数越多，运算量的节约就越大，这就是FFT的优越性。
离散傅里叶变换X(k)可看成是z变换在单位圆上的等距离采样值。同样，X(k)也可看作是序列傅氏变换
的采样，采样间隔为ωN=2π/N 。由此看出，离散傅里叶变换实质上是其频谱的离散频域采样，对频率具有选择性(ωk=2πk/N)，在这些点上反映了信号的频谱。
根据采样定律，一个频带有限的信号，可以对它进行时域采样而不丢失任何信息，FFT变换则说明对于时间有限的信号(有限长序列)，也可以对其进行频域采样，而不丢失任何信息。所以只要时间序列足够长，采样足够密，频域采样也就可较好地反映信号的频谱趋势，所以FFT可以用以进行连续信号的频谱分析。
第四章 MATLAB界面设计
两个画图按钮简化成了一个下拉菜单，坐标图一和二分别显示时域和频域波形，计算按钮进行计算并把数值放到下面四个文本框中。判断按钮判断是否故障并通过发送按钮把特征值传输到单片机中。
第五章各模块设计及程序
1、时域和频域分析及画图
%时域波形
% figure(1);
set(gcf,'CurrentAxes',handles.time);
plot(xdata);
xlabel('样本序号 n/个');
ylabel('电压 V/v');
%db10小波进行4层分解
%一维小波分解
[c,l]=wavedec(xdata,4,'db10');
%第一层细节信号的包络谱
y=hilbert(d1);
ydata=abs(y);
y=y-mean(y);
nfft=200;
p=abs(fft(ydata,nfft));
% figure(3);
set(gcf,'CurrentAxes',handles.frequency);
plot((0:nfft/2-1)/nfft*fs,p(1:nfft/2));
xlabel('频率 f/Hz');
ylabel('功率谱 p/w')；
故障信号的时域和频域图如下
正常信号的时域和频域图如下
可见其频域有较大的差别
2、计算部分
这里，我们调用matlab中已有的程序来计算需要得到的指标，所编程序如下
fs=10000;
fid=fopen('bearingout.dat','r');
load Normal.mat
fid=fopen('Normal.mat');
N=1000;
xdata=fread(fid,N,'int16');
xdata=fread(fid,N,'int16');
fclose(fid);
xdata=(xdata-mean(xdata))/std(xdata,1);
A=mean(xdata); % 平均值
set(handles.edit1,'string',num2str(A));
B=max(xdata); % 最大值
C=min(xdata); % 最小值
D=B-C; % 峰峰值
set(handles.edit2,'string',num2str(D));
E=kurtosis(xdata);%峭度
set(handles.edit3,'string',num2str(E));
F=skewness(xdata);%偏度;
set(handles.edit4,'string',num2str(F));
3、判断程序
这里我们根据两种不同信号偏度的不同来判断是否有故障，定义偏度<=-0.05时为正常信号。
global F
global J
global FF
if (FF<=-0.05)
J=0;
else J=1;
end
set(handles.edit6,'string',num2str(J))
Matlab中本实验所得数据计算界面如下：故障数据
正常数据:
第六章 设计心得及总结
刚开始设计时，因为对matlab的使用不是十分地了解，所以有很多的困难，觉得无从下手，经过这么多天的学习以及查阅资料，已经能初步完成一些简单的语言编程，这使我感到非常有成就感，使我体会到了自主学习的乐趣。通过这次课程设计，我对数字信号处理及傅里叶变换的实际应用有了深刻的理解，尤其加深了我对快速傅里叶变换的理解。在这次课程设计中，我还对小波变换在数字信号处理中的应用有了一定的理解，在对信号的细节进行分析时，小波变换是非常常见和方便的一种方法。
总之，经过这次数字信号处理的课程设计，使我收获颇多。不仅巩固了理论知识，增强了自主学习能力，还学会了团队合作和互相帮助，使我受益良多。
参考文献
1）《微型计算机控制系统》 赖寿宏，机械工业出版社（教材）
2）《单片机及应用》李大友，高等教育出版社（教材）
3）《信号处理原理及应用》 谢平 等 机械工业出版社（教材）
4)《Matlab程序设计及其在信号处理中的应用》聂祥飞等 西南交通大学出版社
附录一 matlab总程序
function varargout = keshe0628(varargin)
gui_Singleton = 1;
gui_State = struct('gui_Name', mfilename, ...
'gui_Singleton', gui_Singleton, ...
'gui_OpeningFcn', @keshe0628_OpeningFcn, ...
'gui_OutputFcn', @keshe0628_OutputFcn, ...
'gui_LayoutFcn', [] , ...
'gui_Callback', []);
if nargin && ischar(varargin{1})
gui_State.gui_Callback = str2func(varargin{1});
end
if nargout
[varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
gui_mainfcn(gui_State, varargin{:});
end
function keshe0628_OpeningFcn(hObject, eventdata, handles, varargin)
handles.output = hObject;
guidata(hObject, handles);
function varargout = keshe0628_OutputFcn(hObject, eventdata, handles)
function jisuan_Callback(hObject, eventdata, handles)
global xdata
global J
global F
global FF
A=mean(xdata); % 平均值
set(handles.edit1,'string',num2str(A));
B=max(xdata); % 最大值
C=min(xdata); % 最小值
D=B-C; % 峰峰值
set(handles.edit2,'string',num2str(D));
E=kurtosis(xdata);%峭度
set(handles.edit3,'string',num2str(E));
F=skewness(xdata);%偏度;
FF=F
set(handles.edit4,'string',num2str(F));
function edit1_Callback(hObject, eventdata, handles)
function edit1_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
set(hObject,'BackgroundColor','white');
end
function edit2_Callback(hObject, eventdata, handles)
function edit2_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
set(hObject,'BackgroundColor','white');
end
function edit3_Callback(hObject, eventdata, handles)
function edit3_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
set(hObject,'BackgroundColor','white');
end
function edit4_Callback(hObject, eventdata, handles)
function edit4_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
set(hObject,'BackgroundColor','white');
end
function edit5_Callback(hObject, eventdata, handles)
function edit5_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
set(hObject,'BackgroundColor','white');
end
function pushbutton2_Callback(hObject, eventdata, handles)
global F
global J
global FF
global ss
ss=serial('com1');
ss.BaudRate=2400;
ss.DataBits=8;
ss.Parity='none';
ss.StopBits=1;
ss.TimeOut=60;
ss.DataTerminalReady='off';
ss.RequestToSend='off';
ss.FlowControl='none';
ss.InputBufferSize=1000;
fopen(ss);
fwrite(ss,J);
fclose(ss);
clear ss
function pushbutton3_Callback(hObject, eventdata, handles)
global F
global J
global FF
if (FF<=-0.05)
J=0;
else J=1;
end
set(handles.edit6,'string',num2str(J))
function edit6_Callback(hObject, eventdata, handles)
function edit6_CreateFcn(hObject, eventdata, handles)
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
set(hObject,'BackgroundColor','white');
end
function time_CreateFcn(hObject, eventdata, handles)
function frequency_CreateFcn(hObject, eventdata, handles)
function huatu1_Callback(hObject, eventdata, handles)
global xdata
global F
global FF
global J
fs=10000;
load Abnormal.mat
fid=fopen('Abnormal.mat');
N=1000;
xdata=fread(fid,N,'int16');
fclose(fid);
xdata=(xdata-mean(xdata))/std(xdata,1);
set(gcf,'CurrentAxes',handles.time);
plot(xdata);
xlabel('样本序号 n/个');
ylabel('电压 V/v');
%db10小波进行4层分解
%一维小波分解
[c,l]=wavedec(xdata,4,'db10');
d4=wrcoef('d',c,l,'db10',4);
d3=wrcoef('d',c,l,'db10',3);
d2=wrcoef('d',c,l,'db10',2);
d1=wrcoef('d',c,l,'db10',1);
%第一层细节信号的包络谱
y=hilbert(d1);
ydata=abs(y);
y=y-mean(y);
nfft=200;
p=abs(fft(ydata,nfft));
set(gcf,'CurrentAxes',handles.frequency);
plot((0:nfft/2-1)/nfft*fs,p(1:nfft/2));
xlabel('频率 f/Hz');
ylabel('功率谱 p/w');
function huatu2_Callback(hObject, eventdata, handles)
global xdata
global F
global FF
global J
fs=10000;
load Normal.mat
fid=fopen('Normal.mat');
N=1000;
xdata=fread(fid,N,'int16');
fclose(fid);
xdata=(xdata-mean(xdata))/std(xdata,1);
set(gcf,'CurrentAxes',handles.time);
plot(xdata);
xlabel('样本序号 n/个');
ylabel('电压 V/v');
[c,l]=wavedec(xdata,4,'db10');
d4=wrcoef('d',c,l,'db10',4);
d3=wrcoef('d',c,l,'db10',3);
d2=wrcoef('d',c,l,'db10',2);
d1=wrcoef('d',c,l,'db10',1);
y=hilbert(d1);
ydata=abs(y);
y=y-mean(y);
nfft=200;
p=abs(fft(ydata,nfft));
set(gcf,'CurrentAxes',handles.frequency);
plot((0:nfft/2-1)/nfft*fs,p(1:nfft/2));
xlabel('频率 f/Hz');
ylabel('功率谱 p/w');
function Untitled_2_Callback(hObject, eventdata, handles)
function Untitled_1_Callback(hObject, eventdata, handles)
function Untitled_3_Callback(hObject, eventdata, handles)
function Untitled_4_Callback(hObject, eventdata, handles)
燕山大学专业综合训练评审意见表
