振动工程学报

冲击信号成分	谐波信号成分
冲击类型	幅值A、B	时间间隔T/s	共振频率f/Hz	谐波编号	幅值C	共振频率f/Hz	初相位φ/rad
故障冲击	1	1/29	2000	谐波1	0.025	7	π/6
随机冲击	3	—	5100	谐波2	0.025	14	−π/3

冲击信号成分

谐波信号成分

冲击类型

幅值A、B

时间间隔T/s

共振频率f/Hz

谐波编号

幅值C

共振频率f/Hz

初相位φ/rad

故障冲击

1/29

2000

谐波1

0.025

π/6

随机冲击

—

5100

谐波2

0.025

−π/3

冲击信号成分	谐波信号成分
冲击类型	幅值A、B	时间间隔T/s	共振频率f/Hz	谐波编号	幅值C	共振频率f/Hz	初相位φ/rad
故障冲击	1	1/29	2000	谐波1	0.025	7	π/6
随机冲击	3	—	5100	谐波2	0.025	14	−π/3

冲击信号成分

谐波信号成分

冲击类型

幅值A、B

时间间隔T/s

共振频率f/Hz

谐波编号

幅值C

共振频率f/Hz

初相位φ/rad

故障冲击

1/29

2000

谐波1

0.025

π/6

随机冲击

—

5100

谐波2

0.025

−π/3

输入为信号x，信号长度为N，初始化j=1
步骤1：令时间长度τ_i = j，对信号x进行时域同步平均处理，得到结果y_j；步骤2：计算y_j的基尼指标GI(y_j)，令j=j+1；步骤3：重复步骤1和2，直至 $τ j = N 2$ ，获得基尼指标集合GI_set={GI(y₁),GI(y₂),⋯,GI(y_N/2)}；步骤4：从基尼指标集合中找出最大的基尼指标值，将其对应的时间长度作为故障信号周期的估计： $T = arg ⁡ max τ j ⁡ G I (y j)$ 。

输入为信号x，信号长度为N，初始化j=1

步骤1：令时间长度τ_i = j，对信号x进行时域同步平均处理，得到结果y_j；
步骤2：计算y_j的基尼指标GI(y_j)，令j=j+1；
步骤3：重复步骤1和2，直至

τ j = N 2

，获得基尼指标集合GI_set={GI(y₁),GI(y₂),⋯,GI(y_N/2)}；
步骤4：从基尼指标集合中找出最大的基尼指标值，将其对应的时间长度作为故障信号周期的估计：

T = arg ⁡ max τ j ⁡ G I (y j)

。

输入为信号x，信号长度为N，初始化j=1
步骤1：令时间长度τ_i = j，对信号x进行时域同步平均处理，得到结果y_j；步骤2：计算y_j的基尼指标GI(y_j)，令j=j+1；步骤3：重复步骤1和2，直至 $τ j = N 2$ ，获得基尼指标集合GI_set={GI(y₁),GI(y₂),⋯,GI(y_N/2)}；步骤4：从基尼指标集合中找出最大的基尼指标值，将其对应的时间长度作为故障信号周期的估计： $T = arg ⁡ max τ j ⁡ G I (y j)$ 。

输入为信号x，信号长度为N，初始化j=1

步骤1：令时间长度τ_i = j，对信号x进行时域同步平均处理，得到结果y_j；
步骤2：计算y_j的基尼指标GI(y_j)，令j=j+1；
步骤3：重复步骤1和2，直至

τ j = N 2

T = arg ⁡ max τ j ⁡ G I (y j)

。

部件	故障特征频率/Hz
内圈	5.41f_r
外圈	3.59f_r
保持架	4.71f_r
滚动体	0.4f_r

部件

故障特征频率/Hz

内圈

5.41f_r

外圈

3.59f_r

保持架

4.71f_r

滚动体

0.4f_r

部件	故障特征频率/Hz
内圈	5.41f_r
外圈	3.59f_r
保持架	4.71f_r
滚动体	0.4f_r

部件

故障特征频率/Hz

内圈

5.41f_r

外圈

3.59f_r

保持架

4.71f_r

滚动体

0.4f_r

周期精炼的最大相关峭度解卷积在滚动轴承微弱故障特征提取中的应用

PDF下载

苗永浩 ¹^,² , 石惠芳 ¹ , 李晨辉 ¹ , 顾晓辉 ²

振动工程学报 | 2025,38(6): 1317-1325

收起

振动工程学报 | 2025, 38(6): 1317-1325

周期精炼的最大相关峭度解卷积在滚动轴承微弱故障特征提取中的应用

全屏

苗永浩¹^,², 石惠芳¹, 李晨辉¹, 顾晓辉²

作者信息

^1.北京航空航天大学可靠性与系统工程学院，北京 100191

^2.石家庄铁道大学省部共建交通工程结构力学行为与系统安全国家重点实验室，河北石家庄 050043

通讯作者:

苗永浩（1992—），男，博士，副教授。E-mail：miaoyonghao@buaa.edu.cn

Period-refined maximum correlated kurtosis deconvolution method for weak fault feature extraction in rolling bearings

Yonghao MIAO¹^,², Huifang SHI¹, Chenhui LI¹, Xiaohui GU²

Affiliations

^1.School of Reliability and Systems Engineering，Beihang University，Beijing 100191，China

^2.State Key Laboratory of Mechanical Behavior and System Safety of Traffic Engineering Structures，Shijiazhuang Tiedao University，Shijiazhuang 050043，China

出版时间: 2025-06-10 doi: 10.16385/j.cnki.issn.1004-4523.2025.06.020

文章导航

摘要

收起

最大相关峭度解卷积（maximum correlated kurtosis deconvolution, MCKD）以相关峭度作为解卷目标，能有效提取兼具周期性和冲击性的机械故障特征，是当前解决滚动轴承故障诊断问题的常用手段。然而，MCKD性能的发挥严重依赖先验故障周期信息，而且现有的解决方案也只关注于迭代周期的估计，难以在低信噪比条件下奏效。为此，提出了一种周期精炼的最大相关峭度解卷积（period-refined MCKD, PRMCKD）方法，借助时域同步平均（time synchronous averaging, TSA）方法进行迭代周期的精炼以用于解卷，能在强噪声情况下准确提取微弱的轴承故障特征。该方法利用滤波器组进行共振频带的初步定位，以明确解卷方向；以相关峭度为目标函数，基于TSA技术精炼出的周期信息，迭代求解最优滤波器系数；通过滤波信号中的故障特征实现滚动轴承故障定位。仿真与试验分析结果表明，所提PRMCKD方法相比传统解卷积方法在滚动轴承微弱故障特征提取方面更具优势。

关键词

滚动轴承 / 故障诊断 / 时域同步平均 / 滤波器初始化 / 解卷积

Abstract

收起

Maximum correlated kurtosis deconvolution (MCKD), which uses correlated kurtosis as its deconvolution target, effectively extracts both periodic and impulsive features of mechanical faults. This is a widely used method for solving rolling bearing fault diagnosis problems. However, the performance of MCKD heavily relies on accurate prior fault period information. Existing solution often only focus on period estimation during the iterative process, making them ineffective under low signal-to-noise ration (SNR) conditions. To address this limitation, a period-refined maximum corrlated kurtosis deconvolution (PRMCKD) method is proposed. This approach refines the iteration period using time synchronous averaging (TSA) for reconolution, enabling accurate extraction of subtle bearing fault features even in strong noise environments. The method operates by first utilizing a filter bank for preliminary localization of the resonance frequency band, thus defining the correct deconvolution direction. With correlated kurtosis as the objective function, and leveraging the period information refined by TSA technology, the optimal filter coefficients are iteratively solved. Rolling bearing fault localization is achieved through the fault features present in the filtered signal. Simulation and experimental analysis results demonstrate that the proposed PRMCKD method offers significant advantages over traditional deconvolution methods for extracting weak fault features in rolling bearings.

Key words

rolling bearing / fault diagnosis / time synchronous averaging / initialization filter / deconvolution

引用本文

苗永浩, 石惠芳, 李晨辉, 顾晓辉. 周期精炼的最大相关峭度解卷积在滚动轴承微弱故障特征提取中的应用. 振动工程学报, 2025 , 38 (6) : 1317 -1325 . DOI: 10.16385/j.cnki.issn.1004-4523.2025.06.020

Yonghao MIAO, Huifang SHI, Chenhui LI, Xiaohui GU. Period-refined maximum correlated kurtosis deconvolution method for weak fault feature extraction in rolling bearings[J]. Journal of Vibration Engineering, 2025 , 38 (6) : 1317 -1325 . DOI: 10.16385/j.cnki.issn.1004-4523.2025.06.020

正文

收起

滚动轴承是旋转机械中最为常见的机械部件，其健康状态与装备的服役安全密切相关。严酷恶劣的运行环境使得滚动轴承成为旋转机械的薄弱环节^[1]。因此，开展滚动轴承故障诊断对保障设备安全运行、提高生产效益来说至关重要。

振动分析是滚动轴承故障诊断的重要手段。当滚动轴承发生局部故障时，其振动信号中会出现一系列的瞬态冲击，这也是诊断轴承故障的重要依据^[2]。事实上，在滚动轴承故障萌发的早期阶段，故障引发的冲击特征能量十分微弱，而复杂传递路径产生的衰减影响和信号测试过程中的各种干扰成分会使得故障特征被淹没。因此，从测试信号中提取出故障相关的特征是实现滚动轴承故障诊断的关键步骤^[3-4]。

解卷积方法通过设置目标函数、迭代选择逆滤波器系数构造滤波器，能自适应提取最优目标特征的信号成分，是故障诊断的常用手段。ENDO等^[5]在机械故障诊断领域引入了最小熵解卷积（minimum entropy deconvolution, MED）方法，使用峭度指标作为目标函数，用于提取振动信号中的故障冲击。然而，峭度指标对单冲击噪声更为敏感，使得MED方法倾向于增强信号中的随机冲击干扰，而非周期性的故障冲击导致该方法广为诟病。

为此，学者们在经典MED方法的基础上开展了大量工作，构造或引入新的指标代替峭度，从而提出多种改进方法。许多性能更加优秀的稀疏性指标被相继提出。例如，PEETERS等^[6]设计了基于L2/L1范数、Hoyer指数和谱负熵三种稀疏性指标的解卷积方法，并通过瑞利熵迭代优化滤波器系数。然而，基于稀疏性指标的方法未能全面评价机械故障的典型特征。为兼顾故障冲击的周期性和冲击性，MCDONALD等^[7]将周期嵌入峭度指标，提出了全新的相关峭度指标，能准确衡量信号的周期性冲击强度，由此提出了最大相关峭度解卷积（maximum correlated kurtosis deconvolution, MCKD）。类似地，BUZZONI等^[8]则选择循环平稳指标作为目标函数，提出了最大二阶循环平稳盲解卷积（maximum second-order cyclostationarity blind deconvolution, CYCBD）。上述兼顾周期的解卷积方法都成功克服了传统MED对于随机冲击噪声敏感的缺点，但它们的优异性能严重依赖精确的先验故障周期，在一定程度上限制了其应用^[9]。

在工业实际中，由于转速信息难以精确测量以及轴承存在随机滑动等原因，故障周期的精确估计或计算变得十分困难。因此，无需先验周期的解决方案被相继提出。MIAO等^[10]利用迭代过程中的滤波信号的自相关函数进行周期估计，成功改进了MCKD，实现了无需先验周期的最大相关峭度解卷积（improved MCKD, IMCKD）的应用。随后，谐波乘积谱被用于识别周期信息，从而改进传统CYCBD方法^[11]。此外，周期性调制强度指标、相关峭度以及多点峭度等周期评估指标也被证实可用于旋转机械的故障周期识别^[12]。然而，上述周期估计手段局限于从给定的信号中识别周期，当信号中的干扰成分能量过高时，这类方法极容易失效，进而难以准确提取出微弱的故障特征。

为了解决现有方法在低信噪比情况下难以奏效的问题，本文提出了一种基于时域同步平均（time synchronous averaging, TSA）的周期精炼方法。该方法首先通过滤波器组初步定位故障频带，随后利用细粒度划分的采样长度对信号进行TSA处理，并依据TSA结果的基尼指标（Gini index, GI）选取最合适的采样长度作为迭代周期。得益于采样长度的细粒度划分和TSA的梳状滤波特性，与已有方法相比，该方法可以在得到准确故障特征周期信息的同时增强信号中的周期成分^[13]。基于所提周期精炼方法和滤波器组初始化策略，本文构建周期精炼的最大相关峭度解卷积（period-refined MCKD, PRMCKD），通过仿真数据与试验数据验证所提PRMCKD方法在滚动轴承微弱故障特征提取方面的有效性和实用性。

1　最大相关峭度解卷积方法

收起

最大相关峭度解卷积方法以相关峭度指标（correlated kurtosis, CK）作为目标函数，能兼顾故障信号成分的周期性和冲击性特征。MCKD的滤波器系数求解目标是最大化滤波信号的CK，该过程可以看作下式的最优化问题：

max f ⁡ C K (y, T s) ， s . t . y = X T f ， C K (y, T s) = ∑ n = 1 N (y n y n − T s) 2 (∑ n = 1 N y n 2) 2

(1)

式中，f为滤波器系数向量；N为原信号长度；T_s为需预先提供的故障信号周期对应的采样点数；CK(y, T_s)为滤波信号y的相关峭度；y_n表示信号y的第n个元素；X为根据原信号x构造的Toeplitz矩阵，其表达式如下：

X = [x 1 x 2 ⋯ x N 0 x 1 ⋯ x N − 1 ⋮ ⋮ ⋱ ⋮ 00 ⋯ x N − L + 1] L × N

(2)

式中，L为滤波器长度。

为了求解最优滤波器f，令目标函数对滤波器的偏导为0：

∂ C K (y, T s) ∂ f = 0

(3)

进而得到最优滤波器系数的表达式：

f = ‖ y ‖ 22 2 ‖ β ‖ 22 (X X T) − 1 (X α 0 + X T α 1)

(4)

式中，

‖ ⋅ ‖ 2

表示L2范数。其中：

α 0 = [y 1 y 1 − T s 2 y 2 y 2 − T s 2 ⋮ y N y N − T s 2], α 1 = [y 12 y 1 − T s y 22 y 2 − T s ⋮ y N 2 y N − T s], β = [y 1 y 1 − T s y 2 y 2 − T s ⋮ y N y N − T s], X T = [x 1 − T s x 2 − T s ⋯ x N − T s 0 x 1 − T s ⋯ x N − T s − 1 ⋮ ⋮ ⋮ 00 ⋯ x N − T s − L + 1] L × N 。

式(4)是关于f的非线性方程，本文使用固定点迭代法求解其最优解。综上，MCKD方法的步骤如下：

步骤1：输入原信号x，选择滤波器长度L和迭代次数M，根据轴承的尺寸和转速计算各部件的故障特征频率，换算出故障信号周期对应的采样点数T；

步骤2：设置解卷积初始滤波器为

f 0

=[0,1,0,

⋯

,0]，令迭代次数i=0；

步骤3：使用滤波器

f i

对原信号x进行滤波，获得滤波信号

y i

；

步骤4：利用原信号x、滤波信号

y i

与先验周期信息T，根据式(4)获得更新后的滤波器系数

f i + 1

；令i=i+1；

步骤5：如果i<M，重复步骤3和4；否则结束算法，得到最优滤波器系数

f M

，对原信号进行滤波得到滤波信号。

2　周期精炼的最大相关峭度解卷积

收起

2.1　滤波器初始化策略

滤波器初始化对于解卷积结果具有显著影响。HE等^[14]证明覆盖故障频带且包含更多故障信息的滤波器初始化方式可以提高解卷过程的效率和解卷结果的准确性。因此，本文提出采用频带均分的初始化策略，构建覆盖全频带的滤波器组，并从中优选出故障含量最为丰富的带通滤波器作为初始滤波器。

为了比较不同初始化策略对于解卷结果的影响，本文首先进行了仿真测试。为了模拟滚动轴承外圈故障，信号x(t)由外圈故障振动冲击d(t)、谐波成分h(t)、随机冲击成分r(t)和背景噪声n(t)组成，其表达式如下：

x (t) = d (t) + h (t) + r (t) + n (t)

(5)

其中：

d (t) = ∑ i = 1 M 0 A i s i (t − i T a − θ i)

(6)

h (t) = ∑ k = 1 2 C k cos (2 π f k t + φ k)

(7)

r (t) = ∑ j = 1 2 B j s j (t − T j)

(8)

s i (t) = e − α i t cos ⁡ (2 π f i t)

(9)

式中，A_i为故障脉冲的幅值；s_i(t)为冲击响应函数；α_i和f_i分别表示系统的阻尼系数（取值为300）和共振频率（取值为2000 Hz）；M₀为故障冲击的数量；T_a为故障冲击的时间间隔；θ_i表示时间波动，用于模拟轴承中滚动体随机滑移对故障冲击时间间隔造成的波动影响，取值为m倍的T_a，其中m服从均匀分布U(0.01, 0.02)；C_k、f_k和φ_k分别为第k个谐波成分的振幅、频率和初相位；B_j和T_j分别为随机冲击信号的幅值和时刻。

另外，仿真信号的采样频率设置为20000 Hz，信噪比设置为−18 dB，其他参数的取值情况如表1所示。图1为仿真信号的时域波形，由图可见，其背景噪声强烈，完全淹没了滚动轴承的故障特征。

传统MCKD方法的解卷性能极大程度上依赖于先验周期的准确性^[9]。IMCKD方法在此基础上引入自相关函数进行故障周期估计，从而摆脱了先验周期的限制，其被用于测试两种滤波器初始化策略的性能。

首先，设置滤波器长度为100，分别实施传统初始化策略（滤波器系数设置为[0,1,0,

⋯

,0]）与带通滤波器策略，即通过MATLAB的fir1函数自动生成汉宁窗滤波器，其截止频率设定为[1500, 2500] Hz，以覆盖故障频带。随后，基于得到的两种初始滤波器，利用IMCKD处理信号。

图2显示了两种滤波器初始化策略下的IMCKD迭代过程中滤波信号CK值的变化。可以清楚地看到，在基于传统初始化策略的IMCKD的迭代过程中，滤波信号的CK值存在巨大的波动且在40 次迭代后仍未收敛。相比之下，带通滤波器策略使得IMCKD仅仅经历10次迭代就收敛至理想状态。这是由于初始值对优化求解具有重要影响。基于两种初始化策略的解卷积的求解效率有所不同：传统初始化滤波器在求解初期需要经过多次迭代才能寻找到最优解的收敛方向，因此需要更多的迭代次数；与之相比，带通滤波器对信号进行一次滤波即可滤除大量的噪声，在求解初期即可迅速锁定目标函数最优解。该仿真结果证实了带通滤波器初始化策略可以有效提升解卷积的收敛速度和故障特征提取的准确性。

然而，在实际应用中，由于滚动轴承早期故障冲击能量微弱，故障频带充满噪声，难以直接定位。因此，本文拟采用基尼指标作为度量准则，初步锁定故障频带位置。

基尼指标能衡量信号的稀疏性，且与峭度指标相比，其对随机冲击的鲁棒性更强，被广泛应用于评价复杂信号中的故障信息^[15-16]。对于长度为N的序列x，其基尼指标的定义如下：

G I (x) = 1 − 2 ∑ k = 1 N | x (k) | ‖ x ‖ 1 (N − k + 0.5 N)

(10)

式中，

| x (1) | ⩽ | x (2) | ⩽ ⋯ ⩽ | x (N) |

为信号中各元素按绝对值升序排列的结果；

‖ x ‖ 1

为x的L1范数。

结合GI和带通滤波器组策略，本文提出的滤波器初始化策略可描述如下：

步骤1：将信号全频带等分为10个相邻的区间，构造每个频率区间对应的汉宁窗滤波器；

步骤2：计算上述滤波器组得到的滤波信号的基尼指标；

步骤3：选取基尼指标最大值对应的带通滤波器作为解卷积的初始滤波器。

2.2　基于时域同步平均的周期精炼方法

尽管传统IMCKD等方法克服了对于先验周期的依赖，但是其优异性能仍难以在低信噪比条件下凸显。这些方法仅仅利用迭代过程中的滤波信号进行周期估计，其估计精度随原信号信噪比的降低而降低，从而难以奏效。鉴于此，本文提出利用时域同步平均方法对迭代滤波信号中的周期特征实现进一步增强，以提高迭代周期的估计精度，达到提取低信噪比下微弱故障特征的目的。

时域同步平均方法通过同一信号中相同长度信号片段的叠加增强特定周期的信号成分，因此对信号进行TSA处理相当于对该信号进行梳状滤波。其定义如下：

y (n) = 1 K ∑ r = 0 K − 1 x (n − r τ)

(11)

式中，x为待处理信号；K为截断的相同长度信号片段的数量；τ为信号片段的长度；y为处理后的信号。

那么，通过TSA处理，信号y中周期为τ的成分被保留，而其他成分被抑制。而且，TSA物理意义明确，操作简单，计算效率高。因此，TSA被用于MCKD迭代过程中迭代滤波信号的周期增强与准确周期估计。

本文依据细粒度划分的时间长度对迭代滤波信号进行TSA处理，使用基尼指标作为衡量准则，选择最大基尼指标对应的时间长度作为估计周期。其具体操作流程如表2所示。

为了比较基于TSA的周期精炼方法与传统基于自相关函数的周期估计方法在不同信噪比下的周期估计性能，本节进行了仿真测试。首先，仿真信号仍使用2.1节中的信号模型，其他参数保持不变，仅调节仿真信号的信噪比，将其从10 dB逐渐降低至−18 dB。利用不同信噪比下的仿真信号对这两种方法进行测试，其周期估计误差随信噪比变化的情况如图3所示。信噪比较高时，两种方法都可以准确估计出故障信号的周期，估计误差几乎为0；而当信噪比降低到−12 dB时，基于自相关函数的估计开始出现较为明显的误差，且随着信噪比的降低，估计误差逐渐变大。反观提出的基于TSA的周期精炼方法，在极低信噪比下也能获得准确、稳定的周期估计值。

2.3　周期精炼的最大相关峭度解卷积

基于以上介绍，本文提出的PRMCKD方法流程如图4所示。

具体步骤描述如下：

步骤1：输入原信号x，选择滤波器长度L和迭代次数M；

步骤2：构造汉宁窗滤波器组，以滤波信号的基尼指标作为度量准则，确定解卷积初始滤波器

f 0

，设定迭代次数i=0；

步骤3：使用滤波器

f i

对原信号x进行滤波，获得滤波信号

y i

；

步骤4：使用所提的基于TSA的周期精炼方法从滤波信号

y i

中得到故障周期的估计值

T i

，并据此构造矩阵X_T；

步骤5：利用滤波信号

y i

、故障周期估计值

T i

和矩阵X_T，根据式(4)获得更新后的滤波器系数

f i + 1

，令i=i+1；

步骤6：如果i<M，返回步骤3；否则进入步骤7；

步骤7：使用滤波器

f i

对原信号进行滤波，对滤波结果进行时域分析和包络谱分析。

需要注意的是，所提PRMCKD方法并不能直接用于非平稳转速下的滚动轴承故障诊断。但是，利用阶次跟踪技术借助转速信息将非平稳信号转换为平稳的角度域信号，可以将所提方法的应用范围扩展到交变转速场合。

3　仿真分析

收起

为了验证所提PRMCKD方法在滚动轴承故障诊断中的有效性，本节进行了仿真分析。并且使用应用广泛的IMCKD和MED两种解卷积方法作为对比参照。

3.1　强背景噪声情况下的故障特征提取

为分析PRMCKD方法面对强背景噪声时的故障特征提取效果，使用2.1节中的轴承故障信号模型和表1的参数配置生成滚动轴承外圈故障仿真信号，信噪比设置为−16 dB。

图5分别展示了仿真信号的时域波形及其包络谱，包络谱中的红色线条表示外圈故障特征频率（BPFO=29 Hz）及其整数倍谐波在谱图中所在的位置。由于强烈的背景噪声，从仿真信号的波形中难以发现故障冲击，其包络谱中含有大量干扰信息，说明该仿真信号具有极低的信噪比。

为了公平比较，对MED、IMCKD方法与所提PRMCKD方法设置相同的参数，滤波器长度为50，迭代次数为30。

图6为3种方法处理相同信号得到的结果。其中，在MED和IMCKD方法的解卷结果中，信号的波形和包络谱含有强烈的噪声，无法识别出周期性的冲击或突出的外圈故障特征频率成分。而PRMCKD方法的滤波信号中周期性的故障冲击清晰可见，其包络谱中也存在极为明显的外圈故障特征频率及其倍频成分，这表明滚动轴承出现了外圈故障，与实际情况相符。

通过以上分析可以得出，相比MED和IMCKD方法，所提PRMCKD方法能够抵抗强背景噪声的影响，从信号中提取微弱的滚动轴承故障特征。

3.2　含有随机冲击噪声情况下的故障特征提取

进一步，为了评价提出的方法对于随机冲击干扰的鲁棒性，设置故障冲击的时间间隔为1/37.7 s，令随机冲击信号的幅值为4，信噪比为−12 dB，其他参数与表1相同。此外，使用频率为7 Hz的谐波对故障信号进行幅值调制，生成如图7所示的滚动轴承内圈故障仿真信号的时域波形及其包络谱，包络谱中的红色线条、黑色线条分别表示内圈故障特征频率（BPFI=37.7 Hz）及其倍频、转轴转频（f_r=7 Hz）以及由转频调制产生的边频在谱图中的位置。由于随机冲击成分能量极高，信号中的故障特征变得十分微弱，在时域波形和包络谱中都难以识别出内圈故障的存在。

对三种解卷积方法使用相同的参数设置，得到如图8所示的解卷积结果。单从三者的滤波信号时域可以看出，由于随机冲击能量的强大，MED与IMCKD方法都倾向于增强两个单一的冲击成分。相比之下，在所提PRMCKD方法得到的滤波信号中可以看到极具辨识度的内圈故障脉冲序列，而包络谱的能量集中在内圈故障特征频率以及倍频处，这些频率周围存在着间隔为7 Hz的边频带，这些特征都暗示着内圈故障的存在。

另外，基于本数值仿真案例，比较了自相关函数法与TSA方法对周期信息估计的性能。图9展示了本仿真分析中，IMCKD和PRMCKD迭代过程中的故障特征周期估计值及其误差。可以发现，IMCKD方法估计的故障周期值收敛于真实值的8.88倍左右，显然与真实故障周期相差甚远，这也是该方法难以提取故障信息的根源。与之相比，所提方法得益于TSA对信号周期信息的精准精炼，能准确估计故障周期并用于解卷滤波，得到理想的诊断结果。

PRMCKD和IMCKD方法的算法复杂度处于同一量级。与此同时，本文还重复多次数值仿真用于评价3种解卷积方法的运行效率。测试结果显示，在运行环境Intel Core (TM) i7-10700 （2.90 GHz）CPU，16 GB RAM下运用MATLAB 2018a处理数据长度为20000采样点的信号，MED、IMCKD和PRMCKD方法的运行时间分别在0.7、1.0和1.8 s左右波动。可见，所提PRMCKD方法在提升性能的同时并没有牺牲大量的计算效率。

该仿真案例不仅证明了PRMCKD方法在故障信息提取方面的优势，也证实了PRMCKD方法采用的周期精炼方法比经典的自相关函数周期估计法具有更出色的周期估计能力。

4　试验分析

收起

为了验证所提方法在实际应用中的有效性，本节基于风力发电齿轮箱故障模拟试验平台的实测数据开展验证工作。图10为风力发电机滚动轴承试验台的实物图。其工作原理为：驱动电机提供动力，通过行星齿轮箱和平行轴齿轮箱进行扭矩传递，磁粉制动器用于加载和制动。其中，滚动轴承故障试件安装于轴承座中，其外圈固定，内圈随转轴匀速旋转。

根据滚动轴承的尺寸参数计算出各部件的故障特征频率并列于表3中，其中f_r表示转轴转频。滚动轴承的振动信号由轴承座上安装的加速度传感器测试得到，采用的加速度传感器为PCB 356A32型三轴加速度传感器，采样频率为20000 Hz。

4.1　轴承外圈故障试验分析

图11为利用线切割加工得到的含外圈裂纹故障的滚动轴承。将该轴承替换试验台中的正常轴承用于试验，其中，转轴转频设置为20 Hz，截取1 s长的测试信号，其时域波形及其包络谱如图12所示，包络谱中的红色线条表示外圈故障特征频率（BPFO=71.8 Hz）及其整数倍谐波在谱图中所在的位置。由于损伤微小，故障能量非常微弱，信号的波形和包络谱中充斥着噪声成分，导致信噪比极低。

同样对3种方法设置相同参数，即滤波器长度为50，迭代次数为30，得到如图13所示的结果。可以发现，MED和IMCKD方法都不能从原信号中提取出故障特征，意味着这两种方法受到了强背景噪声的影响。因此，所提PRMCKD方法不仅能提取出具有高信噪比的冲击成分，其包络谱也证实了它的有效性。

4.2　轴承内圈故障试验分析

使用如图14所示的内圈磨损故障轴承进行测试。试验中转轴的转频为10 Hz，截取1 s长的测试信号，其时域波形及其包络谱如图15所示。包络谱中的红色线条、黑色线条分别表示内圈故障特征频率（BPFI=54.1 Hz）及其倍频、转轴转频（f_r=10 Hz）以及由于转频调制产生的边频在谱图中的位置。从其时域和包络谱中都难以发现故障信息，而且该测试信号中有能量突出的随机冲击干扰。

使用相同的解卷积参数设置，图16展示了3种解卷积方法处理测试信号的结果。从图16(a)~(d)中可以看出，MED和IMCKD方法受原信号中的随机冲击影响明显，尽管包络谱中能看到故障特征频率，但是滤波信号的能量完全集中于随机冲击。而所提PRMCKD方法可以不受随机冲击的影响，有效提取故障冲击成分，并且转频调制现象在时域和包络谱中都能清晰呈现。

以上试验结果可以证实，相比于传统MED和IMCKD方法，所提PRMCKD方法对随机冲击干扰具有极强的鲁棒性，进一步凸显了其有效性与优势。

5　结论

收起

（1）本文提出基于GI和带通滤波器组的滤波器初始化策略，并通过仿真证明了该初始化策略可以明显提高解卷积方法的收敛速度和全局最优性。

（2）本文提出的基于TSA和GI的周期精炼方法可以对迭代滤波信号进行特征增强，提高周期估计的准确度，保障所提方法在极低信噪比下的有效性。

（3）仿真与试验分析表明，相比传统MED和IMCKD方法，所提PRMCKD方法具有更强的抑制强背景噪声和随机冲击的能力，因此PRMCKD方法更适用于滚动轴承的微弱故障特征提取，有助于实现滚动轴承的早期故障诊断。

基金

收起

国家重点研发计划资助项目(2021YFB2500604)

参考文献

收起

文献

收起

参考文献引证文献

排序方式：

[1]

ZHANG

, ZHAO

B Y

, LIN

. Machine learning based bearing fault diagnosis using the case western reserve university data: a review[J]. IEEE Access, 2021, 9: 155598-155608.

[2]

MCFADDEN

P D

, SMITH

J D

. Model for the vibration produced by a single point defect in a rolling element bearing[J]. Journal of Sound and Vibration, 1984, 96(1): 69-82.

[3]

陈仁祥，黄鑫, 杨黎霞，等.基于卷积神经网络和离散小波变换的滚动轴承故障诊断[J]. 振动工程学报, 2018, 31(5):883-891.

CHEN

Renxiang

, HUANG

Xin

, YANG

Lixia

, et al.Rolling bearing fault identification based on convolution neural network and discrete wavelet transform[J]. Journal of Vibration Engineering, 2018, 31(5): 883-891.

[4]

吴超，崔玲丽, 张建宇，等.改进MBCV法在滚动轴承故障诊断中的应用[J]. 振动工程学报, 2022, 35(4):942-948.

Chao

, CUI

Lingli

, ZHANG

Jianyu

, et al.Fault diagnosis method of rolling bearing based on improved MBCV method[J]. Journal of Vibration Engineering, 2022, 35(4): 942-948.

[5]

ENDO

, RANDALL

R B

. Enhancement of autoregressive model based gear tooth fault detection technique by the use of minimum entropy deconvolution filter[J]. Mechanical Systems and Signal Processing, 2007, 21(2): 906-919.

[6]

PEETERS

, ANTONI

, HELSEN

. Blind filters based on envelope spectrum sparsity indicators for bearing and gear vibration-based condition monitoring[J]. Mechanical Systems and Signal Processing, 2020, 138: 106556.

[7]

MCDONALD

G L

, ZHAO

, ZUO

M J

. Maximum correlated Kurtosis deconvolution and application on gear tooth chip fault detection[J]. Mechanical Systems and Signal Processing, 2012, 33: 237-255.

[8]

BUZZONI

, ANTONI

, D’ELIA

. Blind deconvolution based on cyclostationarity maximization and its application to fault identification[J]. Journal of Sound and Vibration, 2018, 432: 569-601.

[9]

黄包裕，张永祥, 赵磊. 基于布谷鸟搜索算法和最大二阶循环平稳盲解卷积的滚动轴承故障诊断方法[J]. 机械工程学报, 2021, 57(9):99-107.

HUANG

Baoyu

, ZHANG

Yongxiang

, ZHAO

Lei

. Research on fault diagnosis method of rolling bearings based on cuckoo search algorithm and maximum second order cyclostationary blind deconvolution[J]. Journal of Mechanical Engineering, 2021, 57(9): 99-107.

[10]

MIAO

Y H

, ZHAO

, LIN

, et al.Application of an improved maximum correlated kurtosis deconvolution method for fault diagnosis of rolling element bearings[J]. Mechanical Systems and Signal Processing, 2017, 92: 173-195.

[11]

ZHANG

B Y

, MIAO

Y H

, LIN

, et al.Adaptive maximum second-order cyclostationarity blind deconvolution and its application for locomotive bearing fault diagnosis[J]. Mechanical Systems and Signal Processing, 2021, 158: 107736.

[12]

CHEN

B Y

, ZHANG

W H

, SONG

D L

, et al.Blind deconvolution assisted with periodicity detection techniques and its application to bearing fault feature enhancement[J]. Measurement, 2020, 159: 107804.

[13]

BRAUN

. The synchronous(time domain) average revisited[J]. Mechanical Systems and Signal Processing,2011,25(4): 1087-1102.

[14]

, YI

, WANG

, et al.Optimized minimum generalized Lp/Lq deconvolution for recovering repetitive impacts from a vibration mixture[J]. Measurement, 2021, 168: 108329.

[15]

HOU

B C

, WANG

, XIA

T B

, et al.Generalized Gini indices: complementary sparsity measures to Box-Cox sparsity measures for machine condition monitoring[J]. Mechanical Systems and Signal Processing, 2022, 169: 108751.

[16]

WANG

, PENG

Z K

, XI

L F

. The sum of weighted normalized square envelope: a unified framework for kurtosis, negative entropy, Gini index and smoothness index for machine health monitoring[J]. Mechanical Systems and Signal Processing, 2020, 140: 106725.

2025年第38卷第6期

PDF下载

引用本文

BibTeX

文章信息

doi: 10.16385/j.cnki.issn.1004-4523.2025.06.020

接收时间：2024-01-26
首发时间：2026-02-12
出版时间：2025-06-10

补充材料

相关文章

文章信息

作者

出版历史

收稿日期：2024-01-26
修回日期：2024-05-30

基金

国家重点研发计划资助项目(2021YFB2500604)

作者信息

^1.北京航空航天大学可靠性与系统工程学院，北京 100191

^2.石家庄铁道大学省部共建交通工程结构力学行为与系统安全国家重点实验室，河北石家庄 050043

通讯作者:

苗永浩（1992—），男，博士，副教授。E-mail：miaoyonghao@buaa.edu.cn

参考文献

分享链接

https://castjournals.cast.org.cn/joweb/zdgcxb/CN/10.16385/j.cnki.issn.1004-4523.2025.06.020

分享至

全文二维码

扫描看全文

引用本文

BibTeX

本文的引用情况

2种不同金属材料的力学参数

科 Family	属数 Number of genus	种数 Number of species	占总种数比例 Percentage of total species (%)	属 Genus	种数 Number of species	占总种数比例 Percentage of total species (%)
鹅膏菌科Amanitaceae	2	11	5.26	鹅膏菌属 Amanita	10	4.78
小菇科 Mycenaceae	2	12	5.74	丝盖伞属 Inocybe	5	2.39
多孔菌科 Polyporaceae	8	14	6.70	蜡蘑属 Laccaria	5	2.39
红菇科 Russulaceae	3	23	11.00	小皮伞属 Marasmius	6	2.87
				小菇属 Mycena	11	5.26
				光柄菇属 Pluteus	5	2.39
				红菇属 Russula	17	8.13
				栓菌属 Trametes	5	2.39

关闭全屏

BibTeX
EndNote
RefWorks
TxT

冲击信号成分	谐波信号成分
冲击类型	幅值A、B	时间间隔T/s	共振频率f/Hz	谐波编号	幅值C	共振频率f/Hz	初相位φ/rad
故障冲击	1	1/29	2000	谐波1	0.025	7	π/6
随机冲击	3	—	5100	谐波2	0.025	14	−π/3

冲击信号成分

谐波信号成分

冲击类型

幅值A、B

时间间隔T/s

共振频率f/Hz

谐波编号

幅值C

共振频率f/Hz

初相位φ/rad

故障冲击

1/29

2000

谐波1

0.025

π/6

随机冲击

—

5100

谐波2

0.025

−π/3

输入为信号x，信号长度为N，初始化j=1
步骤1：令时间长度τ_i = j，对信号x进行时域同步平均处理，得到结果y_j；步骤2：计算y_j的基尼指标GI(y_j)，令j=j+1；步骤3：重复步骤1和2，直至 $τ j = N 2$ ，获得基尼指标集合GI_set={GI(y₁),GI(y₂),⋯,GI(y_N/2)}；步骤4：从基尼指标集合中找出最大的基尼指标值，将其对应的时间长度作为故障信号周期的估计： $T = arg ⁡ max τ j ⁡ G I (y j)$ 。

输入为信号x，信号长度为N，初始化j=1

步骤1：令时间长度τ_i = j，对信号x进行时域同步平均处理，得到结果y_j；
步骤2：计算y_j的基尼指标GI(y_j)，令j=j+1；
步骤3：重复步骤1和2，直至

τ j = N 2

T = arg ⁡ max τ j ⁡ G I (y j)

。

部件	故障特征频率/Hz
内圈	5.41f_r
外圈	3.59f_r
保持架	4.71f_r
滚动体	0.4f_r

部件

故障特征频率/Hz

内圈

5.41f_r

外圈

3.59f_r

保持架

4.71f_r

滚动体

0.4f_r