数字信号-机械振动信号分析及故障报警.doc

燕山大学课 程 设 计 说 明 书题目：机械振动信号分析及故障报警学院（系）： 电气工程学院 年级专业： 10级仪表3班 学 号： 1001030201 学生姓名： 指导教师： 教师职称： 电气工程学院课程设计任务书课程名称： “单片机原理及应用数字信号处理”课程设计 院（系）：电气工程学院 基层教学单位： 自动化仪表系学号1001030201学生姓名（专业）班级自动化仪表3班设计题目机械振动信号分析及故障报警设计技术参数1、了解加速度传感器工作原理。2、掌握单片机串口通信协议。3、连接硬件实验线路，实现所要求的功能。4、基于matlab信号处理工具箱实现功能5、完成设计，提交课程设计报告。设计要求1、利用加速度传感器采集轴承正常与故障状态下的振动信号，基于matlab对振动信号进行时域和频域(包括fft，功率谱，倒谱)分析，提取时域波形指标如均值、峰峰值、峭度、偏度等。设计相应的信号分析及显示界面。2、基于振动信号分析结果，通过串口发送命令给单片机系统，根据振提取的特征指标判断轴承是否故障，若有故障驱动蜂鸣器报警。3、扩展：也控制相应的数码管或led发光显示相应振动信号参数(幅度、频率参数等)。工作量软件编程与硬件调试相结合，绘制设计流程图，编制相应软件界面，实现单片机控制与信号处理任务的综合应用参考资料1）微型计算机控制系统 赖寿宏，机械工业出版社（教材）2）单片机及应用李大友，高等教育出版社（教材）3）信号处理原理及应用 谢平 等 机械工业出版社（教材）4)matlab程序设计及其在信号处理中的应用聂祥飞等 西南交通大学出版社周次第一周第二周应完成内容熟悉伟福单片机编程环境，调试单片机各基本功能模块；熟悉matlab信号处理工具箱，信号处理系统基本功能模块学习和调试单片机系统与信号处理系统综合进行硬件调试，撰写课程设计报告指导教师签字基层教学单位主任签字说明：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变换大约就需要n2次运算。当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=2k/n)，在这些点上反映了信号的频谱。 根据采样定律，一个频带有限的信号，可以对它进行时域采样而不丢失任何信息，fft变换则说明对于时间有限的信号(有限长序列)，也可以对其进行频域采样，而不丢失任何信息。所以只要时间序列足够长，采样足够密，频域采样也就可较好地反映信号的频谱趋势，所以fft可以用以进行连续信号的频谱分析。第4章 matlab界面设计两个画图按钮简化成了一个下拉菜单，坐标图一和二分别显示时域和频域波形，计算按钮进行计算并把数值放到下面四个文本框中。判断按钮判断是否故障并通过发送按钮把特征值传输到单片机中。第5章 各模块设计及程序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.matfid=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 fglobal jglobal ff if (ff=-0.05) j=0;else j=1;endset(handles.edit6,string,num2str(j)matlab中本实验所得数据计算界面如下：故障数据正常数据:第6章 设计心得及总结刚开始设计时，因为对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(varargin1) gui_state.gui_callback = str2func(varargin1);endif nargout varargout1:nargout = gui_mainfcn(gui_state, varargin:);else gui_mainfcn(gui_state, varargin:);endfunction 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 xdataglobal jglobal fglobal ffa=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=fset(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);endfunction edit2_callback(hobject, eventdata, handles)function edit2_createfcn(hobject, eventdata, handles)if ispc & isequal(get(hobject,backgroundcolor), get(0,defaultuicontrolbackgroundcolor) set(hobject,backgroundcolor,white);endfunction edit3_callback(hobject, eventdata, handles)function edit3_createfcn(hobject, eventdata, handles)if ispc & isequal(get(hobject,backgroundcolor), get(0,defaultuicontrolbackgroundcolor) set(hobject,backgroundcolor,white);endfunction edit4_callback(hobject, eventdata, handles)function edit4_createfcn(hobject, eventdata, handles)if ispc & isequal(get(hobject,backgroundcolor), get(0,defaultuicontrolbackgroundcolor) set(hobject,backgroundcolor,white);endfunction edit5_callback(hobject, eventdata, handles)function edit5_createfcn(hobject, eventdata, handles)if ispc & isequal(get(hobject,backgroundcolor), get(0,defaultuicontrolbackgroundcolor) set(hobject,backgroundcolor,white);endfunction pushbutton2_callback(hobject, eventdata, handles)global fglobal jglobal ffglobal ssss=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 ssfunction pushbutton3_callback(hobject, eventdata, handles)global fglobal jglobal ffif (ff=-0.05) j=0;else j=1;endset(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);endfunction time_createfcn(hobject, eventdata, handles)function frequency_createfcn(hobject, eventdata, handles)function huatu1_callback(hobject, eventdata, handles)global xdataglobal fglobal ffglobal jfs=10000;load abnormal.matfid=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);ylab