CN104776908A - 基于emd广义能量的轮轨振动信号故障特征提取方法 - Google Patents

Abstract

本发明公开了一种基于EMD广义能量的轮轨振动信号故障特征提取方法。该方法首先采集实时运行列车的振动加速度信号并对车速积分确定车轮旋转一圈对应的起止时刻，用该起止时刻去截取对应时间历程的加速度信号；其次对采集到的振动加速度信号进行小波分解、各层小波系数的阈值处理和小波重构，实现小波去噪；而后对获取的轴箱振动加速度信号进行经验模态分解，得到一系列的本征模态函数；最后结合由车辆轨道耦合动力学模型仿真得到的存在故障激励下的振动加速度信号确定能量权系数，计算经验模态分解广义能量，根据该值确定故障特征。本发明具有成本低、特征提取分辨率高、实时性强等优点。

Description

基于EMD广义能量的轮轨振动信号故障特征提取方法

技术领域

本发明涉及城市轨道列车在途监测与安全预警技术领域，特别是一种基于经验模态分解(EMD)广义能量的轮轨振动信号故障特征提取方法。

背景技术

根据中国城市轨道交通协会统计，“十二五”期间城市轨道交通将要建成投运2500公里左右，年均500公里左右，到2020年末，全国建成城市轨道交通总里程将达7000公里左右。城市轨道交通向着高速度、高密度、技术***构成复杂、业务***联动性高等方面发展，对车辆行车控制、故障诊断与检测提出了越来越高的要求。

城市轨道交通车辆的走行系是车辆的关键***之一，它不但承受车体重量，还传递钢轨与车辆间的驱动力与制动力。作为走行系的核心零部件之一，城市轨道交通车辆轮对是列车与轨道的耦合部位，承载着整个列车的重量并保证列车在轨道上的运行。平轮故障是城市轨道交通车辆重大行车事故隐患之一，车辆的平轮故障在车辆行驶时会产生周期性的噪声和轮轨冲击力，这种噪声不仅严重影响了乘客乘车舒适度，给线路周边环境带来干扰；而且轮轨的冲击能量由道床向路基内层传递，是造成混凝土轨枕和钢轨断裂的重要原因之一，并且这类冲击也是车辆轮轴疲劳冷切，轴承损坏的重要原因，在车辆速度高、负载重和擦伤加深的情况下，这种高强度冲击会造成更严重的破坏。

目前城市轨道交通运营部门虽然对列车轮对制定了各类不同周期的检修计划，但是这种静态检测既占用车辆服役时间又带来大量人工成本，无法实时监测车辆轮对的平轮故障，越来越不能满足城市轨道交通运营安全系数和运营服务质量要求。

发明内容

本发明的目的在于提供一种成本低、工程实施性好的基于EMD广义能量的轮轨振动信号故障特征提取方法，通过采集轴箱振动、构架振动和车体振动等信号，并计算EMD广义能量，实时在线监测轮对平轮故障。

实现本发明目的的技术解决方案是：一种基于EMD广义能量的轮轨振动信号故障特征提取方法，包括以下步骤：

第1步，实时采集运行列车的振动加速度信号并对车速积分确定车轮旋转一圈对应的起止时刻，用该起止时刻去截取对应时间历程的加速度信号；

第2步，对第1步截取得到的振动加速度信号依次进行小波分解、各层小波系数的阈值处理和小波重构，实现小波去噪；

第3步，对第2步小波去噪后的振动加速度信号进行经验模态分解，得到一系列本征模态函数IMF_n；

第4步，确定本征模态函数IMF_n的能量权系数，并计算EMD广义能量，根据该EMD广义能量确定故障特征。

本发明与现有技术相比，其显著优点在于：(1)成本低，该方法所涉及到的硬件***只包含振动传感器和数据处理主机；(2)从信号的固有模态函数出发，引入基于信号局部特性的IMF，在时间和频率以高的分辨率来处理非线性、非平稳的轮轨平轮振动信号，求取能量特征；(3)利用故障激励的动力学响应与无故障激励的响应的差值确定能量权系数；(4)EMD广义能量计算量小，可实现在线实时监测，通过EMD广义能量的阈值判断及时发现突发的和长期累积的平轮状态变化，从而提供及时的维护预警。

附图说明

图1为本发明基于EMD广义能量的轮轨振动信号故障特征提取方法的流程图。

图2为车辆轨道耦合动力学仿真模型图。

图3为本发明平轮故障为7mm的轴箱振动信号进行小波消噪前后信号的时域图，其中(a)为小波消噪前信号的时域图，(b)为小波消噪后信号的时域图。

图4为本发明平轮故障为7mm的轴箱振动信号进行小波消噪前后信号的频域图，其中(a)为小波消噪前信号的频域图，(b)为小波消噪后信号的频域图。

图5为本发明无故障的轴箱振动信号进行EMD分解得到的8个IMF分量和1个剩余分量。

图6为本发明平轮故障为7mm的轴箱振动信号进行EMD分解得到的8个IMF分量和1个剩余分量。

图7为平轮故障为7mm的轴箱振动信号进行EMD分解后各分量能量与正常信号EMD分解后各分量能量图，其中(a)为无平轮故障能量特征向量，(b)为7mm平轮能量特征向量。

具体实施方式

下面结合附图及具体实施例对本发明作进一步详细说明。

结合图1，本发明基于EMD广义能量的轮轨振动信号故障特征提取方法，首先采集实时运行列车的振动加速度信号并对车速积分确定车轮旋转一圈对应的起止时刻，用该起止时刻去截取对应时间历程的加速度信号；其次对采集到的振动加速度信号进行小波分解、各层小波系数的阈值处理和小波重构，实现小波去噪；而后对获取的轴箱振动加速度信号进行经验模态分解，得到一系列的本征模态函数IMF_n；最后结合由车辆轨道耦合动力学模型仿真得到的存在平轮故障激励下的振动加速度信号确定能量权系数，计算EMD广义能量，根据该值确定故障特征。包括以下步骤：

第1步，实时采集运行列车的振动加速度信号并对车速积分确定车轮旋转一圈对应的起止时刻，用该起止时刻去截取对应时间历程的加速度信号；

所述列车的振动加速度信号记为x_i(t)，i＝1,2，…，m，下标i表示对应的轴箱、构架和车体等走行系部件，m表示走行系部件总数。

第2步，对第1步截取得到的振动加速度信号依次进行小波分解、各层小波系数的阈值处理和小波重构，实现小波去噪；

对采集到的振动加速度信号进行小波去噪。一维含噪声信号可表示为

x(t)＝f(t)+σe(t),t＝0,1,…，n-1

式中，f(t)表示真实信号，e(t)表示噪声，x(t)表示含有噪声的信号。在实际工程应用中，有用信号和噪声的小波变换系数特性不同，含噪声信号进行小波分解时，得到不同尺度下的小波系数，通过保留真实信号的小波系数和衰减噪声的小波系数，然后经信号重构即可实现信号降噪。具体步骤为：

(2.1)信号的小波分解：选择合适的小基函数波并确定小波分解层数N，对含噪信号进行N层小波分解。

(2.2)各层小波系数的阈值处理：依据阈值处理函数(硬阈值或者软阈值等)对信号分解得到的各层小波系数进行阈值处理，得到新的小波系数。

(2.3)信号重构：对处理完成的小波系数进行重构，获得去噪后的信号。

第3步，对第2步小波去噪后的振动加速度信号进行经验模态分解，得到一系列本征模态函数IMF_n；具体步骤如下：

(3.1)设原始信号为x(t)，确定信号整个时域上所有的局部极值点；

(3.2)采用三次样条函数分别对信号x(t)所有的极大值点和极小值点进行插值拟合得到信号x(t)的上包络u₀(t)和下包络v₀(t)，则由u₀(t)和v₀(t)得信号的均值包络m₀(t)：

$m_0\left(t\right)=\frac{1}{2}(u_0\left(t\right)+v_0\left(t\right))$

(3.3)用原始信号x(t)减去均值包络m₀(t)，得到剩余部分，记为h₁(t)：

h₁(t)＝x(t)-m₀(t)

步骤(3.1)到步骤(3.3)可认为是一次“筛”的过程，原始信号x(t)经过一次“筛”的过程变为h₁(t)，通常情况下h₁(t)不是满足条件的固有模态函数，需要进一步分解。此时将h₁(t)当作新的原始信号，重复步骤(3.1)～步骤(3.3)，假设重复k次，直至得到h_k(t)满足固有模态函数的条件，过程如下：

$\begin{array}{c}{m}_1\left(t\right)=\frac{1}{2}(u_1\left(t\right)+v_1\left(t\right))\\ h_2\left(t\right)=h_1\left(y\right)-m_1\left(t\right)\\ .\\ .\\ .\\ m_{k-1}\left(t\right)=\frac{1}{2}(u_{k-1}\left(t\right)+v_{k-1}\left(t\right))\\ h_k\left(t\right)=h_{k-1}\left(y\right)-m_{k-1}\left(t\right)\end{array}$

此时h_k(t)是原始信号x(t)的第一个固有模态函数，记为c₁(t)，它包含原始信号x(t)的频率最高的分量。

(3.4)用原始信号x(t)减去c₁(t)得到剩余部分记为r₁(t)：

r₁(t)＝x(t)-c₁(t)

对剩余信号r₁(t)重复上述步骤(3.1)～步骤(3.4)，则可得到第二个固有模态函数c₂(t)，重复n次分解得到各个固有模态函数，其过程如下：

$\begin{array}{c}{r}_1\left(t\right)-c_2\left(t\right)=r_2\left(t\right)\\ .\\ .\\ .\\ r_{n-1}\left(t\right)-c_n\left(t\right)=r_n\left(t\right)\end{array}$

(3.5)当剩余分量r_n(t)满足给定的终止准则时，则停止整个信号分解过程，得到最后一个固有模态函数c_n(t)和剩余分量r_n(t)，此时原始信号x(t)表示为多个固有模态函数和剩余分量之和：

$x\left(t\right)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{c}_i\left(t\right){+r}_n\left(t\right)$

式中c_i(t)表示第i个固有模态分量，c_i(t)代表原始信号不同频率段的分量，r_n(t)代表信号中平稳成分的剩余分量。

第4步，确定本征模态函数IMF_n的能量权系数，并计算EMD广义能量，根据该EMD广义能量确定故障特征；

EMD广义能量的表达式为：

$Q=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\sqrt{\frac{1}{n}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\left(f_j\right)E_i\left({\mathrm{\η}}_j\right)}$

式中：Q为EMD广义能量指数；i＝1,2，…，m，i表示走行系部件，m表示走行系部件总数；j＝1,2，…，n，j为EMD分解后得到各个本征模函数的标号，n为本征模函数总数；η_j为第j项本征模函数提取的信号，长度为车轮旋转一圈对应的起止时刻截取对应时间历程的振动加速度；λ_i(f_j)为第i项振动信号对应的第j项本征模函数的能量权系数；E_i(η_j)为第i项振动信号对应的第j项本征模函数的能量，能量E_i用下式表示：

$E_i={\∫}_0^x{G}_i\left(f\right)\mathrm{dx}={\mathrm{\ψ}}_i^2$

式中：G_i为轨道功率谱密度函数；f为频率；ψ_i ²为均方值，而均方值ψ_i ²计算式为

${\mathrm{\ψ}}_i^2=\frac{1}{x}{\∫}_0^x{\mathrm{\η}}^2\left(x\right)\mathrm{dx}\≈\frac{1}{m}\underset{\mathrm{\μ}=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\η}}_{\mathrm{\μ}}^2$

能量系数λ(f_j)需满足如下量纲规一化条件：

$\underset{j=1}{\overset{n}{\mathrm{\Σ}}}\mathrm{\λ}\left(f_j\right)=1$

不同的本征模函数有着不同的能量系数，能量系数越大表示与之对应的频率段的振动信号IMF对轨道质量影响越大。

采用特征能量法计算能量系数，首先建立车辆轨道耦合动力学模型，仿真获取无故障运行和存在故障激励下的振动加速度信号x_r(t)和x_f(t)；其次对x_r(t)和x_f(t)进行第3步的EMD分解，得到各个波段下的IMF模态函数；再次计算x_r(t)和x_f(t)的各个波段下的本征模函数的能量E_rj和E_fj；而后将x_f(t)各个IMF能量E_fj减去x_r(t)各个IMF能量E_rj得到E_rfj；最后将各IMF能量E_rfi与总能量E的比值作为能量系数的特征向量，即 $\mathrm{\λ}\left(f_j\right)=E_{\mathrm{rfj}}/\underset{j=}{\overset{n}{\mathrm{\Σ}}}{E}_{\mathrm{rfj}}.$ 具体步骤为：

(4.1)采用翟婉明方法建立车辆轨道耦合动力学模型。车辆-轨道垂向耦合模型如图2，包括车辆子***、轨道子***和轮轨接触模型。车辆子***为两个构架、四组轮对组成的多自由度刚体***，共有10个自由度，其中车体和构架考虑沉浮和点头两个方向自由度，轮对只考虑沉浮方向自由度；轨道子***简化为钢轨、轨枕、道床和路基，而钢轨、轨枕和道床只考虑沉浮方向自由度；轮轨接触模型应用Hertz接触理论。

车辆轨道耦合***动力学方程可根据D’Alembert原则写成矩阵形式：

$M\stackrel{..}{x}+C\stackrel{.}{x}+\mathrm{Kx}=P$

式中，M，C，K分别表示质量、阻尼、刚度矩阵；P表示载荷矢量；和x分别表示相关自由度的位移矢量、速度矢量、加速矢量。在轨道不平顺激励作用下，应用快速显式积分法可进行数值积分，求解动力学模型的响应，具体见文献(翟婉明著.车辆-轨道耦合动力学(第三版).科学出版社.2007)。

(4.2)获取无故障运行和存在车轮扁巴激励下的振动加速度信号x_r(t)和x_f(t)，并加入白噪声。对x_r(t)和x_f(t)重复执行第3步，获取得到一系列的本征模态函数IMF_n。

(4.3)按照式：

$E_j={\mathrm{\ψ}}_j^2=\frac{1}{x}{\∫}_0^x{\mathrm{\η}}^2\left(x\right)\mathrm{dx}$

计算x_r(t)和x_f(t)的各个波段下的本征模函数的能量E_rj和E_fj；

将x_f(t)各个IMF能量E_fj减去x_r(t)各个IMF能量E_rj得到E_rfj，即：

E_rfj＝E_fj-E_rj

上式会因为E_fj小于E_rj在特定f波段为负数的情况，在实际计算中，若出现负数，则直接定义E_rfj为零。它的物理意义表示故障信号中几乎没有对应该波长的故障成分，***的输出中也没有对应该波长的能量。

将各IMF能量E_rfi与总能量的比值作为能量系数的特征向量，即：

$\mathrm{\λ}\left(f_j\right)=E_{\mathrm{rfj}}/\underset{j=}{\overset{n}{\mathrm{\Σ}}}{E}_{\mathrm{rfj}}$

实施例1

本实施例基于动力学仿真模型获取在车轮平轮故障和轨道不平顺叠加激励下的轮对振动信号响应。仿真模型的车辆***的参数为高速机车，轨道***的参数为我国高速线路HST60基本动力参数。仿真过程中分别设定列车运行速度v＝20m/s，迭代空间步长Δs＝1mm。车轮平轮故障模型表示为式中x表示擦伤长度方向的位移，z(x)表示擦伤对应的轨道垂向不平顺，L表示擦伤长度，D_f＝L²/16R表示有效的擦伤深度。定义L＝0，1，2，…，7mm长度的车轮扁疤为激励数据，应用快速显式积分法可进行数值积分，求解动力学模型的响应。实验采集轴箱振动加速度信号，采样频率为20000Hz，持续时间为0.1s。

结合图3～4，选取平轮故障为7mm的轴箱振动信号进行小波消噪前后信号的时域和频域分析。小波基函数选择为sym8小波，尺度序列长度为256，可见轮轨冲击信号主要频率在1000Hz以内，冲击持续时间约为16毫秒，在时域图上冲击形成若干个区域峰值。通过比较去噪前后的时频图与时频幅值图，可见小波阈值去噪除去了信号中大部分噪声，而冲击信号的时频特征被良好地保留，验证了轮轨振动信号小波阈值去噪的有效性。

结合图5～6，选取无故障和平轮故障为7mm的轴箱振动信号进行EMD分解，轮轨振动信号分解得到8个IMF分量和1个剩余分量，图中signal为原始信号，imf1-imf8为8个IMF分量，res为剩余分量。从图5中可以看出，每一个模态分量上的信号均无明显的冲击线性，信号呈平稳性，因此无平轮故障特征。从图6中可以看出，imf5、imf6、imf7和imf8中无明显轮对冲击信号，应为原始信号中高频噪声，imf2和imf3为幅值最大的分量，是平轮故障特征成分的主要集中区域，因此选取前5个分量进行能量分析。

结合图7，平轮故障为7mm的轴箱振动信号进行EMD分解后各分量能量与正常信号EMD分解后各分量能量。可见，相对正常信号的各分量能量，imf1和imf2为幅值最大的分量，是平轮故障特征成分的主要集中区域。平轮故障为7mm的振动信号的各分量能量与正常信号各分量能量相减得到E_rfj，即：E_rfj＝E_fj-E_rj。进一步根据确定能量权系数为，λ_(1～5)＝[0.16030.68000.12980.02820.0016]。

对其余平轮故障等级采集的振动信号根据进行EMD广义能量计算，同时计算λ_(1～5)＝1/5情况下的广义总能量，与EMD广义能量进行对比。实验仅采集轴箱振动信号，即m＝1；选取了EMD分解后的前5个分量进行能量分析，即n＝5。

结合表1，广义总能量值与平轮故障等级呈现线性递增关系，在区分第3级和第4级平轮故障时，由于轨道不平顺叠加激励的干扰，振动信号的广义总能量分别为4.25×10³和3.93×10³，无法具体区分故障等级。EMD广义能量幅值较总能量值有一定幅度的增大，同样与平轮故障等级呈现线性关系。由于引入了各个波段的能量系数，使得由平轮故障引起的特定IMF波段的能量对EMD广义能量的贡献更大，因此避免了在区分第3级和第4级平轮故障时，由于轨道不平顺叠加激励的干扰引起的难以区分故障等级。

表1总能量和EMD广义能量对比

综上所述，本发明与现有技术相比，能够实时、准确的提取平轮故障引起的轮轨振动信号特征，实现平轮故障在线实时监测。

Claims (5)

1.一种基于EMD广义能量的轮轨振动信号故障特征提取方法，其特征在于，包括以下步骤：

第1步，实时采集运行列车的振动加速度信号并对车速积分确定车轮旋转一圈对应的起止时刻，用该起止时刻去截取对应时间历程的加速度信号；

第2步，对第1步截取得到的振动加速度信号依次进行小波分解、各层小波系数的阈值处理和小波重构，实现小波去噪；

第3步，对第2步小波去噪后的振动加速度信号进行经验模态分解，得到一系列本征模态函数IMF_n；

第4步，确定本征模态函数IMF_n的能量权系数，并计算EMD广义能量，根据该EMD广义能量确定故障特征。

2.根据权利要求1所述的基于EMD广义能量的轮轨振动信号故障特征提取方法，其特征在于，第1步所述列车的振动加速度信号记为x_i(t)，i＝1,2，…，m，下标i表示对应的走行系部件，m表示走行系部件总数。

3.根据权利要求1所述的基于EMD广义能量的轮轨振动信号故障特征提取方法，其特征在于，第2步所述小波去噪的具体步骤为：

(2.1)信号的小波分解：选择小基函数波并确定小波分解层数N，对含噪信号进行N层小波分解；

(2.2)各层小波系数的阈值处理：依据阈值处理函数对信号分解得到的各层小波系数进行阈值处理，得到新的小波系数；

(2.3)小波重构：对新的小波系数进行重构，获得去噪后的信号。

4.根据权利要求1所述的基于EMD广义能量的轮轨振动信号故障特征提取方法，其特征在于，第3步所述对小波去噪后的振动加速度信号进行经验模态分解，得到一系列本征模态函数IMF_n，具体步骤如下：

(3.1)设原始信号为x(t)，确定信号整个时域上所有的局部极值点；

(3.2)采用三次样条函数分别对信号x(t)所有的极大值点和极小值点进行插值拟合得到信号x(t)的上包络u₀(t)和下包络v₀(t)，则由u₀(t)和v₀(t)得信号的均值包络m₀(t)：

$m_0\left(t\right)=\frac{1}{2}(u_0\left(t\right)+v_0\left(t\right))$

(3.3)用原始信号x(t)减去均值包络m₀(t)，得到剩余部分，记为h₁(t)：

h₁(t)＝x(t)-m₀(t)

进一步分解，将h₁(t)当作新的原始信号，重复步骤(3.1)～步骤(3.3)，假设重复k次，直至得到h_k(t)满足固有模态函数的条件，过程如下：

$\begin{array}{c}{m}_1\left(t\right)=\frac{1}{2}(u_1\left(t\right)+v_1\left(t\right))\\ h_2\left(t\right)=h_1\left(y\right)-m_1\left(t\right)\\ .\\ .\\ .\\ m_{k-1}\left(t\right)=\frac{1}{2}(u_{k-1}\left(t\right)+v_{k-1}\left(t\right))\\ h_k\left(t\right)=h_{k-1}\left(y\right)-m_{k-1}\left(t\right)\end{array}$

此时h_k(t)是原始信号x(t)的第一个固有模态函数，记为c₁(t)，c₁(t)包含原始信号x(t)的频率最高的分量；

(3.4)用原始信号x(t)减去c₁(t)得到剩余信号记为r₁(t)：

r₁(t)＝x(t)-c₁(t)

对剩余信号r₁(t)重复上述步骤(3.1)～(3.4)，则得到第二个固有模态函数c₂(t)，重复n次分解得到各个固有模态函数，其过程如下：

$\begin{array}{c}{r}_1\left(t\right)-c_2\left(t\right)=r_2\left(t\right)\\ .\\ .\\ .\\ r_{n-1}\left(t\right)-c_n\left(t\right)=r_n\left(t\right)\end{array}$

(3.5)当剩余分量r_n(t)满足给定的终止准则时，则停止整个信号分解过程，得到最后一个固有模态函数c_n(t)和剩余分量r_n(t)，此时原始信号x(t)表示为多个固有模态函数和剩余分量之和：

$x\left(t\right)=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{c}_i\left(t\right)+r_n\left(t\right)$

式中c_i(t)表示第i个固有模态分量，c_i(t)代表原始信号不同频率段的分量，r_n(t)代表信号中平稳成分的剩余分量。

5.根据权利要求1所述的基于EMD广义能量的轮轨振动信号故障特征提取方法，其特征在于，第4步所述确定本征模态函数IMF_n的能量权系数，并计算EMD广义能量，根据该EMD广义能量确定故障特征，具体如下：

EMD广义能量的表达式为：

$Q=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\sqrt{\frac{1}{n}\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\left(f_j\right)E_i\left({\mathrm{\η}}_j\right)}$

式中：Q为EMD广义能量指数；i＝1,2，…，m，i表示走行系部件，m表示走行系部件总数；j＝1,2，…，n，j为EMD分解后得到各个本征模函数的标号，n为本征模函数总数；η_j为第j项本征模函数提取的信号，长度为车轮旋转一圈对应的起止时刻截取对应时间历程的振动加速度；λ_i(f_j)为第i项振动信号对应的第j项本征模函数的能量权系数；E_i(η_j)为第i项振动信号对应的第j项本征模函数的能量，能量E_i用下式表示：

$E_i={\∫}_0^x{G}_i\left(f\right)\mathrm{dx}={\mathrm{\ψ}}_i^2$

式中：G_i为轨道功率谱密度函数；f为频率；为均方值，而均方值计算式为

${\mathrm{\ψ}}_i^2=\frac{1}{x}{\∫}_0^x{\mathrm{\η}}^2\left(x\right)\mathrm{dx}\≈\frac{1}{m}\underset{\mathrm{\μ}=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\η}}_{\mathrm{\μ}}^2$

能量系数λ(f_j)需满足如下量纲规一化条件：

$\underset{j=1}{\overset{n}{\mathrm{\Σ}}}\mathrm{\λ}\left(f_j\right)=1$

不同的本征模函数有着不同的能量系数，能量系数越大表示与之对应的频率段的振动信号IMF对轨道质量影响越大；

采用特征能量法确定本征模态函数IMF_n的能量权系数，具体步骤为：

(4.1)采用翟婉明方法建立车辆轨道耦合动力学模型，车辆轨道耦合***动力学方程根据D’Alembert原则写成矩阵形式：

$M\stackrel{\·\·}{x}+C\stackrel{\·}{x}+\mathrm{Kx}=P$

式中，M，C，K分别表示质量、阻尼、刚度矩阵；P表示载荷矢量；和x分别表示相关自由度的位移矢量、速度矢量、加速矢量；

(4.2)获取无故障运行的振动加速度信号x_r(t)和存在车轮扁巴激励下的振动加速度信号x_f(t)，并加入白噪声，对x_r(t)和x_f(t)重复执行第3步，得到一系列的本征模态函数IMF_n；

(4.3)按照下式计算x_r(t)各个波段下的本征模函数的能量E_rj和x_f(t)各个波段下的本征模函数的能量E_fj：

$E_j={\mathrm{\ψ}}_j^2=\frac{1}{x}{\∫}_0^x{\mathrm{\η}}^2\left(x\right)\mathrm{dx}$

将x_f(t)各个IMF能量Efj减去xr(t)各个IMF能量Erj得到Erfj，即：

E_rfj＝E_fj-E_rj

将各IMF能量Erfi与总能量的比值作为能量系数的特征向量λ(f_j)，即：

$\mathrm{\λ}\left(f_j\right)=E_{\mathrm{rfj}}/\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{E}_{\mathrm{rfj}}.$