CN112665858A - 一种用于对轴承故障信号进行nsst时频分析的方法 - Google Patents

Abstract

本发明提出了一种用于对轴承故障信号进行NSST时频分析的方法，包括以下步骤：对非平稳信号进行非线性稀疏分解，得到分解后的各个信号稀疏分量；对分解后的各个信号稀疏分量依次进行同步压缩变换SST，得到压缩后的稀疏分量；将压缩后的稀疏分量按时频重排，并将重排后的结果进行累加整合，得到最终处理完成的非平稳信号时频分布。本发明的有益效果是：能够较好的提高轴承故障信号时频分布的聚集度，更为精准地获取由故障所产生的频率曲线，在轴承故障的诊断中得到良好的应用。

Description

一种用于对轴承故障信号进行NSST时频分析的方法

技术领域

本发明涉及非平稳信号处理领域，尤其涉及一种用于对轴承故障信号进行NSST时频分析的方法。

背景技术

滚动轴承是机械设备中的重要零部件，运行状态往往直接影响整台机器的性能，但存在故障率高、易损坏等问题，若能及时地有效提取出故障特征信息，实现轴承运行状态的准确判断，并对轴承进行及时的更换或修复，可以有效避免连锁故障的发生，对于减少经济损失意义重大。目前对于滚动轴承故障信号的诊断，主流的方案是根据谱峭度进行故障的诊断。而在日益复杂的工业背景下，传统的方案，对轴承故障数据难以做到准确提取，且即使能够提取，也会存在时频分布聚集度不高的问题。

发明内容

针对以上技术问题，本发明提出一种用于对轴承故障信号进行NSST时频分析的方法，该方法是基于非线性稀疏算法(NMP:Nonlinear Matching Pursuit)的同步压缩变换(SST:Synchrosqueezing Transform)的改进方法。稀疏实现过程中设定不同的稀疏分量，得到稀疏分解信号；将SST应用于稀疏分解时频谱的优化，以提高时频分布的聚集度，且很大程度上减少噪声干扰对信号时频表示的影响。

本发明提出的一种用于对轴承故障信号进行NSST时频分析的方法，具体包括以下步骤：

S101：对非平稳信号进行非线性稀疏分解，得到分解后的各个信号稀疏分量；

S102：对分解后的各个信号稀疏分量依次进行同步压缩变换SST，得到压缩后的稀疏分量；

S103：将压缩后的稀疏分量按时频重排，并将重排后的结果进行累加整合，得到最终处理完成的非平稳信号时频分布。

进一步地，步骤S101具体为：

S201：首先，信号f(t)的初始化残差为r⁰＝f(t)；

S202：其次，构建字典集D＝{a(t)cosθ(t):θ'(t)≥0,a(t)≥0；a(t),θ'(t)∈V(θ)}；其中，θ(t)为待求分解系数；θ′(t)为θ(t)的转置；a(t)为待稀疏分解信号；V(θ)为超完备傅里叶基张成的线性空间；同时，将V(θ)进一步构造为一组超完备傅里叶基张成的线性空间V(θ,λ)；V(θ,λ)是所有比cosθ(t)平滑的函数的集合；λ是控制V(θ)平滑性的参数；可以得到非平稳信号f(t)在字典集D上的非线性稀疏分解和V(θ,λ)空间的表示；

S203：利用初始化残差，求解正则化非线性最小二乘问题，得到更新的残差；

S204：若更新的残差小于且等于预设的阈值，则非线性稀疏分解结束，得到分解后的各个信号稀疏分量；否则，返回步骤S203。

进一步地，步骤S202中，

其中

进一步地，步骤S203中，求解正则化非线性最小二乘问题的公式如下：

P:

s.t.θ'_k≥0,a_k(t)∈V(θ,λ)

式中，

是a_k在V(θ,λ)空间的表示；k表示迭代次数。

进一步地，步骤S203中，更新的残差如式：

进一步地，步骤S102具体为：

S301：分解后的各个信号稀疏分量作为SST的输入s(t)，根据式：

得到小波系数W_s(a,b)；其中a为小数位数，b为时移，Ψ(·)表示母小波；

S302：根据所述小波系数计算s(t)的瞬时频率估计值：

S303：将瞬时频率估计值作为同步压缩算法SST的重排规则，并依据该规则对所述小波系数进行压缩，得到压缩后的小波系数谱T_s(ω,b)，即压缩后的稀疏分量。

步骤S303中对所述小波系数进行压缩，压缩的具体公式如下：

l表示连续频段[ω_l-Δω/2,ω_l+Δω/2]和Δω＝ω_l-ω_l-1的中心。

本发明提供的有益效果是：能够较好的提高轴承故障信号时频分布的聚集度，更为精准地获取由故障所产生的频率曲线，在轴承故障的诊断中得到良好的应用。

附图说明

图1是一种用于对轴承故障信号进行NSST时频分析的方法的流程图；

图2是在-1dB噪声下，仿真信号x(t)的时频分布(a)HHT；(b)STFT；(c)WVD；(d)ST；(e)STST；(f)NSST；

图3是在不同噪声背景下x(t)的SST与NSST的时频分布(a)15dB(b)0dB(c)-12dB；

图4是不同SNR下的CM值比较；

图5是SST、NSST处理不同风扇端轴承故障数据的时频分布(a)278.mat；(b)282.mat；

图6是时频局部放大图(a)SST；(b)NSST；

图7是SST、NSST处理不同驱动端故障信号的时频分布(a)105.mat；(b)156.mat；

图8是时频的局部放大图(a)SST；(b)NSST。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方式作进一步地描述。

在具体阐述本发明前，先对相关专业术语统一翻译如下：

(1)同步压缩变换(SST：Synchrosqueezing Transform)，SST是基于STFT(或者WT)的一种后处理工具，从频率方向对STFT的结果进行压缩，从而提升时频谱的能量聚集性；

(2)非线性稀疏(NMP:Nonlinear Matching Pursuit)；

请参考图1，本发明提供一种用于对轴承故障信号进行NSST时频分析的方法，具体包括以下步骤：

S101：对非平稳信号进行非线性稀疏分解，得到分解后的各个信号稀疏分量；

步骤S101具体为：

S201：首先，信号f(t)的初始化残差为r⁰＝f(t)；

S202：其次，构建字典集D＝{a(t)cosθ(t):θ'(t)≥0,a(t)≥0；a(t),θ'(t)∈V(θ)}；其中，θ(t)为待求分解系数；θ′(t)为θ(t)的转置；a(t)为待稀疏分解信号；V(θ)为超完备傅里叶基张成的线性空间，

其中

同时，将V(θ)进一步构造为一组超完备傅里叶基张成的线性空间V(θ,λ)；V(θ,λ)是所有比cosθ(t)平滑的函数的集合；λ是控制V(θ)平滑性的参数；可以得到非平稳信号f(t)在字典集D上的非线性稀疏分解和V(θ,λ)空间的表示；

S203：利用初始化残差，求解正则化非线性最小二乘问题，得到更新的残差；求解正则化非线性最小二乘问题的公式如下：

P:

s.t.θ'_k≥0,a_k(t)∈V(θ,λ)

式中，

是a_k在V(θ,λ)空间的表示；k表示迭代次数。

S204：若更新的残差小于且等于预设的阈值，则非线性稀疏分解结束，得到分解后的各个信号稀疏分量；否则，返回步骤S203。

S102：对分解后的各个信号稀疏分量依次进行同步压缩变换SST，得到压缩后的稀疏分量；

步骤S102具体为：

S301：分解后的各个信号稀疏分量作为SST的输入s(t)，根据式：

得到小波系数W_s(a,b)；其中a为小数位数，b为时移，Ψ(·)表示母小波；

S302：根据所述小波系数计算s(t)的瞬时频率估计值：

S303：将瞬时频率估计值作为同步压缩算法SST的重排规则，并依据该规则对所述小波系数进行压缩，得到压缩后的小波系数谱T_s(ω,b)，即压缩后的稀疏分量。压缩的具体公式如下：

l表示连续频段[ω_l-Δω/2,ω_l+Δω/2]和Δω＝ω_l-ω_l-1的中心；

S103：将压缩后的稀疏分量按时频重排，并将重排后的结果进行累加整合，得到最终处理完成的非平稳信号时频分布。

实施例一：

本发明对最终处理完成的非平稳信号时频分布进行仿真验证；

时频聚集度是指在由时频分析方法获得的时频分布中，估计的频率分量在真实频率位置的聚集程度。时频聚集度越高就表明时频分布对于信号频率分量的局部定位和分辨越准确，从而说明该时频分析方法的处理效果越好。本发明将CM值作为衡量时频分布聚集度的性能评价指标，其计算公式为：

其中，n为时间窗长度，ω为信号在某点的频率。CM的值越大，表示时频分布的聚集度越高。

NSST方法可以绘制信号的时-频分布图，能有效的提高时频聚集度，且抗噪声性能比较优越。通过下面的实验验证了该方法的优点，设计了多组仿真实验分别加以证明，可作为对后面实际信号分析可靠性的支持，通过对不同信号的处理结果进行分析，证明了NSST方法应用的广泛性以及有效性。

仿真信号x(t)：x(t)＝cos(8t)+cos(t²+t+cos(t))

(1)对比经典时频分析方法

信号时间区间为[0s,15s]，采样频率为800。该仿真信号由恒频信号与线性调频信号组成。在SNR＝-1dB高斯噪声环境下。通过处理仿真信号x(t)，来对NSST与HHT、STFT、WVD、ST、SST、等方法来进行对比，从而来验证NSST方法的优越性。

图2为通过以上六种时频分析方法处理混合信号得到的结果。对于含噪声的混合信号，图2(a)的SST和图2(d)的ST算法对噪声处理较差，存在严重的噪声干扰，时频曲线非常模糊；图2(b)的STFT可以大致看出两信号分量，但是时频聚集度较差，频率曲线呈带状分布；图2(c)的WVD也能够对两分量有基本的提取，但是受到交叉项的影响导致时频表达效果较差；图2(f)是本改进算法的整个过程展现，可以看出图2(e)的SST时频图聚集度较差，两信号分量频率附近存在较为明显的噪声干扰，而NSST能有有效的抑制噪声的影响，准确提取信号分量，同时保证时频聚集度的进一步提高，达到更为理想的效果；

从时频分布的视觉观察上，SST和NSST都能表现信号的幅度和频率调制，为了进一步比较两者的时频聚集度，表1列出了不同时频分布图的CM值。CM的值越大表示时频聚集度越高，NSST的CM值最大，对应在图2中WVD和ST时频图的聚集度都比较差，而SST和NSST的CM较大对应两者的时频聚集度较高，也说明了CM值作为评价时频聚集度的性能指标的可靠性。NSST方法的CM值为高于SST方法，即NSST的时频聚集度较高，同样也证明本改进方法在提高聚集度方面有效。

表1信号x(t)不同时频分布的CM(10^-5)

HHT	STFT	WVD	ST	SST	NSST
						51.340	7.8043	3.5758	5.5369	38.833	57.595

综合上述结果可以看出其他的时频分析方法在时频表现上比较差，因此后面单独与SST进行详细比较。

(2)不同噪声背景下与SST进行详细对比

时间区间取为[0s,15s]，采样频率为800。同样的，为了表现出一般性，仿真信号采用混合信号(恒频信号与线性调频信号)。图4是基于不同分析方法得出的信号的时频分布结果，对比可以看出NSST算法的优越性。

图3为不同噪声背景下SST方法和NSST方法处理相同仿真信号的时频结果图。从图3(a)到图3(c)，随着噪声水平的增加，SST方法的时频表示越来越差，当SNR＝-12dB时，频率分量已经被噪声完全覆盖，与之对应的是无论噪声水平如何，NSST都能够以更高的精度提供清晰的时频曲线以及时频聚集度，表明NSST方法可以有效抑制背景噪声，尽管时频曲线随着噪声的增强而略有失真。

表2信号x(t)不同噪声背景的时频分布的CM(10^-4)

为进一步详细证明这一结论，表2中列出了SST和NSST方法在不同噪声SNR下的CM值，图4进行了比较分析。结果表明NSST方法的CM始终高于SST，并且当SNR>3dB时，CM几乎不变，这表明，针对不同的噪声水平，本发明提出的方法有着较强的鲁棒性，同时在信号包含强噪声的情况下，即使噪声能量高于信号能量，NSST方法可以得出较高时频聚集度的时频表示，也能够更为准确的反映时频信号，适合信号的精准分析.

实施例二：

本发明算法可以较好的得到该信号的时间-频率的联合分布情况，有助于观察故障信号的时频分布变化，为故障识别诊断提供依据和帮助。本发明选取凯斯西储备大学轴承数据中心网站数据，测试台由2hp电动机、扭矩传感器/编码器、测力计和控制电子设备组成。当前网站数据电机转速为1797r/min。轴承故障数据主要分为风扇端轴承故障和驱动端轴承故障两类，其中每类分别取不同的两组数据进行处理对比，每组数据分为DE、FE、BA三部分，分别代表了驱动端加速数据、风扇端加速数据和基础加速度数据，根据不同的故障类型选择不同的数据进行处理，根据官方提供的采样频率、轴承参数和转速，可以计算出驱动端轴承故障和风扇端轴承故障的故障基本特征频率分别为162.2Hz和148.2Hz，因此后面的实验也主要是提取该范围的信号特征。对不同的轴承故障信号进行时频分析，并与SST方法进行比较和分析，可以得出的不同方法的时频结果。

首先，选择风扇端轴承故障数据(文件名：278.mat、282.mat)，故障直径0.007’，无电动机负载，电动机速度为1787rpm，故障特征频率为148.2Hz。图5(a)分别是SST和NSST处理278.mat数据集的结果，图5(b)是SST和NSST处理282.mat的结果。由于有实际环境噪声的干扰，可以看出SST处理后得出的时频图不清晰，噪声干扰严重，且时频分布较为分散，频率波动范围也大，故障信号的频率分量无法直接识别并定位，对比之下NSST的时频表示明显要优于SST，首先能够有效的分离噪声和故障信号，提取有效的故障信号信息，虽然图形上方存在少量的噪声干扰，但对于轴承故障信号的识别与分析没有什么影响，可得出初步结论，NSST方法能够有效的处理实际轴承故障信号，分离干扰噪声识别故障信号信息，并提供高时频聚集度的时频图。

图6是图5(a)中SST和NSST时频分布的局部放大对比，可以看出SST方法处理结果的时频表示很模糊，尤其是噪声对于轴承故障信号的干扰较为严重，时频分布也较分散，导致无法准确有效的识别出有效信息；而图6(b)中NSST方法处理后的时频图明显地清晰，干扰的噪声也能有效的滤除，对于故障信号频率分布的细节描述准确清晰，可看出频率集中在150Hz的位置附近，也符合故障信号的特征频率，可诊断为风扇端轴承故障。

图7为SST和NSST处理驱动端轴承故障数据(文件名：105.mat、156.mat)的结果,基本参数为：故障直径为0.007’，无电动机负载，电动机速度为1797rpm。经计算故障信号的频率为162.2Hz。同样的，SST处理的结果具有时频图不清晰，噪声干扰严重，且存在能量分布分散，频率波动范围的问题，因此很难去识别故障信号以及其具有的信息特征。相对应的，NSST能够有效的分离噪声和故障信号，处理得到的时频结果明显地更加清晰、且易于识别故障信号，同时故障信号分布较为集中，表明时频分布的时频聚集度大。图8为图7(a)中的局部放大对比，SST方法处理的结果还是存在时频表示模糊、噪声干扰严重等问题，而NSST虽然在0.3S左右有一定的失真，但总体时频表现明显优于SST方法，且很好地分离出轴承故障信号，频率集中在160Hz附近，也符合驱动端轴承故障的特征频率，可诊断为驱动端轴承故障。

从与SST时频分析方法的处理结果比较来看，传统的SST方法无法有效的分离噪声和目标，尤其是故障信号频率附近受噪声干扰严重，导致时频图信息混乱，同时处理得到的时频图时频分布分散，表明能量分布较为分散、时频聚集度较差，时频曲线细节刻画较为粗糙。本发明所提出的NSST方法可以有效地分离故障信号和噪声，得出的时频分布显示故障信号频率分布较为集中，能够对轴承故障信号的准确识别和定位，从而实现故障类型的诊断。

本发明解决了实际复杂工业背景下轴承故障数据处理中有效信息不能准确提取以及时频分布聚集度不高的问题。通过大量函数仿真验证了NSST时频分析方法可以较为理想的处理所研究的数据，可以更加准确的估计信号瞬时频率，得到高聚集度的时间-频率分布，更为精准地获取由故障所产生的频率曲线，区分和辨别故障类型。

本发明提供的有益效果是：能够较好的提高轴承故障信号时频分布的聚集度，更为精准地获取由故障所产生的频率曲线，在轴承故障的诊断中得到良好的应用。

以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种用于对轴承故障信号进行NSST时频分析的方法，其特征在于：具体包括以下步骤：

S101：对非平稳信号进行非线性稀疏分解，得到分解后的各个信号稀疏分量；

S102：对分解后的各个信号稀疏分量依次进行同步压缩变换SST，得到压缩后的稀疏分量；

S103：将压缩后的稀疏分量按时频重排，并将重排后的结果进行累加整合，得到最终处理完成的非平稳信号时频分布。

2.如权利要求1所述的一种用于对轴承故障信号进行NSST时频分析的方法，其特征在于：

步骤S101具体为：

S201：首先，信号f(t)的初始化残差为r⁰＝f(t)；

S202：其次，构建字典集D＝{a(t)cosθ(t):θ'(t)≥0,a(t)≥0；a(t),θ'(t)∈V(θ)}；其中，θ(t)为待求分解系数；θ′(t)为θ(t)的转置；a(t)为待稀疏分解信号；V(θ)为超完备傅里叶基张成的线性空间；同时，将V(θ)进一步构造为一组超完备傅里叶基张成的线性空间V(θ,λ)；V(θ,λ)是所有比cosθ(t)平滑的函数的集合；λ是控制V(θ)平滑性的参数；可以得到非平稳信号f(t)在字典集D上的非线性稀疏分解和V(θ,λ)空间的表示；

S203：利用初始化残差，求解正则化非线性最小二乘问题，得到更新的残差；

S204：若更新的残差小于且等于预设的阈值，则非线性稀疏分解结束，得到分解后的各个信号稀疏分量；否则，返回步骤S203。

3.如权利要求2所述的一种用于对轴承故障信号进行NSST时频分析的方法，其特征在于：步骤S202中，

其中

4.如权利要求2所述的一种用于对轴承故障信号进行NSST时频分析的方法，其特征在于：步骤S203中，求解正则化非线性最小二乘问题的公式如下：

s.t.θ'_k≥0,a_k(t)∈V(θ,λ)

式中，

是a_k在V(θ,λ)空间的表示；k表示迭代次数。

5.如权利要求2所述的一种用于对轴承故障信号进行NSST时频分析的方法，其特征在于：步骤S203中，更新的残差如式：

6.如权利要求1所述的一种用于对轴承故障信号进行NSST时频分析的方法，其特征在于：步骤S102具体为：

S301：分解后的各个信号稀疏分量作为SST的输入s(t)，根据式：

得到小波系数W_s(a,b)；其中a为小数位数，b为时移，Ψ(·)表示母小波；

S302：根据所述小波系数计算s(t)的瞬时频率估计值：

S303：将瞬时频率估计值作为同步压缩算法SST的重排规则，并依据该规则对所述小波系数进行压缩，得到压缩后的小波系数谱T_s(ω,b)，即压缩后的稀疏分量。

7.如权利要求6所述的一种用于对轴承故障信号进行NSST时频分析的方法，其特征在于：步骤S303中对所述小波系数进行压缩，压缩的具体公式如下：

l表示连续频段[ω_l-Δω/2,ω_l+Δω/2]和Δω＝ω_l-ω_l-1的中心。

Images