自适应多小波基函数构造与机械故障诊断应用研究

何正嘉，孙海亮，訾艳阳

(西安交通大学机械制造系统工程国家重点实验室，西安 710049)

1 前言

机械、运载、能源、冶金、石化、国防等国民经济重要行业的关键机械设备，长期在重载、疲劳、腐蚀、高温等复杂恶劣的工况下运行，设备中的核心零部件和重要机械结构不可避免地发生不同程度的故障。机械设备一旦出现事故，将带来巨大的经济损失和人员伤亡。在工程实践中，为了对事故的发生做到防患于未然，因此必须开展对关键机械设备早期故障和复合故障诊断与预示的研究。为有效开展故障诊断，获取机械设备运行过程中的状态信息，提取信号特征是实现故障诊断的重要手段，其本质是信号变换和特征提取。

近年来工程中广泛应用的小波变换和多小波变换等非平稳信号处理方法，都是基于Hilbert空间的信号与不同“基函数”的内积变换。这些信号处理方法的本质都是探求信号中包含与“基函数”最相似或相关的分量，其关键在于构造和选择与故障特征波形相匹配、且具有优良性质的合适基函数，使所构造的自适应基函数与动态信号物理特征达到最佳匹配，获得不同物理意义并符合工程实际的故障特征信息，实现科学、正确的状态监测与故障诊断［1］。

基于故障既是状态(设备性能和状况)又是过程(故障萌生和扩展)这一本质属性，笔者深入研究机械故障动态信号处理的基函数选择与内积变换原理，集中分析小波变换和多小波变换的内积变换表述，揭示它们的内积变换共性并探讨其物理本质和关键问题。为提取运行状态下微弱动态信号的特征，以解决复合故障诊断难题，基于两尺度相似变换、频域提升变换、对称提升变换和时域提升变换自适应构造多小波基函数，阐述机械故障诊断的内积变换原理及其工程应用效果。设备出现早期损伤时，其微弱故障特征被强背景噪声所淹没，考虑小波变换系数间的相关性，深入研究平移不变多小波相邻系数降噪、多小波滑动窗局部阈值降噪和多小波改进相邻系数降噪方法，阐述运行状态下微弱动态信号的特征提取及工程应用效果。

2 内积变换原理与基函数构造

在平方可积空间L2中的函数x(t)和y(t)，它们的内积定义为［2］

式(1)中，y*(t)是y(t)的共轭。设一个Ψ域中的信号x可以表示为以下展开形式:

其中，{ψn}n∈Z是 Ψ 域中的基本函数集。如果{ ψn}n∈Z是Ψ空间的完备序列，即Ψ域内的信号均可由式(2)表示，存在一个对偶函数集，使得其展开系数可由内积函数计算:

从式(3)可以看出，内积结果an越大，则表明x与对偶函数越接近。在上述内积变换中，对偶函数可视为一种“基函数”，则内积变换即为信号x与“基函数”关系紧密度或相似性的一种度量，其关键在于基函数的合理构造，使构造的基函数与机械设备运行过程的动态信号物理特征达到最佳匹配，获得不同物理意义并符合工程实际的故障特征信息，实现科学、正确的状态监测与故障诊断。

2．1 小波变换的内积表述

对信号x(t)进行小波变换，关键是选择小波基函数ψa，b(t)，这里a是反映频率变量的尺度因子，b是对信号进行扫描的时移因子。信号x(t)的小波变换 WTx(a，b)为

这一内积运算旨在探求信号x(t)中包含与小波基函数ψa，b(t)最相关或最相似的分量。因此，构造出一个小波基函数ψa，b(t)，就能够进行一种小波变换。不同类型的机械故障在动态信号中往往表征为不同的特征波形，如果采用了不恰当的小波基函数，则会冲淡特征信息，给特征提取与故障诊断造成困难。同时，小波基函数的正交性、正则性、消失矩、紧支性、对称性、相似性、冗余度等性质与信号特征提取直接有关。因此，小波变换的关键取决于构造和选择与故障特征波形相匹配、且具有优良性质的小波基函数。

2．2 多小波变换的内积表述

多小波变换是小波理论的新发展，它是由两个或两个以上的函数作为尺度函数生成的小波。多小波分析也是建立在内积变换的基础上，多小波两尺度关系可以表述为

利用式(5)对矢量输入信号x进行内积运算，得到多小波分解表达:

对于 Φj，k与 Φj－1，n和 Ψj－1，n有基函数分解关系

这里符号*表示共轭转置。

与单小波不同，多小波是矢量运算，滤波器均为多维向量。由于多小波可以兼备单小波所不能同时满足的优良性质［3］，有助于基函数在内积匹配过程中更加准确地提取微弱故障特征，在早期故障诊断方面具备优势;多小波拥有多个小波基函数，可以同时匹配信号中的多个故障特征，从而实现复合故障耦合特征的一次性分离与识别。

3 自适应多小波基函数构造方法

为了实现基函数的按需构造，使得多小波拥有自适应能力，研究基于两尺度相似变换、频域提升变换、对称提升变换和时域提升变换的自适应多小波基函数构造，以及针对不同应用对象的基函数优选准则，增强在强噪声背景下的早期、微弱和复合故障特征提取和识别能力。

3．1 基于两尺度相似变换的自适应多小波基函数构造

逼近阶是刻画尺度函数对原始函数逼近性能或分析精度的重要特征之一。V Strela博士引入两尺度相似变换［4］，可以改变多尺度函数的逼近阶，从而构造具有优良特性的多小波基函数。随后，F Keinert研究了m重的两尺度相似变换理论［5］，并给出了基于两尺度相似变换的双正交多小波构造方法，完善了两尺度相似变换理论。

在众多多小波中，GHM多小波是最常用的2重多小波［6］，它具有紧支性、对称性、正交性及2阶逼近阶。文章基于两尺度相似变换的转换矩阵M1(ω)对GHM多小波进行改造，构造具有阶逼近阶的双正交多小波，且保持多尺度函数的紧支性和对称性。

式(9)中，参数a和b为非零常数;e为自然常数。

在多小波基函数逼近阶提高为3时，对偶多小波基函数的逼近阶降低为1，为改善其逼近阶特性，用式(10)的两尺度相似变换转换矩阵M2(ω)，将其逼近阶提升为 2［7］。式(10)中，参数c，d和e均为非零常数。

某炼油厂重油催化裂化装置由烟机、风机、齿轮箱、电机组成，结构如图 1所示。机组转速为5850 r/min(转频为97．656 Hz);齿轮箱传动比为5745/1484，低速轴转频为25．22 Hz。监测系统采用电涡流传感器对烟机1#瓦、烟机2#瓦、风机3#瓦、风机4#瓦和齿轮箱5#瓦的垂直和水平方向共10个测点进行实时监测。采样频率为2 kHz。

采集齿轮箱5#瓦的振动信号，利用两尺度相似变换自适应构造多小波基函数，分析结果如图2所示。图2中，每两个冲击I2和一个幅值较大冲击的I1交替出现。冲击I2之间的平均间隔约为0．01 s，与机组高速轴转频接近。而冲击I1的时间间隔约等于冲击I2的3倍，刚好对应于高速轴转频的1/3。该1/3次谐波成分说明机组存在碰摩故障。经现场验证，该碰摩是由于压缩机与齿轮箱之间的膜片式离合器因补偿对中不良而引起的错动所致［7］。

图1 催化裂化装置结构图Fig．1 The structural scheme of a heavy oil catalytic cracking unit

图2 齿轮箱5#瓦振动信号的自适应基函数分析结果Fig．2 The result decomposed by multiwavelet constructed for the signal of 5#bearing bush in the gearbox

3．2 基于频域提升变换的自适应多小波基函数构造

提升方法［8］(lifting scheme)是由贝尔实验室W Sweldens提出的一种构造小波及其滤波器组的强有力工具，基本思想是在双正交小波和完全重构滤波器组的理论基础上，通过设计不同的提升算子改变原有小波滤波器的特性，得到不同性质的双正交小波。研究人员尝试将提升方法引入到多小波的构造中，并成功构造出了多小波的实例。

定理［9］:给定一个初始多小波滤波器组，则新的多小波滤波器组为:

其中，S(z)和T(z)均为有限阶，T(z)的行列式为单项式。它与单小波的提升框架式是类似的，不同之处在于这里的滤波器均为矩阵形式，且在Gnew(z)的表达中多了一项T(z2)。

利用多小波两尺度关系及式(8)可得:

提升框架可以改造现有多小波，使多小波函数具有任意阶消失矩，从而构造出具有理想特性的自适应多小波基函数。提升框架实现的关键是矩阵S(z)和 T(z)的设计［10］。

某炼钢厂连铸连轧机组的精轧机机架的主传动系统为单级减速器，齿数比为Z22/Z65。机组运行时该机架主减速器高速轴靠电机侧有明显异响，但在其例行检测中，机组振动值处于正常范围。该侧轴承为双列短圆柱轴承，型号为SKF316077A，和以往轴承滚动体为实心圆柱体不同，该轴承为柱销式焊接保持架结构，滚动体中空，套装在柱销轴上，轴承参数见表1。在主减速器外侧安装加速度传感器对其进行监测，采样频率为2560 Hz。

表1 高速轴输入侧轴承参数Table 1 The bearing parameters of input side of high speed shaft

基于三次Hermite多小波，利用局部故障域谱熵最小原则得到具有7阶消失矩的自适应多小波基函数。分解结果中，最敏感频带是三层多小波包分解的第四频带第一分支，其平方包络及其频谱见图3。

图3 主减速器振动信号的自适应基函数分析结果Fig．3 The result decomposed by customized multiwavelets for signals of the main reducer

由图 3(b)可以看出，谱线 1．875 Hz、11．88 Hz及25．63 Hz在包络谱中占优，分别接近于保持架故障特征频率2．05 Hz(即保持架自转频率)、滚动体通过频率半频及滚动体通过频率26．34 Hz(即滚动体自转频率)。由于轧钢过程的非平稳性，这三个明显谱线并不是严格对应于相应特征频率。鉴于滚动体自转频率的半频分量为11．88 Hz，而半频等次谐波成分是典型的碰摩故障特征，经现场验证，如图4所示，该柱销式焊接保持架滚动轴承发生保持架脱焊失效，导致了套装在柱销上的滚动体在运行中出现碰摩故障［2］。

图4 轴承保持架结构及损坏状况图Fig．4 The structural scheme of rolling bearing cage and its defects

3．3 基于对称提升变换的自适应多小波基函数构造

根据第3．2节中多小波提升变换基础理论，本节介绍一种基于对称提升变换的自适应多小波基函数构造方法。对称性可以保证滤波器具有线性相位或广义线性相位，能够有效地避免重构误差;其次对称性有利于多小波运算时的边界处理，利用对称扩展的方式可以减少边界畸变。如果初始的多尺度函数与多小波函数均为对称或反对称，若要构造出的新多小波函数也是对称的，则要求初始多尺度函数与多小波函数的平移必须满足对称条件。为确保提升过程的对称性，使用提升过程的“对称选择”方法，来选择用于修正多小波的其他函数的平移量［11］。

某空气分离压缩机(简称空分机)大修后开机，发现齿轮箱振动剧烈，并伴随尖叫声。需要通过对其进行详细地检测与诊断来查清故障。使用便携式诊断仪对该机组进行了全面测试。轴承座振动用加速度传感器测量，采样频率为15 kHz。

利用对称提升方法自适应构造多小波基函数，分解结果如图5所示，第一个分量(图5(a))清晰地表示出与4．7 ms对应的周期性冲击特征，对应于空分机齿轮箱高速轴的转频(213 Hz)。第二个分量(图5(b))中的高频振荡则较为突出，这种特征同样具有周期性，且为4．7 ms。对图5中的信号取其能量，如图6所示，该图中的周期性冲击更加明显。在图5中所观察到的多小波输出的两个分支的结果有所不同，是因为多小波的两个小波基函数各不相同，同时提取出冲击和调制的故障特征。诊断为由于大修过程中安装不善，导致齿轮箱止推夹板和大齿轮端面接触部位由于相对运动产生撞击和摩擦。经开箱验证，诊断结论正确［11］。

图5 空分机信号多小波分解结果Fig．5 The result decomposed by multiwavelet constructed for the signals of an air compressor

图6 空分机信号多小波分解结果的能量Fig．6 The power of the result decomposed by multiwavelet constructed for the signals of an air compressor

3．4 基于时域提升变换的自适应多小波基函数构造

与第二代小波变换相类似，多小波的时域提升框架(为了区别于前面的频域提升原理，将其命名为时域提升框架)主要包含剖分、预测和更新三个步骤，所不同的是多小波提升框架都是矩阵运算，同时在多小波整体提升框架变换中还有前处理和后处理过程。基于时域提升框架的2重多小波分解过程和重构过程如图7所示。

图7 时域提升变换的多小波分解和重构过程Fig．7 The decomposition and reconstruction of multiwavelet constructed using lifting scheme in time domain

信息熵是描述信号概率分布均匀性的有力工具。在此选用最小熵原则来选取与给定信号最佳匹配的多小波基函数。在此过程中，选择对第一层分解的细节信号进行最小熵约束，由此可以得到时域提升框架中最优的自由参数及其相应的自适应多小波基函数。

该方法和滑动窗降噪相结合应用于精轧机齿轮箱复合故障诊断。

4 多小波降噪方法研究

在机械设备运行过程中，所采集的动态信号往往受到背景噪声干扰。因此，如何提高信号信噪比，突出并提取有用的微弱故障特征信息是故障诊断的一个关键环节。

设含噪信号x(t)是所关心的特征信号f(t)和噪声信号r(t)的叠加，通常表示为

其中，r(t)为单位方差、零均值的高斯白噪声;σ为噪声信号的方差。对信号进行降噪的目的就是抑制信号中的噪声部分，从含噪信号x(t)中恢复真实特征信号f(t)。

传统的小波阈值降噪方法存在以下不足:a．由于小波基函数是由母小波函数经伸缩和平移得到，该平移是非一致性取样，引起信号在奇异点(突变、尖点等不规则的瞬变结构)附近产生的急剧振荡的Gibbs现象;b．采用“逐点比较”方式，忽略了小波系数之间的相关性;c．采用全局阈值设定方式，将过扼杀同一尺度下相对较小的信号幅值，导致微弱故障特征信息的丢失。

4．1 平移不变多小波相邻系数降噪方法

结合多小波和平移不变方法的优势，平移不变多小波降噪方法［12］利用时域平移运算获得与输入信号具有一定相位差的新信号，改变原始数据结构中奇异点的位置，从而降低或消除因奇异点位置特殊性所引起的Gibbs现象。对于一个给定信号，通过选择最佳平移量可以实现由Gibbs现象引起的异常振幅的最小化。为了解决最佳平移量问题，采用对一定范围的平移量做循环平移运算，再平均所获得的结果，可以更好地保持信号的光滑性且具有优越的降噪性能。相邻系数降噪方法将若干紧邻的多小波系数组合，进行局部阈值设定，当多小波系数之间存在较强的相关性时，相邻系数降噪方式明显优于逐点比较的传统阈值处理方式，它能在滤去噪声的同时更加有效地保留信号的局部特征信息。

为了验证机车轴承轻微剥落故障，在型号为552732QT的滚动轴承外圈上加工如图8所示的长形小凹槽损伤。将该故障轴承装配到某型客运电力机车，并在机车走行部轴承箱安装加速度传感器对其进行顶轮测试，如图9所示。轴承相关参数见表2。采样频率为 12．8 kHz，机车车轴的转速为650 r/min，即10．833 Hz。根据滚动轴承故障特征频率计算得到，机车滚动轴承外圈损伤的特征频率fo=78．169 Hz。

图8 机车滚动轴承外圈损伤Fig．8 The defect in outer ring of the locomotive rolling bearing

图9 机车滚动轴承顶轮试验Fig．9 The test rig of the locomotive rolling bearing

表2 机车滚动轴承参数Table 2 The parameters of an locomotive rolling bearing

降噪结果中可以看出平移不变多小波相邻系数降噪方法可以全面、准确地提取出所有周期性冲击特征，其周期约为12．8 ms，即78．125 Hz(见图10)，与滚动轴承外圈损伤特征频率相符［13］。

图10 顶轮测试信号的平移不变多小波相邻系数降噪结果Fig．10 The denoising result for signals of the test rig using translation invariant multiwavelet denoising with improved neighboring coefficients

4．2 多小波滑动窗局部阈值降噪方法

滑动窗截断技术如图11所示，可见如果采用传统的全局阈值处理，第二个和第四个微弱冲击信号很容易被过扼杀，造成微弱故障信息的丢失。滑动窗截断技术处理保证每个故障特征被独立、完整地保留到相应子信号中，采用局部阈值设定，以高阶统计量［14］作为优化目标，对降噪结果进行盲解卷积处理［15］。

图11 采用滑动窗截断技术的数据分割Fig．11 The data partition using sliding window truncation

某炼钢厂连铸连轧机组的精轧机机架的主传动系统为单级减速器，齿数比为Z22/Z65。对该精轧机机架进行监测的过程中，发现其主传动系统振动发生异常。在齿轮箱外部安装速度传感器，以采集其运行动态信号，其采样频率为5120 Hz，高速轴小齿轮的转频为4．54 Hz，低速轴大齿轮的转频为1．54 Hz。

采用多小波滑动窗局部阈值降噪方法提取信号特征。由于信号中除了间隔0．22 s的冲击外，其他振动信号没有明显的奇异点，并且信号点间相关性较强，因此局部阈值处理方法采用的是相邻系数局部降噪技术。选择GHM多小波基函数，分解层数选为5层，重复采样的前处理方式，滑动窗窗口宽度由小齿轮转频fr=4．54 Hz确定。多小波滑动窗相邻系数降噪结果如图12所示。图12中交替出现了两种强弱不等的周期冲击I1和I2。在多小波滑动窗相邻系数降噪后，采用基于峭度指标的最优盲解卷积技术以充分突显故障特征。滑动窗相邻系数降噪信号的最优盲解卷积分析结果如图13所示，更加清晰地显示出周期性交替出现两类冲击I1和I2，故可以确定该齿轮箱小齿轮存在两处损伤程度不同的局部故障。并且，大冲击I1反映严重损伤和小冲击I2反映较严重损伤的两处故障的时间间隔约为小齿轮旋转周期的1/3，说明两处局部故障相距约为小齿轮圆周的1/3，即相隔约8个齿［16］。

图12 轧机振动信号的多小波滑动窗相邻系数降噪结果Fig．12 The denoising result for signals of the rolling mills using multiwavelet sliding window

图13 轧机振动信号多小波滑动窗相邻系数降噪和最优盲解卷积分析结果Fig．13 The denoising result for signals of rolling mills using multiwavelet sliding window and optimal blind deconvolution

停机检修后发现，该齿轮箱的高速轴小齿轮存在两处损伤，损伤照片见图14，高速轴小齿轮存在两处局部磨损和胶合故障，且两处故障大概相距1/3个圆周，图14(a)为在齿宽范围内的局部故障，(b)为在整个齿宽范围内的故障。小齿轮的故障位置和损伤状况与以上分析结论十分吻合。

4．3 改进的相邻系数降噪方法

图14 精轧机齿轮箱损伤照片Fig．14 The defects on the pinion of the fine rolling mills

传统多小波相邻系数降噪方法采用的是固定的邻域范围，即对于不同小波分解尺度，其邻域范围是固定不变的，这与多小波分解的特点不能准确吻合，影响了多小波相邻系数的降噪效果。为此，提出一种改进的相邻系数降噪方法［17］，在考虑相邻的多小波系数相关性的基础上，根据多小波分解的特点，构造了随分解层变化的邻域格式，对多小波分解的细节信号进行阈值处理，并进行重构以提纯故障信号。

行星齿轮箱结构尺寸小、输出扭矩大、传动比大、效率高、性能安全可靠，作为核心零部件广泛应用于风力发电、重载起重机、直升机和航天测量船等国民经济重要行业的关键设备中。针对船载卫星通信地球站机械传动机构，测试图15所示的行星齿轮箱(减速器)振动信号，采样频率为12．8 kHz，驱动电机为永磁同步交流电机，转速为255 r/min。

行星减速器第一级传动的啮合频率可由以下公式计算得到:

式(14)中，f啮合为行星减速器第一级传动的啮合频率;N输入为太阳轮输入轴转速;N输出为行星轮系杆输出轴转速;Z中心轮为中心轮齿数。根据结构参数，可以计算出行星减速器第一级传动的啮合频率为太阳轮输入轴转频的10．5倍。

图15 行星减速器第一级传动机构Fig．15 The first transmission device of the planetary gearbox

图16 行星减速器齿面损伤振动信号Fig．16 The collected signals of the planetary gearbox with defects

图16(a)和(b)分别是采集到的振动信号时域波形及其频谱，图16(c)是多小波改进相邻系数的降噪结果，太阳轮旋转一个周期，交替出现三组类似冲击波形A、B和C，表明太阳轮轮齿存在损伤。每组均包含两个邻近的冲击波形，其平均间距约为0．023 s，即43 Hz，约为第一级行星传动太阳轮转频(4．25 Hz)的10．5倍，和行星减速器第一级行星传动的啮合频率接近，表明太阳轮上两个邻近的轮齿出现损伤而引发冲击波形。经开箱验证，情况属实，见图17。

图17 行星减速器太阳轮齿损伤Fig．17 The defect in the sun gear of the planetary gearbox

5 结语

深入研究机械故障动态信号处理的基函数选择与内积变换原理，集中分析小波变换和多小波变换的内积变换表述，揭示它们的内积变换共性和物理本质，指出构造与故障特征相匹配的基函数是机械故障成功诊断的关键。

基于两尺度相似变换、频域提升变换、对称提升变换和时域提升变换自适应构造多小波基函数，可准确提取运行状态下微弱动态信号的特征，成功诊断出重油催化裂化装置、连铸连轧机组、空分机等装备的碰摩、保持架脱焊、撞击和摩擦等多类故障。

所研究的平移不变多小波相邻系数降噪、滑动窗局部阈值降噪和多小波改进相邻系数降噪等方法，准确提取出电力机车、连铸连轧机组和船载卫通站传动系统在运行状态下，滚动轴承外圈轻微剥落、齿轮局部磨损和胶合、行星减速器齿轮早期损伤等多种故障，实现了微弱动态信号特征增强和复合故障特征提取。

致谢:感谢袁静、王晓冬、李臻、陈雪峰、张周锁、李兵和曹宏瑞等为本文做出的贡献。

［1］何正嘉，訾艳阳，陈雪峰，等．内积变换原理与机械故障诊断［J］．振动工程学报，2007，20(5):528－533．

［2］袁 静．机械故障诊断的内积变换原理与多小波特征提取方法研究［D］．西安交通大学，2011．

［3］Daubechies I．Ten Lectures on Wavelets［M］．CBMS－ Conference Lecture Notes．Philadephia:SIMA，1992．

［4］Strela V．Multiwavelets:theory and application［D］．Massachusetts Institute of Technology，1996．

［5］Keinert F．Wavelets and Multiwavelets［M］．USA:A CRC Press Company，2004．

［6］Geronimo J S，Hardin D P，Massopust P R．Fractal functions and wavelet expansions based on several scaling functions［J］．Journal of Approximation Theory，1994(78):373 －401．

［7］Yuan J，He Z J，Zi Y Y，et al．Adaptive multiwavelets via two －scale similarity transforms for rotating machinery fault diagnosis［J］．Mechanical Systems and Signal Processing，2009，23(5):1490－1508．

［8］Sweldens W．The lifting scheme:a construction of second generation wavelet［J］．Siam Journal on Mathematical Analysis，1998，29(2):511－546．

［9］Davis G，Strela V，Turcajova R．Multiwavelet construction via the lifting scheme［J］．Wavelet Analysis and Multiresolution Methods，Dekker，New York，2000，57 －59．

［10］Keinert F．Raising multiwavelet approximation order through lifting［J］．Siam．Journal on Mathematical Analysis，2001(32):1032－1049．

［11］Wang X D，Zi Y Y ，He Z J．Multiwavelet construction via an adaptive symmetric lifting scheme and its application for rotating machinery fault diagnosis［J］．Measurement Science ＆ Technology，2009，20(2):1－17．

［12］Bui T D ，Chen G Y．Translation-invariant de-noising using multiwavelets［J］．IEEE Transactions on Signal Processing，1998，46(12):414－420．

［13］袁 静，何正嘉，王晓东，等．平移不变多小波相邻系数降噪方法及其在监测诊断中的应用［J］．机械工程学报，2009，45(4):155－160．

［14］Sawalhi N，Randall R B，Endo H．The enhancement of fault detection and diagnosis in rolling element bearings using minimum entropy deconvolution combined with spectral kurtosis［J］．Mechanical Systems and Signal Processing，2007，21(6):2616 －2633．

［15］Lee J Y，Nandi A K．Extraction of impacting signals using blind deconvolution［J］．Journal of Sound and Vibration，2000，232(5):945－962．

［16］Yuan J，He Z J，Zi Y Y，et al．Gearbox fault diagnosis of rolling mills using multiwavelet sliding window neighboring coefficient denoising and optimal blind deconvolution［J］．Science in China Series E:Technological Sciences，2009，52(10):2801－2809．

［17］Wang X D，Zi Y Y，He Z J．Multiwavelet denoising with improved neighboring coefficients for application on rolling bearing［J］．Mechanical Systems and Signal Processing，2011，25(1):285－304．