基于模糊C均值聚類小波包重構(gòu)的車輛頻域非平穩(wěn)振動譜模擬方法研究

叢 楠， 任焱晞， 陳俊達， 王彬星

(1.工程裝備系統(tǒng)工程研究所，北京 100093； 2.北京市標準化研究院，北京 100013)

采集車輛在實際道路上行駛時的振動響應(yīng)，是準確開展道路模擬試驗的重要前提。由于受到采集場地與設(shè)備的限制，在有限里程內(nèi)采集的實車響應(yīng)譜往往無法滿足在試驗臺架上或虛擬仿真環(huán)境下模擬數(shù)百乃至數(shù)十萬公里行程的需求。因此，在開展模擬試驗時，實際輸入振動臺架作為激勵的往往是根據(jù)采集譜的某些特征進行“重構(gòu)”處理后所得到的重構(gòu)譜。

按照重構(gòu)目標特征的不同，對于車輛振動譜的重構(gòu)方法可分為時域循環(huán)、頻域(功率譜)重構(gòu)和幅值概率密度重構(gòu)三個大類[1]。其中，以功率譜為目標的重構(gòu)方法應(yīng)用最為廣泛，并形成了以諧波疊加法、線性濾波法、ARMA法為代表的一系列成熟方法[2-5]。然而，該類方法均以采集譜的總體功率譜為重構(gòu)目標，因而無法重構(gòu)出非平穩(wěn)振動譜。

隨著車輛設(shè)計與試驗水平的不斷提高，對于如何準確模擬車輛的各種非平穩(wěn)隨機振動得到了越來越多的關(guān)注。王巖松等[6]由時變狀態(tài)方程建立一種在車速變化條件下的非平穩(wěn)振動輸入模型；張立軍等[7]基于等效協(xié)方差法建立了一種通用的車輛非平穩(wěn)振動輸入模型；萇亮等[8]由轉(zhuǎn)移函數(shù)法與狀態(tài)空間法對車輛對路面的非平穩(wěn)動載荷進行了模擬。上述方法均須對車輛振動系統(tǒng)建模，并對該系統(tǒng)做線性假設(shè)，僅考慮由于車速變化造成路面激勵關(guān)于時間的非平穩(wěn)性(假定路面激勵關(guān)于空間是平穩(wěn)的)所導(dǎo)致的車輛振動的在振動強度上的非平穩(wěn)性。而事實上，由于車輛振動系統(tǒng)是個復(fù)雜的非線性系統(tǒng)，車輛行駛時，即使路面激勵輸入是(關(guān)于時間)平穩(wěn)的，其上某一部位的振動也會受到各種結(jié)構(gòu)非線性因素以及振動傳遞路線多樣性的影響，在頻域上呈現(xiàn)出一定的非平穩(wěn)特征，即振動的瞬時功率譜具有時變性。在時變瞬時功率譜信號的模擬方面，包括小波(包)、經(jīng)驗?zāi)Ｊ椒纸?、希爾伯特譜分析在內(nèi)多種時頻分析方法，已被應(yīng)用于包括地震、風(fēng)載、水輪機以及車輛振動的模擬，但這些方法大多要求功率譜的變化須滿足特定的隨機分布[9-11]，因此上述方法得到的重構(gòu)譜并不能認為是嚴格意義上頻域非平穩(wěn)的，不能很好的與實際車輛振動瞬時功率譜的變化特點相適應(yīng)。

本文通過對實車采集譜進行小波包變換，將一維振動譜的重構(gòu)問題轉(zhuǎn)化為二維小波包系數(shù)矩陣(圖像)的重構(gòu)問題?；趯π〔ò禂?shù)歸一化列向量的模糊C均值(Fuzzy C-Means，F(xiàn)CM)聚類，提取包括振動譜所包含的主要振動類型、各類型發(fā)生的總概率與相鄰條件概率、各類型中圍繞聚類中心的分散程度等采集譜在頻域中的非平穩(wěn)特征，并依據(jù)該特征實現(xiàn)對采集譜的頻域非平穩(wěn)化模擬。通過將該方法應(yīng)用于某試驗道路響應(yīng)譜的重構(gòu)工作，所得到的重構(gòu)譜相比于原采集譜，在具有相似的頻域非平穩(wěn)特性的同時，其在時域波形連續(xù)性與平滑性、總體功率譜以及幅值概率密度等方面均具有良好的相似程度。本文為車輛振動譜的重構(gòu)與模擬工作，提供了一種新的思路與方法。

1 振動譜的小波包系數(shù)矩陣及其行列向量

1.1 小波包變換與小波包系數(shù)矩陣

小波變換是目前廣泛應(yīng)用于非平穩(wěn)信號分析與處理的時頻分析方法。在信號合成與重構(gòu)方面，其已在包括地震波、脈動風(fēng)載、軌道不平度等諸多領(lǐng)域取得了成功的應(yīng)用[12-14]。

小波包變換由小波變換發(fā)展而來，其基本思想是在對小波分解后的高頻部分與低頻部分等同對待，均做進一步分解，從而在全頻帶上具有一致的時頻分辨率。相比小波變換，更適合于對已經(jīng)過低通濾波等前處理的在全頻帶上均須做出準確分析的車輛振動譜進行分解與合成。

對于長度為N的離散信號，工程中通常使用Mallat快速(遞推)算法[15]來進行小波包分解

(1)

相應(yīng)地，小波包重構(gòu)(反變換)的遞推算法為

(2)

(3)

通常，依據(jù)上述系數(shù)矩陣中系數(shù)取值的不同，以圖像的形式繪制出該矩陣可以更清晰的顯示出信號的時頻特征(由此又被稱為時頻圖像)。根據(jù)小波包變換性質(zhì)，對于采樣頻率Fs的信號進行j層小波包分解，其時頻分辨率分別為

(4)

即隨著層數(shù)的增加，小波包分解的時間分辨率下降而頻率分辨率上升。

圖1(a)所示為經(jīng)過預(yù)處理的某型車駛過試驗場塊石路時的單一輪軸振動譜，采樣頻率為128 Hz，使用db6小波包對其進行6級分解的系數(shù)圖像如圖1(b)所示，由式(4)可知，該圖像的時間/頻率分辨率分別為0.5 s/1 Hz。

(a) 振動譜

(b) 小波包圖像圖1 某試驗場塊石路的車輛振動譜及其小波包系數(shù)圖像Fig.1 Vehicle vibration spectrum and its wavelet packet coefficients image of a proving ground gravel road

1.2 小波包變換與小波包系數(shù)矩陣

在式(3)所示的小波系數(shù)矩陣中，其第i行系數(shù)組成的行向量

ri={d[i,1],d[i,2],…，d[i,n]}

(5)

即為某一小波包系數(shù)，根據(jù)小波包分解的性質(zhì)，其代表了在頻帶[(i-1)Δf,iΔf]內(nèi)的信號能量隨時間的變化，且其上各系數(shù)的平方和

(6)

代表了信號在該頻帶內(nèi)的廣義能量。由帕賽瓦爾定理，其與信號的功率譜之間有如下對應(yīng)關(guān)系

(7)

式中：Cφ為一個與小波包基的類型有關(guān)的常數(shù)；P(f)為原信號的功率譜密度函數(shù)。在全頻帶上，亦有

(8)

式(7)、式(8)說明了各小波包系數(shù)的平方和直接反映了該頻帶內(nèi)信號能量的強弱。同理，對于系數(shù)矩陣中的第k列系數(shù)向量及其平方和

ck={d[1,k],d[2,k],…，d[2j,k]}T

(9)

(10)

其分別代表了在[(k-1)Δt,kΔt]時段內(nèi)，信號的廣義瞬時功率譜與廣義瞬時能量。

此外，為便于進行平穩(wěn)性分析，對于如式(5)、式(9)所示的小波包系數(shù)行/列系數(shù)向量，分別定義其歸一化版本

(11)

圖2給出了如圖1所示小波包系數(shù)進行歸一化后的行向量與列向量圖像。

(a) 歸一化行向量

(b) 歸一化列向量圖2 小波包歸一化系數(shù)向量圖像Fig.2 Normalized wavelet packet coefficients image

2 基于FCM聚類小波包系數(shù)矩陣重構(gòu)的頻域非平穩(wěn)振動模擬

正如前文所述，當(dāng)前討論關(guān)于車輛振動信號的平穩(wěn)性的討論，大多單指其振動強度關(guān)于時間的平穩(wěn)性。事實上，由圖2(a)中所示的歸一化行向量圖像即可看出，當(dāng)路面條件單一時，車輛響應(yīng)在所有頻帶內(nèi)的強度可認為是平穩(wěn)的(此時振動譜可順利通過平穩(wěn)性檢驗)。而圖2(b)所示的歸一化列向量圖像則表現(xiàn)出明顯的非平穩(wěn)特性，由歸一化列向量的含義，說明該振動譜的瞬時功率譜具有明顯的時變性。

2.1 歸一化列向量的FCM聚類

(12)

步驟2計算隸屬度

(13)

步驟3根據(jù)隸屬度修正聚類中心

(14)

由于FCM聚類方法的靈活性，其無需預(yù)先設(shè)定聚類中心的確切數(shù)量及數(shù)值(即瞬時功率譜的樣式)，僅需任意指定一個較大的聚類數(shù)量M(該值由受試車輛的主要激振源的數(shù)量決定，一般取10～20即可滿足要求)作為初始值，即可滿足對車輛瞬時功率譜的聚類要求。

2.2 基于FCM聚類結(jié)果的歸一化列向量非平穩(wěn)重構(gòu)

經(jīng)過FCM聚類，原有n組列向量將被歸入m(m≤M)個類別里，按照類別內(nèi)所包含列向量數(shù)量由多到少依次排序，各類別包含的數(shù)量為nci|i=1,2,…,m，其對應(yīng)的聚類中心(曲線)ccenter-i|i=1,2,…,m。為了進行歸一化列向量重構(gòu)，須對聚類結(jié)果做進一步統(tǒng)計：

(2) 統(tǒng)計每一類別中所屬的各歸一化列向量相對聚類中心的擴散程度。由于各列向量及聚類中心均可表示為曲線，且類別中所包含的原列向量數(shù)量往往比較有限，因此，在不對擴散的分布特征做前提假設(shè)時，本文使用百分位對曲線上的各點做單獨統(tǒng)計，即對于第i|1,2,…,m個分類上所有歸一化列向量上的第k|1,2,…,2j個值(對應(yīng)于相應(yīng)的列歸一化小波包系數(shù))，統(tǒng)計出若干個百分位點：{Qik-q1,Qik-q2,…,Qik-qn}|0

在完成上述對聚類結(jié)果的統(tǒng)計與分析后，即可對歸一化列向量進行重構(gòu)。

重構(gòu)過程分為兩步：

(1) 指定各重構(gòu)歸一化列向量的類別序列，并使相鄰的概率為pij。具體可由如下遞推過程實現(xiàn)：

(i) 隨機指定一個類別i∈{1,2,…,m}；

(ii) 產(chǎn)生一個在[0,1]范圍內(nèi)的均勻分布的隨機數(shù)u；

(iii) 令i+1←min(j|pcum-ij≥u)；

(iv) 令i←i+1并重復(fù)第(ii)步，直至確定所有重構(gòu)列向量的分類序列。

(2) 按照既定的分類序列生成重構(gòu)歸一化列向量曲線，并使各重構(gòu)曲線滿足與原聚類結(jié)果中各類別中列向量圍繞其聚類中心的散布樣式相一致。對于基于百分位統(tǒng)計出的原分散結(jié)果，可由如下過程生成各重構(gòu)曲線：

(ii) 對應(yīng)于同一個分類中的第k|1,2,…,2j個值，在其各百分位點qn+1形成的qn-1個分段中按均勻分布隨機插值，同時產(chǎn)生Nci|i=1,2,…,m個隨機數(shù)；

(iii) 將所有2j組隨機數(shù)合并，則構(gòu)成了第i類的所有Nci|i=1,2，…，m個重構(gòu)歸一化列向量；

(iv) 令i←i+1并重復(fù)第(ii)步，直至生成所有類別所需的所有重構(gòu)列向量；

(v) 對所有最大值不為1的重構(gòu)列向量進行校正，使其歸一化；

(vi) 按照已生成的重構(gòu)列向量的分類序號，從所生成的各類重構(gòu)歸一化列向量中依次抽取，生成重構(gòu)歸一化列向量矩陣。

2.3 小波包系數(shù)矩陣重構(gòu)

對于已生成的重構(gòu)歸一化列向量陣，僅需補充表示強度的信息，即可將其恢復(fù)為重構(gòu)小波包系數(shù)矩陣。由于由式(10)所示的廣義瞬時能量被認為是平穩(wěn)的，因此利用信號小波包系數(shù)矩陣的該廣義能量信息，即可完成重構(gòu)小波包系數(shù)矩陣的工作，具體步驟如下：

步驟2將重構(gòu)廣義瞬時能量賦予各重構(gòu)的歸一化列向量

(15)

3 塊石路響應(yīng)譜重構(gòu)實例

對于如圖1所示的塊石路采集譜及其小波包系數(shù)圖像，將圖2(b)所示的歸一化列向量圖像用曲線集合的形式表示，如圖3所示。由小波包分解的含義，圖3即表示了原信號的時間分辨率為0.5 s、頻率分辨率為1 Hz的歸一化瞬時功率譜譜線的集合。

圖3 歸一化小波包列系數(shù)向量(局部)Fig.3 Normalized wavelet packet column coefficient vectors(partial)

通過FCM對如圖3所示的歸一化列向量曲線集合進行聚類，形成了9個有效分類(有少量僅包含1條列向量曲線的分類合并另作為一個類別，在重構(gòu)時可根據(jù)需要從中隨機抽取或忽略)。按照聚類中心峰值頻率、峰值頻率變化范圍以及中高頻成分強弱的不同，該9個有效分類又進一步可劃分為4個大類，如圖4所示。

(a) 第I大類(第2、第3有效分類)

(b) 第II大類(第4、第6有效分類)

(c) 第III大類(第1、第7有效分類)

(d) 第IV大類(第5、第8、第9有效分類)圖4 聚類中心Fig.4 Centers of clustering

由圖4可知：第I類中心線具有尖銳的峰值，代表了較為理想的車輛振動響應(yīng)，此時各種外部非線性因素影響較小，車輛運行平穩(wěn)，中心曲線峰值頻率的輕微變化表明了車輛振動系統(tǒng)結(jié)構(gòu)動力學(xué)參數(shù)具有輕微的非線性變化；第II類中心曲線具有一強一弱兩處峰值，代表了由采集處振動受到某一外部(此處為對側(cè)車輪)激勵的相對強烈的影響；第III類中心曲線具有較寬的主峰，其所包含的曲線的峰值頻率具有較大的變化范圍，表明采集處振動受到多處外部激勵的影響以及振動系統(tǒng)結(jié)構(gòu)動力學(xué)參數(shù)的較強的非線性變化；第IV類中心曲線在中高頻段上的具有多處峰值，表明振動中包含的間歇性中高頻或沖擊性激勵成分。

對各有效分類進行統(tǒng)計，可得到如圖5所示的所屬曲線集合圍繞中心線的散布情況(以第3有效類為例，并以箱線圖的形式表示)，以及如表1所示的各有效類別之間的一階相鄰條件概率(以第1、第3、第5、第7個有效分類為例)。

圖5 第3有效分類中所屬各曲線圍繞聚類中心的散布情況Fig.5 Scatter of curves in valid category 3

表1 有效分類間的一階相鄰條件概率(部分)Tab.1 First-order adjacent conditional probabilities among valid categories (partial)

由表1可知，各有效分類之間并非均勻分布的隨機出現(xiàn)，而是呈較明顯的成組及連續(xù)出現(xiàn)的趨勢，從而顯現(xiàn)出明顯的頻域非平穩(wěn)特性。根據(jù)上文所述的重構(gòu)方法，首先由上述一階相鄰概率依次指定重構(gòu)列向量的類別，由全概率公式，當(dāng)重構(gòu)列向量數(shù)量足夠大時，重構(gòu)列向量的各類別的總體概率將與原信號保持一致。圖6展示了使用本文方法重構(gòu)出10倍于原信號長度的重構(gòu)信號時，重構(gòu)信號及原信號中歸一化列向量各有效分類所占的比例。圖7對比了原信號的歸一化列向量與重構(gòu)歸一化列向量圖像。由圖6、圖7可知，對于歸一化列向量的按類別及散布特征重構(gòu)，可對車輛振動中的頻域非平穩(wěn)特征進行有效模擬。

圖6 原信號與重構(gòu)信號中歸一化列向量各有效分類的概率Fig.6 Probabilities of valid clustering categories of original vs reconstructed signal

圖7 原信號與重構(gòu)信號的歸一化列向量圖像Fig.7 Normalized column vector images of original vs reconstructed signal

原始采集譜及重構(gòu)的廣義能量序列如圖8(a)所示。當(dāng)原信號為強度平穩(wěn)時，可認為該廣義能量序列為滿足某一分布(此處為對數(shù)正態(tài)分布)的隨機序列，根據(jù)其生成的重構(gòu)廣義能量序列及兩信號廣義能量序列的概率密度分布分別如圖8(b)、圖9所示。將該重構(gòu)廣義能量序列賦予圖7所示的重構(gòu)歸一化向量，得到的重構(gòu)小波包系數(shù)圖像如圖10所示。

(a) 原信號(b) 重構(gòu)信號圖8 廣義能量序列Fig.8 Generalized energy sequence圖9 廣義能量序列的幅值概率密度Fig.9 Probability density of generalized energy se-quences 圖10 重構(gòu)小波包系數(shù)圖像Fig.10 Wavelet packet coeffi-cients image of recon-structed signal

基于式(2)所示的小波包重構(gòu)，將重構(gòu)小波系數(shù)矩陣進行小波包重構(gòu)，即可得到重構(gòu)振動譜。該重構(gòu)譜與原始譜在時域波形、總體功率譜以及幅值概率密度的對比情況分別如圖11～圖13所示。

(a) 總體波形

(b) 局部放大圖11 不同尺度下的時域波形Fig.11 Time domain waveform in different scales

圖12 總體功率譜密度Fig.12 Gross PSDs

圖13 幅值概率密度分布Fig.13 Amplitude PDFs

通過圖11～圖13可知，采用本文方法得到的重構(gòu)譜在具有與原采集譜相似的非平穩(wěn)特性的基礎(chǔ)上，同時在時域波形、總體功率譜及幅值概率密度上與原始采集譜具有較高的一致性。并且，根據(jù)dB小波包自身性質(zhì)，通過合理選擇小波包消失矩的階數(shù)，可使得重構(gòu)譜在保持原采集譜的平滑性連續(xù)性與功率譜一致性之間，進行靈活的取舍——對于包含更多沖擊樣高幅值激勵且平滑性差的振動譜，可選用消失矩階數(shù)較低的小波包；相反，對于較為平緩的振動譜，應(yīng)選用消失矩較高的小波包。

相比于傳統(tǒng)振動譜的重構(gòu)方法，本文方法得到的重構(gòu)振動譜能夠更真實地反映出車輛的實際振動情況，為開展各種車輛振動模擬試驗提供更為精確的振動輸入。

4 結(jié) 論

受到地面-車輛系統(tǒng)中多種因素的影響，實測車輛振動譜在頻域上往往是非平穩(wěn)的；該非平穩(wěn)特征體現(xiàn)為若干種類別振動模式的間歇性交替出現(xiàn)并相互交疊；通過對振動譜小波包系數(shù)歸一化列向量進行FCM聚類，可有效地對各振動模式(及其主要交疊模式)進行提取與分析；通過重構(gòu)與聚類結(jié)果相一致的歸一化列向量進而生成重構(gòu)系數(shù)矩陣，經(jīng)由小波包逆變換，可得到與原采集譜具有相似振動強度、總體功率譜及幅值概率密度，且同時擁有相似頻域非平穩(wěn)特性的重構(gòu)振動譜；通過選取合適的小波包類型及分解層數(shù)，原采集譜的平滑性與連續(xù)性也可得到有效的模擬。

通過將該方法應(yīng)用于某試驗場塊石路采集譜的重構(gòu)工作，驗證了本文方法在車輛振動譜模擬方面的正確性和有效性。此外，由于小波包分析及FCM聚類方法本身對振動譜分析的普遍適用性，本文所述方法亦可應(yīng)用于其它非平穩(wěn)振動譜的模擬工作，為各種振動模擬試驗提供與實際更為吻合的振動輸入。