CN114757230A - 一种用于轴承故障信号的时频分析方法及系统 - Google Patents

一种用于轴承故障信号的时频分析方法及系统 Download PDF

Info

Publication number
CN114757230A
CN114757230A CN202210395577.1A CN202210395577A CN114757230A CN 114757230 A CN114757230 A CN 114757230A CN 202210395577 A CN202210395577 A CN 202210395577A CN 114757230 A CN114757230 A CN 114757230A
Authority
CN
China
Prior art keywords
time
frequency
group delay
initial
delay estimation
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.)
Pending
Application number
CN202210395577.1A
Other languages
English (en)
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.)
University of Jinan
Original Assignee
University of Jinan
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 University of Jinan filed Critical University of Jinan
Priority to CN202210395577.1A priority Critical patent/CN114757230A/zh
Publication of CN114757230A publication Critical patent/CN114757230A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • 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
    • G06F2218/06—Denoising by applying a scale-space analysis, e.g. using wavelet analysis
    • 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
    • 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

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Signal Processing (AREA)
  • General Engineering & Computer Science (AREA)
  • Complex Calculations (AREA)

Abstract

本公开提供了一种用于轴承故障信号的时频分析方法及系统，所述方案包括：获取待分析的轴承强频变调制信号；基于小波变换对所述强频变调制信号进行处理，获得初始时频分布结果；基于获得的初始时频分布结果，获得信号的初始群延迟估计算子；基于所述初始群延迟估计算子及时间重排同步压缩变换，对获得初始时频分布结果中的时频能量沿时间方向进行重排，获得初始重分配结果；基于不动点迭代法获得迭代群延迟估计算子，并基于所述迭代群延迟估计算子对获得的初始重分配结果中的时频能量沿时间方向进行二次优化重排，获得优化后的时频分布结果；基于所述优化后的时频分布结果实现对所述强频变调制信号的分析。

Description

一种用于轴承故障信号的时频分析方法及系统
技术领域
本公开属于非平稳旋转机械故障信号时频分析技术领域，尤其涉及一种用于轴承故障信号的时频分析方法及系统。
背景技术
本部分的陈述仅仅是提供了与本公开相关的背景技术信息，不必然构成在先技术。
旋转机械故障诊断是工程研究的一个重要领域，其目的是确保安全、减少故障、降低维护成本。旋转机械的故障通常由滚动轴承的损坏引起。轴承故障信号的动态变化规律为旋转机械故障检测提供了基础。基于轴承振动信号的状态监测是一种有效的非侵入式技术，在工业中得到了广泛应用。轴承振动信号是一种工业中易于获取的信号，但是其通常运行在非平稳运行的情况下，且具有明显的瞬态特征。此外，检测到的信号还包含来自其他信号源和背景噪声的污染，使得对故障特征进行准确描述变得十分困难。信号处理技术的进一步发展是提高故障诊断准确性的最重要途径之一。
通过传统的小波变换和短时傅里叶变换通过将观察窗口内信号的瞬时频率或群延迟近似为准线性，以实现对非平稳故障信号的描述，但是，受海森堡不确定性原则的限制，其时频表示会在时频平面上形成沿着脊线的扩散，这将严重影响信号的时频分辨率。
针对群延迟恒为定值的弱频变故障信号，基于短时傅里叶的时间重排同步压缩变换可以实现对其群延迟的准确估计，并利用沿时间方向的同步压缩算子将时频能量聚集到真实群延迟脊线上，但是，短时傅里叶变换的观察窗一经选定便无法改变，缺乏自适应的能力，一旦窗口选择不合适，其压缩能力将大大削弱。现有改进小波变换拥有了在不同频率范围自适应调整合适的观察窗的能力，结合时间重排同步压缩变换能更有效地重排时频系数，但其无法应对在短时间内群延迟变化规律更复杂且特征微弱的轴承故障信号。发明人发现，在弱频变调制信号模型框架下的改进小波的时间重排同步压缩变换，无法准确地将强频变调制信号的时频系数挤压到真实群延迟的位置。利用改进小波的时间重排同步压缩变换处理强频变调制信号的时频分布结果与理想状态下强频变调制信号的时频分布结果差距过大。
发明内容
本公开为了解决上述问题，提供了一种用于轴承故障信号的时频分析方法及系统，所述方案利用不动点迭代的思想，减小真实的群延迟与估计的群延迟算子之间误差项的影响，利用最终的迭代群延迟算子替换改进小波的时间重排同步压缩变换中的原始群延迟估计算子进行再次压缩，使获得的时频能量分布结果更逼近强频变调制信号的理想时频分布，在旋转机械嘈杂的工作环境中，该方法可以更准确地提取到轴承的故障特征，有效提高信号时频分析的准确性。
根据本公开实施例的第一个方面，提供了一种用于轴承故障信号的时频分析方法，包括：
获取待分析的轴承强频变调制信号；
基于小波变换对所述强频变调制信号进行处理，获得初始时频分布结果；
基于获得的初始时频分布结果，获得信号的初始群延迟估计算子；基于所述初始群延迟估计算子及时间重排同步压缩变换，对获得初始时频分布结果中的时频能量沿时间方向进行重排，获得初始重分配结果；
基于不动点迭代法获得迭代群延迟估计算子，并基于所述迭代群延迟估计算子对获得的初始重分配结果中的时频能量沿时间方向进行二次优化重排，获得优化后的时频分布结果；
基于所述优化后的时频分布结果实现对所述强频变调制信号的分析。
进一步的，所述小波变换采用改进小波变换公式，具体表示如下：
Figure RE-RE-GDA0003690429240000031
其中，所述改进小波变换的小波基采用解析小波，小波函数定义为ψ(t)＝g((t-b)/a)·exp(itω0/a)/a，其中a为局部尺度变量，b为局部时间变量，t为全局时间变量，ω0为小波函数初始频点，g()为高斯窗函数，*代表取复共轭。
进一步的，所述基于获得的初始时频分布结果，获得信号的初始群延迟估计算子，具体为：利用改进小波变换与t·g(t)的改进小波变换间的比值关系，实现群延迟估计算子的构造，所述初始群延迟估计算子具体表示如下：
Figure RE-RE-GDA0003690429240000032
其中，
Figure RE-RE-GDA0003690429240000033
表示窗口为t·g(t)的改进小波变换结果，g()为高斯窗函数。
进一步的，所述基于不动点迭代法获得迭代群延迟估计算子，具体为：
在所述初始群延迟估计算子中引入群延迟的二阶泰勒展开的信号模型，获得高阶信号的具体群延迟估计算子，具体表示如下：
Figure RE-RE-GDA0003690429240000034
其中，
Figure RE-RE-GDA0003690429240000035
为信号的实际群延迟，
Figure RE-RE-GDA0003690429240000036
为实际群延迟的变化率，σ为窗函数的支撑因子，ω为频率因子，a为局部尺度变量，b为局部时间变量；
将群延迟估计算子
Figure RE-RE-GDA0003690429240000037
中局部时间变量b不断用
Figure RE-RE-GDA0003690429240000038
替代，获得迭代群延迟估计算子。
进一步的，所述优化后的时频分布结果与改进小波变换具有相同的重构特性，经时间多次重排同步压缩变换处理的强频变调制信号通过下式进行恢复：
Figure RE-RE-GDA0003690429240000039
其中，
Figure RE-RE-GDA00036904292400000310
为窗函数g(t)的频域表示
Figure RE-RE-GDA00036904292400000311
在ω＝0处的值，a为局部尺度变量，u为局部时间变量，ω0为小波函数初始频点，i为虚数单位，t为全局时间变量。
根据本公开实施例的第二个方面，提供了一种用于轴承故障信号的时频分析系统，包括：
数据获取单元，其用于获取待分析的轴承强频变调制信号；
时频结果获取单元，其用于基于小波变换对所述强频变调制信号进行处理，获得初始时频分布结果；基于获得的初始时频分布结果，获得信号的初始群延迟估计算子；基于所述初始群延迟估计算子及时间重排同步压缩变换，对获得初始时频分布结果中的时频能量沿时间方向进行重排，获得初始重分配结果；
时频结果优化单元，其用于基于不动点迭代法获得迭代群延迟估计算子，并基于所述迭代群延迟估计算子对获得的初始重分配结果中的时频能量沿时间方向进行二次优化重排，获得优化后的时频分布结果；
分析单元，其用于基于所述优化后的时频分布结果实现对所述强频变调制信号的分析。
根据本公开实施例的第三方面，提供了一种计算机可读存储介质，其上存储有程序，该程序被处理器执行时实现上所述的一种用于轴承故障信号的时频分析方法。
根据本公开实施例的第四方面，提供了一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的程序，所述处理器执行所述程序时实现如上所述的一种用于轴承故障信号的时频分析方法。
与现有技术相比，本公开的有益效果是：
本公开所述方案提出一种用于轴承故障信号的时频分析方法及系统，所述方案针对改进小波的时间重排同步压缩变换中估计的强频变调制信号的群延迟误差项，利用不动点迭代的思想，减小真实的群延迟与估计的群延迟算子之间误差项的影响，利用最终的迭代群延迟算子替换改进小波的时间重排同步压缩变换中的原始群延迟估计算子进行再次压缩，使获得的时频能量分布结果更逼近强频变调制信号的理想时频分布，特别是在旋转机械嘈杂的工作环境中，该方法可以更准确地提取到轴承的故障特征，有效提高信号时频分析的准确性。
本公开附加方面的优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本公开的实践了解到。
附图说明
构成本公开的一部分的说明书附图用来提供对本公开的进一步理解，本公开的示意性实施例及其说明用于解释本公开，并不构成对本公开的不当限定。
图1为本公开实施例中所述的基于时间多次重排压缩变换的信号时频分析方法与相关方法的演变过程示意图；
图2(a)和图2(b)分别为本公开实施例中所述的强频变调制信号时域波形和理想时频分布示意图；
图3(a)和图3(b)分别为本公开实施例中所述的基于改进小波变换处理结果处理强频变调制信号时频分布结果以及基于改进小波的时间重排同步压缩变换处理结果示意图；
图4为本公开实施例中采用基于时间多次重排压缩变换的信号时频分析方法处理强频变调制信号的结果示意图。
图5(a)至图5(c)为本公开实施例中所述的基于改进小波变换方法处理强频变调制信号时频分布及不同频率处沿时间方向的时频切片，其中，图5(a) 为f＝44Hz处频率切片，图5(b)为f＝68Hz处频率切片，图5(c)为f＝82Hz 处频率切片；
图6(a)至图6(c)为本公开实施例中所述的基于改进小波的时间重排同步压缩变换方法处理强频变调制信号时频分布及不同频率处沿时间方向的时频切片(迭代次数N＝1)，其中，图6(a)为f＝44Hz处频率切片，图6(b)为f＝68Hz 处频率切片，图6(c)为f＝82Hz处频率切片。
图7为本公开实施例中所述的基于时间多次重排压缩变换的信号时频分析方法处理强频变调制信号时频分布及不同频率处沿时间方向的时频切片(迭代次数N＝10)，其中，图7(a)为f＝44Hz处频率切片，图7(b)为f＝68Hz处频率切片图7(c)为f＝82Hz处频率切片；
图8为本公开实施例中所述的三种时频分析方法得到的时频分布结果的 Rényi熵对比结果示意图；
图9(a)和图9(b)为本公开实施例中所述的基于时间多次重排压缩变换的信号时频分析方法的重构性能结果示意图，其中，图9(a)为重构结果比较，图9(b)为重构误差；
图10为本公开实施例中所述的基于时间多次重排压缩变换的信号时频分析方法流程图。
具体实施方式
下面结合附图与实施例对本公开做进一步说明。
应该指出，以下详细说明都是示例性的，旨在对本公开提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本公开所属技术领域的普通技术人员通常理解的相同含义。
需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本公开的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，当在本说明书中使用术语“包含”和/或“包括”时，其指明存在特征、步骤、操作、器件、组件和/或它们的组合。
在不冲突的情况下，本公开中的实施例及实施例中的特征可以相互组合。
实施例一：
本实施例的目的是提供一种用于轴承故障信号的时频分析方法。
首先，基于Dirac信号作为典型的弱频变调制信号，具有恒为定值的群延迟和理想的时间定位能力，下面以Dirac信号为例，对改进小波的时间重排同步压缩变换的工作原理进行说明，Dirac信号的时域表达式为：
xδ(t)＝A·δ(t-t0) (1)
该信号仅在t＝t0处幅值为A，其余时刻的幅值为0。对于改进的小波变换其表达式为
Figure RE-RE-GDA0003690429240000071
其中，小波基为解析小波，即小波函数可以定义为一个被特定频率ω0调制的窗函数：ψ(t)＝g((t-b)/a)·exp(itω0/a)/a，ψ*为小波函数并取复共轭，i为虚数单位，b为局部时间变量，ω0为小波函数初始频点。
对Dirac信号进行改进小波变换，其表达式为：
Figure RE-RE-GDA0003690429240000072
其中，g(t)＝exp(-t2/(2·σ))为高斯窗函数，频率和尺度存在反比关系ω＝ω0/a， A为信号幅度，δ()为狄拉克函数。对式(3)取模可以得到Dirac信号的改进小波变换的时频分布：
Figure RE-RE-GDA0003690429240000073
可以看出窗函数的紧支撑域的存在使能量无法集中在信号的群延迟脊线上，无法体现出Dirac信号理想的时间定位能力。
改进小波的时间重排同步压缩变换利用改进小波变换与t·g的改进小波变换间的比值关系构造群延迟估计算子。由式(3)可推得，对t·g的改进小波变换的时频结果如下：
Figure RE-RE-GDA0003690429240000074
由此可知，当
Figure RE-RE-GDA0003690429240000075
时，信号xδ(t)的群延迟
Figure RE-RE-GDA0003690429240000076
在时频空间中可以表示为：
Figure RE-RE-GDA0003690429240000081
在此基础上，改进小波的时间重排同步压缩变换通过构造沿时间方向的同步压缩算子，将扩散的时频能量压缩至真实群延迟脊线处，其表达式为：
Figure RE-RE-GDA0003690429240000082
改进小波的时间重排同步压缩变换通过对时频平面上改进小波变换的系数进行重排，从而实现了能量聚集的目的。
改进小波的时间重排同步压缩变换处理后的时频分布结果仍具有重构特性，其处理后的时频分布结果沿时间方向积分，可以得到如下结果：
Figure RE-RE-GDA0003690429240000083
由式(8)可知，改进小波的时间重排同步压缩变换拥有与改进小波变换相同的重构特性。通过对改进小波变换后的结果沿时间方法进行积分，得到如下结果：
Figure RE-RE-GDA0003690429240000084
其中：
Figure RE-RE-GDA0003690429240000085
为窗函数g(t)的频域表示
Figure RE-RE-GDA0003690429240000086
在ω＝0处的值，
Figure RE-RE-GDA0003690429240000087
为弱频变信号xδ(t)的频域表示。
根据式(9)可知，弱频变信号xδ(t)的频域表示
Figure RE-RE-GDA0003690429240000088
可通过下式表示：
Figure RE-RE-GDA0003690429240000089
当已知频域表示
Figure RE-RE-GDA00036904292400000810
利用下式可恢复其对应的时域信号：
Figure RE-RE-GDA0003690429240000091
综合式(8)～式(11)可知，当时间重排同步压缩变换的时频分布结果已知时，可以通过下式获得其对应的时域信号：
Figure RE-RE-GDA0003690429240000092
然而，在强频变调制信号模型框架下的改进小波的时间重排同步压缩变换，无法准确地将强频变调制信号的时频系数挤压到真实群延迟的位置。利用改进小波的时间重排同步压缩变换处理强频变调制信号的时频分布结果与理想状态下强频变调制信号的时频分布结果差距过大。
尽管利用改进小波的时间重排同步压缩变换处理弱频变调制信号能够获得强于其他方法时频结果的能量聚集性。但是在面对强频变调制信号时改进小波的时间重排同步压缩变换仍旧会存在严重的能量扩散。为此，提出了改进小波的时间多次重排同步压缩变换，依据不动点迭代的思想，对改进小波的时间重排同步压缩变换中的群延迟估计算子进行反复迭代精炼，进一步提高时间重排同步压缩变换对强频变调制信号时频结果的能量聚集性。
利用改进小波的时间重排同步压缩变换处理强频变调制信号时，选择高斯窗构造小波基函数可以给出群延迟为二阶泰勒展开下的信号模型的群延迟估计算子。但群延迟估计算子与真实群延迟之间存在过大的差距。误差项可以看作是一个小于1的加权系数与群延迟之间的误差乘积。为此，针对群延迟估计算子中的误差项，依据不动点迭代的思想，消减误差项的影响使群延迟估计算子逼近真实群延迟，为时间方向的重排操作确定一个更准确的位置指向。
具体的，基于上述问题，如图10所示，本公开提供了一种用于轴承故障信号的时频分析方法，包括：
获取待分析的轴承强频变调制信号；
基于小波变换对所述强频变调制信号进行处理，获得初始时频分布结果；
基于获得的初始时频分布结果，获得信号的初始群延迟估计算子；基于所述初始群延迟估计算子及时间重排同步压缩变换，对获得初始时频分布结果中的时频能量沿时间方向进行重排，获得初始重分配结果；
基于不动点迭代法获得迭代群延迟估计算子，并基于所述迭代群延迟估计算子对获得的初始重分配结果中的时频能量沿时间方向进行二次优化重排，获得优化后的时频分布结果；
基于所述优化后的时频分布结果实现对所述强频变调制信号的分析。
进一步的，所述小波变换采用改进小波变换公式，具体表示如下：
Figure RE-RE-GDA0003690429240000101
其中，所述改进小波变换的小波基采用解析小波，小波函数定义为ψ(t)＝g((t-b)/a)·exp(itω0/a)/a，其中a为局部尺度变量，b为局部时间变量，t为全局时间变量，ω0为小波函数初始频点，g()为高斯窗函数，*代表取复共轭。
进一步的，所述基于获得的初始时频分布结果，获得信号的初始群延迟估计算子，具体为：利用改进小波变换与t·g(t)的改进小波变换间的比值关系，实现群延迟估计算子的构造，所述初始群延迟估计算子具体表示如下：
Figure RE-RE-GDA0003690429240000102
其中，
Figure RE-RE-GDA0003690429240000103
表示窗口为t·g(t)的改进小波变换结果，g()为高斯窗函数。
进一步的，所述基于不动点迭代法获得迭代群延迟估计算子，具体为：
在所述初始群延迟估计算子中引入群延迟的二阶泰勒展开的信号模型，获得高阶信号的具体群延迟估计算子，具体表示如下：
Figure RE-RE-GDA0003690429240000104
其中，
Figure RE-RE-GDA0003690429240000105
为信号的实际群延迟，
Figure RE-RE-GDA0003690429240000106
为实际群延迟的变化率，σ为窗函数的支撑因子，ω为频率因子，a为局部尺度变量，b为局部时间变量；
将群延迟估计算子
Figure RE-RE-GDA0003690429240000111
中局部时间变量b不断用
Figure RE-RE-GDA0003690429240000112
替代，获得迭代群延迟估计算子。
具体的，利用所述初始群延迟估计算子，可以得到改进小波变换的时频系数的初始重排结果，具体表达式为：
Figure RE-RE-GDA0003690429240000113
其中，δ()为狄拉克函数，u为局部时间变量。
通过引入群延迟的二阶泰勒展开的信号模型：
Figure RE-RE-GDA0003690429240000114
其中
Figure RE-RE-GDA0003690429240000115
为信号的频域形式，A(ω)为信号幅度，
Figure RE-RE-GDA0003690429240000116
为信号相位，
Figure RE-RE-GDA0003690429240000117
为信号群延迟，ξ为全局频率变量，ω为局部频率变量。
将公式(17)代入(14)，可以得到高阶信号的具体群延迟估计算子，即公式(15)；
进一步的，通过将群延迟估计算子
Figure RE-RE-GDA0003690429240000118
中局部时间变量b不断用
Figure RE-RE-GDA0003690429240000119
替代，可以得到：
Figure RE-RE-GDA00036904292400001110
Figure RE-RE-GDA00036904292400001111
Figure RE-RE-GDA00036904292400001112
Figure RE-RE-GDA00036904292400001113
其中，
Figure RE-RE-GDA00036904292400001114
为时间重排同步压缩变换中的
Figure RE-RE-GDA00036904292400001115
其中N为迭代次数，由于不同信噪比下改进小波的多次时间重排同步压缩变换拥有几乎相同的收敛速度，并且为了减小计算负担，N 通常取值在10以内。可以看出，由于上式中加号右边的分式项是小于1且在逐渐减小的，因此增加迭代次数可以使群延迟估计算子
Figure RE-RE-GDA0003690429240000121
逼近信号真实的群延迟
Figure RE-RE-GDA0003690429240000122
当对改进小波变换的时间重排同步压缩变换的结果进行N＝2的二次优化重排时可以得到：
Figure RE-RE-GDA0003690429240000123
其中，
Figure RE-RE-GDA0003690429240000124
δ()为狄拉克函数。τ,u为局部时间变量。由上式的结果可以看出，如果对上式迭代N次，则有
Figure RE-RE-GDA0003690429240000125
进一步的，所述优化后的时频分布结果与改进小波变换具有相同的重构特性，经时间多次重排同步压缩变换处理的强频变调制信号通过下式进行恢复：
Figure RE-RE-GDA0003690429240000126
其中，
Figure RE-RE-GDA0003690429240000127
为窗函数g(t)的频域表示
Figure RE-RE-GDA0003690429240000128
在ω＝0处的值，a为局部尺度变量，u为局部时间变量，ω0为小波函数初始频点，i为虚数单位，t为全局时间变量。
进一步的，为了证明本公开所述方案的可行性，本实施例中进行了实验证明：
首先，构造的强频变调制信号频域表示如下：
Figure RE-RE-GDA0003690429240000129
该信号的采样频率为200Hz，采样时间为10s。该信号对应的群延迟函数为：
GD2(ω)＝0.1π(sin(0.05ω)+0.05ωcos(0.05ω))+5 (23)
该强频变调制信号的群延迟在t＝5s附近快速变化，其对应的时域波形和理想状态下的时频分布结果如图2(a)至图2(b)所示。该信号的时域波形近似描述为t＝5s附近的一簇类脉冲信号。理想时频分布结果中，该强频变调制信号的群延迟快速变化，群延迟位置处的时频能量幅值为一指数函数。
利用改进小波变换和改进小波的时间重排同步压缩变换处理该强频变调制信号得到的时频表示如图3(a)和图3(b)所示。可以看出，改进小波变换处理后的时频表示仍然存在较大的能量扩散现象；利用改进小波的时间重排同步压缩变换处理后的时频结果表明在信号低频处，局部群延迟近似为线性时能较好地聚集信号的能量。但在高频处该信号的群延迟在短时间内变化过快，改进小波的时间重排同步压缩变换处理后的时频结果仍然与理想状态下时频分布结果存在较大差距。为了验证所提时频分析方法的有效性，利用所提方法处理式 (17)所示强频变调制信号。图4给出了所提方法的时频结果，可以看出所提方法拥有与理想时频分布一致的时频脊线。图5(a)至图5(c)，图6(a)至图 6(c)和图7(a)至图7(c)分别给出了改进小波变换，改进小波的时间重排同步压缩变换以及改进小波的时间多次重排同步压缩变换获得的时频结果在 f＝7Hz、f＝75Hz和f＝85Hz处沿时间方向的时频切片。可以看出，由于强频变调制信号的群延迟在局部变换率较大，改进小波的时间重排同步压缩变换对原始改进小波变换的能量分布的聚集效果不理想，但是所提方法能够将强频变调制信号的能量分布接近理想地汇聚至群延迟脊线位置处的一根直线。
Rényi熵可用于定量评价利用不同时频分析方法获得的时频分布的能量聚集性。时频能量聚集性越高，对应的Rényi熵值越小。图8给出了利用不同时频分析方法处理该强频变调制信号的时频分布结果对应的Rényi熵，可以看出改进小波的时间多次重排同步压缩变换处理该强频变调制信号所得的时频结果具有更集中的能量。随着信噪比的减小，改进小波变换和改进小波的时间重排同步压缩变换处理该强频变调制信号所得的时频分布结果对应的Rényi熵值均在增大，这意味着噪声的增加会严重破坏这两种方法的时频结果的能量聚集性。然而，从整体上看，所提方法对应的Rényi熵相对于改进小波的时间重排同步压缩变换波动较小且始终拥有最小的熵值。由此证明，所提方法相对于改进小波的时间重排同步压缩变换能进一步提高时频结果的能量聚集性。
重构特性将允许所提方法拥有更广泛的应用途径，图9(a)至图9(b)给出了所提方法恢复的信号与原始信号之间的比较，可以看出所提方法可以做到以较小的误差重构出原始信号。由此证明，所提方法在提高时频结果的能量聚集性的同时允许恢复信号的幅值和相位信息。
实施例二：
本实施例的目的是提供一种用于轴承故障信号的时频分析系统。
一种用于轴承故障信号的时频分析系统，包括：
数据获取单元，其用于获取待分析的轴承强频变调制信号；
时频结果获取单元，其用于基于小波变换对所述强频变调制信号进行处理，获得初始时频分布结果；基于获得的初始时频分布结果，获得信号的初始群延迟估计算子；基于所述初始群延迟估计算子及时间重排同步压缩变换，对获得初始时频分布结果中的时频能量沿时间方向进行重排，获得初始重分配结果；
时频结果优化单元，其用于基于不动点迭代法获得迭代群延迟估计算子，并基于所述迭代群延迟估计算子对获得的初始重分配结果中的时频能量沿时间方向进行二次优化重排，获得优化后的时频分布结果；
分析单元，其用于基于所述优化后的时频分布结果实现对所述强频变调制信号的分析。
在更多实施例中，还提供：
一种电子设备，包括存储器和处理器以及存储在存储器上并在处理器上运行的计算机指令，所述计算机指令被处理器运行时，完成实施例一中所述的方法。为了简洁，在此不再赘述。
应理解，本实施例中，处理器可以是中央处理单元CPU，处理器还可以是其他通用处理器、数字信号处理器DSP、专用集成电路ASIC，现成可编程门阵列FPGA或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等。
存储器可以包括只读存储器和随机存取存储器，并向处理器提供指令和数据、存储器的一部分还可以包括非易失性随机存储器。例如，存储器还可以存储设备类型的信息。
一种计算机可读存储介质，用于存储计算机指令，所述计算机指令被处理器执行时，完成实施例一中所述的方法。
实施例一中的方法可以直接体现为硬件处理器执行完成，或者用处理器中的硬件及软件模块组合执行完成。软件模块可以位于随机存储器、闪存、只读存储器、可编程只读存储器或者电可擦写可编程存储器、寄存器等本领域成熟的存储介质中。该存储介质位于存储器，处理器读取存储器中的信息，结合其硬件完成上述方法的步骤。为避免重复，这里不再详细描述。
本领域普通技术人员可以意识到，结合本实施例描述的各示例的单元即算法步骤，能够以电子硬件或者计算机软件和电子硬件的结合来实现。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本公开的范围。
上述实施例提供的一种用于轴承故障信号的时频分析方法及系统可以实现，具有广阔的应用前景。
以上所述仅为本公开的优选实施例而已，并不用于限制本公开，对于本领域的技术人员来说，本公开可以有各种更改和变化。凡在本公开的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本公开的保护范围之内。

Claims (10)

1.一种用于轴承故障信号的时频分析方法，其特征在于，包括：
获取待分析的轴承强频变调制信号；
基于小波变换对所述强频变调制信号进行处理，获得初始时频分布结果；
基于获得的初始时频分布结果，获得信号的初始群延迟估计算子；基于所述初始群延迟估计算子及时间重排同步压缩变换，对获得初始时频分布结果中的时频能量沿时间方向进行重排，获得初始重分配结果；
基于不动点迭代法获得迭代群延迟估计算子，并基于所述迭代群延迟估计算子对获得的初始重分配结果中的时频能量沿时间方向进行二次优化重排，获得优化后的时频分布结果；
基于所述优化后的时频分布结果实现对所述强频变调制信号的分析。
2.如权利要求1所述的一种用于轴承故障信号的时频分析方法，其特征在于，所述小波变换采用改进小波变换公式，具体表示如下：
Figure RE-FDA0003690429230000011
其中，所述改进小波变换的小波基采用解析小波，小波函数定义为ψ(t)＝g((t-b)/a)·exp(itω0/a)/a，其中a为局部尺度变量，b为局部时间变量，t为全局时间变量，ω0为小波函数初始频点，g()为高斯窗函数，*表示取复共轭。
3.如权利要求1所述的一种用于轴承故障信号的时频分析方法，其特征在于，所述基于获得的初始时频分布结果，获得信号的初始群延迟估计算子，具体为：利用改进小波变换与t·g(t)的改进小波变换间的比值关系，实现群延迟估计算子的构造，所述初始群延迟估计算子具体表示如下：
Figure RE-FDA0003690429230000012
其中，
Figure RE-FDA0003690429230000021
表示窗口为t·g(t)的改进小波变换结果，g()为高斯窗函数。
4.如权利要求1所述的一种用于轴承故障信号的时频分析方法，其特征在于，所述基于不动点迭代法获得迭代群延迟估计算子，具体为：
在所述初始群延迟估计算子中引入群延迟的二阶泰勒展开的信号模型，获得高阶信号的具体群延迟估计算子，具体表示如下：
Figure RE-FDA0003690429230000022
其中，
Figure RE-FDA0003690429230000023
为信号的实际群延迟，
Figure RE-FDA0003690429230000024
为实际群延迟的变化率，σ为窗函数的支撑因子，ω为频率因子，a为局部尺度变量，b为局部时间变量；
将群延迟估计算子
Figure RE-FDA0003690429230000025
中局部时间变量b不断用
Figure RE-FDA0003690429230000026
替代，获得迭代群延迟估计算子。
5.如权利要求4所述的一种用于轴承故障信号的时频分析方法，其特征在于，基于所述迭代群延迟估计算子，获得优化重排后的时频分布结果，具体表示：
Figure RE-FDA0003690429230000027
其中，W(a,b)为改进小波变换，δ()为狄拉克函数，u为局部时间变量，
Figure RE-FDA0003690429230000028
为N次迭代群延迟估计算子，a为局部尺度变量，b为局部时间变量。
6.如权利要求1所述的一种用于轴承故障信号的时频分析方法，其特征在于，所述优化后的时频分布结果与改进小波变换具有相同的重构特性，经时间多次重排同步压缩变换处理的强频变调制信号通过下式进行恢复：
Figure RE-FDA0003690429230000029
其中，
Figure RE-FDA00036904292300000210
为窗函数g(t)的频域表示
Figure RE-FDA00036904292300000211
在ω＝0处的值，a为局部尺度变量，u为局部时间变量，ω0为小波函数初始频点，i为虚数单位，t为全局时间变量。
7.一种用于轴承故障信号的时频分析系统，其特征在于，包括：
数据获取单元，其用于获取待分析的轴承强频变调制信号；
时频结果获取单元，其用于基于小波变换对所述强频变调制信号进行处理，获得初始时频分布结果；基于获得的初始时频分布结果，获得信号的初始群延迟估计算子；基于所述初始群延迟估计算子及时间重排同步压缩变换，对获得初始时频分布结果中的时频能量沿时间方向进行重排，获得初始重分配结果；
时频结果优化单元，其用于基于不动点迭代法获得迭代群延迟估计算子，并基于所述迭代群延迟估计算子对获得的初始重分配结果中的时频能量沿时间方向进行二次优化重排，获得优化后的时频分布结果；
分析单元，其用于基于所述优化后的时频分布结果实现对所述强频变调制信号的分析。
8.如权利要求7所述的一种用于轴承故障信号的时频分析系统，其特征在于，所述基于不动点迭代法获得迭代群延迟估计算子，具体为：
在所述初始群延迟估计算子中引入群延迟的二阶泰勒展开的信号模型，获得高阶信号的具体群延迟估计算子，具体表示如下：
Figure RE-FDA0003690429230000031
其中，
Figure RE-FDA0003690429230000032
为信号的实际群延迟，
Figure RE-FDA0003690429230000033
为实际群延迟的变化率，σ为窗函数的支撑因子，ω为频率因子，a为局部尺度变量，b为局部时间变量；
将群延迟估计算子
Figure RE-FDA0003690429230000034
中局部时间变量b不断用
Figure RE-FDA0003690429230000035
替代，获得迭代群延迟估计算子。
9.一种计算机可读存储介质，其上存储有程序，其特征在于，该程序被处理器执行时实现如权利要求1-6任一项所述的一种用于轴承故障信号的时频分析方法。
10.一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的程序，其特征在于，所述处理器执行所述程序时实现如权利要求1-6任一项所述的一种用于轴承故障信号的时频分析方法。
CN202210395577.1A 2022-04-15 2022-04-15 一种用于轴承故障信号的时频分析方法及系统 Pending CN114757230A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210395577.1A CN114757230A (zh) 2022-04-15 2022-04-15 一种用于轴承故障信号的时频分析方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210395577.1A CN114757230A (zh) 2022-04-15 2022-04-15 一种用于轴承故障信号的时频分析方法及系统

Publications (1)

Publication Number Publication Date
CN114757230A true CN114757230A (zh) 2022-07-15

Family

ID=82330774

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210395577.1A Pending CN114757230A (zh) 2022-04-15 2022-04-15 一种用于轴承故障信号的时频分析方法及系统

Country Status (1)

Country Link
CN (1) CN114757230A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115728612A (zh) * 2022-12-01 2023-03-03 广州广电计量检测股份有限公司 变压器放电故障诊断方法及装置
CN117740377A (zh) * 2023-12-19 2024-03-22 济南大学 一种基于幅值谱二阶瞬态提取变换的轴承故障诊断方法

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115728612A (zh) * 2022-12-01 2023-03-03 广州广电计量检测股份有限公司 变压器放电故障诊断方法及装置
CN117740377A (zh) * 2023-12-19 2024-03-22 济南大学 一种基于幅值谱二阶瞬态提取变换的轴承故障诊断方法

Similar Documents

Publication Publication Date Title
Li et al. Fault diagnosis for rolling bearing based on VMD-FRFT
CN114757230A (zh) 一种用于轴承故障信号的时频分析方法及系统
Yu A multisynchrosqueezing-based high-resolution time-frequency analysis tool for the analysis of non-stationary signals
CN101561314B (zh) 随机共振-混沌微弱信号检测方法
Ming et al. Weak fault feature extraction of rolling bearing based on cyclic Wiener filter and envelope spectrum
CN107783938B (zh) 一种旋转设备瞬时转速估计方法
Wang et al. Maximum envelope-based autogram and symplectic geometry mode decomposition based gear fault diagnosis method
CN113405823B (zh) 一种基于迭代扩展本征模态分解的旋转机械故障诊断方法
WO2020057066A1 (zh) 基于增强调制双谱分析的滚动轴承故障诊断方法
CN105527077A (zh) 一种基于振动信号的通用旋转机械故障诊断与检测的方法
CN111414893A (zh) 基于vmd精细复合多尺度散布熵的转子故障特征提取方法
Liu et al. Asymmetric penalty sparse model based cepstrum analysis for bearing fault detections
CN107831013A (zh) 一种利用概率主分量分析增强循环双谱的轴承故障诊断方法
Yan et al. Adaptive synchroextracting transform and its application in bearing fault diagnosis
Antoni et al. Cyclostationarity in condition monitoring: 10 years after
Dan et al. An automatic real-time cable modal frequency identification and tracking algorithm by combining recursive band-pass filter and recursive Hilbert transform
CN116361733A (zh) 一种故障诊断方法、装置、系统以及存储介质
Li et al. An iterative reassignment based energy-concentrated TFA post-processing tool and application to bearing fault diagnosis
Dong et al. Time–frequency-multisqueezing transform
Lv et al. Longitudinal synchroextracting transform: A useful tool for characterizing signals with strong frequency modulation and application to machine fault diagnosis
CN116067657A (zh) 一种滚动轴承故障诊断方法及系统
Yi et al. High-order synchrosqueezing superlets transform and its application to mechanical fault diagnosis
Yi et al. Time-varying fault feature extraction of rolling bearing via time–frequency sparsity
CN107941511B (zh) 一种基于信号时频分解的频率—峭度图的实现方法
CN111353415B (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