第八章 基于EMD的时频分析方法及其应用在机械动态分析、设备状态监测与故障诊断过程中，存在着大量的非平稳信号，短时傅里叶变换(STFT)、Wigner-Ville分布、小波和小波包分析等方法不同程度地对此类信号的时变性给予了恰当的描述，在工程实际中获得了广泛的应用[1]。
对非平稳、非线性信号比较直观的分析方法是使用具有局域性的基本量和基本函数，如瞬时频率。
1998年，美籍华人Norden E. Huang等人在对瞬时频率的概念进行了深入研究之后，创造性地提出了本征模式函数(Intrinsic Mode Function, IMF)的概念以及将任意信号分解为本征模式函数组成的新方法——经验模式分解(Empirical Mode Decomposition, EMD)[2]，从而赋予了瞬时频率合理的定义和有物理意义的求法，初步建立了以瞬时频率表征信号交变的基本量，以本征模式分量为时域基本信号的新的时频分析方法体系，并迅速地在水波研究[3]、地震学[4]、合成孔径雷达图像滤波[5]及机械设备故障诊断[6～8]等领域得到应用。
8.第八章 基于EMD的时频分析方法及其应用 (201)8.1 EMD的基本理论和算法 (201)1 EMD的基本理论和算法 (202)8.1.1 基本概念 (202)8.1.2 EMD的基本原理 (204)8.1.3 EMD方法的完备性和正交性 (207)8.1.4基于EMD的Hilbert变换(HHT)的基本原理和算法 (209)8.2 EMD实用化技术研究 (211)8.2.1局部均值的求解 (211)8.2.2 端点效应处理方法 (213)8.3 基于EMD的Laplace小波结构模态参数识别方法研究 (214)8.3.1 基于EMD的Laplace小波模态参数识别方法 (215)8.3.2 应用实例 (218)8.4 EMD方法在机械设备故障诊断中的应用 (220)8.4.1机车轮对轴承损伤定量识别方法 (221)8.4.2 烟气轮机摩擦故障诊断 (223)1 EMD 的基本理论和算法8.1.1 基本概念在讨论基于EMD 的时频分析方法之前，必须先建立两个基本概念：一个是瞬时频率的概念，信号的瞬时能量与瞬时包络的概念已被广泛接受，而瞬时频率的概念在Hilbert 变换方法产生之前，却一直具有争议性；另一个是基本模式分量的概念，相对于原信号的Hilbert 变换的结果，只有对基本模式分量进行Hilbert 变换出来的时频谱才具有具体的物理意义。
1） 瞬时频率在Hilbert 变换方法产生之前，有两个主要原因使得接受瞬时频率的概念较为困难：一是受到傅里叶变换分析的影响；二是瞬时频率没有唯一的定义。
当可以使离散数据解析化的Hilbert 变换方法产生以后，瞬时频率的概念得到了统一[9]。
对任意时间序列)(t x ，可得到它的Hilbert 变换)(t y 为： τττπd t x t y ∫∞∞−−=)(1)( (8.1.1) 构造解析函数)(t z ()()()()()i t z t x t iy t a t e Φ=+= (8.1.2)其中幅值函数22)()()(t y t x t a += (8.1.3)相位函数()()arctan()y t t x t Φ= (8.1.4) 而相位函数的导数即为瞬时频率 ()()d t t dt Φω=(8.1.5) 或 1()()2d t f t dtΦπ= (8.1.6) 然而按上述定义求解的瞬时频率在某些情况下是有问题的，可能会出现没有意义的负频率。
考虑如下信号)(2121)()()()(21t j t j t j e t A e A e A t x t x t x ϕωω=+=+= (8.1.7)为了简单起见，假设信号幅值1A 和2A 是恒定的，1ω和2ω是正的。
信号)(t x 的频谱应由两个在1ω和2ω的δ函数组成，即)()()(2211ωωδωωδω−+−=A A X (8.1.8)既然认为1ω和2ω是正的，所以这个信号是解析的，按式(8.1.3)和(8.1.4)可以求解其相位和幅值，得到11221122sin sin ()arctancos cos A t A t t A t A t ωωΦωω+=+ (8.1.9) t A A A A t A )cos(2)(122122212ωω−++= (8.1.10)取相位的导数，得到其瞬时频率，有222121212()11()()()22()A A d t t dt A t Φωωωωω−==−+− (8.1.11) 当两个正弦频率取101=ωHz ，202=ωHz 两个频率时，幅值的取值不同，其瞬时频率也有很大的不同。
如图8.1.1(a)所示，2.01=A ，12=A 时，其瞬时频率是连续的。
而在图8.1.1(b)中，2.11−=A ，12=A ，虽然信号是解析的，瞬时频率却出现了负值，而已知信号的频率是离散的和正的。
可见，对任一信号做简单的Hilbert 变换可能会出现无法解释的、缺乏实际物理意义的频率成分。
Norden E. Huang 等人对瞬时频率进行深入研究后发现，只有满足一定条件的信号才能求得具有物理意义的瞬时频率，并将此类信号称之为本征模式函数(Intrinsic Mode Function, IMF) 或基本模式分量。
具体的推导过程见文献[2]。
(a)2.01=A ，12=A （b ）2.11−=A ，12=A图8.1.1 两个正弦波叠加的瞬时频率2） 基本模式分量基本模式分量的概念是为了得到有意义的瞬时频率而提出的。
基本模式分量)(t f 需要满足的两个条件为：(1) 在整个数据序列中，极值点的数量e N(包括极大值点和极小值点)与过零点的数量z N 必须相等，或最多相差不多于一个，即)1()1(+≤≤−z e z N N N (8.1.12)(2) 在任一时间点i t 上，信号局部极大值确定的上包络线)(max t f 和局部极小值确定的下包络线)(min t f 的均值为零。
即02/)]()([min max =+i i t f t f ],[b a i t t t ∈ (8.1.13)其中[,]a b t t 为一段时间区间。
第一个限定条件非常明显，类似于传统平稳高斯过程的分布。
第二个条件是创新的地方，它把传统的全局性的限定变为局域性的限定。
这种限定是必须的，可以去除由于波形不对称而造成的瞬时频率的波动。
第二个限定条件的实质是要求信号的局部均值为零。
而对于非平稳信号而言，“局部均值”又涉及到用于计算局部均值的“局部时间”，这是很难定义的。
因而用局部极大值和极小值的包络作为代替和近似，强迫信号局部对称。
钟佑明等人在对基本模式分量的数学模型进行分析之后，论证了局部对称性的必要性和用极值点拟合包络线的合理性[10]。
满足以上两个条件的基本模式分量，其连续两个过零点之间只有一个极值点，即只包括一个基本模式的振荡，没有复杂的叠加波存在。
需要注意的是，如此定义的基本模式分量并不被限定为窄带信号，可以是具有一定带宽的非平稳信号，例如纯粹的频率和幅度调制函数。
一个典型的基本模式分量如图8.1.2所示。
图8.1.2 一个典型的基本模式分量8.1.2 EMD 的基本原理对满足基本模式分量两个限定条件的信号可以通过Hilbert 变换求出其瞬时频率。
但不幸的是，大多数信号或数据并不是基本模式分量，任何时刻，信号中可能包括多个振荡模式，这就是为什么简单的Hilbert 变换不能给出一般信号的完全的频率内容的原因。
我们必须把复杂的非平稳信号按一定的规则提取出所包含的基本模式分量。
基于此，Norden E. Huang 等人创造性地提出了如下假设：任何信号都是由一些不同的基本模式分量组成的；每个模式可以是线性的，也可以是非线性的，满足IMF 的两个基本条件；任何时候，一个信号可以包含多个基本模式分量；如果模式之间相互重叠，便形成复合信号。
在此基础上，Huang 进一步指出，可以用EMD 方法将信号的基本模式提取出来，然后再对其进行分析。
该分解算法也称为筛选过程。
这种方法的本质是通过数据的特征时间尺度来获得基本模式分量，然后分解数据[2]。
基于基本模式分量的定义，我们可以提出信号的模式分解原理，信号模式分解的目的就是要得到使瞬时频率有意义的时间序列－基本模式分量。
而基本模式分量必须满足两个条件，即式（8.1.12）和（8.1.13）。
因而，其分解原理如下：（1） 把原始信号)(t x 作为待处理信号，确定该信号的所有局部极值点（包括极大值和极小值点），然后将所有极大值点和所有极小值点分别用三次样条曲线连接起来，得到)(t x 的上、下包络线，使信号的所有数据点都处于这两条包络线之间。
取上、下包络线均值组成的序列为)(t m 。
如图8.1.3所示，N 表示数据点数，A 表示幅值，实线为原始信号)(t x ，“○”和“*”分别表示了原始信号中的极大值和极小值，两条虚线表示用这些极大极小值拟合的上、下包络线，点划线表示均值序列)(t m 。
（2） 从待处理信号()x t 中减去其上、下包络线均值)(t m ，得到：)()()(1t m t x t h −= （8.1.14）检测)(1t h 是否满足基本模式分量的两个条件。
如果不满足，则把)(1t h 作为待处理信号，重复上述操作，直至)(1t h 是一个基本模式分量，记)()(11t h t c = （8.1.15）图8.1.3 信号)(t x 的上、下包络线及均值)(t m（3） 从原始信号)(t x 中分解出第一个基本模式分量)(1t c 之后，从)(t x 中减去)(1t c，N得到剩余值序列)(1t r)()()(11t c t x t r −= （8.1.16）（4） 把)(1t r 作为新的“原始”信号重复上述操作，依次可得第二、第三直至第n 个基本模式分量，记为)(,),(),(21t c t c t c n "，这个处理过程在满足预先设定的停止准则后即可停止，最后剩下原始信号的余项)(t r n 。
这样就将原始信号)(t x 分解为若干基本模式分量和一个余项的和：1()()()ni n i x t c t r t ==+∑ （8.1.17）上述第（4）步中的停止条件被称为分解过程的停止准则，它可以是如下两种条件之一：①当最后一个基本模式分量)(t c n 或剩余分量)(t r n ，变得比预期值小时便停止；②当剩余分量)(t r n 变成单调函数，从而从中不能再筛选出基本模式分量为止。
基本模式分量的两个限定条件只是一种理论上的要求，在实际的筛选过程中，很难保证信号的局部均值绝对为零。