CN112577746A - 一种转速波动下滚动轴承包络阶次谱故障特征的提取方法 - Google Patents

一种转速波动下滚动轴承包络阶次谱故障特征的提取方法 Download PDF

Info

Publication number
CN112577746A
CN112577746A CN202011427884.0A CN202011427884A CN112577746A CN 112577746 A CN112577746 A CN 112577746A CN 202011427884 A CN202011427884 A CN 202011427884A CN 112577746 A CN112577746 A CN 112577746A
Authority
CN
China
Prior art keywords
time
signal
fault characteristic
spectrum
fault
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202011427884.0A
Other languages
English (en)
Other versions
CN112577746B (zh
Inventor
胡建中
方波
许飞云
贾民平
黄鹏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Southeast University
Original Assignee
Southeast University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Southeast University filed Critical Southeast University
Priority to CN202011427884.0A priority Critical patent/CN112577746B/zh
Publication of CN112577746A publication Critical patent/CN112577746A/zh
Application granted granted Critical
Publication of CN112577746B publication Critical patent/CN112577746B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
    • G01M13/00—Testing of machine parts
    • G01M13/04—Bearings
    • G01M13/045—Acoustic or vibration analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/02—Preprocessing
    • G06F2218/04—Denoising
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/08—Feature extraction

Landscapes

  • Physics & Mathematics (AREA)
  • Acoustics & Sound (AREA)
  • General Physics & Mathematics (AREA)
  • Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
  • Testing Of Devices, Machine Parts, Or Other Structures Thereof (AREA)

Abstract

本发明公开了一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，解决了由于滚动轴承转速估计不准确导致无法精确获取故障特征阶次的技术问题，其技术方案要点是先对振动信号滤波降噪，然后在时频域构建时变带通滤波器对包络信号进行滤波得到滤波后的包络信号，通过计算相邻显著性极大值时间间隔的倒数得到瞬时故障特征频率序列，对瞬时故障特征频率序列进行拟合得到故障特征频率曲线，消除由测量、滚动体打滑和计算等原因造成的误差。该方案能够精确估计故障特征频率曲线，而不依赖转速测量传感器。获取的包络信号阶次谱精度高，故障特征阶次清晰明显，提高了转速波动下滚动轴承的故障诊断的精度与可靠度。

Description

一种转速波动下滚动轴承包络阶次谱故障特征的提取方法
技术领域
本公开涉及滚动轴承特征提取技术领域，尤其涉及一种转速波动下滚动轴承包络阶次谱故障特征的提取方法。
背景技术
滚动轴承是旋转机械的重要部件，其失效是造成机械故障的最常见原因之一。为了保证设备运行的可靠性，减少安全隐患，需要对滚动轴承的运行状态进行检测，包络谱(或者平方包络谱)能够提取信号的故障特征频率，是滚动轴承故障诊断领域最常用的技术之一。但是包络谱分析只适用于恒定转速的场合，在转速波动的工况下，故障特征频率发生混叠或者完全消失，此时包络谱无法提供有用的诊断信息。针对转速波动工况下滚动轴承的故障诊断问题，广泛使用的方法是利用阶次跟踪算法将信号从时域变换到与转轴相锁的角度域，时域信号的故障特征频率对应角域信号的故障特征阶次。由于故障特征阶次与轴转速无关，因此角域信号的包络阶次谱能够实现变转速滚动轴承故障诊断。
阶次跟踪算法的关键是获取参考轴的转速随时间变化的曲线。一种思路是通过转速表或编码器提供实时转速信息，这不仅增加了测量成本，也给安装和调整带来不便，在某些情况下，甚至不可能安装那些传感器。另一种思路是从振动信号中提取转速信息，现有方法是从振动信号的包络时频谱中提取瞬时故障特征频率(Instantaneous FaultCharacteristic Frequency, IFCF)或者瞬时转频，但转速精度一方面受时频分析方法时间分辨率和频率分辨率的制约，另一方面对时频分析的参数敏感；因此这种方法估计的转速精度不高，信号从时间域变换到角度域误差较大，导致包络阶次谱中特征阶次谱线不能与其他成分区分开来，提取的阶次谱特征不明显，影响诊断精度的同时，制约了诊断的可靠度。
发明内容
本公开提供了一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，其技术目的是精确估计转速曲线，使得包络信号阶次谱中特征阶次谱线清晰明显，从而提高转速波动下滚动轴承故障特征提取的精度与可靠度。
本公开的上述技术目的是通过以下技术方案得以实现的：
一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，包括：
S1：以等时间间隔Ts采集转速波动下故障滚动轴承的振动信号，记录额定转速v，依据轴承型号确定故障特征阶次o；
S2：采用快速谱峭度图算法对所述振动信号进行滤波，得到滤波后的第一信号，使用最小熵反卷积算法对所述第一信号进行滤波，得到降噪后的第二信号；
S3：通过希尔伯特变换获取所述第二信号的包络信号，通过短时傅里叶变换获取所述包络信号的第一时频谱，通过峰值搜索算法从所述第一时频谱中提取第一瞬时故障特征频率序列，对所述第一瞬时故障特征频率序列进行拟合，得到第一故障特征频率曲线；设置时变带通滤波器的带宽参数和时变中心频率，所述时变中心频率为所述第一故障特征频率曲线，在所述第一时频谱中提取所述时变带通滤波器范围内的时频谱，得到所述包络信号的第二时频谱，对所述第二时频谱进行短时傅里叶逆变换得到第三信号；
S4：求解所述第三信号中的极大值，根据相邻极大值之间的时间间隔估计第二瞬时故障特征频率序列，对所述第二瞬时故障特征频率序列进行拟合，得到第二故障特征频率曲线，结合所述第二故障特征频率曲线和所述故障特征阶次o，得到转速曲线；
S5：依据所述转速曲线结合三次样条插值技术得到等相位间隔的采样时间序列，依据所述采样时间序列对所述包络信号进行三次样条插值得到角度域包络信号，对所述角度域包络信号进行离散傅里叶变换获得所述角度域包络信号阶次谱，提取所述角度域包络信号阶次谱中的故障特征。
进一步地，所述步骤S2具体包括：
构建所述振动信号的快速谱峭度图，选择所述快速谱峭度图中峭度最大的通频带，根据所述通频带构建带通滤波器对所述振动信号进行滤波，得到滤波后的第一信号x(t)；
构建滤波器F＝{f1,f2,...,fl,...,fL}，L表示所述滤波器的长度，fl表示所述滤波器的滤波器系数，则使用最小熵反卷积算法对所述第一信号x(t)进行滤波的过程表示为：
Figure BDA0002819812700000031
y(t)即为所述第一信号x(t)滤波后得到的第二信号；
其中，对所述第一信号x(t)进行滤波的迭代过程包括：随机初始化所述滤波器的滤波器系数得到滤波器F1＝{f11,f12,...,f1l,...,f1L}，再由式(1)得到y1(t)，将y1(t)的峭度作为目标函数则有
Figure 1
式(2)中OF1为y1(t) 的峭度，N1表示的是信号y1(t)的长度，当
Figure BDA0002819812700000042
时，获取OF1的最大值并得到对应的滤波器的滤波器系数，即F2＝{f21,f22,...,f2l,...,f2L}，再根据式(1)、 (2)、(3)反复进行迭代得到第二信号y(t)。
进一步地，所述步骤S3具体包括：
S31：通过希尔伯特变换获取所述第二信号y(t)的解析信号，再根据所述解析信号获取所述第二信号y(t)的包络信号，即：sa(t)＝y(t)+iH[y(t)] (4)， sA(t)＝|sa(t)|(5)；式中H[·]表示希尔伯特变换，i为虚数单位，sa(t)表示所述解析信号，sA(t)表示所述包络信号；
S32：通过短时傅里叶变换获取所述包络信号sA(t)的第一时频谱，即：
Figure BDA0002819812700000043
其中，j表示虚数单位，g(t-τ)表示中心位于τ时刻的时间窗，STFTs(τ,f)表示所述包络信号sA(t)的第一时频谱；
S33：根据所述额定转速v和所述故障特征阶次o估计故障特征频率范围，获取所述STFTs(τ,f)在所述故障特征频率范围内所有τ时刻对应的频率幅值最大处的频率值，得到第一瞬时故障特征频率序列；
S34：使用加权最小二乘法和二次多项式局部回归模型对所述第一瞬时故障特征频率序列进行拟合，得到第一故障特征频率曲线；
在拟合的过程中，去除拟合残差大于第一预设阈值的数据点，重复拟合，直至拟合残差大于第一预设阈值的数据点的数量为0，得到最终的第一故障特征频率曲线；
S35：设置时变带通滤波器的带宽参数为ΔF、时变中心频率为所述第一故障特征频率曲线，在所述包络信号sA(t)的所述第一时频谱STFTs(τ,f)中提取所述时变带通滤波器范围内的时频谱，得到所述包络信号sA(t)的第二时频谱，对所述第二时频谱进行短时傅里叶逆变换，得到时变带通滤波后的第三信号 sb(t)。
进一步地，所述步骤S4具体包括：
S41：标记所述第三信号sb(t)的所有极值，获取所有极大值与其相邻极小值之间的两个高度差，将两个高度差中较小的高度差定义为所述极大值的显著性水平，则所述显著性水平与所述极值中最大值的比值表示相对显著性水平，使用Xn表示所述极值的序列，n∈[1,Z]，Z表示极值的总数，若Xk为Xn中的一个极大值，则Xk-1、Xk+1为Xk的相邻极小值，则所述Xk的显著性水平为P＝min(Xk-Xk-1,Xk-Xk+1) (7)，所述Xk的相对显著性水平为
Figure BDA0002819812700000051
Figure BDA0002819812700000052
其中，max(·)表示求最大值，k∈n；
S42：提取所述P'中大于第二预设阈值的极大值序列所对应的时间序列 Tm，则Tm＝{t1,t2,...,tM}，M表示所述时间序列Tm的长度，有tq,tq+1∈Tm，则时刻
Figure BDA0002819812700000061
的瞬时故障特征频率为
Figure BDA0002819812700000062
根据式(9)计算所述时间序列Tm所有时间中点的瞬时故障特征频率，得到第二瞬时故障特征频率序列；其中，所述第二预设阈值∈(0.05,0.4)；
S43：使用加权最小二乘法和二次多项式局部回归模型对所述第二瞬时故障特征序列进行拟合，得到第二故障特征频率曲线；
在拟合过程中，去除拟合残差大于第三预设阈值的数据点，重复拟合，直至拟合残差大于第三预设阈值的数据点的数量为0，得到最终的第二故障特征频率曲线fc(t)；
S44：结合所述故障特征阶次o，得到转速曲线
Figure BDA0002819812700000063
进一步地，所述步骤S5具体包括：
S51：根据所述fr(t)得到相位-时间映射
Figure BDA0002819812700000064
Figure BDA0002819812700000065
则角度域采样的相位间隔
Figure BDA0002819812700000066
其中fs表示采样频率，tend表示
Figure BDA0002819812700000067
的信号时间长度；
S52：通过对相位-时间映射
Figure BDA0002819812700000068
进行三次样条插值得到所有等相位间隔序列
Figure BDA0002819812700000069
对应的时间序列{t1,t2,...,tend}；
S53：将时间序列{t1,t2,...,tend}作为插值点，对所述包络信号sA(t)进行三次样条插值，得到等角度采样的角度域包络信号sB(t)；
S54：通过离散傅里叶变换求解所述角度域包络信号sB(t)阶次谱；
S55：提取所述角度域包络信号sB(t)阶次谱中的故障特征。
本公开的有益效果在于：时变带通滤波信号显著性极大值的横坐标与滚动轴承产生故障冲击的时刻相对应，本发明首先通过快速谱峭度图算法和最小熵反卷积算法依次对振动信号进行滤波，然后在时频域构建时变带通滤波器对包络信号进行滤波得到滤波后的包络信号，通过计算相邻显著性极大值时间间隔的倒数得到瞬时故障特征频率序列，对瞬时故障特征频率序列进行拟合得到故障特征频率曲线，消除由测量、滚动体打滑和计算等原因造成的误差，从而能够精确估计故障特征频率曲线，而不依赖转速测量传感器。获取的包络信号阶次谱精度高，故障特征阶次清晰明显，提高了转速波动下滚动轴承的故障诊断的精度与可靠度。
附图说明
图1为本公开方法流程图；
图2为滚动轴承振动信号及其局部放大图；
图3为振动信号的快速谱峭度图；
图4为降噪后的第二信号及其局部放大图；
图5为第二信号的局部包络信号示意图；
图6为包络信号的第一时频谱局部图；
图7为第一瞬时故障特征频率、第二故障特征频率曲线及其异常点示意图；
图8为时变带通滤波之后的第三信号的局部示意图；
图9为第二瞬时故障特征频率、第二故障特征频率曲线及其异常点示意图；
图10为角度域包络信号示意图；
图11为角度域包络信号阶次谱示意图。
具体实施方式
下面将结合附图对本公开技术方案进行详细说明。在本公开的描述中，需要理解地是，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量，仅用来区分不同的组成部分。
图1为本公开方法流程图，具体包括：S1：以等时间间隔Ts采集转速波动下故障滚动轴承的振动信号，记录额定转速v，依据轴承型号确定故障特征阶次o。S2：采用快速谱峭度图算法对所述振动信号进行滤波，得到滤波后的第一信号，使用最小熵反卷积算法对所述第一信号进行滤波，得到降噪后的第二信号。S3：通过希尔伯特变换获取所述第二信号的包络信号，通过短时傅里叶变换获取所述包络信号的第一时频谱，通过峰值搜索算法在所述第一时频谱中提取第一瞬时故障特征频率序列，对所述第一瞬时故障特征频率序列进行拟合，得到第一故障特征频率曲线；设置时变带通滤波器的带宽参数和时变中心频率，所述时变中心频率为所述第一故障特征频率曲线，在所述第一时频谱中提取所述时变带通滤波器范围内的时频谱，得到所述包络信号的第二时频谱，对所述第二时频谱进行短时傅里叶逆变换得到第三信号。S4：求解所述第三信号中的极大值，根据相邻极大值之间的时间间隔估计第二瞬时故障特征频率序列，对所述第二瞬时故障特征频率序列进行拟合，得到第二故障特征频率曲线，结合所述第二故障特征频率曲线和所述故障特征阶次o，得到转速曲线。S5：依据所述转速曲线结合三次样条插值技术得到等相位间隔的采样时间序列，依据所述采样时间序列对所述包络信号进行三次样条插值得到角度域包络信号，对所述角度域包络信号进行离散傅里叶变换获得所述角度域包络信号阶次谱，提取所述角度域包络信号阶次谱中的故障特征。
作为具体实施例地，S1：本发明通过加速度传感器采集滚动轴承振动信号，本实施例中轴承信息及采样信息如表1所述，采集的振动信号如图2所示，图2右上角为矩形框对应的局部放大图。
轴承信号 SKF6205
故障类型 内圈故障
故障特征阶次 5.415
采样频率 25600Hz
采样时间 5s
额定转速 1632r/min
表1
S2：采用快速谱峭度图算法获取振动信号的快速谱峭度图如图3所示，振动信号的最优共振带为[6400，9600]Hz。构建相应带通的有限长单位冲激响应(Finite ImpulseResponse，FIR)滤波器对振动信号进行滤波处理得到第一信号，采用最小熵反卷积算法对上述第一信号进行滤波降噪。经过带通滤波和最小熵反卷积降噪后的第二信号如图4所示，由于通频带不是能量集中带，经过带通滤波信号的幅值大幅衰减，但是信号的冲击特征得到明显增强。
S3：时变带通滤波：希尔伯特变换获取第二信号的包络信号后，通过短时傅里叶变换求包络信号的第一时频谱，然后沿时间轴拟合故障特征频率的峰值曲线，以此确定时变带通滤波的时变中心频率和带宽，计算带通范围内的第二时频谱，进行短时傅里叶逆变换。具体包括以下步骤：
S31：通过希尔伯特变换获得第二信号的包络信号sA(t)，0.5秒到0.6 秒的局部包络信号如图5所示。
S32：通过短时傅里叶变换获取包络信号sA(t)的第一时频谱STFTs(τ,f)，已知额定转速v＝1632r/min，即27.2Hz，瞬时转速在v附近波动，故障特征阶次为5.415，因此可大致估计故障特征频率的范围。本实施例选择第一时频谱中100-200Hz范围内的数据分析，包络信号的第一时频谱局部图如图6 所示。
S33：求解STFTs(τ,f)所有τ时刻对应的频率轴中频率的最大值，即通过峰值算法搜索从第一时频谱中提取第一瞬时故障特征频率，即图7中细实线表示的IFCF。
S34：使用加权最小二乘法和二次多项式局部回归模型对第一瞬时故障特征频率序列进行拟合，得到平滑的第一故障特征频率曲线。
在拟合的过程中，去除拟合残差大于第一预设阈值的数据点(异常点)，重复拟合，直至拟合残差大于第一预设阈值的数据点的数量为0，得到最终的第一故障特征频率曲线。第一预设阈值可根据实际情况而定，例如可以是标准差的3倍。所有的异常点在图7中标出，去除异常点后的拟合曲线为图 7中连续的曲线。
S35：本实施例设置带宽参数ΔF＝30Hz，提取第一时频谱中以第一故障特征频率曲线为中心、宽度为ΔF的第二时频谱，对所述第二时频谱进行短时傅里叶逆变换，得到时变带通滤波之后的第三信号sb(t)，如图8所示。
S4:计算转速曲线：求第三信号sb(t)的极大值点的相对显著性水平(即显著性极大值)，根据相邻极大值之间的时间间隔计算故障特征频率曲线，剔除异常值后进行平滑处理，之后结合额定转速v计算转速曲线。具体包括以下步骤：
S41:标记第三信号sb(t)的所有极值(极大值与极小值)，求解极大值与相邻极小值之间的两个高度差，将其中较小的高度差定义为该极大值的显著性水平，显著性水平与极值中最大值的比值定义为相对显著性水平，设极值序列为Xn，n∈[1,Z]，Z代表极值的总数量，设Xk为Xn中的某个极大值：
P＝min(Xk-Xk-1,Xk-Xk+1)；
Figure BDA0002819812700000111
式中Xk-1，Xk+1为Xk的相邻极小值，P代表Xk的显著性水平，P'代表Xk的相对显著性水平；
S42:提取P'中大于第二预设阈值的极大值序列所对应的时间序列Tm， Tm＝{t1,t2,...,tM}，M表示所述时间序列Tm的长度；设tq,tq+1∈Tm，则
Figure BDA0002819812700000121
时刻的瞬时故障特征频率为：
Figure BDA0002819812700000122
依据式(9)计算时间序列Tm所有时间中点的瞬时故障特征频率，得到第二瞬时故障特征频率序列；其中，第二预设阈值∈(0.05,0.4)；
S43:使用加权最小二乘法和二次多项式局部回归模型对第二瞬时故障特征序列进行拟合，得到平滑的第二故障特征频率曲线；
在拟合的过程中，去除拟合残差大于第三预设阈值的数据点(异常点)，重复拟合，直至拟合残差大于第三预设阈值的数据点的数量为0，得到最终的第二故障特征频率曲线，结果如图9所示。第三预设阈值可根据实际情况而定，例如可以是标准差的3倍。
S44:设拟合后的第二故障特征频率曲线为fc(t)，结合上述获取的故障特征阶次o，则转频曲线fr(t)为：
Figure BDA0002819812700000123
S5：包络信号阶次谱分析：通过等角度采样将包络信号变换到角度域，通过离散傅里叶变换获得包络阶次谱，提取故障特征。具体包括以下步骤：
S51：根据fr(t)得到相位-时间函数
Figure BDA0002819812700000131
Figure BDA0002819812700000132
为了保证角度域采样数据量与原包络信号数据量一致，确定角度域采样的相位间隔
Figure BDA0002819812700000133
其中fs表示采样频率，tend表示
Figure BDA0002819812700000134
的信号时间长度；实施例中fs为 25600Hz，tend为5s，可求得
Figure BDA0002819812700000135
为0.0066rad。
S52：通过对相位-时间映射
Figure BDA0002819812700000136
进行三次样条插值求得所有等相位间隔序列
Figure BDA0002819812700000137
对应的时间序列{t1,t2,...,tend}。
S53：将时间序列{t1,t2,...,tend}作为插值点，对包络信号sA(t)进行三次样条插值，得到等角度采样的角度域包络信号sB(t)，sB(t)曲线及局部放大图如图 10所示。
S54：通过离散傅里叶变换求解角度域包络信号sB(t)阶次谱，该阶次谱如图11所示。
S55：从图11中可以看出，滚动轴承的转轴特征阶次、内圈故障特征阶次等成分非常突出。
以上为本公开示范性实施例，本公开的保护范围由权利要求书及其等效物限定。

Claims (5)

1.一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，其特征在于，包括：
S1：以等时间间隔Ts采集转速波动下故障滚动轴承的振动信号，记录额定转速v，依据轴承型号确定故障特征阶次o；
S2：采用快速谱峭度图算法对所述振动信号进行滤波，得到滤波后的第一信号，使用最小熵反卷积算法对所述第一信号进行滤波，得到降噪后的第二信号；
S3：通过希尔伯特变换获取所述第二信号的包络信号，通过短时傅里叶变换获取所述包络信号的第一时频谱，通过峰值搜索算法从所述第一时频谱中提取第一瞬时故障特征频率序列，对所述第一瞬时故障特征频率序列进行拟合，得到第一故障特征频率曲线；设置时变带通滤波器的带宽参数和时变中心频率，所述时变中心频率为所述第一故障特征频率曲线，在所述第一时频谱中提取所述时变带通滤波器范围内的时频谱，得到所述包络信号的第二时频谱，对所述第二时频谱进行短时傅里叶逆变换得到第三信号；
S4：求解所述第三信号中的极大值，根据相邻极大值之间的时间间隔估计第二瞬时故障特征频率序列，对所述第二瞬时故障特征频率序列进行拟合，得到第二故障特征频率曲线，结合所述第二故障特征频率曲线和所述故障特征阶次o，得到转速曲线；
S5：依据所述转速曲线结合三次样条插值技术得到等相位间隔的采样时间序列，依据所述采样时间序列对所述包络信号进行三次样条插值得到角度域包络信号，对所述角度域包络信号进行离散傅里叶变换获得所述角度域包络信号阶次谱，提取所述角度域包络信号阶次谱中的故障特征。
2.如权利要求1所述的一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，其特征在于，所述步骤S2具体包括：
构建所述振动信号的快速谱峭度图，选择所述快速谱峭度图中峭度最大的通频带，根据所述通频带构建带通滤波器对所述振动信号进行滤波，得到滤波后的第一信号x(t)；
构建滤波器F＝{f1,f2,...,fl,...,fL}，L表示所述滤波器的长度，fl表示所述滤波器的滤波器系数，则使用最小熵反卷积算法对所述第一信号x(t)进行滤波的过程表示为：
Figure FDA0002819812690000021
y(t)即为所述第一信号x(t)滤波后得到的第二信号；
其中，对所述第一信号x(t)进行滤波的迭代过程包括：随机初始化所述滤波器的滤波器系数得到滤波器F1＝{f11,f12,...,f1l,...,f1L}，再由式(1)得到y1(t)，将y1(t)的峭度作为目标函数则有
Figure FDA0002819812690000022
式(2)中OF1为y1(t)的峭度，N1表示的是信号y1(t)的长度，当
Figure FDA0002819812690000023
时，获取OF1的最大值并得到对应的滤波器的滤波器系数，即F2＝{f21,f22,...,f2l,...,f2L}，再根据式(1)、(2)、(3)反复进行迭代得到第二信号y(t)。
3.如权利要求2所述的一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，其特征在于，所述步骤S3具体包括：
S31：通过希尔伯特变换获取所述第二信号y(t)的解析信号，再根据所述解析信号获取所述第二信号y(t)的包络信号，即：sa(t)＝y(t)+iH[y(t)] (4)，sA(t)＝|sa(t)|(5)；式中H[·]表示希尔伯特变换，i为虚数单位，sa(t)表示所述解析信号，sA(t)表示所述包络信号；
S32：通过短时傅里叶变换获取所述包络信号sA(t)的第一时频谱，即：
Figure FDA0002819812690000031
其中，j表示虚数单位，g(t-τ)表示中心位于τ时刻的时间窗，STFTs(τ,f)表示所述包络信号sA(t)的第一时频谱；
S33：根据所述额定转速v和所述故障特征阶次o估计故障特征频率范围，获取所述STFTs(τ,f)在所述故障特征频率范围内所有τ时刻对应的频率幅值最大处的频率值，得到第一瞬时故障特征频率序列；
S34：使用加权最小二乘法和二次多项式局部回归模型对所述第一瞬时故障特征频率序列进行拟合，得到第一故障特征频率曲线；
在拟合的过程中，去除拟合残差大于第一预设阈值的数据点，重复拟合，直至拟合残差大于第一预设阈值的数据点的数量为0，得到最终的第一故障特征频率曲线；
S35：设置时变带通滤波器的带宽参数为ΔF、时变中心频率为所述第一故障特征频率曲线，在所述包络信号sA(t)的所述第一时频谱STFTs(τ,f)中提取所述时变带通滤波器范围内的时频谱，得到所述包络信号sA(t)的第二时频谱，对所述第二时频谱进行短时傅里叶逆变换，得到时变带通滤波后的第三信号sb(t)。
4.如权利要求3所述的一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，其特征在于，所述步骤S4具体包括：
S41：标记所述第三信号sb(t)的所有极值，获取所有极大值与其相邻极小值之间的两个高度差，将两个高度差中较小的高度差定义为所述极大值的显著性水平，使用Xn表示所述极值的序列，n∈[1,Z]，Z表示极值的总数，若Xk为Xn中的一个极大值，则Xk-1、Xk+1为Xk的相邻极小值，则所述Xk的显著性水平为P＝min(Xk-Xk-1,Xk-Xk+1) (7)，定义Xk的相对显著性水平为
Figure FDA0002819812690000041
其中，max(·)表示求最大值，k∈n；
S42：提取所述P'中大于第二预设阈值的极大值序列所对应的时间序列Tm，则Tm＝{t1,t2,...,tM}，M表示所述时间序列Tm的长度，有tq,tq+1∈Tm，则时刻
Figure FDA0002819812690000042
的瞬时故障特征频率为
Figure FDA0002819812690000043
根据式(9)计算所述时间序列Tm所有时间中点的瞬时故障特征频率，得到第二瞬时故障特征频率序列；其中，所述第二预设阈值∈(0.05,0.4)；
S43：使用加权最小二乘法和二次多项式局部回归模型对所述第二瞬时故障特征序列进行拟合，得到第二故障特征频率曲线；
在拟合过程中，去除拟合残差大于第三预设阈值的数据点，重复拟合，直至拟合残差大于第三预设阈值的数据点的数量为0，得到最终的第二故障特征频率曲线fc(t)；
S44：结合所述故障特征阶次o，得到转速曲线
Figure FDA0002819812690000051
5.如权利要求4所述的一种转速波动下滚动轴承包络阶次谱故障特征的提取方法，其特征在于，所述步骤S5具体包括：
S51：根据所述fr(t)得到相位-时间映射
Figure FDA0002819812690000052
Figure FDA0002819812690000053
则角度域采样的相位间隔
Figure FDA0002819812690000054
其中fs表示采样频率，tend表示
Figure FDA0002819812690000055
的信号时间长度；
S52：通过对相位-时间映射
Figure FDA0002819812690000056
进行三次样条插值得到所有等相位间隔序列
Figure FDA0002819812690000057
对应的时间序列{t1,t2,...,tend}；
S53：将时间序列{t1,t2,...,tend}作为插值点，对所述包络信号sA(t)进行三次样条插值，得到等角度采样的角度域包络信号sB(t)；
S54：通过离散傅里叶变换求解所述角度域包络信号sB(t)阶次谱；
S55：提取所述角度域包络信号sB(t)阶次谱中的故障特征。
CN202011427884.0A 2020-12-07 2020-12-07 一种转速波动下滚动轴承包络阶次谱故障特征的提取方法 Active CN112577746B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011427884.0A CN112577746B (zh) 2020-12-07 2020-12-07 一种转速波动下滚动轴承包络阶次谱故障特征的提取方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011427884.0A CN112577746B (zh) 2020-12-07 2020-12-07 一种转速波动下滚动轴承包络阶次谱故障特征的提取方法

Publications (2)

Publication Number Publication Date
CN112577746A true CN112577746A (zh) 2021-03-30
CN112577746B CN112577746B (zh) 2022-07-22

Family

ID=75127932

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011427884.0A Active CN112577746B (zh) 2020-12-07 2020-12-07 一种转速波动下滚动轴承包络阶次谱故障特征的提取方法

Country Status (1)

Country Link
CN (1) CN112577746B (zh)

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113191317A (zh) * 2021-05-21 2021-07-30 江西理工大学 一种基于极点构造低通滤波器的信号包络提取方法和装置
CN113326782A (zh) * 2021-06-01 2021-08-31 西安交通大学 基于包络谱形态拟合的滚动轴承故障特征自动提取方法
CN113567127A (zh) * 2021-07-23 2021-10-29 西安交通大学 一种基于时频特征分离的滚动轴承退化指标提取方法
CN114486252A (zh) * 2022-01-28 2022-05-13 西北工业大学 一种矢量模极大值包络的滚动轴承故障诊断方法
CN114778114A (zh) * 2022-04-01 2022-07-22 西南交通大学 一种基于信号冲击性和周期性的轴承健康指标构建方法
CN116380467A (zh) * 2023-05-24 2023-07-04 成都工业职业技术学院 一种基于多时频脊线提取的滚动轴承故障诊断方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106769033A (zh) * 2016-11-30 2017-05-31 西安交通大学 基于阶次包络时频能量谱的变转速滚动轴承故障识别方法
CN110617964A (zh) * 2019-07-29 2019-12-27 中国铁道科学研究院集团有限公司城市轨道交通中心 用于滚动轴承故障诊断的同步压缩变换阶比分析法
CN111307460A (zh) * 2020-03-14 2020-06-19 中国石化销售股份有限公司华南分公司 基于计算阶次跟踪与谱峭度的滚动轴承故障诊断方法
CN111665051A (zh) * 2020-07-01 2020-09-15 天津大学 基于能量权重法的强噪声变转速条件下轴承故障诊断方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106769033A (zh) * 2016-11-30 2017-05-31 西安交通大学 基于阶次包络时频能量谱的变转速滚动轴承故障识别方法
CN110617964A (zh) * 2019-07-29 2019-12-27 中国铁道科学研究院集团有限公司城市轨道交通中心 用于滚动轴承故障诊断的同步压缩变换阶比分析法
CN111307460A (zh) * 2020-03-14 2020-06-19 中国石化销售股份有限公司华南分公司 基于计算阶次跟踪与谱峭度的滚动轴承故障诊断方法
CN111665051A (zh) * 2020-07-01 2020-09-15 天津大学 基于能量权重法的强噪声变转速条件下轴承故障诊断方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
林京等: "变转速下机械设备动态信号分析方法的回顾与展望", 《中国科学:技术科学》 *

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113191317A (zh) * 2021-05-21 2021-07-30 江西理工大学 一种基于极点构造低通滤波器的信号包络提取方法和装置
CN113191317B (zh) * 2021-05-21 2022-09-27 江西理工大学 一种基于极点构造低通滤波器的信号包络提取方法和装置
CN113326782A (zh) * 2021-06-01 2021-08-31 西安交通大学 基于包络谱形态拟合的滚动轴承故障特征自动提取方法
CN113567127A (zh) * 2021-07-23 2021-10-29 西安交通大学 一种基于时频特征分离的滚动轴承退化指标提取方法
CN114486252A (zh) * 2022-01-28 2022-05-13 西北工业大学 一种矢量模极大值包络的滚动轴承故障诊断方法
CN114486252B (zh) * 2022-01-28 2023-06-30 西北工业大学 一种矢量模极大值包络的滚动轴承故障诊断方法
CN114778114A (zh) * 2022-04-01 2022-07-22 西南交通大学 一种基于信号冲击性和周期性的轴承健康指标构建方法
CN114778114B (zh) * 2022-04-01 2022-11-22 西南交通大学 一种基于信号冲击性和周期性的轴承健康指标构建方法
CN116380467A (zh) * 2023-05-24 2023-07-04 成都工业职业技术学院 一种基于多时频脊线提取的滚动轴承故障诊断方法
CN116380467B (zh) * 2023-05-24 2024-01-23 成都工业职业技术学院 一种基于多时频脊线提取的滚动轴承故障诊断方法

Also Published As

Publication number Publication date
CN112577746B (zh) 2022-07-22

Similar Documents

Publication Publication Date Title
CN112577746B (zh) 一种转速波动下滚动轴承包络阶次谱故障特征的提取方法
US7133801B2 (en) System and methodology for vibration analysis and condition monitoring
Wang et al. Rolling element bearing fault diagnosis via fault characteristic order (FCO) analysis
US7031873B2 (en) Virtual RPM sensor
CN110617964A (zh) 用于滚动轴承故障诊断的同步压缩变换阶比分析法
US20070032966A1 (en) System and methodology for vibration analysis and conditon monitoring
CN112507769A (zh) 一种基于仿真传感器谐振增强特征的轴承故障诊断方法
CN110763462B (zh) 一种基于同步压缩算子的时变振动信号故障诊断方法
CN111665051A (zh) 基于能量权重法的强噪声变转速条件下轴承故障诊断方法
CN112105907B (zh) 用于监测齿轮系统的方法和设备
KR101966270B1 (ko) 전자 진동 성분의 제거 방법, 회전 기계 진단 방법, 및 회전 기계 진단 장치
Ding et al. Sparsity-based algorithm for condition assessment of rotating machinery using internal encoder data
JP2012042338A (ja) 転がり軸受の異常診断装置および歯車の異常診断装置
CN110987438A (zh) 水轮发电机变转速过程周期性振动冲击信号检测的方法
CN110108467A (zh) 基于便携式移动设备的主动发声测速方法
CN108398260B (zh) 基于混合概率方法的齿轮箱瞬时角速度的快速评估方法
Wu et al. A modified tacho-less order tracking method for the surveillance and diagnosis of machine under sharp speed variation
CN117836599A (zh) 检测旋转系统中轴承缺陷的方法以及实施该方法的监控系统
JP6430234B2 (ja) 回転機械の振動解析装置およびプログラム
CN112465068A (zh) 一种基于多传感器数据融合的旋转设备故障特征提取方法
DK2630453T3 (en) A method for monitoring a rotating member belonging to a mechanical transmission of a wind turbine
JP6550476B2 (ja) 信号を分析するための方法およびその方法を実行するための装置
Jiang et al. A wavelet cluster-based band-pass filtering and envelope demodulation approach with application to fault diagnosis in a dry vacuum pump
CN109948487B (zh) 基于时频谱相关性分析的旋转机械故障特征提取方法
CN114486252B (zh) 一种矢量模极大值包络的滚动轴承故障诊断方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
EE01 Entry into force of recordation of patent licensing contract
EE01 Entry into force of recordation of patent licensing contract

Application publication date: 20210330

Assignee: NANJING HIGH ACCURATE MARINE EQUIPMENT Co.,Ltd.

Assignor: SOUTHEAST University

Contract record no.: X2023990000339

Denomination of invention: A Method for Extracting Envelope Order Spectral Fault Characteristics of Rolling Bearing under Rotational Speed Fluctuation

Granted publication date: 20220722

License type: Exclusive License

Record date: 20230322