代 煜，張建勛

(南開大學(xué) 機(jī)器人與信息自動化研究所，天津 300071)

小波變換是一種能對觀測信號進(jìn)行時(shí)頻聯(lián)合分析的方法，系統(tǒng)的自由衰減響應(yīng)信號在連續(xù)小波變換域上的系數(shù)分布呈現(xiàn)出脊的特征，由此發(fā)展出多種根據(jù)脊上小波變換系數(shù)的相位、幅值信息使系統(tǒng)自動解耦，并提取線性系統(tǒng)的模態(tài)參數(shù)［1-3］或弱非線性阻尼及剛度［4］的方法。

在對系統(tǒng)的自由衰減響應(yīng)信號進(jìn)行小波變換時(shí)，小波變換系數(shù)必然會受到邊界效應(yīng)的影響。邊界效應(yīng)產(chǎn)生的根本原因是小波變換中的卷積運(yùn)算，實(shí)際采集的信號都是在時(shí)域離散而且有限長的，對這樣的序列作卷積相當(dāng)于加矩形窗，因此在小波變換系數(shù)的兩個(gè)邊界處會發(fā)生畸變，這必然影響到模態(tài)參數(shù)的識別。對于模態(tài)參數(shù)識別中遇到的邊界效應(yīng)問題，目前主要有以下幾種解決方案:Boltezar和Slavic［5］提出用反射窗法、等窗面積法和自適應(yīng)小波函數(shù)法來改變小波函數(shù)以消除邊界效應(yīng);伊廷華等［6］提出運(yùn)用自回歸滑動平均模型進(jìn)行信號預(yù)測以消除邊界效應(yīng)的影響;徐增丙等［7］提出運(yùn)用最小二乘支持向量機(jī)對小波骨架進(jìn)行預(yù)測延拓的方法;Kijewski和Kareem［8］提出采用鏡像延拓技術(shù)以消除邊界效應(yīng)對小波變換系數(shù)的影響;Yan等［9］定量計(jì)算了不受邊界效應(yīng)影響的小波變換系數(shù)的范圍，指出只有這部分系數(shù)才能用于估計(jì)模態(tài)參數(shù)，其余受到邊界效應(yīng)影響的部分應(yīng)當(dāng)剔除。這幾種方法可以分為兩類:一類是通過調(diào)節(jié)小波函數(shù)中的參數(shù)實(shí)現(xiàn)對邊界效應(yīng)的抑制，一類是對信號進(jìn)行延拓以獲得較多不受邊界效應(yīng)影響的小波變換系數(shù)。受到種種條件的限制，小波函數(shù)的改變程度有限，因此研究信號的延拓方法是克服邊界效應(yīng)的最重要途徑。在眾多的延拓方法中，類似鏡像法、補(bǔ)零法的一些方法不具有對待分析信號的自適應(yīng)性，這將對模態(tài)參數(shù)的估計(jì)精度帶來一定程度的影響。

考慮到人工神經(jīng)網(wǎng)絡(luò)具有通過自學(xué)習(xí)逼近任意映射的能力，突破了需要建立參數(shù)模型的傳統(tǒng)方法，本文將利用反向傳播網(wǎng)絡(luò)對采樣時(shí)長較短的自由衰減響應(yīng)信號進(jìn)行雙向延拓，然后對信號進(jìn)行連續(xù)小波變換，從而提取各階模態(tài)的阻尼比和無阻尼自振頻率。

1 基于連續(xù)小波變換的模態(tài)參數(shù)識別方法的理論基礎(chǔ)

小波變換是一種能對信號進(jìn)行多尺度(分辨率)分析的線性變換。設(shè)ψ(t)是L2(R)(L2(R)表示平方可積的函數(shù)空間)內(nèi)的實(shí)值或復(fù)值函數(shù)，當(dāng)且僅當(dāng)ψ(t)的傅里葉變換ψ^(ω)滿足以下不等式時(shí)函數(shù)ψ(t)可以成為母小波:

函數(shù)x(t)∈L2(R)在尺度因子為a，時(shí)移因子為b時(shí)的連續(xù)小波變換Wx(a，b)表示為:

式中，a，b和t均是連續(xù)變量，且a＞0，ψ*(·)表示復(fù)共軛，〈·〉表示卷積。本文采用最常用的Morlet復(fù)小波作為母小波，定義如下式所示:

式中fc為小波中心頻率，fb為帶寬參數(shù)，j是虛數(shù)單位。二者都可因分析對象的不同而改變，這就提供了極大的靈活性。

線性多自由度系統(tǒng)的自由衰減響應(yīng)為［9］:

式中Ai表示第i階模態(tài)的振幅，φi表示初相角，fi是第i階模態(tài)的無阻尼自然頻率是有阻尼自然頻率，ξi表示阻尼比，n表示模態(tài)的階數(shù)。根據(jù)漸進(jìn)信號的相關(guān)理論［10］，x(t)的Morlet小波變換系數(shù)的近似表達(dá)式為:

由式(5)可見，對于每一個(gè)時(shí)移因子b，在尺度ai=取得局部極大值，因此可以提取到小波脊，并且認(rèn)為在這個(gè)尺度下僅有第i階模態(tài)的響應(yīng)決定著小波變換系數(shù)，其余(n-1)個(gè)模態(tài)的影響可以忽略不計(jì)，此時(shí)多模態(tài)被解耦。因此，尺度為ai時(shí)的小波變換系數(shù)可以表示為［9］:

小波變換系數(shù)的模和輻角可表示為［9］:

對式(7)和式(8)應(yīng)用對數(shù)和微分運(yùn)算，第i階模態(tài)的無阻尼自然頻率fi和阻尼比ξi由下式估計(jì)［9］:

式中l(wèi)n［·］表示取自然對數(shù)，Arg［·］表示求輻角。

為了分離信號中兩個(gè)較接近的頻率成份fi和fi+1，可以通過改變fb和fc提高頻率分辨率，但同時(shí)時(shí)間分辨率會下降，時(shí)間分辨率下降帶來的后果是邊界效應(yīng)更加嚴(yán)重。Yan等［9］研究了時(shí)頻分辨率對模態(tài)解耦和邊界效應(yīng)的影響，計(jì)算了小波參數(shù)的選擇范圍，如式(11)所示:

式中，α、β和γ均是用于調(diào)節(jié)范圍的常數(shù):α=2、β=4和γ=0.4;T表示自由衰減響應(yīng)信號的采樣時(shí)長。此外，為滿足小波函數(shù)的容許性條件，還應(yīng)當(dāng)使。Yan 等［9］還給出了邊界效應(yīng)的影響范圍:

式中，f表示小波分析窗的中心頻率。值得注意的是，當(dāng)待分析信號的采樣時(shí)長較短或者模態(tài)較密集時(shí)，式(11)將會無解，因此必須考慮增加信號的采樣時(shí)長。

2 利用反向傳播網(wǎng)絡(luò)對信號的雙向延拓

利用反向傳播神經(jīng)網(wǎng)絡(luò)對自由衰減響應(yīng)信號的延拓建立在一定數(shù)量的信號樣本之上。神經(jīng)網(wǎng)絡(luò)的實(shí)質(zhì)就是依據(jù)已有的樣本數(shù)據(jù)，通過訓(xùn)練獲得并存儲樣本所隱含的特征關(guān)系，因此它具有很強(qiáng)的學(xué)習(xí)適應(yīng)能力。

本文提出的用于信號延拓的三層反向傳播網(wǎng)絡(luò)模型如圖1所示，由輸入層，隱藏層和輸出層組成，同層神經(jīng)元間無聯(lián)系，前后層間全連接。圖1中輸入層樣本的個(gè)數(shù)p表示前p個(gè)信號的離散采樣點(diǎn)對下一個(gè)采樣點(diǎn)有直接影響，而最新觀測值將比早期觀測值包含更多信息。網(wǎng)絡(luò)中隱藏層的傳遞函數(shù)都采用S型正切函數(shù):

輸出層用線性傳遞函數(shù)。對于在0到T時(shí)刻等時(shí)間采樣的信號序列{x(1)，x(2)，…，x(q)}，利用反向傳播網(wǎng)絡(luò)對信號進(jìn)行雙向延拓的步驟如下:

(1)離散時(shí)間序列的歸一化處理。為了避免神經(jīng)網(wǎng)絡(luò)的不穩(wěn)定性及收斂速度慢等問題，必須先對原始數(shù)據(jù)進(jìn)行歸一化處理。本文通過式(14)將原序列變換到［0，1］區(qū)間:

式中，max{·}和min{·}分別表示求序列的最大值和最小值。

表1 反向傳播網(wǎng)絡(luò)的訓(xùn)練樣本(1)Tab.1 Training samples of back propagation network(1)

表2 反向傳播網(wǎng)絡(luò)的訓(xùn)練樣本(2)Tab.2 Training samples of back propagation network(2)

(2)對反向傳播網(wǎng)絡(luò)的訓(xùn)練。從上一步驟得到的序列{y(1)，y(2)，…，y(q)}中獲取訓(xùn)練樣本，如表1和表2所示，表1中的p個(gè)訓(xùn)練樣本用于預(yù)測在T時(shí)刻之后的信號序列，而表2中的p個(gè)訓(xùn)練樣本用于預(yù)測在0時(shí)刻之前的信號序列。訓(xùn)練樣本由輸入層經(jīng)隱藏層傳遞到輸出層，輸出層將訓(xùn)練結(jié)果與期望值進(jìn)行比較，按照減少誤差的方向，從輸出層反向經(jīng)隱藏層修正各連接權(quán)值，最后回到輸入層，這個(gè)過程反復(fù)進(jìn)行，直到反向傳播網(wǎng)絡(luò)的輸出誤差小于設(shè)定值時(shí)才認(rèn)為訓(xùn)練完成。

(3)對信號序列的雙向延拓。為預(yù)測T時(shí)刻之后的h點(diǎn)信號序列，對上一步驟中訓(xùn)練完的網(wǎng)絡(luò)依次輸入{y(q-p+1)，y(q-p+2)，…，y(q)}，{y(q-p+2)，y(q-p+3)，…，y(q+1)}，…，{y(q-p+h)，y(q-p+h+1)，…，y(q+h-1)};為了預(yù)測0時(shí)刻之前的h點(diǎn)信號序列，對上一步驟中獲得的另一個(gè)反向傳播網(wǎng)絡(luò)依次輸入{y(p)，y(p-1)，…，y(1)}，{y(p-1)，y(p-2)，…，y(0)}，…，{y(p-h+1)，y(p-h)，…，y(2-h)}。雙向延拓的好處在于，它能夠使得原始信號序列的小波變換系數(shù)不受邊界效應(yīng)的影響，而受到邊界效應(yīng)影響的小波變換系數(shù)均來自于對延拓信號序列的連續(xù)小波變換，這就能夠使得原始信號的小波變換系數(shù)對最終模態(tài)參數(shù)識別結(jié)果的貢獻(xiàn)最大。

(4)通過式(15)將神經(jīng)網(wǎng)絡(luò)輸出序列的幅度還原:

圖1 用于信號雙向延拓的反向傳播網(wǎng)絡(luò)Fig.1 Back propagation network developed for signal extension in both directions

3 數(shù)值仿真

兩自由度有阻尼系統(tǒng)的自由衰減響應(yīng)如式(4)所示，令式中i=2;f1=10 Hz;ξ1=0.005;φ1=π/6;f2=12 Hz;ξ2=0.01;φ2= π/4;Ai=1(i=1，2);采樣頻率fs=200 Hz以及采樣時(shí)長T=1.5 s。將這些數(shù)據(jù)代入式(11)，無解。對于圖1所示的反向傳播網(wǎng)絡(luò)，輸入層的神經(jīng)元個(gè)數(shù)設(shè)為40，隱藏層的神經(jīng)元個(gè)數(shù)為5，利用表1和表2中的樣本對神經(jīng)網(wǎng)絡(luò)進(jìn)行訓(xùn)練，取p=200。初始學(xué)習(xí)速率參數(shù)設(shè)置為0.05，動量系數(shù)設(shè)為 0.9，目標(biāo)誤差設(shè)為 1 ×10-4［11］。圖 2 給出了利用反向傳播網(wǎng)絡(luò)對信號進(jìn)行雙向延拓的結(jié)果(延拓序列的長度h=100)，圖中兩條豎線間的區(qū)域是原始信號，實(shí)線是延拓的信號，虛線表示該自由衰減響應(yīng)信號的理論幅值。從圖2中可見只是在信號幅度的突變處延拓的信號與其理論值有少量差異，其余時(shí)刻二條曲線基本吻合。將信號的采樣時(shí)長改為2.5 s，代入式(11)解得3.5≤≤5。根據(jù)式(12)計(jì)算得到邊界效應(yīng)的影響范圍，將這部分受到邊界效應(yīng)影響的小波變換系數(shù)剔除，而后利用式(9)和式(10)辨識各階模態(tài)的阻尼比和無阻尼自振頻率。

考慮到實(shí)際采集的自由衰減響應(yīng)信號常常被噪聲污染，表3給出了不同噪聲水平下模態(tài)參數(shù)識別的相對誤差，所加的噪聲為零均值的高斯噪聲。信號的噪聲水平由信噪比(SNR)來衡量，SNR定義為信號標(biāo)準(zhǔn)差與噪聲標(biāo)準(zhǔn)差的比值［12］。設(shè)r表示某個(gè)量的精確值，r0表示它的測量值或者估計(jì)值，相對誤差定義為:

從表3可見，基于反向傳播網(wǎng)絡(luò)的雙向延拓方法削弱了邊界效應(yīng)對模態(tài)參數(shù)識別的影響，這種方法對于信號采樣時(shí)長較短或者模態(tài)較密集的場合有效。值得注意的是，由于在信號幅度的突變處延拓信號的質(zhì)量較差，這給阻尼比的估計(jì)帶來一定影響。從提取模態(tài)參數(shù)消耗的時(shí)間來看，SNR越低，訓(xùn)練神經(jīng)網(wǎng)絡(luò)所需要的時(shí)間越長。

圖2 利用反向傳播網(wǎng)絡(luò)對自由衰減響應(yīng)信號的延拓Fig.2 Extension of free decay response signal using back propagation network

如果僅采用基于連續(xù)小波變換的模態(tài)參數(shù)識別方法，由于不能計(jì)算得到合適的，也就是邊界效應(yīng)將影響到所有的小波變換系數(shù)，這將導(dǎo)致模態(tài)參數(shù)識別結(jié)果的完全錯(cuò)誤。如果采用直接設(shè)計(jì)神經(jīng)網(wǎng)絡(luò)辨識模態(tài)參數(shù)的方法，就必須要求一定數(shù)量的訓(xùn)練樣本，如果訓(xùn)練樣本的SNR和待辨識信號的SNR差距較大就不能得到較好的辨識結(jié)果，而小波變換能夠在一定程度上克服噪聲的影響，這就避免了單純使用基于神經(jīng)網(wǎng)絡(luò)的辨識方法帶來的不確定性。

4 實(shí)際問題的仿真分析

本文提出的模態(tài)參數(shù)識別方法在對用于微創(chuàng)外科手術(shù)機(jī)器人的力傳感器的動態(tài)性能測試實(shí)驗(yàn)中得到了應(yīng)用。具備力感知功能的四自由度微創(chuàng)外科手術(shù)機(jī)器手如圖3所示，機(jī)器手通過右端的快換接口安裝在微創(chuàng)外科手術(shù)機(jī)器人的機(jī)械臂末端，機(jī)械手中部是外徑10 mm，長約500 mm的空心金屬圓桿，手術(shù)開始時(shí)機(jī)器手經(jīng)開在患者腹腔上的通道插入體內(nèi)，由機(jī)器手左端鉗形工具夾持手術(shù)針或其他器械完成整個(gè)手術(shù)。為了測量手術(shù)過程中鉗形工具與患者內(nèi)臟組織的接觸力，將應(yīng)變片陣列布置在圓桿表面上靠近鉗形工具的位置，并對圓桿進(jìn)行了必要的切削加工以滿足應(yīng)力集中和輸入輸出引線的需要。

表3 不同信噪比下模態(tài)參數(shù)識別的相對誤差(fb=4，fc=2)Tab.3 Relative errors of identified modal parameters with different SNR(fb=4，fc=2)

力傳感器的動態(tài)性能測試采用如下方法:將粘貼了應(yīng)變片陣列的圓桿一端裝夾在傳感器的標(biāo)定裝置上，此時(shí)的圓桿可視為一個(gè)變截面懸臂梁，圓桿的另一端拴一根細(xì)銅絲，細(xì)銅絲繞過安裝在標(biāo)定裝置上的滑輪吊起砝碼，通過剪斷銅絲，突然釋放砝碼，使得圓桿做彎曲自由振動，采集應(yīng)變信號并從中提取各階應(yīng)變模態(tài)的無阻尼自振頻率和阻尼比。由于應(yīng)變信號很微弱且易被噪聲污染，經(jīng)過調(diào)理(包括調(diào)零、約3萬倍的放大和截止頻率為318 Hz的低通濾波)之后僅能獲得0.15 s的自由衰減響應(yīng)信號(如圖4所示)，信號的采樣頻率為2 000 Hz。由于信號的采樣時(shí)長過短，不能直接使用基于連續(xù)小波變換的模態(tài)參數(shù)識別方法。

圖3 具備力感知功能的微創(chuàng)外科手術(shù)機(jī)器手Fig.3 Minimally invasive surgery manipulator with force sensing ability

圖4 調(diào)理后的應(yīng)變橋輸出信號Fig.4 Output of strain gauge bridge after signal conditioning

表4給出了對信號雙向延拓后，應(yīng)用連續(xù)小波變換對力傳感器的應(yīng)變模態(tài)參數(shù)進(jìn)行辨識的結(jié)果(所用的反向傳播神經(jīng)網(wǎng)絡(luò)如圖1所示，參數(shù)與第3節(jié)數(shù)值仿真中所用的參數(shù)相同)，由于低通濾波器的限制，僅能辨識到前兩階應(yīng)變模態(tài)。表4給出的辨識結(jié)果表征了力傳感器的動態(tài)特性，為傳感器的建模以及進(jìn)一步設(shè)計(jì)補(bǔ)償器提高傳感器的動態(tài)響應(yīng)速度提供了實(shí)驗(yàn)依據(jù)。當(dāng)阻尼比較大時(shí)，自由衰減響應(yīng)信號的衰減速度很快，而在實(shí)際信號中所包含噪聲的幅度卻不會以同樣的規(guī)律衰減，這就導(dǎo)致有用信號逐漸被噪聲淹沒，本文提出的方法為這樣一類問題提供了解決方案:將這部分受噪聲污染較重的信號直接剔除，而后對信號進(jìn)行雙向延拓，延拓的信號比原信號具有較高的信噪比，從而能夠更加準(zhǔn)確地利用基于小波變換的方法估計(jì)模態(tài)參數(shù)。

表4 力傳感器應(yīng)變模態(tài)參數(shù)的辨識結(jié)果(fb=1，fc=2)Tab.4 Identified strain modal parameters of force sensor(fb=1，fc=2)

5 結(jié)論

考慮到人工神經(jīng)網(wǎng)絡(luò)具有自學(xué)習(xí)和自適應(yīng)的能力，本文提出了一種基于連續(xù)小波變換和反向傳播神經(jīng)網(wǎng)絡(luò)的模態(tài)參數(shù)辨識方法。利用反向傳播神經(jīng)網(wǎng)絡(luò)對自由衰減響應(yīng)信號的雙向延拓，解決了小波變換的邊界效應(yīng)問題，擴(kuò)展了傳統(tǒng)模態(tài)參數(shù)辨識方法的應(yīng)用領(lǐng)域。數(shù)值仿真以及對力傳感器進(jìn)行的動態(tài)特性測試實(shí)驗(yàn)表明，本文提出的方法在信號的采樣時(shí)長很短時(shí)也能夠準(zhǔn)確地辨識各階模態(tài)的阻尼比和無阻尼自振頻率，實(shí)用性較強(qiáng)，并體現(xiàn)出一定的抗噪聲能力。

［1］ Le T P，Argoul P.Continuous wavelet transform for modal identification using free decay response［J］.Journal of Sound and Vibration，2004，277(1-2):73-100.

［2］滕 軍，朱焰煌，周 峰，等.基于復(fù)Morlet小波變換的大跨空間結(jié)構(gòu)模態(tài)參數(shù)識別研究［J］.振動與沖擊，2009，28(8):25-29.

［3］徐曉霞，任偉新，韓建剛.基于響應(yīng)協(xié)方差小波變換和SVD的結(jié)構(gòu)工作模態(tài)參數(shù)識別［J］.振動工程學(xué)報(bào)，2010，23(2):194-199.

［4］代 煜，孫和義，李慧鵬，等.基于小波變換的弱非線性阻尼和剛度辨識方法［J］.振動與沖擊，2009，28(2):51-55.

［5］ Boltezar M，Slavic J.Enhancements to the continuous wavelet transform for damping identifications on short signals［J］.Mechanical Systems and Signal Processing，2004，18(5):1065-1076.

［6］伊廷華，李宏男，王國新.基于小波變換的結(jié)構(gòu)模態(tài)參數(shù)識別［J］.振動工程學(xué)報(bào)，2006，19(1):51-56.

［7］徐增丙，軒建平，史鐵林，等.LS-SVM在基于小波變換的模態(tài)分析中端部效應(yīng)的應(yīng)用［J］.中國機(jī)械工程，2008，19(13):1614-1620.

［8］Kijewski T，Kareem A.Wavelettransformsforsystem identification in civil engineering［J］.Computer-Aided Civil and Infrastructure Engineering，2003，18(5):339-355.

［9］Yan B F，Miyamoto A，Bruhwiler E.Wavelet transform-based modal parameter identification considering uncertainty［J］.Journal of Sound and Vibration，2006，291(1-2):285-301.

［10］ Delprat N，Escudie B，Guillemain P，et al.Asymptotic wavelet and gaboranalysis:extraction ofinstantaneous frequencies［J］.IEEE Transactions on Information Theory，1992，38(2):644-664.

［11］ Fredric M H，Ivica K.(著)，葉世偉，王海娟(譯).神經(jīng)計(jì)算原理［M］.北京:機(jī)械工業(yè)出版社，2007.

［12］ Yu K，Ye J，Zou J，et al.Missile flutter experiment and data analysis using wavelet transform ［J］.Journal of Sound and Vibration，2004，269(3-5):899-912.