用加窗插值法抑制振動信號工頻噪聲

張曉濤，李偉光

（1.廣東機電職業技術學院汽車學院, 廣州 510515；2.華南理工大學機械與汽車工程學院, 廣州 510640）

可傾瓦軸承因各油楔的壓力總是指向轉軸中心而形成最佳油膜，具有能承受各方向徑向載荷、高速旋轉穩定性好、承載能力大等優點，被廣泛用于大型汽輪機、燃氣輪機、風機、泵等旋轉機械。采用電渦流位移傳感器采集可傾瓦軸承的振動信號，耦合到振動信號的工頻噪聲對振動特征分析產生影響，因而需要研究有效的工頻噪聲抑制方法。

工頻噪聲抑制是重要研究課題。工頻噪聲實際上是來自電力系統的工頻基波及其諧波，工頻基波不是恒定不變的，而是在50 Hz 上下有小幅波動。為抑制工頻噪聲國內外學者提出了諸多方法，按基本原理主要可分為六類[1]，基于濾波器理論的方法[2]、基于盲源分離的方法[3]、基于神經網絡的方法[4]、基于參數擬合的方法[5]、基于信號分解的方法[6]、基于奇異值分解的方法[7]等。這些方法在特定條件下效果明顯，但又存在一定的局限性，如基于濾波器理論的方法要求工頻頻率為恒定的50 Hz，顯然這與實際情況不符，因而該方法既抑制了工頻噪聲同時也會損害了有用信號。

工頻噪聲來自電力系統，可將電力信號的分析方法擴展到振動信號領域，用于振動信號的工頻噪聲抑制，其中加窗插值法[8]是常見的電力信號分析方法。電力信號主要包含基波、二次諧波、三次諧波等高次諧波及間諧波，加窗插值法用于計算各次諧波的參數（幅值、頻率和相位）。通常電力信號的高次諧波比較微弱，故可以忽略。電力信號各成分的頻譜間距比較大，而受工頻噪聲干擾的振動信號的成分則復雜得多，振動信號中既包含有用信號成分也包含工頻噪聲，有用信號成分頻譜與工頻噪聲頻譜的間距不確定，甚至會重疊在一起。加窗插值法應用時既要考慮窗函數特性也要考慮被分析信號特點[9]。因信號特征不同，用于電力信號的加窗插值法不能直接照搬到振動信號，加窗插值法用于振動信號分析雖具有可行性，但相關研究報道并不多見。張營等[6]用譜插值法消除工頻成分，譜插值法就是將時域插值方法應用于頻域分析，具體做法是先對信號作快速傅里葉變換（Fast Fourier Transform，FFT），并將對應于工頻50 Hz 的單一譜線去掉，用相鄰的兩個頻譜值對其作線性內插以得到新值，并用新值替換原值，最后作逆傅里葉變換（Ⅰnverse Fast Fourier Transform，ⅠFFT），得到消除工頻噪聲的信號。但筆者認為該做法值得商榷，因為將工頻噪聲簡化為一根譜線是不準確的。

本文先仿真分析信號加不同窗函數時的工頻插值計算結果，進而對實測振動信號進行加窗插值分析，最終選用合適的窗函數以打我的最準確的插值計算結果，即使得工頻噪聲的抑制效果最好。

1 加窗插值法分析

所謂加窗插值法就是先對信號加窗函數，再對相鄰的若干條頻譜線應用插值公式，計算出信號的頻率、幅值和相位。按照所用譜線數的不同，插值方法可分為雙譜線插值[8]、三譜線插值[10]、四譜線插值[11]等，本文采用雙譜線插值。加窗插值法用于振動信號分析，首先要選取合適的窗函數。窗函數種類較多，大體可分為冪窗、三角函數窗和指數窗三類，每一類又含有多種窗函數。各窗函數的主要區別是主瓣寬度、旁瓣峰值和旁瓣衰減速度，窗函數的選擇就是依據這三項指標。矩形窗是最簡單形式的冪窗，通常不加窗就是加了矩形窗，而漢寧窗和六項余弦窗都是余弦窗。矩形窗、漢寧窗和六項余弦窗三種窗函數的項數和系數不同[12]，使得各窗函數的頻譜特征差別很大，如圖1 所示，矩形窗主瓣窄、旁瓣多且旁瓣衰減慢，六項余弦窗主瓣寬、旁瓣窄且旁瓣衰減快，而漢寧窗的主瓣及旁瓣特征位于矩形窗和六項余弦窗之間。矩形窗、漢寧窗和六項余弦窗的主瓣寬度、旁瓣峰值、旁瓣衰減速度的具體數值分dB/octave；24 π/N、-88 dB、-66 dB/octave。三種窗函數具有特征完全不同的主瓣和旁瓣，因而能夠比較全面地反映窗函數這一類對象。

周期信號加窗函數并作FFT 后，可由主瓣上相鄰的兩個譜峰經插值公式[13-14]計算出其參數（幅值、頻率、相位）。若周期信號用x0（n）（n=0,1,···,N-1）表示，采樣頻率用fs表示，采樣間隔用T表示，則。用窗函數w（n）對x0（n）作N點截斷，得離散信號序列x（n）：

對x（n）作FFT，分別用hm、fm、phase（hm）表示主瓣上左側譜線的幅值、頻率、相位，分別用hm+1和fm+1表示右側譜線的幅值、頻率。待計算的頻率、幅值、相位分別用f、A、θ表示。

當信號加矩形窗時，令插值系數δ為：

插值公式為：

當信號加漢寧窗時，令插值系數δ為：

插值公式為：

當信號加六項余弦窗時，令插值系數δ為：

插值公式為：

應用加窗插值抑制法，通過上述插值公式可計算工頻噪聲參數的精確值，從而將其從信號中消除，達到工頻噪聲抑制的目的。

2 數值仿真試驗分析

可傾瓦軸承振動試驗需要在不同的轉速下進行，經振動傳感器采集的振動信號的特征頻率也隨之不同，而耦合到振動信號中的工頻噪聲的頻率則基本穩定。應用數值仿真方法對此工況進行模擬，并針對仿真信號研究合適的窗函數，以使經插值公式計算得到的工頻參數最準確。雖然在實際試驗過程中振動信號的振動分量和工頻噪聲分量的參數都是變化的，但此處為便于仿真分析，將工頻噪聲分量的參數設為某一定值，振動分量的幅值和相位也設為定值，而只改變振動分量的頻率。

仿真信號用x（t）表示，如公式（8）所示，包括表征工頻噪聲的分量x1（t），表征振動特征的分量x2（t），及高斯白噪聲ξ（t），信噪比為30 dB。

設x1（t）分量參數為：α1=2，f1=50.2 Hz，?1=π/18。設x2（t）分量的幅值和相位為：α2=1，?2=7π/9，每次試驗時頻率f2取不同的值，f2值取為：30、33、35、38、40、41、42、43、44、45、46、47、48、49、51、52、53、54、55、56、57、58、59、60、62、65、68、70，共28 個值，單位為Hz。當f2=30 Hz時，仿真信號x（t）的波形如圖2所示。

圖2 仿真信號波形（當f2 為30 Hz時）

按以上方法可生成28個仿真信號，每個仿真信號分別加矩形窗、漢寧窗、六項余弦窗，再用插值公式（3）、式（5）、式（7）計算x1（t）的參數（幅值、頻率、相位）。當圖2 所示仿真信號，加矩形窗時計算的幅值、頻率、相位分別為：2.001 296、50.502 924、9.732 154，誤差分別為：0.064 814%、0.005 791%、2.678 465%；加漢寧窗時計算的幅值、頻率、相位分別為：2.002 086、50.500 167、10.155 155，誤差分別為：0.104 282%、0.000 331%、1.551 547%；加六項余弦窗時計算的幅值、頻率、相位分別為：1.997 974、50.503 487、9.528 856，誤差分別為：0.101 276 %、0.006 906%、4.711 439%。無論加哪種窗函數，由加窗插值法計算的工頻參數的誤差都很小，參數值非常準確。進一步地，將28 個仿真信號的工頻分量x1（t）的幅值誤差繪制為圖3（a），頻率誤差繪制為圖3（b），將相位誤差繪制為圖3（c）。由圖3 可見，在48 Hz～53 Hz頻率范圍內，信號加矩形窗時x1（t）參數的誤差最小，加六項余弦窗時x1（t）參數的誤差最大，而加漢寧窗時的x1（t）參數的誤差則位于兩者之間。

圖3 加不同窗函數時工頻參數誤差

圖3 中當超出48 Hz～53 Hz 頻率范圍時，因為加不同窗函數時的x1（t）參數的誤差值都比較小，三者完全重疊在一起而難以分辨。進一步地，將縱軸數值最低的部分區域放大，圖3（a）放大為圖4（a），圖3（b）放大為圖4（b），圖3（c）放大為圖4（c）。由圖4可更清晰地區分開加不同窗函數時x1（t）參數的計算誤差。在圖4（a）中，信號分別加三個窗函數時幅值誤差都有較大波動，三者間沒有明顯區分；在圖4（b）中，信號加矩形窗時頻率誤差最大，加六項余弦窗時頻率誤差次之，而加漢寧窗時頻率誤差最小；在圖4（c）中，信號加矩形窗和六項余弦窗時相位誤差波動都比較大，而加漢寧窗時相位誤差波動最小。就綜合性能而言，信號加漢寧窗時的工頻參數計算誤差最優，顯然此時將工頻噪聲刪除所導致的誤差也最小。

圖4 加不同窗函數時的工頻參數誤差（放大）

數值仿真試驗結論：當振動信號的特征頻率位于48 Hz～53 Hz 范圍時，信號加矩形窗時工頻噪聲的抑制效果最好；當振動信號的特征頻率超出此范圍時，加漢寧窗時工頻噪聲的抑制效果最好。

3 實測振動信號分析

可傾瓦軸承性能優良、應用廣泛，為開發特殊應用需求的新型可傾瓦軸承需開展相關振動試驗，以使所開發的可傾瓦軸承具有最優的減振性能。前期項目已為此搭建了可傾瓦軸承轉子試驗系統[7]，該系統示意圖如圖5 所示，系統包括基座、伺服電機、傳動裝置、轉子、軸承座、轉子兩端推力軸承、稀油站、測試系統、控制系統等。在兩個軸承座的內側固定四個支架，安裝四個Kaman 電渦流位移傳感器進行振動量監測，傳感器型號為KD2306，精度為0.1 μm，頻響速度為50 kHz，安裝現場如圖6所示，圖中圓圈內的就是電渦流位移傳感器。

圖5 可傾瓦軸承轉子試驗系統示意圖

圖6 電渦流位移傳感器安裝現場

試驗中發現有電網工頻噪聲耦合到振動信號，工頻噪聲對振動信號的特征分析造成干擾，本文擬采用加窗插值法加以抑制。當可傾瓦軸承轉子系統工作在轉速3 900轉/分時，即轉頻為65 Hz時系統工作于高速下，此時是特別需要關注的工況。若此工況下的振動信號用x（t）表示，x（t）時域波形及頻譜如圖7所示，圖7（b）中工頻噪聲已被標出。轉子工作頻率為65 Hz，與工頻噪聲頻率50 Hz 之間相差15 Hz，根據前面仿真試驗結論，振動信號加漢寧窗時工頻噪聲的抑制效果最好，為進一步驗證此結論，這里同時給出了加三種窗函數時的計算結果。

圖7 振動信號波形及頻譜

信號x（t）分別加矩形窗、漢寧窗、六項余弦窗，并分別用插值公式（3）、式（5）、式（7）計算出工頻噪聲的參數，將計算結果填入表1。顯然僅從表1中的數值無法得知加哪個窗函數更合適。

若表1 中的幅值、頻率和相位分別用字母A、F和P表示，則工頻噪聲可用式（9）表示，即：

表1 振動信號工頻噪聲參數

由x（t）減x1（t），可得工頻噪聲被抑制后的振動信號x0（t），再對x0（t）作FFT，生成如圖8 所示的頻譜圖。分別將圖8（a）和圖8（b）與圖7（b）相比，發現工頻噪聲被有效地抑制了，而將圖8（c）與圖7（b）相比，發現工頻噪聲并沒有被有效抑制，這說明加矩形窗和漢寧窗后用插值公式計算的工頻參數值比較準確，而加六項余弦窗的插值計算結果最差。

圖8 加不同窗函數時工頻噪聲被抑制的振動信號

工頻噪聲的頻率雖然存在一定波動，但波動范圍很小，因為工頻噪聲來源于電網電源，而電網頻率是比較穩定的。如果多次重復采集振動信號進行試驗，由加窗插值法計算的頻率參數理應比較穩定?；诖艘幝?，在相同的試驗工況和采樣率條件下另外截取七段長度相同的信號，分別加矩形窗、漢寧窗、六項余弦窗，再用插值公式（3）、式（5）、式（7）計算工頻參數，并將8次的試驗結果匯總到圖9。

由圖9（a）可見，雖然加三種窗函數時幅值參數的波動都比較大，而其中加漢寧窗時的幅值參數的波動略小，但三者變化趨勢基本一致，出現幅值波動的原因可能主要是受現場工作環境及電網負荷的影響；由圖9（b）可見，加六項余弦窗時頻率參數波動最大，加漢寧窗和矩形窗時的頻率參數則比較穩定，其中加漢寧窗時的頻率參數的波動更??；由圖9（c）可見，加六項余弦窗時的相位參數的波動較大，而加漢寧窗、矩形窗的相位值雖然有波動，但二者趨勢基本一致，由于每次試驗的信號都是在不同時刻所采集的，因而存在相位波動也是正常的。

圖9 加窗插值法所提取的工頻參數

分析結果表明：加窗插值法用于實測振動信號的工頻噪聲抑制具有可行性。雖然文中是以轉速3 900 轉/分的工況為例來分析的，但同樣也可以分析其他工況，這里受篇幅所限而不再贅述。

4 與陷波器方法對比

陷波器是常用的工頻噪聲抑制方法，由數字濾波器理論[15]可設計50 Hz陷波器，其計算式為：

ω=2 πf/fs=0.306 6，f=50 Hz，fs為采樣率。α值確定陷波器的寬度和深度，其值越大則凹陷越深、寬度越窄。為得到最佳濾波效果，該值取為0.98。設計的陷波器頻譜如圖10所示。

圖10 陷波器頻譜

振動信號經陷波器濾波后的結果如圖11所示，將其與圖8（b）相比，發現工頻噪聲被濾掉，同時與工頻噪聲相鄰的部分有用信號也被濾掉，從而對后續的振動特征的量化分析產生影響。這表明加窗插值法同陷波器法相比具有明顯優勢。

圖11 陷波器濾波后振動信號頻譜

5 結語

為將加窗插值法應用于振動信號的工頻噪聲抑制，本文選取矩形窗、漢寧窗和六項余弦窗等三個典型的窗函數，并基于仿真信號和實測可傾瓦軸承振動信號對比了加不同窗函數時的工頻插值計算結果，當振動信號加漢寧窗時工頻噪聲抑制效果最好，表明加窗插值法用于振動信號工頻噪聲抑制是可行的，能較準確地抑制工頻噪聲。加窗插值法從信號自身特征出發，應用精確計算公式抑制工頻噪聲，對有用信號損害小，同傳統陷波器法相比具有顯著優勢。

實際上對于理想無隨機噪聲的信號，當其加六項余弦窗時工頻噪聲的插值計算結果最準確，下一步將研究如何降低隨機噪聲以獲得更準確的工頻參數，進而使得抑制效果更好。窗函數眾多，如何選取更適合工頻噪聲抑制的窗函數也是后續需要研究的內容。