CN110646199B

CN110646199B - 一种基于加权导数动态时间规整的齿轮微弱故障诊断方法

Info

Publication number: CN110646199B
Application number: CN201910878474.9A
Authority: CN
Inventors: 尚志武; 高茂生; 李万祥; 刘霞; 俞燕; 周士琦
Original assignee: Tianjin Polytechnic University
Current assignee: Tianjin Polytechnic University
Priority date: 2019-09-06
Filing date: 2019-09-06
Publication date: 2022-01-07
Anticipated expiration: 2039-09-06
Also published as: CN110646199A

Abstract

本发明公开了一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，包括以下步骤：基于集合经验模态分解(EEMD)算法和频率相关性选择法的信号预处理、基于加权导数动态时间规整(WDDTW)算法的残差向量信号提取、基于希尔伯特变换方法的包络谱故障特征分析。本发明能够更加高效的提取齿轮微弱局部故障特征，能够及时对微弱故障做出诊断，更加适用于实际工业生产环境。

Description

一种基于加权导数动态时间规整的齿轮微弱故障诊断方法

技术领域

本发明涉及振动故障诊断领域，特别是一种基于加权导数动态时间规整的齿轮微弱故障诊断方法。

背景技术

齿轮箱目前被广泛应用于工业生产领域，在其中扮演十分重要的角色，发挥着至关重要的作用，一旦发生故障，将严重影响工业生产，造成经济损失或者人身伤害。对齿轮箱进行故障诊断显得至关重要，其中齿轮作为动力传递的关键部件长时间、高负荷工作容易发生故障。大部分齿轮故障都是从微弱状态逐渐恶化，最终演变成严重的故障状态，所以有效的对齿轮的微弱故障进行及时检测和诊断，便可以有效的预防事故的发生。齿轮在制造和安装过程中都会存在合理的误差，这些误差将有可能作为振动的激励源，引起齿轮箱振动，在故障诊断时所采集到的齿轮箱振动信号会包含这些振动信号，微弱故障发生时由于故障水平微弱，其所引起的振动水平也比较微弱，微弱故障信息容易淹没在振动信号中。因此，对齿轮进行准确的故障诊断必须要解决以下关键问题：

(1)提高信号的信噪比，使信号尽可能包含更多的原始故障信息。

(2)排除齿轮自身存在的加工制造缺陷或者安装误差的干扰。

(3)能够对微弱故障做出及时诊断，对微弱故障敏感。

(4)能够准确处理非线性非平稳信号，具有较强的鲁棒性和快速性。

目前大多数的研究都集中在滚动轴承的故障诊断方面，齿轮故障诊断也取得许多成果，但是对于齿轮微弱故障的诊断效果不佳。时域分析方法的指标较多，不同的应用场合需要使用不同的时域指标，小波包分析方法依赖于小波基的选取，深度学习方法计算量较大、计算速度慢。

发明内容

本发明的目的是为了解决上述问题，设计了一种基于加权导数动态时间规整的齿轮微弱故障诊断方法。

实现上述目的本发明的技术方案为，一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，包括以下步骤：

步骤1：基于集合经验模态分解(EEMD)算法和频率相关性选择法的信号预处理；

步骤2：基于加权导数动态时间规整(WDDTW)算法的残差向量信号提取；

步骤3：基于希尔伯特变换方法的包络谱故障特征分析。

进一步地，所述的步骤1中，基于集合经验模态分解(EEMD)算法和频率相关性选择法的信号预处理，其具体步骤为：

步骤1.1：使用EEMD算法分解原始健康振动信号OH(t)和原始故障振动信号OF(t)；

步骤1.2：通过频率相关性选择法分别选择所需的本征模态函数(IMF)分量。

进一步地，所述的步骤2中，基于加权导数动态时间规整(WDDTW)算法的残差向量信号提取，其具体步骤为：

步骤2.1：构建WDDTW算法的最优规整匹配路径WR；

步骤2.2：通过求残差运算获得原始残差向量信号OR；

步骤2.3：使用重采样技术将序列恢复原长。

进一步地，所述的步骤3中，基于希尔伯特变换方法的包络谱故障特征分析，其具体步骤为：

步骤3.1：对残差向量信号R进行希尔伯特变换；

步骤3.2：包络谱故障特征分析。

进一步地，所述的步骤1.1中，使用EEMD算法分解原始健康振动信号OH(t)和原始故障振动信号OF(t)，其具体步骤为：

步骤1.1.1：齿轮箱在无故障状态下测得的振动信号作为原始健康振动信号OH(t)，齿轮箱中的齿轮发生微弱故障状态下测得的振动信号作为原始故障振动信号OF(t)；

步骤1.1.2：EEMD方法分别分解原始健康振动信号OH(t)和原始故障振动信号OF(t)，获得多个本征模态函数(IMF)分量。

进一步地，所述的步骤1.2中，通过频率相关性选择法分别选择所需的IMF分量，其具体步骤为：

步骤1.2.1：将原始健康振动信号OH(t)及其经EEMD算法分解所得到的所有IMF分量、原始故障振动信号OF(t)及其经EEMD算法分解所得到的IMF分量进行快速傅里叶变换；

步骤1.2.2：计算原始健康振动信号OH(t)和原始故障振动信号OF(t)的频率相关值；

步骤1.2.3：分别计算原始健康振动信号OH(t)和原始故障振动信号OF(t)经EEMD算法分解后得到的所有IMF分量的频率相关值；

步骤1.2.4：计算各个IMF分量的频率相关性；

步骤1.2.5：按照从小到的原则排列各个IMF分量的频率相关性的大小，选择频率相关性最大的IMF分量作为后续处理信号。

进一步地，所述的步骤2.1中，构建WDDTW算法的最优规整匹配路径W_R，其具体步骤为：

步骤2.1.1：通过频率相关性选择法选择的原始健康振动信号称为健康信号H(t)，作为WDDTW算法的参考信号A，通过频率相关性选择法选择的原始故障振动信号称为故障信号F(t)，作为WDDTW算法的测试信号B；

步骤2.1.2：将参考信号A与测试信号B同时输入到WDDTW算法，获得最优匹配规整路径W_R。

进一步地，所述的步骤2.2中，通过求残差运算获得原始残差向量信号OR，其具体步骤为：

步骤2.2.1：根据最优规整匹配路径W_R，计算新的健康振动信号和新的故障振动信号；

步骤2.2.2：将新的故障振动信号和新的健康振动信号作差，获得原始残差向量信号OR。

进一步地，所述的步骤3.2中，包络谱故障特征分析，其具体步骤为：

步骤3.2.1：根据齿轮箱中齿轮的基本参数计算故障特征频率；

步骤3.2.2：从包络谱中读取特征频率，与理论计算的故障特征频率进行对比，实现故障诊断及定位。

进一步地，所述的步骤1.2.4中，计算各个IMF分量的频率相关性，其具体计算方式如下：

其中，N表示IMF分量经快速傅里叶变换后的频率值的总数量，原始健康振动信号OH(t)和原始故障振动信号OF(t)统称为原始信号，M表示原始信号经快速傅里叶变换后的频率值的总数量，F_m表示原始信号的第m个频率值，A_m表示原始信号的第m个频率值的幅值，F_in表示第i个IMF分量的第n个频率值，A_m表示第i个IMF分量的第n个频率值的幅值。S_iN表示第i个IMF分量的频率相关值，S_M表示原始信号的频率相关值。

利用本发明的技术方案制作的一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，通过该方法可以首先通过EEMD方法对原始健康振动信号OH(t)和原始故障振动信号OF(t)进行分解，能够部分分离信号的噪声，利用频率相关性选择法选择合适的IMF分量，可以科学准确的选择所需要的IMF分量，实现降噪，提高了信号的信噪比，再利用WDDTW算法将健康信号与故障信号规整对齐，通过求残差运算获得原始残差向量信号OR，此时便排除了齿轮自身存在的合理加工与安装误差对微弱故障诊断的影响，使用重采样技术重采样原始残差向量信号OR，得到恢复原长的残差向量信号R，最后利用希尔伯特变换方法的强大的分解解调能力分解处理残差向量信号，对残差向量信号进行包络谱故障特征分析，与理论计算的故障特征频率进行对比，实现故障定位与故障诊断。该方法可以有效的实现齿轮微弱故障的诊断，通过信号的预处理减少了处理信号的数量，具有更快的计算速度，WDDTW算法具有强大的规整对齐能力，能够高效地找到最优规整路径，EEMD方法能够有效的分解和分离信号中的噪声，提高信号的信噪比。因此，该发明具有实用价值，能够有效的为设备检测与维修人员提供有价值的设备状态信息，帮助减少经济损失，满足故障检测与诊断的实际需要。

附图说明

图1是本发明所述的一种基于加权导数动态时间规整的齿轮微弱故障诊断方法的流程图。

图2是本实施例实验所用的定轴齿轮箱的结构简图。

图3是本发明所述的基于集合经验模态分解(EEMD)算法和频率相关性选择法的信号预处理的流程图。

图4实验所采集到的健康振动信号、轮毂裂纹故障振动信号和轮齿部分缺损故障振动信号的时域图。

图5是水平方向(Y向)健康信号的EEMD分解结果图。

图6是轮齿部分缺损故障信号的EEMD分解结果图。

图7 Y向健康振动信号各个IMF分量的频谱图。

图8 Y向轮齿部分缺损故障信号各个IMF分量的频谱图。

图9是本发明所述的基于加权导数动态时间规整(WDDTW)算法的残差向量信号提取的流程图。

图10是原始残差向量信号的时域图。

图11是残差向量信号的时域图。

图12是轮毂裂纹故障残差信号的包络谱图。

图13是轮齿部分缺损故障残差信号的包络谱图。

具体实施方式

下面结合附图对本发明进行具体描述，如图1所示，一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，包括以下步骤：

步骤3：基于希尔伯特变换方法的包络谱故障特征分析。

下面通过一个实施例对本发明及其效果进行具体说明。

本实施例的实验数据来自风力涡轮机动力传动故障诊断综合试验台(WTDS)的定轴齿轮箱，通过加速度检测传感器采集齿轮箱的振动信号，该型加速度检测传感器为三向传感器，可以同时测量三个方向的振动信号，即水平方向、竖直方向和轴向。该型定轴齿轮箱的结构简图如图2所示。

图2中，1代表高速级小齿轮，2代表高速级大齿轮，3代表低速级大齿轮，4代表低速级小齿轮，5表示高速级传动整体，6表示低速级传动整体，7表示高速级大齿轮与低速级小齿轮的支承轴承，8表示高速级大齿轮与低速级小齿轮的轴，9表示低速级大齿轮轴，10表示高速级小齿轮与低速级大齿轮的支承轴承，11表示平行定轴齿轮箱箱体，12表示驱动电机轴与平行定轴齿轮箱输入轴相连接的联轴器、编码器和扭矩传感器。本实施例将故障模拟在高速级小齿轮上，三向加速度检测传感器放置在靠近12的箱体上，由于所模拟故障为齿轮微弱局部故障，该位置较靠近故障位置，所以放置在该位置可以尽可能地采集到故障信息。

本实施例所采用的模拟故障类型为轮毂裂纹故障和轮齿部分缺损故障。固定轴齿轮箱的基本参数如表2所示，轮毂裂纹故障和轮齿部分缺损故障的基本参数如表3所示。

表2齿轮主要参数

表3故障基本参数

步骤1：基于集合经验模态分解(EEMD)算法和频率相关性选择法的信号预处理；步骤1的具体过程如图3所示。

实验所采集到的健康振动信号、轮毂裂纹故障振动信号和轮齿部分缺损故障振动信号的时域图如图4所示。每一个振动信号的时长均为4s，传感器采样频率为3000Hz。驱动电机的转频设置为40Hz。

步骤1.1：使用EEMD算法分解原始健康振动信号和原始故障振动信号；

步骤1.1.2：EEMD方法分别分解原始健康振动信号OH(t)和原始故障振动信号OF(t)，获得多个本征模态函数(IMF)分量。EEMD方法的具体步骤如下：

(1)输入OH(t)，作为待分解信号x(t)；

(2)向x(t)中添加正态分布白噪声n(t)得到新的待分解信号y(t)：

y(t)＝x(t)+n(t)

(3)对于y(t)，确定其所有极值点；

(4)利用三次样条插值函数拟合得到y(t)的上下包络线；

(5)用原始数据序列与上下包络线的均值m₁(t)做差，得到新数据序列h₁(t)；

(6)计算h₁(t)是否小于筛选阈值St，St一般取值0.2-0.3，其计算公式如下式所示：

若满足条件，则得到一个IMF分量，若不满足条件，则需要返回(2)继续重复上述过程，i＝i+1，直到得到一个满足条件的IMF分量，得到一个IMF分量后，此时的信号可以表示成如下式所示：

表示分解所得到的第1个IMF分量之和，r₁(t)表示分解得到第1个IMF分量之和后的残余信号分量。

(7)在上一步中得到的IMF分量需要在计算集合均值之后才能作为最终的第1个IMF分量。下式为集合均值的计算公式。

(8)向得到的残余分量中继续添加正态分布白噪声，重复(1)、(2)和(3)步骤，直到新的待分解信号y(t)被完全分解。最终分解结果可表示为如下式所示。

(9)输出分解结果矩阵MH，OH(t)被完全分解完成后再输入OF(t)继续进行上述分解步骤，得到分解结果矩阵MF。分解结果矩阵的每一个列向量为一个IMF分量，分解结果矩阵按照列向量表示法表示成如下形式：

图5为Y向健康信号的EEMD分解结果时域图，图6为Y向轮齿部分缺损故障信号的EEMD分解结果时域图。其他各路信号均需要按照上述的步骤进行分解。每一个原始信号分解后均得到12个IMF分量，这样可以使得分解的结果具有更大的选择空间。

步骤1.2：通过频率相关性选择法分别选择所需的IMF分量。

步骤1.2.1：将原始健康振动信号OH(t)及其被EEMD分解得到的IMF分量、原始故障振动信号OF(t)及其被EEMD分解得到的IMF分量进行快速傅里叶变换；

步骤1.2.4：计算各个IMF分量的频率相关性，其具体计算方式如下：

其中，N表示IMF分量经快速傅里叶变换后的频率值的总数量，原始健康振动信号OH(t)和原始故障振动信号OF(t)统称为原始信号，M表示原始信号经快速傅里叶变换后的频率值的总数量，F_m表示原始信号的第m个频率值，A_m表示原始信号的第m个频率值的幅值，F_in表示第i个IMF分量的第n个频率值，A_in表示第i个IMF分量的第n个频率值的幅值。S_iN表示第i个IMF分量的频率相关值，S_M表示原始信号的频率相关值。

步骤1.2.5：按照从小到大的原则排列各个IMF分量的频率相关性的大小，选择频率相关性最大的IMF分量作为后续处理信号。实验信号的计算结果如表4所示。Y向健康振动信号的IMF分量的频谱图如图7，Y向轮齿部分缺损故障振动信号的频谱图如图8所示。

根据表4的计算结果，每一个信号分解所得到的所有IMF分量中，IMF1分量的相关性最大，并且是远大于其他分量的相关性，因此仅选取每一个信号的IMF1分量作为后续处理的信号。从图7和图8中也可以发现，IMF1分量的频率分布是最丰富的，并且幅值也是最大的。

步骤2：基于加权导数动态时间规整(WDDTW)算法的残差向量信号提取；其具体流程如图9所示。

步骤2.1：构建WDDTW算法的最优规整匹配路径W_R；

步骤2.1.1：在MH的列向量中，通过频率相关性选择法选择的IMF分量称为健康信号H(t)，作为WDDTW算法的参考信号A，在MF的列向量中，通过频率相关性选择法选择的IMF分量称为故障信号F(t)，作为WDDTW算法的测试信号B；

步骤2.1.2：将参考信号A与测试信号B同时输入到WDDTW算法，获得最优规整匹配路径W_R。WDDTW的具体计算过程如下：

(1)输入参考信号A和测试信号B；

(2)计算逻辑修正权重函数(MLWF)，消除参考信号A与测试信号B的相位差，MLWF的计算公式如下所示：

表4健康信号与故障信号的IMF的S值及其相关性

n为测试信号B的长度，n＝1，2，3，……，n_c为测试信号B的中点，w_max为权重w_m的期望上限，g为惩罚因子，可以根据参考信号A与测试信号B的相位差的不同水平来设置不同g值。(3)分别计算参考信号A和测试信号B的转化信号序列。将参考信号A和测试信号B中的每个元素用平均斜率来代替，得到参考信号A的转化序列D_A和测试信号B的转化序列D_B，参考信号中第a个元素的平均斜率的计算公式如下所示：

A_a表示参考信号序列中的第a个元素，D_Aa表示参考信号转化序列的第a个元素，K表示参考信号与测试信号中元素的个数。测试信号B的计算方式与参考信号A的计算方式和表示方法相同。D_A和D_B表示为如下形式：

(4)定义D_A和D_B的起始元素和终止元素。从平均斜率的计算公式中可以发现，得到的D_A和D_B没有计算D_A1、D_B1、D_AK和D_RK，所以现在对其进行如下定义：

(5)动态时间规整计算。动态时间规整的计算过程就是将D_A和D_B中的元素进行匹配，顺序的以D_A中每一个元素作为计算欧氏距离的标准元素，寻找D_B中与之欧式距离最小的元素，该最小元素便于标准元素匹配成功，对应的两点形成最优规整匹配路径W_R的一个坐标，所有坐标按照顺序排列形成完整的最优规整匹配路径W_R，第i对匹配点表示为：w_i＝(D_Aa，，D_Bb)，第i对匹配点的最优规整匹配路径W_R的坐标表示为：w_Ri＝(a，b)，最优规整匹配结果表示为：W＝w₁，w₂，w₃，……，w_i，……，w_j，最优规整匹配路径表示为W_R＝w_R1，w_R2，w_R3，……，w_Ri，……，w_Rj，j表示匹配元素对的总数。最优规整匹配结果需要根据最优规整路径W_R来计算获得。因此动态时间规整必须满足以下3个条件：

条件1(边界条件)：必须满足w₁＝(D_A1，D_B1)，w_m＝(D_AK，D_BK)；

条件2(连续性条件)：给定w_i＝(D_Aa，D_Bb)，且w_i-1＝(D_Aa′，D_Bb′)，则a-a′≤1，且b-b′≤1；

条件3(单调性条件)：给定w_i＝(D_Aa，D_Bb)，且w_i-1＝(D_Aa′，D_Bb′)，则a-a′≥0，且b-b′≥0。

步骤2.2：通过求残差运算获得原始残差向量信号OR；

步骤2.2.1：根据最优规整匹配路径W_R，计算新的健康振动信号和新的故障振动信号，最优规整匹配路径W_R的坐标W_Ri中的数值为对应信号中元素的顺序号，按照最优规整匹配路径W_R的坐标W_Ri中包含的顺序号抽取相应的信号中的对应元素，分别形成新的故障振动信号和新的健康振动信号；

步骤2.2.2：将新的故障振动信号和新的健康振动信号作差，获得原始残差向量信号OR，OR＝OR₁，OR₂，……，OR_i，……，OR_j，原始残差向量信号的时域图如图10所示。

步骤2.3：使用重采样技术将序列恢复原长；

重采样技术的重采样过程步骤如下：

(1)输入原始残差向量信号OR和最优规整匹配路径W_R，令w_Ri＝(a，b)，w_R(i+1)＝(a′，b′)，m为从1到K的整数，i和m的初始值均为1；

(2)判断b′是否与b相等。若二者相等，则令R_m取OR_i和OR_i+1中的较大值；若二者不相等，则令R_m＝OR_i；

(3)令i＝i+1；

(4)判断i是否等于j，若不相等，令m＝m+1，返回(2)继续循环，若相等则停止循环；

(5)输出残差向量信号R＝R₁，R₂，……，R_m，……，R_K；最终得到的残差向量信号的时域图如图11所示。从图10和图11中可以发现，重采样技术只是将信号进行横向压缩。

步骤3：基于希尔伯特变换方法的包络谱故障特征分析。

步骤3.1：对残差向量信号进行希尔伯特变换；

步骤3.2：包络谱故障特征分析。

步骤3.2.2：从包络谱中读取特征频率，与理论计算的故障特征频率进行对比，实现故障诊断及定位。轮毂裂纹故障残差信号的包络谱图如图12所示，轮齿部分缺损故障残差信号的包络谱图如图13所示。

从图12和图13中可以发现，这两种故障在包络谱图中均可以找到故障特征频率40Hz，只是轮齿部分缺损故障的特征频率要比轮毂裂纹故障的特征频率更加明显，这是由于故障差异造成。轮毂在齿轮旋转过程中是与轴之间存在相互作用力，轮毂裂纹故障所产生的振动是间接的传递到传感器，所产生的振动强度较小。轮齿部分缺损直接影响到齿轮啮合刚度，引起轮齿啮合刚度下降，啮合时故障齿的形变增加，产生的振动水平增加，振动信号直接传递到传感器，所以轮毂裂纹故障产生的振动强度要小于轮齿部分缺损故障。

由于在齿轮运行时轮齿不断发生啮合，一些外力因素或者轮齿长时间的工作使得轮齿容易发生部分缺损的故障。轮齿部分缺损故障是齿轮故障中常见的齿轮微弱局部故障。齿轮运行时轮毂不断的与轴之间有力的相互作用，长时间的工作也容易使得轮毂发生裂纹故障，因此轮毂裂纹故障与轮齿部分缺损故障都是齿轮箱中常见的齿轮微弱局部故障，所以本实施例中采用这两种故障类型进行实际应用的说明。

图12与图13中，同一种故障不同方向的振动信号的结果略有差异，X向与Y向包含的故障信息无大差别，Z向信号的包含的故障信息要少于X向与Y向。齿轮在旋转过程中，X向与Y向的动力学特性是相同的。

综上，本发明所述的一种基于加权导数动态时间规整的齿轮微弱故障诊断方法可以有效地诊断齿轮微弱故障，将EEMD方法与WDDTW算法结合使用，可以有效地提取故障特征信息，准确定位故障。

上述技术方案仅体现了本发明技术方案的优选技术方案，本技术领域的技术人员对其中某些部分所可能做出的一些变动均体现了本发明的原理，属于本发明的保护范围之内。

Claims

1.一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，其特征在于，包括以下步骤：

步骤1：基于集合经验模态分解EEMD算法和频率相关性选择法的信号预处理；

步骤2：基于加权导数动态时间规整WDDTW算法的残差向量信号提取；

所述步骤1中基于集合经验模态分解EEMD算法，向得到的残余分量中继续添加正态分布白噪声，直到新的待分解信号y被完全分解；

所述的步骤2中，基于加权导数动态时间规整WDDTW算法的残差向量信号提取，其具体步骤为：

步骤2.1：通过WDDTW算法和求残差运算获得原始残差向量信号；

步骤2.2：使用重采样技术将原始残差向量信号恢复原长，获得残差向量信号；

步骤3：基于希尔伯特变换方法的包络谱故障特征分析；

所述步骤3中希尔伯特变换方法的包络谱故障特征分析：包括对残差向量信号R进行希尔伯特变换和包络谱故障特征分析，根据齿轮箱中齿轮的基本参数计算故障特征频率，从包络谱中读取特征频率，与理论计算的故障特征频率进行对比，实现故障诊断及定位。

2.根据权利要求1所述的一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，其特征在于，所述的步骤1中，基于集合经验模态分解EEMD算法和频率相关性选择法的信号预处理，其具体为，通过频率相关性选择法选择EEMD算法处理后所需的本征模态函数IMF分量。

3.根据权利要求2所述的一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，其特征在于，通过频率相关性选择法选择EEMD算法处理后所需的本征模态函数IMF分量，其具体步骤为：

步骤1.1：计算原始健康振动信号、原始故障振动信号及二者经EEMD算法分解后得到的所有IMF分量的频率相关值；

步骤1.2：选择频率相关性最大的IMF分量作为后续处理信号。

4.根据权利要求3所述的一种基于加权导数动态时间规整的齿轮微弱故障诊断方法，其特征在于，所述的步骤1.1中，各个IMF分量的频率相关性的具体计算方式如下：

其中，N表示IMF分量经快速傅里叶变换后的频率值的总数量，原始健康振动信号和原始故障振动信号统称为原始信号，M表示原始信号经快速傅里叶变换后的频率值的总数量，F_m表示原始信号的第m个频率值，A_m表示原始信号的第m个频率值的幅值，F_in表示第i个IMF分量的第n个频率值，A_in表示第i个IMF分量的第n个频率值的幅值，S_iN表示第i个IMF分量的频率相关值，S_M表示原始信号的频率相关值。