下载文档原格式
1．数据导入matlab1.1启动Matlab软件1.2点击载入故障数据中的G2015，Workspace窗口出现：1.3取第一组数据G201，命令窗口输入：G201=G2015(1:1:20000);2. 数据预处理在测试中由数据采集所得的原始信号，在分析前需要进行预处理，以提高数据的可靠性和真实性，并检查信号的随机性，以便正确地选择分析处理方法。
预处理工作主要包括三个方面：一是除去信号中的外界干扰信号和剔除异常数据，如趋势项和异点；二是对原始数据进行适当的平滑或拟合；三是对原始信号的特性进行检验。
当然这些处理工作不是全部必需的，可以选—项或两项容，当认为原始信号获取工作十分可靠或原始数据简单可以直接判断的情况下，也可以不进行这些预处理工作。
以下所做数据预处理，故障轴承以G201为例，正常轴承以Z201为例，观察原始数据经过不同方法做处理前后的变化。
1.1零均值化处理（原理公式见报告P8）命令窗口输入：G201l=G201-sum(G201)/20000;%G201l为零均值处理后的数据。
“20000”为采样点数。
sum为求和语句subplot(2,1,1),plot(G201);subplot(2,1,2),plot(G201l);%显示G201与G201l 得到下面图形：从时域图形上看，是波形整体在Y轴的平移。
再看看频域变化，命令窗口输入：N=20000; %采样点数fs=10000; %采样频率f=(0:N-1)'*fs/N; %进行对应的频率转换G201p=abs(fft(G201)); %进行fft变换，G201p为G201进行fft变换后结果G201lp=abs(fft(G201l)); %进行fft变换，G201lp为G201l进行fft变换后结果subplot(2,1,1),plot(f(1:N/2),G201p(1:N/2));subplot(2,1,2),plot(f(1:N/2),G201lp(1:N/2 )); %显示G201与G201p的频谱图得到下面图形：ω处出现一个由直流分量产生的大谱峰（将从频域图可以明显看出，零均值后消除0=近达到45.4⨯），处理后避免了其对周围小峰值产生的负面影响，便于频域分析。
实用标准课程名称：机械故障诊断设计题目：基于FFT的轴承故障诊断学院：机械工程系班级：学号：姓名：指导老师：李奕璠2017年12月23日摘要滚动轴承是旋转机械中重要的零件，以往的动检工作对滚动轴承强烈振动原因分析不足，不能满足设备维修工作的需要。
所以要定期对旋转机械进行动态监测，根据所测数据做出诊断分析，及时发现滚动轴承强烈震动情况。
傅里叶变换在故障诊断技术中是重要的工具，但傅里叶变换及其逆变换都不适合数字计算机计算，要进行数字计算机处理，必须将连续性信号离散化，无限长数据有限化，再进行采样和截断。
这种算法称为有限离散傅里叶变换（DFT），为了提高效率，在DFT的基础上，运用快速傅里叶变换（FFT）对滚动轴承进行故障诊断。
通过FFT方法分析轴承的信号图，对滚动轴承振动的产生原因进行深入分析，不断总结经验，提高故障分析能力，掌握造成滚动轴承强烈振动的原因，及时消除振动，为设备安全提供可行性措施。
关键词：滚动轴承；故障诊断； FFT23第1章 绪论1.1 滚动轴承概述滚动轴承（rolling bearing ）是将运转的轴与轴座之间的滑动摩擦变为滚动摩擦，从而减少摩擦损失的一种精密的机械元件。
滚动轴承一般由内圈、外圈、滚动体和保持架四部分组成，内圈的作用是与轴相配合并与轴一起旋转；外圈作用是与轴承座相配合，起支撑作用；滚动体是借助于保持架均匀的将滚动体分布在内圈和外圈之间，其形状大小和数量直接影响着滚动轴承的使用性能和寿命；保持架能使滚动体均匀分布，引导滚动体旋转起润滑作用。
图1 滚动轴承结构滚动轴承是各类旋转机械中最常用的通用零件之一，也是旋转机械易损件之一。
据统计，旋转机械的故障越有30%是由轴承故障引起的，它的好坏对机械的工作状况影响很大。
轴承的缺陷会导致机器剧烈振动和产生噪声，甚至会引起设备的损坏。
因此，对重要用途的轴承进行工况检测与故障诊断是非常必要的。
1.2 本次任务本次总共给出了4组通过现场测试得到的滚动轴承运行数据，包括1组正常轴承数据，1组内圈故障数据，1组外圈故障数据，1组滚动体故障数据。
基于支持向量机（SVM）的齿轮箱轴承故障识别一、轴承故障诊断1、概述轴承是旋转设备的一个重要部件，它提供重要的负载承受能力，以支撑转子系统抵抗静态的和动态的外力。
轴承构件，由于它的使用寿命长、负载能力高、能量损失低而被广泛应用于工业和公用设施，是大型机械装备(包括动力机械、机车车辆、泵与风机等)中的关键部件。
高速运转的大型机械装备，其轴承的载荷重且为交变载荷，而且工作环境恶劣，经常发生轴承性能劣化和损坏，影响整个装置的安全可靠性，一旦出现故障将导致严重的损失，有必要对轴承工作状态进行模式识别与诊断。
轴承根据工作的摩擦性质不同可分为滑动摩擦轴承（简称滑动轴承）和滚动摩擦轴承（简称滚动轴承）两大类。
本文所测得的数据来自实验室齿轮箱的滑动轴承，滑动轴承的特点有：(1)在高速重载下能正常工作，寿命长。
(2)精度高。
(3)滑动轴承可做成剖分式的，能满足特殊结构的需要。
(4)液体摩擦轴承具有很好的缓冲和阻尼作用，可以吸收震动，缓和冲击。
(5)滑动轴承的径向尺寸比滚动轴承的小。
(6)起动摩擦阻力较大。
通过对轴承进行故障诊断有以下优势：(1)早期预报、防止事故发生，降低事故发生率;(2)预知性维修，提高设备管理水平，降低维修费用，减少维修时间，增加运行时间;(3)提高设备的设计、制造水平，改进产品质量;(4)确定复杂机器的最佳工作参数，提高效率;(5)降低噪声，泄露等污染，保护环境。
2、滑动轴承失效形式（1）磨粒磨损进入轴承间隙的硬颗粒（如灰尘、砂粒等），在起动、停车或轴颈与轴承发生边缘接触时，都将加剧轴承磨损，导致几何形状改变、精度丧失，轴承间隙加大，使轴承性能在预期寿命前急剧恶化。
（2）刮伤进入轴承间隙中的硬颗粒或轴颈表面粗糙的轮廓峰顶，在轴承上划出线状伤痕，导致轴承因刮伤失效。
（3）咬合（胶合）当轴承温升过高，载荷过大，油膜破裂时，或在润滑油供应不足条件下，轴颈和轴承的相对运动表面材料发生粘附和迁移，从而造成轴承损坏。
滚动轴承故障诊断滚动轴承是应用最为广泛的机械零件质疑，同时，它也是机器中最容易损坏的元件之一。
许多旋转机械的故障都与滚动轴承的状态有关。
据统计，在使用滚动轴承的旋转机械中，大约有30%的机械故障都是由于轴承而引起的。
可见，轴承的好坏对机器工作状态影响极大。
通常，由于轴承的缺陷会导致机器产生振动和噪声，甚至会引起机器的损坏。
而在精密机械中（如精密机床主轴、陀螺等），对轴承的要求就更高，哪怕是在轴承上有微米级的缺陷，都会导致整个机器系统的精度遭到破坏。
最早使用的轴承诊断方法是将听音棒接触轴承部位，依靠听觉来判断轴承有无故障。
这种方法至今仍在使用，不过已经逐步使用电子听诊器来替代听音棒以提高灵敏度。
后来逐步采用各式测振仪器、仪表并利用位移、速度或加速度的均方根值或峰峰值来判断轴承有无故障。
这可以减少对设备检修人员的经验的依赖，但仍然很难发现早期故障。
随着对滚动轴承运动学、动力学的深化研究，对轴承振动信号中频率成分和轴承零件的几何尺寸及缺陷类型的关系有了比较清楚的了解，FFT级数的发展也使得利用频率域分析和检测轴承故障成为一种有效的途径。
也是目前滚动轴承监测诊断的基础。
从发展的历程看，滚动轴承故障检测诊断技术大致经历了以下阶段：1961年，W.F.Stokey完成了轴承圈自由共振频率公式的推导，并发表；1964年，O.G.Gustafsson研究了滚动轴承振动和缺陷、尺寸不均匀及磨损之间的关系，这与目前诊断滚动轴承故障的方法是基本一致的；1969年，H.L.Balderston根据滚动轴承的运动分析得出了滚动轴承的滚动体在内外滚道上的通过频率和滚动体及保持架的旋转频率的计算公式。
至此，有关滚动轴承监测诊断的理论体系已经基本完成；1976年，日本新日铁株式会社研制了MCV-021A机器检测仪，其方法是通过检测低频、中频和高频段轴承的信号特征来判断轴承的工作状态；1976～1983年之间，日本精工公司也积极在滚动轴承检测仪器方面做工作，相继推出了NB系列轴承检测仪，利用1～15kHz范围内的轴承振动信号的有效值（rms）和峰峰值（p－p）来诊断轴承的故障；1980年代至今，以改良频率分析的方法来精密诊断滚动轴承的故障、确定故障位置，一直是精密诊断采取的必备方法，其中包括细化谱分析、倒频谱分析、共振解调技术、包络分析技术等。
1．数据导入matlab载入故障数据中的G2015，Workspace窗口出现：1.3取第一组数据G201，命令窗口输入：G201=G2015(1:1:20000);2. 数据预处理在测试中由数据采集所得的原始信号，在分析前需要进展预处理，以提高数据的可靠性和真实性，并检查信号的随机性，以便正确地选择分析处理方法。
预处理工作主要包括三个方面：一是除去信号中的外界干扰信号和剔除异常数据，如趋势项和异点；二是对原始数据进展适当的平滑或拟合；三是对原始信号的特性进展检验。
当然这些处理工作不是全部必需的，可以选—项或两项容，当认为原始信号获取工作十分可靠或原始数据简单可以直接判断的情况下，也可以不进展这些预处理工作。
以下所做数据预处理，故障轴承以G201为例，正常轴承以Z201为例，观察原始数据经过不同方法做处理前后的变化。
1.1零均值化处理〔原理公式见报告P8〕命令窗口输入：G201l=G201-sum(G201)/20000;%G201l为零均值处理后的数据。
“20000〞为采样点数。
sum为求和语句subplot(2,1,1),plot(G201);subplot(2,1,2),plot(G201l);%显示G201与G201l 得到下面图形：从时域图形上看，是波形整体在Y轴的平移。
再看看频域变化，命令窗口输入：N=20000; %采样点数fs=10000; %采样频率f=(0:N-1)'*fs/N; %进展对应的频率转换G201p=abs(fft(G201)); %进展fft变换，G201p为G201进展fft变换后结果G201lp=abs(fft(G201l)); %进展fft变换，G201lp为G201l进展fft变换后结果subplot(2,1,1),plot(f(1:N/2),G201p(1:N/2));subplot(2,1,2),plot(f(1:N/2),G201lp(1:N/2 )); %显示G201与G201p的频谱图得到下面图形：ω处出现一个由直流分量产生的大谱峰〔将从频域图可以明显看出，零均值后消除0=近达到4105.4⨯〕，处理后防止了其对周围小峰值产生的负面影响，便于频域分析。
第二组实验轴承故障数据：数据打开后应采用X105_DE_time作为分析数据，其他可作为参考，转速1797rpm轴承型号：6205-2RS JEM SKF, 深沟球轴承采样频率：12k Hz1、确定轴承各项参数并计算各部件的故障特征频率通过以上原始数据可知次轴承的参数为：轴承转速r=1797r/min；滚珠个数n=9；滚动体直径d=；轴承节径D=39mm；：滚动体接触角α=0由以上数据计算滚动轴承不同部件故障的特征频率为：外圈故障频率f1=r/60 * 1/2 * n(1-d/D *cosα)=内圈故障频率f2=r/60 * 1/2 * n(1+d/D *cosα)=滚动体故障频率f3=r/60*1/2*D/d*[1-(d/D)^2* cos^2(α)]=保持架外圈故障频率f4=r/60 * 1/2 * (1-d/D *cosα)=2.对轴承故障数据进行时域波形分析将轴承数据导入MATLAB中直接做FFT分析得到时域图如下：并求得时域信号的各项特征：（1）有效值：；（2）峰值：；（3）峰值因子：；（4）峭度：；（5）脉冲因子：；（6）裕度因子：：3.包络谱分析对信号做EMD模态分解，分解得到的每一个IMF信号分别和原信号做相关分析，找出相关系数较大的IMF分量并对此IMF分量进行Hilbert变换。
由图中可以看出经过EMD分解后得到的９个IMF分量和一个残余量。
IMF分量分别和原信号做相关分析后得出相关系数如下：由上表得：IMF1的相关系数明显最大，所以选用IMF1做Hilbert包络谱分析。
所得Hilbert包络谱图如下：对包络谱图中幅值较大区域局部放大得到下图由以上包络图的局部放大图中可以看出包络图中前三个峰值最大也最明显，三个峰值频率由小到大排列分别为、、。
把这三个频率数值和前文计算所得的理论值进行比较可知：频率值最大为和内圈的故障理论计算特征频率f2=相近，说明此轴承的故障发生在轴承的内圈。
滚动轴承故障诊断滚动轴承故障诊断初步1、故障原因滚动轴承的早期故障是滚⼦和滚道剥落、凹痕、破裂、腐蚀和杂物嵌⼊。
即主要故障形式：疲劳剥落、磨损、塑性变形、锈蚀、断裂、胶合、保持架损坏。
产⽣主要原因包括搬运粗⼼、安装不当、不对中、轴承倾斜、轴承选⽤不正确、润滑不⾜或密封失效、负载不合适以及制造缺陷。
2、频谱和波形特征滚动轴承它是由内圈、外圈、滚动体和保持架四部分组成。
当滚动体和滚道接触处遇到⼀个局部缺陷时，就有⼀个冲击信号产⽣。
缺陷在不同的元件上，接触点经过缺陷的频率是不相同的，这个频率就称为滚动轴承的特征频率。
滚动轴承的故障特征频率的数值⼀般在⼏赫兹到⼏百赫兹之间，在频谱图中的1000Hz以内的低频区域轴承故障特征频率如下：1、滚动轴承故障特征频率（外圈静⽌）式中：Z——滚动体个数fr——转频（Hz）D——轴承节径（mm）d——滚动体直径（mm）α——接触⾓（1）滚动轴承内圈故障特征频率（2）滚动轴承外圈故障特征频率（3）滚动轴承滚动体特征频率（4）滚动轴承保持架特征频率2、滚动轴承故障特征频率的计算经验公式：⼆、滚动轴承故障诊断的要素滚动轴承由内圈、外圈、滚动体和保持架四部分组成，每个轴承部件对应⼀个轴承故障特征频率。
滚动轴承的故障频率分布有⼀个明显的特点，往往在低频和⾼频两个频段内都有表现。
所以在频率分析时，可以选择在这两个频段进⾏分析。
根据滚动轴承的故障形式在频域中的表现形式，将整个频域分为三个频段，既⾼频段、中频段和低频段。
l ⾼频阶段指频率范围处于2000－5000Hz 的频段，主要是轴承固有频率，在轴承故障的早期，⾼频段反映⽐较敏感；中频阶段指频率范围处于800－1600Hz 的频段，⼀般是由于轴承润滑不良⽽引起碰磨产⽣的频率范围；l 低频阶段指频率范围处于0－800Hz 的频段，基本覆盖轴承故障特征频率及谐波；在⾼频段和低频段中所体现的频率是否为轴承故障频率，还要通过其他⽅法进⾏印证加以确认。
广州大学学生实验报告开课学院及实验室： 526 室2015年 12月25日学院机械与电气工程学院年级、专业、班机械124班姓名邓庆青学号1207200145实验课程名称机械故障诊断技术成绩实验项目名称滚动轴承频谱分析及故障诊断指导老师郑文一、实验目的1、进一步熟悉常用信号分析仪器的使用；2、了解常规滚动轴承的结构、特征频率及安装；3、掌握滚动轴承的振动测量及分析方法。
通过运用振动分析手段，完成滚动轴承振动信号的测量及分析，从而提高学生进行数据采集、滚动轴承振动分析及状态评估、故障判断等方面的能力。
二、实验设备正常滚动轴承型号为：NTN6201、加速度传感器、Data line数据采集器、ODYSSEY系统、振动试验台轴承故障模块：故障模块中使用的是6024轴承，并利用特殊方法对轴承进行了故障处理。
轴承模块也设计成方便安装的方式（如图所示），可以快速方便的安装在齿轮箱的输入轴上。
在轴承故障模块的顶部有一个英制螺孔（1/4”-28），用来安装传感器。
轴承型号滚动体个数保持架滚动体外环内环6204 8 0.382 1.996 3.054 4.946三、实验要求1、列出主要结构参数，如电机参数、轴承型号、传动比等；2、测出并得到振动试验台测点布置图，说明测量的位置、方向及传感器安装方法等；3、计算各特征频率，如转速，不平衡、对中不良及轴承损坏等的特征频率。
四、实验步骤1、仪器连接；2、测试参数选择，如频率范围（要求能测量滚动轴承的各主要频率成分）；测量单位用速度 mm/s ；3、取轴承模块MD711X一块，套在齿轮轴上。
4、将压紧垫片置于轴承的内环侧面，用紧固螺钉将其压紧，从而固定轴承模块，防止轴承内环在轴承上转动。
5、启动试验台；6、测量各测点的时域波形、频谱，并存储于分析仪中；7、用软件将测量结果传输至计算机，并显示、打印（要求用线性坐标），标出各频谱图中主要频率峰值的频率、幅值；五、实验结果MD7111轴承外环故障模块测试结果MD7112轴承内环故障模块测试结果MD7113轴承滚动体故障模块的测试结果MD7114轴承保持架故障模块的测试结果六、分析1、测量对象的结构简图2、计算滚动轴承的特征频率外环fi=Zn/2（1-dcos fi=Zn/2（1-dcosąß/D）ß--接触角fo--外环频率内环fo=Zn/2（1+dcos ß /D）fb--滚动体频率f--保持架频率保持架f=n/2（1-dcos ß /D）d--滚动体直径D--节园直径滚动体fb=D/2d（1-d²cos² ß /D）Z--滚动体数n--轴频3、综合判断滚动轴承的运行状态及存在的问题。
作业名称：傅里叶分析滚动轴承的故障诊断院系：机械工程系学号：姓名：指导教师：20XX年XX月XXXXXXXXX校区傅里叶分析滚动轴承的故障诊断摘要：简要介绍了快速傅里叶变换（FFT）在滚动轴承故障分析中的应用，滚动轴承在机械设备中使用非常广泛，其工作状态直接影响整个设备的运行品质。
对滚动轴承进行状态监测与故障诊断，能够避免重大事故的发生，获得较大的经济和社会效益。
通过快速傅里叶变换（FFT）对滚动轴承运行时的实时数据信号进行分析，可以实现对滚动轴承的状态监测和故障诊断。
同时，采用对正常轴承和故障轴承信号对比分析、各种故障轴承之间信号的对比分析，加深了快速傅里叶变换（FFT）对轴承实时信号分析的运用和理解，能够更好的对轴承进行状态监测和故障分析。
关键词：快速傅里叶变换（FFT）；滚动轴承；故障诊断；状态监测Abstract:This paper describes a fast Fourier transform (FFT) in the rolling bearing failure analysis applications, bearing in machinery and equipment is widely used, and its working status directly affects the quality of the operation of the entire device. Rolling element bearing condition monitoring and fault diagnosis, able to avoid major accidents and achieve greater economic and social benefits. Through Fast Fourier Transform (FFT) for real-time data bearing signal runtime analysis can be achieved on the rolling bearing condition monitoring and fault diagnosis. Meanwhile, the use of normal bearings and bearing fault signal comparative analysis of various fault signals comparative analysis between the bearings and deepened the fast Fourier transform (FFT) of the bearing using real-time signal analysis and understanding of the bearing can be better condition monitoring and fault analysis.Keywords: fast Fourier transform (FFT); Rolling; fault diagnosis; condition monitoring一、概述通过对快速傅里叶变换（FFT）的原理的理解和学习，利用MATLAB软件编程应用快速傅里叶变换（FFT）的方法，对滚动轴承的1组正常数据和2组故障数据（故障类型不同）进行信号分析和处理，并对正常轴承和故障轴承信号对比分析、各种故障轴承之间信号的对比分析，并得出结论，实现对滚动轴承的状态监测和故障分析。
第二组实验
轴承故障数据：
Test2.mat 数据打开后应采用 X105_DE_time 作为分析数据，其他可作为参考，转速 1797rpm
轴承型号： 6205-2RS JEM SKF, 深沟球轴承
采样频率： 12k Hz
1、确定轴承各项参数并计算各部件的故障特征频率通过以上原始数据可知次轴承的参数为：
轴承转速 r=1797r/min；滚珠个数 n=9；滚动体直径
d=7.938mm；轴承节径 D=39mm；：滚动体接触角α=0 由以上数据计算滚动轴承不同部件故障的特征频率为：外圈故障频率 f1=r/60 * 1/2 * n(1-d/D *cos
α )=107.34Hz 内圈故障频率 f2=r/60 * 1/2 * n(1+d/D *cos α)=162.21Hz 滚动体故障频率
f3=r/60*1/2*D/d*[1-(d/D)^2* cos^2( α)]=70.53Hz 保持架外圈故障频率 f4=r/60 * 1/2 * (1-d/D *cos
α )=11.92Hz
2.对轴承故障数据进行时域波形分析
将轴承数据Test2.mat导入 MATLAB 中直接做 FFT 分析得到时域图如下：
并求得时域信号的各项特征：
1）有效值：0.2909；
3）峰值因子：5.2441；2）峰值： 1.5256；4）峭度： 5.2793；6）裕度因子：
3.包络谱分析
对信号做 EMD 模态分解，分解得到的每一个 IMF 信号分别和原信号做相关分析，找出相关系数较大的 IMF 分量并对此 IMF 分量进行 Hilbert 变换。
Empirical Mode Decomposition
im
由图中可以看出经过 EMD 分解后得到的９个 IMF 分量和一个残余量。
IMF 分量分别和原信号做相关分析后得出相关系数如下：
由上表得：IMF1 的相关系数明显最大，所以选用 IMF1 做 Hilbert 包络谱分析。
所得 Hilbert 包络谱图如下：
对包络谱图中幅值较大区域局部放大得到下图
由以上包络图的局部放大图中可以看出包络图中前三个峰值最大也最明显，三个峰值频率由小到大排列分别为
58.59Hz、 105.5Hz、 164.1Hz。
把这三个频率数值和前文计算所得的理论值进行比较可知：频率值最大为 164.1Hz 和内圈的故障理论计算特征频率 f2=162.21Hz 相近，说明此轴承的故障发生在轴承的内圈。
clc
程序 1：原始信号时域分析及小波去噪处理
clear all
轴承诊断\test2.mat');
x1=z.X105_DE_time(1:4096);
clear z;
N=4096; fs=12000;
n=0:N-1;
t=n/fs;
f=n*fs/N; figure(1);
plot(t,x1); xlabel('t');
ylabel(' 幅值');
title(' 原信号时域图')
%小波去噪[thr,sorh,keepapp]=ddencmp('den','wv',x1);
xd=wdencmp('gbl',x1,'db3',2,thr,sorh,keepapp);
figure(2); plot(t,xd); xlabel('t');
ylabel(' 幅值');
title(' 小波去噪后时域图')
程序 2：EMD分解及 Hilbert 包络
clc clear all
轴承诊断\test2.mat');
x=z.X105_DE_time(1:1024);
N=1024; fs=12000;
n=0:N-1; f=n*fs/N;
lag=N;
n=0:N-1;
t=n/fs;
imf=emd(x);
[m,n]=size(imf); %imf 为一m*n阶矩阵，m是imf 分量，n 为数据点emd_visu(x,1:length(x),imf,m); % 实信号的信号重构及emd结果显示函数for i=1:m
a(i)=kurtosis(imf(i,:));% 峭度b(i)=mean(imf(i,:)); % 均值;
c(i)=var(imf(i,:)); % 方差;
d(i)=std(imf(i,:)); % 均方值
e(i)=std(imf(i,:)).^0.5; % 均方根值
f(i)=skewness(imf(i,:)); % 计算偏度
end
[k,c]=max(a); %k 为峭度最大值， c 为最大元素在数组中的位置
[r,lags]=xcorr(x,lag,'unbiased'); % 计算序列的自相关函数for
i=1:m
[R,lags]=xcorr(imf(i,:),lag,'unbiased'); % 计算序列的自相关函数a=corrcoef(R(1:N/2),r(1:N/2)); % 相关系数矩阵【对称】，主对角元素为 1 xg(i)=abs(a(1,2)); % 相关系数
end
[R,C]=max(xg); %R 为最大值， C 为最大元素在数组中的位置
figure(4);
y = hilbert(imf(C,:)); a = abs(y);% 包络b=fft(a);
mag1=abs(b); mag=mag1*2/N; f1=(0:N-1)*fs/N;
plot(f1(1:N/2),mag(1:N/2));
%set(gca,'xlim',[0,.400]); title(' 包络'); xlabel(' 频率'); ylabel(' 幅值');。
