王良斌 樊 凱 劉 洋

1武漢鋼鐵有限公司煉鋼廠 武漢 430080 2寶鋼股份中央研究院（武鋼有限技術(shù)中心） 武漢 430080

0 引言

在冶金行業(yè)，滾動軸承是鋼鐵冶金機(jī)械設(shè)備的重要組成部分，大多工作在變速、重載的環(huán)境中，其發(fā)生故障的概率隨使用年限的增加而增加[1]。一旦滾動軸承在運(yùn)行過程中發(fā)生故障，將直接導(dǎo)致工業(yè)生產(chǎn)事故的發(fā)生，對工作人員的生命造成極大的威脅，故監(jiān)測和診斷滾動軸承的運(yùn)行狀態(tài)具有重大的工程意義[2]。

目前，滾動軸承故障診斷通常是通過采集振動加速度信號進(jìn)行分析、處理的。然而，實(shí)際生產(chǎn)過程中通過傳感器獲得振動加速度信號常常存在明顯的背景噪聲，而噪聲極大地影響了分析結(jié)果的準(zhǔn)確性，故如何去除信號中的噪聲是分析信號的關(guān)鍵步驟之一[3]。小波閾值去噪[4]、傅里葉濾波去噪[5]等傳統(tǒng)降噪方法的主要思想是在變換域中劃分有用信號和噪聲分量的存在范圍，是通過選擇合適的參數(shù)或基函數(shù)來去除噪聲分量，難以適用于每一次采集到的振動信號。由此，出現(xiàn)了其他的分解方法以分離信號中的噪聲。奇異值分解（Singular Value Decomposition，SVD）[6]是將信號轉(zhuǎn)化為奇異矩陣，并在對奇異矩陣進(jìn)行特征分解時(shí)進(jìn)行本征值選取，從而實(shí)現(xiàn)信號的降噪，但因SVD存在特征值難以選擇的缺點(diǎn)，使得去噪效果存在較大偏差。隨機(jī)奇異值分解及軟閾值（Random Singular Value Decomposition And Soft Threshold，rSVD-ST）去噪算法是奇異值分解、隨機(jī)矩陣降維和軟閾值算法相結(jié)合的去噪方法[7]，它對于有阻尼或無阻尼的諧波信號都有很好效果，且具有處理大數(shù)據(jù)的能力，并在軟閾值處理時(shí)保留了被視為有用信號的頻譜峰值強(qiáng)度，可抑制噪聲區(qū)域，有利于下一步進(jìn)行時(shí)頻分析。

針對變轉(zhuǎn)速工況下機(jī)械設(shè)備振動信號的故障診斷，時(shí)頻分析（Time-Frequency Analysis，TFA）是如今最直觀有效的工具，通過對去噪后的信號進(jìn)行時(shí)頻（Time Frequency，TF）變換，采用時(shí)頻域聯(lián)合分布來描述信號的瞬態(tài)特性，從而獲得滾動軸承故障特征的瞬時(shí)頻率(Instantaneous Frequency，IF)曲線[8]。經(jīng)典時(shí)頻變換算法如短時(shí)傅里葉變換（Short-Time Fourier Transform，STFT）[9]是通過對信號進(jìn)行加窗，并逐一對每個(gè)窗口內(nèi)進(jìn)行時(shí)頻表示，將得到的局部時(shí)頻表示連接起來，最終得到的時(shí)頻圖能量整體比較發(fā)散。Daubechies I等[10]提出的同步壓縮變換（Synchrosqueezing Transform，SST）和于剛等[11]提出的同步提取變換（Synchroextracting Transform，SET）雖在STFT的基礎(chǔ)上進(jìn)一步集中了時(shí)頻圖中的能量，但在強(qiáng)噪聲和多分量的情況下時(shí)頻分析的效果仍然有限。Chirplet 變換（Chirplet Transform，CT）有較高的時(shí)頻能量聚集性，但也有明顯的局限性，即在給定點(diǎn)處的窗口長度中IF軌跡的斜率不變，需要不同的核相位函數(shù)來匹配不同時(shí)刻的斜率；CT的TF基不能同時(shí)匹配任何給定時(shí)刻所有分量的IF軌跡（因?yàn)椴煌至康腎F軌道具有不同的斜率），故CT難以在不同頻率分量下都獲得最佳時(shí)頻分辨率[12]。因此，Li M F等[13]提出尺度基線性調(diào)頻變換（Scaling-Basis Chirplet Transform，SBCT）能解決以上問題，它通過構(gòu)造了一個(gè)新的核相位函數(shù)，能在窗口范圍內(nèi)的時(shí)間中心及其周圍縮放TF基來擴(kuò)展CT，生成隨時(shí)間和頻率變化的線性調(diào)頻波，該線性調(diào)頻波可在整個(gè)窗口長度范圍上匹配多個(gè)IF軌跡的斜率，從而實(shí)現(xiàn)多分量信號的時(shí)頻表達(dá)。

基于上述研究背景，本文提出一種基于rSVD-ST和SBCT的滾動軸承故障診斷方法，該方法利用rSVDST對變轉(zhuǎn)速振動信號進(jìn)行降噪，再對降噪后的信號進(jìn)行SBCT時(shí)頻分析，以獲得清晰地時(shí)頻表達(dá)，并得到變轉(zhuǎn)速工況下的特征頻率，最后通過經(jīng)驗(yàn)公式計(jì)算滾動軸承可能出現(xiàn)的故障特征頻率[14]，通過比較診斷出軸承故障類型。本文通過數(shù)值模擬仿真分析，驗(yàn)證了該方法的有效性，并通過實(shí)驗(yàn)臺數(shù)據(jù)分析進(jìn)行了驗(yàn)證。

1 理論描述

1.1 隨機(jī)奇異值分解-軟閾值去噪

均勻采樣的諧波信號的一般模型可表示為

式中：x（n）為觀測信號，s（n）、w（n）分別為無噪聲信號、噪聲，L為抽樣總數(shù)。

將無噪聲信號s（n）重新排列為M×N的Hankel矩陣H（s），且滿足L≥M+N-1和M，N＞K，即有

式中：s是s（n）的矢量化形式，H(·)表示將向量轉(zhuǎn)換為Hankel矩陣的運(yùn)算符。

在執(zhí)行SVD去噪之前，首先降低SVD在大矩陣上的計(jì)算量，通過應(yīng)用隨機(jī)投影來降低矩陣維數(shù)

(x)的SVD可表示為

提取與前K1個(gè)奇異值對應(yīng)的前K1個(gè)子矩陣，并重構(gòu)信號的Hankel矩陣為

因此，SVD算式（5）等價(jià)于(x)在信號子空間s上的投影，即

式（5）、式（6）的輸出結(jié)果通常不是Hankel矩陣，需將s改為Hankel矩陣，通過求得第n個(gè)反對角線元素的平均值，從s中提取初步去噪后的信號(n)，即

同時(shí)，考慮一個(gè)經(jīng)典的去噪優(yōu)化問題，即

假設(shè)已獲得頻譜中峰值位置的良好估計(jì)（即頻譜權(quán)重向量SW），只需約束頻譜中噪聲區(qū)域的稀疏性，并將式（8）改為

2.2 尺度基線性調(diào)頻波變換

尺度基線性調(diào)頻波變換（SBCT）信號s(u)∈L2(R)的CT可表示為

式中：s（u）為由x（u）希爾伯特變換生成的解析信號；h(u)∈L2(R)為非負(fù)、對稱和歸一化的實(shí)窗函數(shù)，通常是高斯函數(shù)；φ(f，u，tc)為相位函數(shù)；tc∈R、f∈R分別為時(shí)間和頻率中心；C∈R為調(diào)頻率。

對于IF軌跡隨時(shí)間非線性變化的多分量信號，CT有3個(gè)限制。首先，IF軌跡的斜率隨時(shí)間變化，需要不同的調(diào)頻率來匹配不同時(shí)刻的斜率并實(shí)現(xiàn)高能量集中。其次，在時(shí)頻面同一時(shí)刻不能同時(shí)匹配所有分量的IF軌跡，因?yàn)椴煌至康腎F軌跡具有不同的斜率。第3個(gè)限制是對在整個(gè)給定的窗口長度中使用的斜是不改變，即給定時(shí)刻的中心斜率不匹配整個(gè)窗口長度上IF軌跡的時(shí)變斜率，故通過利用一種新的核相位函數(shù)以滿足每個(gè)時(shí)刻都能匹配所有IF軌跡。相位函數(shù)的二階導(dǎo)數(shù)是頻率中心、時(shí)間中心和積分變量的函數(shù)，即

式中：tc和f分別為時(shí)間和頻率中心，且定義所有IF軌跡在tc時(shí)的中心頻率為TF基；u是積分變量；θ是TF基的旋轉(zhuǎn)角。

為突破3個(gè)限制，需構(gòu)造一個(gè)新的相位函數(shù)φs，即

式中：（a1，a2，…，an）為必須確定的參數(shù)。相應(yīng)的第1和第2導(dǎo)數(shù)可寫為

式中：φs'為相位函數(shù)φs的IF軌跡。

以時(shí)間中心tc為特定窗口，其中u服從（tc-L/2，tc+L/2），當(dāng)u等于時(shí)間中心tc時(shí)，式（16）可寫為

式中：θ是TF基在時(shí)間中心tc處的旋轉(zhuǎn)角。

由上式開始，θ值隨頻率中心f而變化，而不是保持恒定。假設(shè)f（t）和f1（t）在時(shí)間中心tc的中心頻率為fc0和fc1，為了同時(shí)匹配這些分量，TF基的相應(yīng)旋轉(zhuǎn)角θ0和θ1應(yīng)等于f（t）和f1（t）在fc0和fc1處的傾斜角。為實(shí)現(xiàn)此目標(biāo)，應(yīng)為a1選擇一個(gè)值，以確保以下方程成立，即

式中：θ0和θ1為TF基的旋轉(zhuǎn)角，α和β為f（t）和f1（t）在時(shí)間中心tc處的傾斜角。

在這種情況下，TF基與時(shí)間中心tc處的IF軌跡與f（t）和f1（t）匹配，表明與傳統(tǒng)CT相比，SBCT為每個(gè)頻率分量提供了令人滿意的最佳能量集中。

在考慮tc+Δu時(shí)，即TF基在單個(gè)窗口長度中旋轉(zhuǎn)角度的情況，通過將時(shí)間中心增加少量Δu∈（-L/2，L/2），對于不同的頻率中心可改寫為

在式（19）和式（20）中，TF基的旋轉(zhuǎn)角隨Δu變化。當(dāng)Δu從–L/2變化到L/2時(shí)，即u值從tc-L/2變化到tc+L/2時(shí)，如果（a1，a2，…，an）值正確分配，則TF基的旋轉(zhuǎn)角等于整個(gè)窗口的IF軌跡的傾斜角。由此，當(dāng)滿足以下條件時(shí)，可以確定參數(shù)（a1，a2，…，an）可表示目標(biāo)信號的IF彼此成比例，即

式中：fcj1和fcj2分別為目標(biāo)信號在中心時(shí)間tc的任意2個(gè)IF的中心頻率。

綜合以上算式，SBCT可表示為

2.3 基于rSVD-ST和SBCT的變轉(zhuǎn)速滾動軸承故障診斷方法

圖1為基于rSVD-ST和SBCT的變轉(zhuǎn)速滾動軸承故障診斷方法的流程圖，其實(shí)現(xiàn)的步驟如下：

圖1 變轉(zhuǎn)速滾動軸承故障診斷方法流程圖

1）利用rSVD-ST對變轉(zhuǎn)速振動信號進(jìn)行降噪，得到除去強(qiáng)背景噪聲的高信噪比信號；

2）利用SBCT 對去噪后的信號進(jìn)行時(shí)頻分析，以獲得時(shí)頻能量集中的時(shí)頻表達(dá)，進(jìn)而得到實(shí)驗(yàn)對象在變轉(zhuǎn)速工況下的故障特征頻率軌跡；

3）通過經(jīng)驗(yàn)公式計(jì)算不同類型故障的理論故障特征頻率，并與SBCT計(jì)算后獲得的故障特征頻率進(jìn)行比較，從而判斷出變轉(zhuǎn)速下滾動軸承的故障類型和故障部位。

3 數(shù)值仿真分析

為了驗(yàn)證該方法的有效性，通過設(shè)計(jì)多分量信號進(jìn)行數(shù)字仿真分析，其中多分量信號的IF軌跡彼此非常接近，這樣更具有說服力。因此，構(gòu)建一個(gè)間隔頻率非常近的模擬信號s（t），即

式中：s1（t）、s2（t）、s3（t）、s4（t）分別為4個(gè)不同的調(diào)頻分量；n（t）為噪聲分量，并設(shè)置信噪比SNR＝-5 dB。

圖2為仿真信號的時(shí)域圖和頻譜圖。由于噪聲的干擾，圖2b中的頻譜結(jié)構(gòu)非常雜亂，無法對定義信號所在頻帶進(jìn)行判斷，故需利用信號分解算法進(jìn)行去噪處理。為了證明rSVD-ST去噪的效果，利用小波閾值去噪、SVD分解去噪和本文所提方法分別對仿真信號進(jìn)行降噪處理。

圖2 仿真信號時(shí)域圖和頻譜

圖3為其降噪后的時(shí)域圖和頻譜圖，從圖中去噪結(jié)果可以看出，rSVD-ST去噪的效果相較于小波閾值去噪和SVD分解去噪的效果更好，并用信噪比衡量降噪的效果。信噪比(SNR)是信號功率與噪聲功率的比值，其表達(dá)式為

圖3 3種算法的去噪結(jié)果

式中：PS為信號功率，PN為噪聲功率。

表1 是不同去噪方法降噪后的信噪比（SNR），信噪比數(shù)值越高，則去噪效果越好。

表1 不同去噪方法降噪后的信噪比

對于時(shí)變信號，傳統(tǒng)頻譜分析難以表征信號的時(shí)變信息，故需利用時(shí)頻分析的方法對去噪后的信號進(jìn)行信號處理。隨后，分別使用STFT、SST、SET和SBCT進(jìn)行時(shí)頻分析，計(jì)算結(jié)果如圖4所示。由圖4可知，SBCT具有最佳處理效果，即使頻率間隔較近的IF脊線也能很好地分隔開，而其他幾種方法均存在不同程度的時(shí)頻能量發(fā)散，在頻率間隔較近時(shí)，相近的IF脊線會出現(xiàn)紊亂，難以直接分辨。

圖4 4種時(shí)頻分析算法的結(jié)果

4 實(shí)驗(yàn)數(shù)據(jù)分析

為進(jìn)一步驗(yàn)證該方法在實(shí)際應(yīng)用中的有效性，利用試驗(yàn)臺數(shù)據(jù)進(jìn)行實(shí)驗(yàn)分析，圖5為試驗(yàn)臺實(shí)物圖及結(jié)構(gòu)簡圖。實(shí)驗(yàn)軸承故障類型為外圈故障，實(shí)驗(yàn)變轉(zhuǎn)速工況為電動機(jī)轉(zhuǎn)速從0 r/min加速到1 900 r/min并持續(xù)運(yùn)行約1 s，隨后從1 900 r/min減速到0 r/min。將振動加速度傳感器徑向放置于外置軸承座上，在與實(shí)驗(yàn)軸承相連的聯(lián)軸器上貼上反光片，并將轉(zhuǎn)速傳感器對準(zhǔn)該聯(lián)軸器進(jìn)行轉(zhuǎn)速測量。此外，試驗(yàn)臺減速器的減速比為2:3，采樣頻率為2 560 Hz，實(shí)驗(yàn)采集的時(shí)長為7 s，采集了加速、勻速、減速的過程。滾動軸承型號為6202，其滾動體個(gè)數(shù)z＝8。采集原始振動信號的時(shí)域圖及頻譜如圖6所示。

圖5 試驗(yàn)臺實(shí)物圖及結(jié)構(gòu)簡圖

圖6 振動信號時(shí)域圖和頻譜

由圖6可知，實(shí)測信號中包含有較多噪聲，為了使信號時(shí)頻分析效果更佳，本文首先應(yīng)用rSVD-ST去噪方法對實(shí)測振動信號進(jìn)行降噪，降低實(shí)驗(yàn)中的大量背景噪聲等無關(guān)因素的影響。rSVD-ST去噪結(jié)果及頻譜如圖7所示。

圖7 rSVD-ST去噪結(jié)果及頻譜

利用rSVD-ST去噪后的信號進(jìn)行SBCT分析，其計(jì)算結(jié)果如圖8所示，由圖8可知，SBCT獲得的時(shí)頻圖中能夠很好地識別故障特征頻率fc及其二倍頻2fc、三倍頻3fc等。根據(jù)轉(zhuǎn)速傳感器測得的轉(zhuǎn)速可計(jì)算出轉(zhuǎn)頻fr，根據(jù)故障診斷的經(jīng)驗(yàn)公式fr＝0.4zfr，可以得到滾動軸承的理論外圈故障特征頻率fo，將SBCT計(jì)算的故障特征頻率fc與理論外圈故障特征頻率fo做出比較，可以發(fā)現(xiàn)fc與fo十分吻合，由此判斷實(shí)驗(yàn)滾動軸承為外圈故障，與實(shí)驗(yàn)預(yù)設(shè)的外圈故障一致。

圖8 SBCT計(jì)算結(jié)果

5 結(jié)論

本文提出了一種基于rSVD-ST和SBCT的強(qiáng)背景噪聲下變轉(zhuǎn)速滾動軸承故障診斷算法，首先利用隨機(jī)投影來降低矩陣維數(shù)，然后通過應(yīng)用SVD提取基本特征向量，并應(yīng)用該特征向量重構(gòu)成平滑信號，再對其進(jìn)行軟閾值處理。實(shí)驗(yàn)結(jié)果表明，rSVD-ST能去除復(fù)雜工況下的強(qiáng)背景噪聲。為了能夠獲得高分辨率的時(shí)頻表達(dá)，對去噪后的信號進(jìn)行SBCT分析，得到了清晰的時(shí)頻圖，最后與理論所計(jì)算出的故障特征頻率相對比，從而實(shí)現(xiàn)滾動軸承的故障診斷。本文首先通過數(shù)值仿真分析驗(yàn)證了方法的可行性，并利用試驗(yàn)臺對滾動軸承外圈故障進(jìn)行特征分析。