郭遠晶　金曉航　魏燕定　楊友東

摘要： 為實現旋轉機械設備故障診斷，提出基于相關性檢測的振動信號改進時域同步平均（Time Synchronous Averaging，TSA）降噪方法，并利用平方包絡譜提取故障特征。在采樣振動信號中任取一段參考子信號，令參考子信號從振動信號初始時間開始，沿時間軸逐點滑動至終了時間，同時計算參考子信號與其在振動信號中遮掩子信號之間的相關系數;利用步長迭代算法獲取相關系數最優(yōu)閾值，選擇與參考子信號相關系數大于最優(yōu)閾值的所有遮掩子信號作為總體同步信號，平均計算后獲得降噪的目標信號;對目標信號進行平方包絡譜分析。仿真信號、齒輪和軸承振動信號的分析結果表明，所提出的方法能夠有效用于旋轉機械設備的振動信號降噪與故障特征提取。

關鍵詞： 故障診斷; 時域同步平均; 相關系數; 步長迭代算法; 平方包絡譜

中圖分類號： TH165+.3; TN911.7 ? ?文獻標志碼： A ? ?文章編號： 1004-4523（2021）02-0402-09

DOI：10.16385/j.cnki.issn.1004-4523.2021.02.021

引 ?言

旋轉機械設備出現故障損傷時，振動信號中通常包含有規(guī)律性的故障特征，其提取對于設備故障診斷具有重要意義。但由于噪聲干擾、傳遞路徑衰減、多振源耦合等因素[1]，往往需要對振動信號進行降噪處理，如濾波器降噪、小波閾值降噪、SVD降噪以及經驗模態(tài)分解降噪等，才能順利提取出故障特征。此外，在旋轉機械設備平穩(wěn)運行工況下，適應于振動信號的周期性變化特征，時域同步平均（Time Synchronous Averaging，TSA）降噪方法也得到了廣泛應用[2?7]。

振動信號TSA降噪的關鍵在于參與平均計算的各段振動信號要保持同步，即無相位差。目前保持振動信號同步的方法主要有按整數周期截取信號和時標脈沖觸發(fā)同步采樣[8]。按整數周期截取信號的方法需要獲取振動信號的準確周期[7]，但這并不是一件容易的事情，原因在于旋轉機械設備的轉速通常不能保持時刻恒定，尤其是在故障損傷狀態(tài)下，轉速的波動使振動信號沒有一個嚴格穩(wěn)定的周期，因此，按照一個固定周期對振動信號進行截斷，會出現相位誤差累積效應[9]，影響降噪效果。時標脈沖觸發(fā)同步采樣則需要在參考軸上安裝轉速計或光電編碼器，輸出時標脈沖觸發(fā)振動信號的同步采樣或者插值重采樣[2，5，10?11]。此類方法能夠保證采樣振動信號的準確同步，但需要額外增加硬件成本，并且在很多實際應用場合中，轉速計或光電編碼器的安裝也并非易事[12]，因此，此類方法的應用也受到一定限制。

為了獲取時域同步振動信號參與TSA降噪，本文從旋轉機械設備上采集的振動信號中選取一段參考子信號，令參考子信號從振動信號的初始時間開始，沿時間軸滑動至終了時間，期間參考子信號不斷地在振動信號中遮掩一段等長度的子信號，同時計算參考子信號與遮掩子信號之間的相關系數;然后利用步長迭代算法獲取相關系數的最優(yōu)閾值;最后選擇與參考子信號相關系數大于最優(yōu)閾值的所有遮掩子信號作為總體同步信號，平均計算后獲得降噪的目標信號。對目標信號進行平方包絡譜分析，提取故障特征，判別旋轉機械設備故障狀態(tài)。

1 TSA降噪原理

旋轉機械設備運行時的振動信號為y（t），采樣時間間隔為Δt，則采樣振動信號為y（kΔt），記作y（k），k=0，1，…，N（N為振動信號采樣數據點數）。若y（k）由周期為NT（NT為一個時間周期內的采樣數據點數）的特征信號s（k）和白噪聲sn（k）組成，即

式中 ?fr為齒輪轉頻，fm為齒輪嚙合頻率，N為沖擊特征個數，An為第n個沖擊特征的最大幅值，[tn ，tn+1）表示第n個沖擊特征的持續(xù)時間范圍。

假設一對嚙合齒輪中發(fā)生局部損傷齒輪的齒數z=20，轉頻fr=25 Hz，嚙合頻率fm=zfr=500 Hz。由于齒輪發(fā)生局部損傷的故障特征頻率與其轉頻fr一致，且考慮諧波分量sh（t）中調制信號的初相位為π/2，經過計算可知，第1個沖擊特征出現的時間t1=0.03 s，持續(xù)時間（t2-t1）=1/fr=0.04 s。沖擊特征的最大幅值序列An由1.5?1.7的均勻分布隨機序列產生。設定信號采樣頻率fs=10 kHz，采樣時間t=12 s，疊加高斯白噪聲，使信噪比為-8 dB，由此產生長度L=120000的仿真振動信號y（k）（k=1，2，…，L），如圖3（a）所示。

此處期望從y（k）中選取的參考子信號包含5個周期的沖擊特征，則參考子信號長度至少為Np×（1/fr）×fs=2000。然后在y（k）中選擇一個時間點k=58930，對應時間t=5.8930 s，以此開始截取一段長度為2000的參考子信號，如圖3（a）中間框選的信號段所示，參考子信號的時域波形如圖3（b）所示，其中的沖擊特征完全被噪聲所淹沒。參考子信號對應的無噪波形如圖3（c）所示。

令參考子信號從仿真振動信號y（k）的初始時間點k=1，即初始時間t=0開始，如圖3（a）左側框選位置所示，沿時間軸逐點滑動至終了時間點k=118001，即滑動終了時間t=11.8001 s，如圖3（a）右側框選位置所示。參考子信號滑動至每個時間點處遮掩一段等長度的子信號，計算參考子信號與各個遮掩子信號之間的相關系數。將相關系數從小到大排列，得到重排的相關系數序列ρ ?_sm （k_t）（kt=1，2，…，118001），舍棄其負值及末尾Nr=10個取值后的曲線如圖3（d）所示。

需要說明的是，ρ_sm^*的選取帶有一定主觀性，選取的ρ_sm^*只是相對最優(yōu)，并不是唯一，但只要ρ_sm^*在一定范圍內，如在本仿真實例中，只需滿足0.10≤ρ_sm^*≤0.15，就可以使改進TSA取得良好的降噪效果。R（r）關于ρth的關系曲線則具有良好的指示意義，借助它可以方便地選取出合適的最優(yōu)閾值ρ_sm^*，用于篩選符合要求的總體同步信號，實現改進TSA降噪。改進TSA降噪后目標信號的平方包絡譜如圖3（g）所示，可見故障特征頻率25 Hz被順利提取出來。

對比圖3（b）和3（f）可見，仿真振動信號中的噪聲得到了較好的抑制。但對比圖3（c）和3（f）可見，降噪后的目標信號幅值存在一定的衰減，原因為在總體同步信號中，以參考子信號為參考，雖然超前的遮掩子信號與滯后的遮掩子信號之和無相位差，但兩者在幅值上存在一定的相互抵消，因此兩者求和之后再求平均，獲得的幅值要比參考子信號小。

為進一步驗證改進TSA降噪方法的有效性，在y（k）中再任選一個時間點k=12689，對應時間t=1.2689 s，以此開始截取一段長度同樣為2000的參考子信號，按前述相同步驟獲取降噪后的目標信號如圖4（a）所示，其平方包絡譜如圖4（b）所示。在y（k）中改選一個時間點k=102371，對應時間t=10.2371 s，按相同方法獲取的降噪后目標信號及其平方包絡譜分別如圖5（a）和5（b）所示。通過對比可以發(fā)現，參考子信號選取的時間點不同，對最終獲取目標信號的降噪效果沒有影響，只是目標信號的相位可能會發(fā)生改變，但這對目標信號平方包絡譜中故障特征的提取并無影響。

5 實際振動信號分析

5.1 齒輪故障振動信號分析

齒輪故障振動信號取自一臺風力發(fā)電機，其齒輪箱中一個小齒輪發(fā)生輪齒斷裂損傷故障，如圖6所示。該小齒輪的齒數z=32，額定轉速n=1770 r/min，轉頻fr=29.5 Hz，故障特征頻率fg=fr=29.5 Hz，嚙合頻率fm=z·fn=960 Hz，振動信號采樣頻率fs=97656 Hz，采樣得到的一段振動信號如圖7（a）所示。

由于該小齒輪的故障特征頻率為29.5 Hz，此處期望截取的參考子信號中至少包含6個故障沖擊特征，那么參考子信號的長度至少為6×（1/29.5）×fs=19531。因此，在采樣的振動信號中任選一個時間點k=50836，對應時間為t=0.5206 s，以此開始截取一段長度為19854的參考子信號，時域波形如圖7（b）所示，其中噪聲較大，故障特征并不顯著。令參考子信號從采樣振動信號初始時間開始，沿時間軸逐點滑動，計算參考子信號與各時間點處遮掩子信號之間的相關系數。滑動結束后，將相關系數從小到大排列，得到重排相關系數序列，舍棄其負值及末尾Nr=10個取值后的曲線如圖7（c）所示。

選擇最優(yōu)閾值ρ_sm^*的搜尋區(qū)間為[0.3ρ_sm^max， ρ_sm^max]，即[0.0997，0.3324]，設定迭代步數M=100，迭代計算后獲得均方根誤差R（r）關于閾值變量ρth的關系曲線，如圖7（d）所示。選取R（r）在低值平穩(wěn)變化過程中最小值處的閾值ρth=0.1928作為最優(yōu)閾值ρ_sm^*，然后選擇與參考子信號相關系數ρsm>ρ_sm^*=0.1928的所有遮掩子信號作為總體同步信號，平均計算后獲得降噪的目標信號，如圖7（e）所示，可見其中顯著的周期性故障沖擊特征。

對降噪后的目標信號進行平方包絡譜分析，如圖7（f）所示，可見小齒輪故障特征頻率fg及其倍頻（2?13）×fg被清晰地提取出來。目標信號的包絡譜如圖7（g）所示，故障特征頻率fg及其倍頻的譜線不如平方包絡譜清晰，因此，相比于包絡譜，平方包絡譜在振動信號分析中具有一定的優(yōu)勢。

作為比較，采用基于周期分段的傳統TSA降噪方法對采樣振動信號進行處理，獲得的降噪信號如圖8（a）所示，其平方包絡譜與包絡譜分別如圖8（b）和8（c）所示。從時域上講，由于齒輪故障損傷會引起轉速波動，振動信號的周期并非嚴格恒定，按周期分段截取的各段振動信號存在一定的相位誤差，且相位誤差具有累計效應，相比于改進TSA降噪方法，傳統TSA降噪方法獲取的目標信號降噪不夠徹底，存在較多的高頻分量，其中的周期性故障沖擊特征不太容易被識別，信號幅值也存在較大程度的衰減。從頻域上講，改進TSA降噪方法與傳統TSA降噪方法對于故障特征頻率的提取效果基本相同，且平方包絡譜的故障特征頻率及其倍頻譜線要比包絡譜更加清晰可辨。

5.2 軸承故障振動信號分析

軸承故障振動信號取自Case Western Reserve University軸承數據中心。選擇試驗臺驅動端的NTN滾動軸承為研究對象，其內滾道上有一直徑為0.7112 mm（0.028"），深度為1.27 mm（0.050"）的局部損傷故障。在試驗臺運行過程中，驅動電機功率為2.2 kW，轉速為1730 r/min，加載電機施加3 N·m載荷，軸承內圈的故障特征頻率fBPFI=123 Hz。振動數據采樣頻率為12 kHz，采樣到的一段振動信號如圖9（a）所示。

在采樣振動信號中任選一個時間點k=20185，對應時間為t=1.6821 s，以此開始截取一段長度為2000的參考子信號，如圖9（b）所示，其中故障沖擊特征并不顯著。令參考子信號從采樣的振動信號初始時間開始，沿時間軸逐點滑動，計算參考子信號與遮掩子信號之間的相關系數?；瑒咏Y束后獲得參考子信號與所有遮掩子信號之間的相關系數，將其從小到大排列后得到重排相關系數序列，舍棄其負值及末尾Nr=10個取值后的曲線如圖9（c）所示。

選擇最優(yōu)閾值ρ_sm^*的搜尋區(qū)間為[0?ρ_sm^max，ρ_sm^max]，即[0，0.5822]，設定迭代步數M=100，迭代計算后獲得均方根誤差R（r）關于閾值變量ρth的關系曲線，如圖9（d）所示。選取R（r）從低值平穩(wěn)變化向急劇增大變化轉變的一個臨界值ρth=0.4366作為最優(yōu)閾值ρ_sm^*，然后選擇與參考子信號相關系數ρsm>ρ_sm^*=0.4366的所有遮掩子信號作為總體同步信號，平均計算后獲得降噪的目標信號如圖9（e）所示，可見噪聲得到了較好的抑制，周期性特征非常顯著。降噪后目標信號的平方包絡譜如圖9（f）所示，明顯可見軸承內圈的故障特征頻率fBPFI及其倍頻（2?5）×f_BPFI。降噪后目標信號的包絡譜如圖9（g）所示，故障特征頻率fBPFI及其倍頻的譜線也基本清晰可辨，但不如平方包絡譜突出。

作為比較，采用傳統TSA降噪方法對軸承的采樣振動信號進行處理，獲得降噪的目標信號如圖10（a）所示。從該目標信號中能夠大致看出信號的周期性，但難以辨別出周期性的故障沖擊特征，原因在于軸承故障損傷引起轉速波動，導致采樣振動信號的周期并非嚴格恒定，按周期分段截取的各段時域同步信號在求和過程中存在較大的相位誤差累計效應。相比于改進TSA降噪方法，傳統TSA降噪方法獲取的目標信號降噪不夠徹底，存在較多的高頻成分以及故障特征頻率的0.5倍、1.5倍、2.5倍等分數倍頻成分，這從傳統TSA降噪目標信號的平方包絡譜（如圖10（b）所示）與包絡譜（如圖10（c）所示）中也可以看出來。雖然傳統TSA降噪結合平方包絡譜或包絡譜也能夠提取出故障特征頻率，但由于故障特征頻率分數倍頻成分的存在，對故障特征頻率識別有一定的干擾，影響故障判別。

6 結 ?論

（1）提出一種基于相關性檢測的改進TSA降噪方法，該方法在旋轉機械設備的采樣振動信號中任選一段參考子信號，令參考子信號在采樣振動信號中逐點滑動，計算參考子信號與遮掩子信號之間的相關系數;然后利用步長迭代算法獲取相關系數的最優(yōu)閾值;最后選擇與參考子信號相關系數大于最優(yōu)閾值的所有遮掩子信號作為總體同步信號，平均計算后可以獲得降噪效果良好的目標信號。

（2）改進TSA降噪方法成功應用于齒輪與軸承故障振動信號的處理，降噪后得到的目標信號經過平方包絡譜分析，能夠有效提取出齒輪與軸承的故障特征頻率，從而實現旋轉機械設備故障診斷。

（3）改進TSA降噪方法對旋轉機械設備轉速的穩(wěn)定性以及振動信號的平穩(wěn)性沒有嚴格要求，無需獲取振動信號的準確周期，避免了振動信號按周期分段而導致的相位誤差累積效應;無需轉速傳感器，硬件成本低、應用場合廣，具有較好的適用性。

參考文獻：

[1] 熊國良，胡俊鋒，陳慧，等. 基于SK-NLM包絡的滾動軸承故障沖擊特征增強[J]. 儀器儀表學報，2016， 37（10）：2176-2184.

Xiong Guoliang， Hu Junfeng， Chen Hui， et al. Rolling bearing fault impact feature enhancement based on spectral kurtosis and non-local means[J]. Chinese Journal of Scientific Instrument， 2016， 37（10）： 2176-2184.

[2] Ahamed N， Pandya Y， Parey A. Spur gear tooth root crack detection using time synchronous averaging under fluctuating speed[J]. Measurement， 2014， 52： 1-11.

[3] Ha J M， Youn B D， Oh H， et al. Autocorrelation-based time synchronous averaging for condition monitoring of planetary gearboxes in wind turbines[J]. Mechanical Systems and Signal Processing， 2015， 70-71： 161-175.

[4] Li H， Xu F， Liu H， et al. Incipient fault information determination for rolling element bearing based on averaging reassigned wavelet scalogram[J]. Measurement， 2015， 65： 1-10.

[5] Roy S K， Mohanty A R， Kumar C S. Fault detection in a multistage gearbox by time synchronous averaging of the instantaneous angular speed[J]. Journal of Vibration and Control， 2016， 22（2）： 468-480.

[6] Camerini V， Coppotelli G， Bendisch S， et al. Impact of pulse time uncertainty on synchronous average： Statistical analysis and relevance to rotating machinery diagnosis[J]. Mechanical Systems and Signal Processing， 2019， 129： 308-336.

[7] 龍海軍，孫燦飛，莫固良. 直升機振動檢測通用算法的研究與實現[J]. 振動、測試與診斷，2016，36（3）：524-528.

Long Haijun， Sun Canfei， Mo Guliang. Research and implementation of general algorithm for helicopter vibration detection[J]. Journal of Vibration， Measurement & Diagnosis， 2016， 36（3）：524-528.

[8] 邵毅敏，周曉君，歐家福. 變周期信號的時域同步平均新算法[J]. 振動工程學報，2009，22（5）：497-502.

Shao Yimin， Zhou Xiaojun， Ou Jiafu. A algorithm of time domain averaging for fluctuant-period signal[J]. Journal of Vibration Engineering， 2009， 22（5）：497-502.

[9] 沈國際，陶利民，徐永成. 時域同步平均的相位誤差累積效應研究[J]. 振動工程學報，2007，20（4）：335-339.

Shen Guoji， Tao Limin， Xu Yongcheng. Research on phase error accumulation effect of time synchronous averaging[J]. Journal of Vibration Engineering， 2007， 20（4）：335-339.

[10] 羅德揚. 時域同步平均原理與應用[J]. 振動、測試與診斷，1999，19（3）：202-207.

Luo Deyang. Principals and applications of time domain synchronous averaging[J]. Journal of Vibration， Measurement & Diagnosis， 1999， 19（3）： 202-207.

[11] 呂琛，宋希庚，鄒積斌. 基于DSP的振動信號階比與時域同步平均分析[J]. 振動與沖擊，2002，21（2）：53-57.

L? Chen， SONG Xigeng， ZOU Jibin. DSP-based order domain and time domain synchronous averaging analysis of vibration signal[J]. Journal of Vibration and Shock， 2002， 21（2）： 53-57.

[12] Combet F， Gelman L. An automated methodology for performing time synchronous averaging of a gearbox signal without speed sensor[J]. Mechanical Systems and Signal Processing， 2007， 21（6）： 2590-2606.

[13] 何正嘉，劉雄，屈梁生. 信號時域平均原理和應用[J]. 信號處理，1986，2（4）：46-53.

He Zhengjia， Liu Xiong， Qu Liangsheng. The principle and application of signal time domain average[J]. Signal Processing， 1986， 2（4）： 46-53.

[14] 郭遠晶，魏燕定，周曉軍. 基于STFT時頻譜系數收縮的信號降噪方法[J]. 振動、測試與診斷，2015，35（6）：1090-1096.

Guo Yuanjing， Wei Yanding， Zhou Xiaojun. Signal denoising method based on STFT time-frequency spectrum coefficients shrinkage[J]. Journal of Vibration， Measurement & Diagnosis， 2015， 35（6）： 1090-1096.

[15] Antoni J. Cyclic spectral analysis of rolling-element bearing signals： Facts and fictions[J]. Journal of Sound and Vibration， 2007， 304（3-5）： 497-529.

[16] Borghesani P， Pennacchi P， Chatterton S. The relationship between kurtosis- and envelope-based indexes for the diagnostic of rolling element bearings[J]. Mechanical Systems and Signal Processing， 2014， 43（1-2）： 25-43.

Fault feature extraction based on improved TSA denoising and squared envelope spectrum

GUO Yuan-jing1， JIN Xiao-hang2， WEI Yan-ding3， YANG You-dong1

（1. Zhijiang College， Zhejiang University of Technology， Shaoxing 312030， China;

2. College of Mechanical Engineering， Zhejiang University of Technology， Hangzhou 310023， China;

3. Key Laboratory of Advanced Manufacturing Technology of Zhejiang Province， Zhejiang University， Hangzhou 310027， China）

Abstract： In order to achieve the fault diagnosis of the rotating machine equipment， a method based on improved time synchronous averaging （TSA） denoising using correlation detection and fault features extracting using squared envelope spectrum is proposed. A reference sub-signal is arbitrarily taken from the sampled vibration signal. The reference sub-signal starts from the initial time of the vibration signal and moves point by point along the time axis to the end time. Meanwhile， the correlation coefficient between the reference sub-signal and its masking sub-signal in the vibration signal is calculated. The correlation coefficient optimal threshold are obtained using the step iterative algorithm. All the masking sub-signals whose correlation coefficients are greater than the optimal threshold are selected as overall synchronous signals whose arithmetic average is calculated for the denoised target signal. The target signal analysis is achieved using squared envelope spectrum. The analysis results of the simulated signal， gear and rolling bearing fault vibration signals have shown that the proposed method is an efficient implementation for rotating machine vibration signal denoising and fault feature extraction.

Key words： fault diagnosis; time synchronous averaging; correlation coefficient; step iterative algorithm; squared envelope spectrum

作者簡介： 郭遠晶（1987-），男，博士，講師。電話：（0575）81112562;E-mail：gyjyn@126.com