基于Welch功率譜的加速度積分改進方法研究

李展銓， 陳太聰,2,3

(1. 華南理工大學(xué) 土木與交通學(xué)院，廣州 510641； 2. 亞熱帶建筑科學(xué)國家重點實驗室，廣州 510641； 3. 琶洲實驗室，廣州 510335)

在工程振動實踐中，結(jié)構(gòu)的加速度響應(yīng)相對于位移和速度響應(yīng)更易測量，常用于結(jié)構(gòu)模態(tài)識別和損傷判別等分析[1-2]。而在結(jié)構(gòu)抗震性能評價、結(jié)構(gòu)阻尼評估等動力效應(yīng)評定中，位移和速度指標更具有關(guān)鍵意義[3-5]。對加速度信號進行時域積分獲得位移和速度是振動信號處理中的常用方法[6]，但由于儀器誤差和環(huán)境影響，測量得到的加速度信號往往帶有干擾噪聲，導(dǎo)致積分得到的速度和位移與真實響應(yīng)有較大差別，甚至完全失真。如何消除干擾噪聲以還原真實的響應(yīng)信號，一直是振動工程領(lǐng)域的重點關(guān)注問題。

不同于傳統(tǒng)的時域積分，徐慶華[7]嘗試頻域積分方法，通過快速傅立葉變換(fast Fourier transform，F(xiàn)FT)實現(xiàn)加速度、速度與位移間的相互轉(zhuǎn)換。Brandt等[8]通過對不同積分方法的對比研究，推薦帶低頻截止的頻域積分方法，參考加速度信號的第一階主頻，將低頻部分進行歸零處理。方新磊等[9]在低頻截止的基礎(chǔ)上引入高頻截止濾波，以同時消除高頻干擾噪聲和低頻趨勢項。以上頻率截止法對截止頻率參數(shù)較為敏感，而相關(guān)參數(shù)的設(shè)置有很大的主觀性。

陳太聰?shù)萚10]針對頻率截止法的局限性提出了有效頻段法，基于加速度信號的FFT譜曲線形態(tài)擬合，進行加速度信號的濾噪積分。該法對于中等噪聲干擾的加速度信號的積分結(jié)果精度較高，有效避免了頻率截止法的參數(shù)敏感性和主觀性等問題。但在高噪聲干擾下，積分精度隨著噪聲水平的增大而迅速變差。其次，該法的實施也需要人為給定主頻近似信息以判定分析頻率范圍，不利于工程應(yīng)用。

事實上，在信號處理領(lǐng)域，除了FFT譜曲線可以反映信號的頻域信息，功率譜曲線也具有類似的功能。傳統(tǒng)的周期圖法直接通過FFT譜的平方求解功率譜曲線，計算方便，但譜線起伏大，弱信號下的頻譜分辨率低[11]；Welch算法引入數(shù)據(jù)分段交疊和加窗函數(shù)的方法，對周期圖法進行了改進，得到修正功率譜，有效降低噪聲影響，所得譜線形狀平滑，易于清晰分辨主頻。如王福杰等[12]通過試驗對比各種功率譜后發(fā)現(xiàn)，應(yīng)用Welch功率譜估計，從寬帶噪聲中檢測出窄帶信號的效果最好。如今Welch功率譜已廣泛應(yīng)用于如脈動風(fēng)譜估計[13]、鐵路軌道不平順度測試[14]、地震監(jiān)測背景噪聲分析[15]等研究和工程實踐中。

針對上述有效頻段法的若干應(yīng)用問題，本文提出引入Welch功率譜的改進有效頻段法，用于加速度信號濾噪積分的自動處理。基于帶噪加速度信號的Welch功率譜，自動實現(xiàn)主頻的識別和分析頻率范圍的定義。繼而分別根據(jù)Welch功率譜曲線和Welch功率譜開方曲線進行形態(tài)擬合，實現(xiàn)有效頻段的自動識別。最后進行有效頻段內(nèi)的頻域積分，得到速度與位移信號。文后通過多頻激勵和隨機激勵下的數(shù)值模擬算例，開展原有效頻段法和本文改進方法的對比分析，以檢驗方法在不同激勵下的積分精度和抗噪性能。

1 有效頻段法的基本原理

本文方法是有效頻段法的改進研究，在此先介紹有效頻段法的基本計算原理。

1.1 分段頻譜曲線形態(tài)假設(shè)

假設(shè)按等時距Δt測量得到的帶噪加速度信號為X(n),n=0,1,2,…,N-1，則根據(jù)離散傅里葉變化的基本原理，離散的時域加速度可轉(zhuǎn)換為頻域內(nèi)若干離散的FFT頻譜信號H(k),k=0,1,2,…,N-1，其相互轉(zhuǎn)化關(guān)系為[16]

(1)

(2)

式中：N為采樣數(shù)據(jù)量；n為時刻點；k為譜線序列點。

設(shè)加速度信號含m個峰值主頻，按從小到大的順序依次為f1,f2,…,fm，則整個頻率范圍按式(3)方式分為用于分析的m段

(3)

式中，fmax為頻譜曲線對應(yīng)的最大頻率。

對m個頻段均進行歸一化處理，即都歸一化為范圍(0,1]。假設(shè)每個頻段內(nèi)的頻譜曲線都符合相應(yīng)的高斯函數(shù)分布[17]

Gi(x)=aie-(x-bi)2/(2c2)

(4)

式中，i=1,2,…,m；x∈(0,1]為每個頻段內(nèi)的歸一化頻率；參數(shù)ai，bi和ci分別為高斯函數(shù)的幅值、中心位置和半徑。

基于式(4)定義的頻譜分布形態(tài)，相應(yīng)的每一頻段內(nèi)的頻譜累計能量分布函數(shù)可由式(5)計算得到

(5)

式中： erf(·)為誤差函數(shù)；Ei(x)在形態(tài)上表現(xiàn)為反Z型的單調(diào)遞增函數(shù)。

對累計能量分布函數(shù)進一步進行歸一化處理，得

(6)

1.2 分段累積信號能量計算

在第i個歸一化頻段內(nèi)，離散的加速度頻譜信號的累積能量由式(7)計算得到

(7)

對信號累計能量函數(shù)進一步做歸一化處理，得

(8)

1.3 參數(shù)擬合與有效頻段確定

基于式(8)計算得到的離散數(shù)據(jù)，用式(6)定義的非線性連續(xù)函數(shù)進行非線性擬合，迭代收斂后得到參數(shù)bi和ci。

根據(jù)概率統(tǒng)計理論的三倍標準差原則[18]，第i階主頻能量分布的歸一化有效頻段取為[bi-3ci,bi+3ci]，相應(yīng)的實際有效頻段為

(9)

最后經(jīng)過有效頻段截斷后得到的加速度、速度和位移分別為[19]

(10)

式中，ωk為加速度信號的頻率，ωk=2πkFs/N。

(11)

有效頻段法用于加速度的濾噪和積分的相關(guān)計算流程，如圖1所示。

圖1 有效頻段法的積分流程圖Fig.1 Integration flowchart of frequency truncation method

2 基于Welch功率譜的有效頻段法

已有研究表明，有效頻段法在處理低于25%噪聲水平的加速度信號積分時有較好的效果，但在更高噪聲干擾下，由于主頻處的FFT譜混入了較大的噪聲，噪聲能量造成了擬合高斯函數(shù)的半徑變大，導(dǎo)致有效頻段的范圍擴大，從而引入更多噪聲影響，積分精度難以保證。此外，有效頻段法需要人為指定主頻的近似取值，不能實現(xiàn)加速度積分的全程自動分析，也不便于實際工程應(yīng)用。

針對有效頻段法的以上問題，本文引入Welch功率譜對原方法進行改進，以實現(xiàn)自動化和高耐噪的加速度積分分析。一方面，Welch譜曲線的形態(tài)平滑，峰值清晰，易于自動識別主頻；另一方面，它作為一種平均功率譜，本身具有良好的抑噪功能，有利于消除高噪聲對擬合高斯函數(shù)參數(shù)的影響。

如圖2為文后算例中多頻激勵情況下的帶噪和不帶噪加速度信號的FFT譜曲線和Welch譜曲線對比，由圖2可知：①與FFT譜相比，Welch譜在信號主頻處的主瓣較寬、旁瓣較低、起伏性小、曲線平滑，在形態(tài)上更加符合高斯曲線的分布；②與FFT譜不同，帶噪Welch譜與不帶噪Welch譜的差別不大，且頻率影響范圍基本一致，表明Welch函數(shù)處理具有一定的消噪功能。

圖2 35%噪聲水平信號時程曲線以及其Welch功率譜圖與FFT頻譜圖Fig.2 Time-history curve, Welch power spectra and FFT spectra of 35% noise level signals

2.1 基于Welch譜的自動選峰和頻域分段

設(shè)實測加速度信號X(n)，對應(yīng)FFT頻譜為H(k)，對應(yīng)Welch功率譜為P(k)。取P(k)max×5%作為識別主頻峰值的下限，得到m階主頻，按從小到大的順序依次為f1,f2,…,fm，并選擇各個主頻間的最低波谷頻率作為頻域分段的依據(jù)，則有

(12)

2.2 分段頻譜曲線形態(tài)擬合假設(shè)

根據(jù)不同類型激勵下的結(jié)構(gòu)振動響應(yīng)特點，本文考察3種不同類型的加速度頻譜曲線：①FFT頻譜曲線，這也是原有效頻段法所采用的頻譜曲線；②Welch功率譜曲線；③Welch功率譜開方曲線，由于功率譜在數(shù)值上是幅值的平方關(guān)系，因此功率譜開方在數(shù)值上與FFT譜具有對應(yīng)關(guān)系。

分別假設(shè)這3種頻譜曲線符合相應(yīng)的高斯曲線分布，即可采用如式(4)所示的帶不同參數(shù)的高斯函數(shù)進行描述。相應(yīng)地，擬合高斯函數(shù)的累計能量分布函數(shù)及其歸一化函數(shù)分別由式(5)和式(6)計算。

2.3 分段累積信號能量計算

對于實測的帶噪加速度信號，在第i個歸一化頻段內(nèi)，離散的加速度頻譜信號的累積能量，根據(jù)3種不同類型的頻譜曲線，分別按式(13)計算得到

(13)

其相應(yīng)的歸一化函數(shù)分別根據(jù)式(8)進行計算。

2.4 參數(shù)擬合與有效頻段確定

基于式(13)和式(8)計算得到的離散數(shù)據(jù)，用式(6)定義的非線性連續(xù)函數(shù)進行非線性擬合，迭代收斂后得到參數(shù)bi和ci。繼而由式(9)確定各有效頻段。最后采用式(10)積分計算得到濾噪后的加速度、速度和位移信號。

最終，改進方法用于加速度的濾噪和積分的相關(guān)計算流程，如圖3所示。由圖3可知，本文改進方法的實現(xiàn)過程與有效頻段法大致相同，不同的是利用了Welch功率譜自動進行頻率分段，并針對基于Welch譜的兩種譜曲線進行高斯函數(shù)形態(tài)擬合，擬合得到的高斯函數(shù)參數(shù)不會因噪聲水平而產(chǎn)生大的偏差，保證有效頻段的精度，避免高噪聲干擾下的積分效果失真情況。

圖3 基于Welch功率譜形態(tài)擬合的積分流程圖Fig.3 Integration flowchart based on morphology fitting of Welch power spectrum

本文方法在MATLAB軟件中進行了實現(xiàn)，相關(guān)技術(shù)開發(fā)應(yīng)用了軟件內(nèi)嵌的傅里葉變換函數(shù)、傅里葉逆變換函數(shù)(inverse fast Fourier transform，IFFT)、Welch功率譜生成函數(shù)和尋峰函數(shù)等若干函數(shù)命令。

3 數(shù)值算例

為了驗證方法效果，本文采用與原有效頻段法一致的數(shù)值算例進行對比分析。

六自由度懸臂梁結(jié)構(gòu)，考慮豎向振動，如圖4所示。其各階模態(tài)阻尼比為0.01，質(zhì)量陣和剛度陣為

圖4 六自由度懸臂梁結(jié)構(gòu)模型Fig.4 Structural model of a 6-dofs cantilever beam

在結(jié)構(gòu)自由端分別施加豎向多頻激勵和白噪聲隨機激勵，其中多頻激勵取為

z(t)=200sin 6πt+100sin 20πt+100sin 40πt

然后按照Newmark-β法(γ=0.5，β=0.25)求出加速度、速度和位移響應(yīng)時程，作為精確參考解。在加速度解中分別加入1%，5%，10%，15%，20%，25%，30%，35%和40%的白噪聲干擾作為觀測加速度。再分別采用原有效頻段法和本文改進方法得到修正的加速度、速度和位移信號。其中，改進方法采用兩種頻譜曲線進行高斯函數(shù)形態(tài)擬合。最后采用如下定義的總體誤差指標來評價積分精度

(14)

3.1 多頻激勵情況

在多頻激勵下，以質(zhì)點3處的豎向響應(yīng)為例，不同噪聲水平下兩種方法的結(jié)果精度對比，如圖5所示。35%噪聲水平下兩種方法濾噪積分后所得的位移時程對比，如圖6所示。

圖5 不同噪聲水平下的結(jié)果誤差Fig.5 Result errors under different noise levels

圖6 35%噪聲下的位移時程Fig.6 Time history of displacement under 35% noise level

由圖5～圖6可知，在多頻激勵作用下：

(1) 基于Welch功率譜曲線擬合的改進方法計算得到的位移結(jié)果，在40%噪聲水平以下均比原方法的精度更高，具有良好的抗噪穩(wěn)定性，在高噪聲情況下仍能保持結(jié)果誤差低于3%。

(2) 基于Welch功率譜開方曲線擬合的改進方法計算得到的位移結(jié)果，在15%噪聲水平以下比原方法的精度更高，在15%～25%噪聲水平間與原方法相近，在25%噪聲水平以上具有一定的抗噪穩(wěn)定性，沒有出現(xiàn)結(jié)果失真。

綜合而言，在多頻激勵的情況下，基于兩種頻譜曲線的改進方法均比原有效頻段法的積分精度更高，具有更好的抗噪穩(wěn)定性。其中，基于Welch功率譜曲線擬合的改進效果最好。其原因是Welch功率譜曲線在數(shù)值上可表示為FFT頻譜曲線和Welch功率譜開方曲線的平方關(guān)系，從而放大各主頻處的真實信號與旁瓣噪聲信號之間的差異，減少了噪聲信號對形態(tài)擬合參數(shù)的影響，避免了高噪聲下有效頻段過寬的問題。因此，針對多頻激勵情況，本文建議使用基于Welch功率譜曲線擬合的改進方法。

3.2 隨機激勵情況

在白噪聲隨機激勵下，同樣以質(zhì)點3處的豎向響應(yīng)為例。不同噪聲水平下兩種方法的結(jié)果精度對比，如圖7所示。35%噪聲水平下兩種方法濾噪積分后所得的位移時程對比，如圖8所示。

圖7 不同噪聲水平下的結(jié)果誤差Fig.7 Result errors under different noise levels

由圖7～圖8可知，在白噪聲隨機激勵作用下：

圖8 35%噪聲下的位移時程Fig.8 Time history of displacement under 35% noise level

(1) 基于Welch功率譜曲線擬合的改進方法計算得到的位移結(jié)果，誤差穩(wěn)定于10%左右，在低噪聲水平時略差于原方法，在高噪聲水平時具有一定的抗噪穩(wěn)定性，沒有出現(xiàn)結(jié)果失真。

(2) 基于Welch功率譜開方曲線擬合的改進方法計算得到的位移結(jié)果，在40%噪聲水平以下均比原方法的精度更高，具有良好的抗噪穩(wěn)定性，在高噪聲情況下仍能保持結(jié)果誤差低于10%。

綜合而言，在隨機激勵的情況下，基于Welch功率譜開方曲線擬合的改進方法比基于Welch功率譜曲線擬合的改進方法的效果更好，也比原有效頻段法的積分精度更高且抗噪穩(wěn)定性更強。其原因是白噪聲作用下的響應(yīng)信號屬于旁瓣較大的信號，Welch功率譜開方曲線相對于Welch功率譜曲線，相對弱化了信號主瓣和旁瓣之間的差異，適當(dāng)保留了旁瓣頻率的貢獻，避免了高噪聲下有效頻段過窄的問題。因此，針對隨機激勵情況，本文建議使用基于Welch功率譜開方曲線擬合的改進方法。

4 結(jié) 論

針對工程振動測試中如何由帶噪加速度信號還原出真實速度和位移的問題，新近發(fā)展的有效頻段法通過對信號FFT頻譜曲線的形態(tài)擬合，實現(xiàn)了中等噪聲干擾下的高精度加速度積分。但方法仍存在高噪聲干擾下積分結(jié)果容易失真，以及需人為給定主頻近似信息等問題，不利于工程應(yīng)用。

本文引入Welch功率譜對有效頻段法進行改進，以實現(xiàn)自動化和高耐噪的加速度積分分析?；趲г爰铀俣刃盘柕腤elch功率譜，自動實現(xiàn)主頻的識別和分析頻率范圍的定義。繼而分別根據(jù)Welch功率譜曲線和Welch功率譜開方曲線進行形態(tài)擬合，實現(xiàn)有效頻段的自動識別。其中，針對多頻激勵情況，推薦使用基于Welch功率譜曲線擬合的改進方法；針對隨機激勵情況，推薦使用基于Welch功率譜開方曲線擬合的改進方法。

還需說明的是，以上方法僅適用于穩(wěn)態(tài)加速度信號的積分處理。針對非穩(wěn)態(tài)加速度信號，本文作者基于有效頻段法，綜合應(yīng)用Welch功率譜和短時傅里葉變換，獲得了良好的積分效果，后續(xù)將另文報告相關(guān)研究成果。