振动信号处理与多功能时频分析系统的研究

孙苗钟，殷 磊，谈炳发，崔世海，徐元利

(天津市轻工与食品工程机械装备集成设计与在线监控重点实验室，天津科技大学机械工程学院，天津 300222)

振动信号处理与分析的研究能优化设计生产，实时监测设备运行状态，及时发现故障和实施故障诊断，保障生产质量，提高设备可靠性等[1]．如何高效地开发出集信号处理与分析方法于一体的振动测试系统，一直是科研人员研究开发的技术难题．不少学者对振动测试系统进行研究与设计，取得了一些成果．如屈威等[2]采用 LabVIEW 开发出多通道的振动测试系统，用于高铁桥梁的振动测试，但该系统不具备时频分析功能；唐进元等[3]采用LabVIEW生产者-消费者模式开发出用于齿轮传动的振动测试系统，分析功能比较单一；唐奕等[4]使用 LabVIEW 编程语言设计了振动信号分析系统，时频分析部分仅实现短时快速傅里叶变换(STFT)和 Helbert-Huang变换(HHT)．以上开发的振动测试系统运用LabVIEW及其他可视化编程工具，如VB、VC++等编程平台，这些平台对各种振动测试信号处理和分析方法编程难度较大[5]．

MATLAB是集数值计算、时频分析、信号处理和图形分析等强大功能于一体的编程平台，提供大量的函数库，具有较强的信号分析和绘图能力，可以生成面向对象编程的组件对象模型(COM)组件；使用MATLAB对信号处理与分析(尤其是时频分析)的算法设计效率高，但不能脱离 MATLAB集成环境工作，运行速度慢，界面开发能力差[6]．Delphi是面向对象的可视化编程语言，适合图形界面的开发，但要实现对各种信号处理与分析方法的编程十分困难．如果将两者混合编程，就可充分利用 MATLAB和Delphi各自的优点，由MATLAB完成信号处理与分析过程和绘图功能，在 Delphi语言开发的程序界面中嵌入 MATLAB COM 组件[7–8]．本文通过MATLAB与 Delphi混合编程，基于 COM 组件技术设计与实现了集32种信号处理与分析功能方法的振动测试系统，可以完成信号数据回放、信号数字滤波处理、平稳信号时域与频域分析和非平稳信号时频联合分析的方法．在时频分析方法中，提出了时频组合方法并实施完成，如经验模态分解(EMD)与快速傅里叶变换(FFT)组合(EMD-FFT)等，这些组合方法能更彻底地提取非平稳信号的特征频率成分，分析更为方便且结果表达更为深刻．该系统可通过扩充 COM组件增加分析功能，且系统软件可以脱离于MATLAB环境独立运行于Windows操作系统．

1 系统编程框图与功能设计

系统的编程框图如图 1所示．“振动信号处理与多功能时频分析”系统由信号分析设置、信号数字滤波处理、数据回放、平稳信号分析和非平稳信号分析组成．

图1 系统的编程框图Fig. 1 Flow chart of system programming

1.1 信号分析设置

信号分析设置分为4类分析：在线采集单通道分析、在线采集双通道分析、非实时采集单通道分析及非实时采集双通道分析．在线采集分析指的是对实时采集后的信号分析．非实时采集分析指的是对以前数据的信号分析(包括分析后保存的数据文件)．信号分析设置中要求输入分析通道数、采样频率、分析点数、输入文件名、输出文件名(要求输出)．

1.2 信号数字滤波处理

采集后的信号可进行选择性频率的数字滤波，增加信号的预处理能力．运用 MATLAB平台中的频域滤波方法和时域滤波方法，能进行高通、低通、带阻和带通 4种方式滤波，时域滤波可选择巴特沃斯、切比雪夫Ⅰ和切比雪夫Ⅱ 3种方法．需要输入不同的截止频率、通带波动系数和阻带波动系数．

1.3 数据回放

要对已存测试信号进行不丢点和准确性的验证，就需要进行数据重新回放，观察信号整个时间历程．图 2中的波形回放主窗口内是采集到的 8个通道波形回放．即是一台双通道 SF2-DDS函数信号发生器分别发生 A路 100Hz、峰值电压 1V的正弦波和B路100Hz、峰值电压1.5V的方波，A路正弦波同时接数据采集卡的 0、2、4、6通道、B 路方波同时接数据采集卡的 1、3、5、7通道．通过回放各个波形完好无缺．下面是数据回放窗口，可检查每个采集点的数据大小．其左边是信号回放的信息和采样信息，包括回放速度、文件长度、通道点数、屏幕位置、量程、通道频率和通道总数．本次选择回放速度点数为20．总数据文件长度达 139420点，每个通道数据为 17427点．

图2 振动信号处理及多功能时频分析系统界面Fig. 2 System interface of vibration signal processing and multi-function time-frequency analysis

1.4 平稳信号分析

平稳信号分析又分为单通道分析和双通道分析，单通道分析的方法有实时信号波形显示、快速傅里叶变换(FFT)、自相关函数、自功率谱密度、倒谱和概率密度函数；双通道分析的方法有互相关函数、互功率谱密度、频率响应函数和相干函数．

1.5 非平稳信号分析

该部分为系统分析的核心，分为变换与分布和组合分析两部分，变换与分布的方法有短时快速傅里叶变换(STFT)、Helbert-Huang变换(HHT)、魏格尔分布(WVD)、伪魏格尔分布(PWVD)、平滑伪魏格尔分布(SPWVD)和 Choi-Williams分布(CWD)．组合分析的方法有小波分解与 FFT、小波降噪与 FFT、小波包分解、小波包降噪与FFT、EMD-FFT、EEMD-FFT、EMD-WVD和EMD-WVD的叠加．

本系统还具有保存或打印分析波形的功能．

2 混合编程过程

2.1 COM组件技术

组件对象模型(COM)为Windows提供了一种面向对象的、可扩充的通信协议，通过此协议，任何编程语言、不同平台和彼此独立的对象都可以实现连通．MATLAB通过调用外部编译器对函数文件编译，使用COM Builder生成COM组件的动态链接库(DLL)文件，Delphi通过安装、调用该 DLL文件，实现MATLAB和Delphi之间的数据通信，从而实现混合编程[8]．

2.2 生成DLL文件并注册

MATLAB平台拥有各种信号处理与分析函数和强大的绘图函数．可编写各种函数 M 文件(子程序)，该文件可以接受和返回参数．在MATLAB命令窗口中输入 mex-setup，选择外部编译器 Microsoft Visual C++ 6.0．再输入命令 comtool，弹出 COM 组件创建窗口，输入组件名字，点击菜单中的 Build COM Object进行编译，编译成功后形成相应的动态链接库(DLL)文件，该文件会在 Windows系统中完成注册[9]．

2.3 安装DLL文件

在 Delphi平台中，点击菜单选项中的 Project Import Type Library，弹出窗会显示上面注册后的各种类型库(DLL)文件，选择所需DLL文件，安装完毕后就会在ActiveX类中出现相应的组件．双击组件会在 Delphi的编程窗口中出现相应的函数类，如图 3为系统主界面窗口调用各函数类情况，在 Delphi的子程序中调用相应的分析方法函数类，实现该方法与MATLAB数据接口的连接与通信．

图3 系统编程主界面窗口调用各函数类Fig. 3 Function classes in the main system programming interface

混合编程的流程框图如图4所示．

图4 混合编程的流程框图Fig. 4 Flow chart of mixed programming

3 时频分析的组合方法实施

在非平稳信号的时频分析中，为了更好地提取非平稳信号的频率特征信息，在实施小波变换(WT)、EMD和EEMD等方法时采用2种或3种方法的组合(组合方法)分析，如 WT-FFT、EMD-WVD 和EEMD-FFT等，这可使结果分析更为深刻和清晰[10].以经验模式分解与魏格尔分布(EMD-WVD)组合方法为例说明其实施过程，该方法采用了EMD、相关系数法和WVD 3种方法进行组合对非平稳信号分析.

经验模式分解由于种种原因存在过分解现象，过分解产生的原因主要有以下几方面：局部均值的数值计算方法的插值误差、边界效应的影响、终止筛选标准的不严格等．过分解使得到的分量个数比原信号组成分量多，把非原信号组成成分的分量称为“伪分量”[10]．用相关系数法来剔除“伪分量”，采用 EMD分解后的各个内禀模式函数(IMF)与原信号进行相关得到各个相关系数，当相关系数很小时，就可判断该分量为伪分量，可能剔除．由于 WVD在非平稳信号分析具有较好的时频聚集性和分辨率[11]，所以再用WVD方法对剩余的IMF进行分析，得到各WVD分布的时频谱图，分析结果更为清晰，信号被分解后的各频率成分更清晰地表达出来．

4 仿真信号验证系统

为了验证系统的分析效果与可靠性，运用系统中多种分析方法对一仿真信号进行分析比较．

设仿真信号的解析表达式为

该信号是由一个线性趋势项 t，两个基频分别为20Hz、200Hz的正弦波和一个基频100Hz并被基频20Hz正弦波调频的余弦波叠加组成的非平稳信号．对其调频部分的频率分析，得角频率：

由式(3)可得其频率变化范围为 80～120Hz，且以余弦波规律变化．

仿真信号的波形和 FFT分析后的频谱图(输出打印得到)如图 5所示．输入采样点数 1000，采样频率 1000Hz．

从频谱图可知 20、100、200Hz基频成分出现分明，这是符合原信号成分的．由于在FFT计算中对信号进行截断，并采用了窗函数(本系统采用汉宁窗)，会产生截断和渗漏效应，明显地在 100Hz处 20Hz正弦波调频成分表达不清楚，在频率为零附近的低频成分是不符合趋势项分析要求的，说明 FFT对非平稳信号分析存在缺陷．

应用时频组合方法(EMD-WVD)对仿真信号EMD分解后各IMF波形以及对各IMF用WVD分析后的时频谱图如图 6所示．在图 6(b)中标出了各IMF与原信号的相关系数，分别为：0.6361、0.5489、0.5544和 0.0327．显然，IMF1、IMF2和 IMF3与原信号相关，IMF4与原信号几乎无相关．从 WVD分析可知，IMF1、IMF3分别是 200Hz、20Hz正弦波，IMF2是基频 100Hz被基频 20Hz正弦波调频的余弦波，IMF4为伪 IMF，应剔除，Res是线性趋势项t．可见采用这种组合方法来分析非平稳仿真信号显示其优点．

图5 仿真信号的波形及FFT分析后的频谱Fig. 5 Simulation signal waveform and its spectrum

图6 仿真信号用EMD-WVD分析后的结果图Fig. 6 Results of simulation signal analyzed with EMD-WVD

5 实例分析

一台 S3SL型砂轮机转速为 2850r/min，即工频f0为 47.5Hz，振动剧烈噪声大，达到 96.2dB．测试系统的硬件组成框图如图7所示．

用传感器提取砂轮机在工频转速下其机壳上的垂直振动信号，选择采样频率为 1000Hz，取数据长度1000点，运用本系统分析．砂轮机原信号、频域方法滤波后的信号及其频谱图如图 8所示．频域方法滤波是利用 FFT快速算法对输入信号采样数据进行离散傅里叶变换分析其频谱，根据滤波要求，将需要滤除的频率成分直接设置成零或加渐变过渡带后再设置成零．滤波时采用低通滤波方式，截止频率取300Hz，由于其频率特性，对频域数据的突然截断造成的谱泄会造成滤波后的时域信号出现失真变形，但频域方法具有较好的频率选择性和灵活性．砂轮机原信号、小波滤波后的信号及其频谱图如图 9．小波分解时，对于不同的信号，应选择不同的小波基．小波变换相当于通过小波的尺度因子和时移因子的变化去观察信号，尺寸因子大可观察信号的总体，尺寸因子小可观察信号的细节．本系统选择广泛使用的Daubechies(db)小波基，适中选择 5层小波，适合大多数信号．分解后取对第一层重构系数 a1(信号的低频成分)，对上面两图中信号滤波后的频谱图进行比较，结果吻合很好，证实了两种滤波方法的可靠性．从滤波后的信号频谱图中反映了砂轮机振动信号低频区域内各倍频的频率特征成分信息，即低频处的1～5倍频峰值十分明显，尤其是1倍频(1f0)．

图7 砂轮机测试系统的硬件组成框图Fig. 7 Hardware framework of the testing system of the grinding machine

图9 砂轮机原信号、小波滤波后信号及其频谱图Fig. 9 Grinding machine original signal，wavelet denoise signal and their spectrograms

使用短时快速傅里叶变换(STFT)和魏格尔分布(WVD)两种时频方法对砂轮机的振动信号进行分析，这两种方法分析结果如图10所示．

图10 砂轮机振动信号STFT、WVD时频谱图Fig. 10 STFT and WVD time-frequency spectrograms of the grinding machine

经STFT分析后，得到1倍频约50Hz频率成分最为突出，而3倍频约150Hz的频率成分也较明显，2、4、5倍频的频率成分隐约可见，在300～400Hz范围内存在比较多频率成分．由 WVD分析结果可见，频率分辨率明显增加，但出现不可知的频率成分交叉项，1倍频的频率成分也最为突出．

上述几种方法的分析得到砂轮机振动信号 1、3倍频的特征频率成分突出，根据旋转机械故障诊断理论[12]可知，砂轮机运行时存在严重偏心而造成运转不平衡，导致振动剧烈．测后更换砂轮，振动与噪声恢复正常．

6 结 语

应用 COM 组件技术，可以将 MATLAB强大的信号处理与分析功能、图形处理能力和 Delphi灵活高效的程序设计能力相结合，实现混合编程．开发集信号处理与多功能分析于一体的振动测试系统可以脱离于 MATLAB环境，以执行文件的形式独立运行于 Windows环境．本系统为更好提取非平稳信号频率特征成分信息提出了多种新的时频组合方法，是对时频分析进行补充和完善．通过仿真信号分析和实例测试，验证了系统的可靠性，可应用于振动信号处理与分析，尤其是非平稳信号的分析．

该系统主要优点：实现 32通道信号同时大容量数据采集显示与存盘；对采集后的数据可以回放分析；根据信号分析要求选择不同的数字滤波处理；根据非平稳信号的特点，选择合适的组合时频分析，使分析结果更为清晰．该系统存在的缺点：缺少边采集信号边实时信号分析功能；对于有关时频分析方法中，输入控制参数不够全面，分析有局限性，如小波分析的基函数、层次数，EMD 分析中的标准差、循环次数等．

本系统具有在线帮助功能，对 MATLAB平台的每个函数应用都有较好的解释和理论指导，对各个分析函数使用时作了详细的说明．

在使用本系统时，强调说明对于不同振动信号进行预处理时，分析要求不同，处理是不同的．如对轴承的振动信号不应进行滤波，要保存其信号的内部频率成分．