一種基于Lamb波的結構損傷檢測波形選擇方法

(1.中國飛行試驗研究院測試技術研究所，陜西 西安 710089；2.中國飛行試驗研究院飛機飛行試驗技術研究所，陜西 西安 710089)

0 引言

近幾十年來，隨著航空航天領域的快速發展，復雜工程結構內部的應力變得越來越高，加之使用環境也越來越惡劣，結構中的微小損傷有可能會導致結構損傷失效事故層。因此對于這些結構需要進行實時在線監控。

由于Lamb波具有傳播距離遠、對結構微小損傷敏感等特點，基于Lamb波的主動監測技術被認為是有效的結構損傷檢測方法，并成為了國際研究熱點之一[1]。在基于Lamb波的結構健康監測中，通過對損傷散射信號的處理分析，能夠實現對結構損傷的監測與評估[2]。Lamb波在結構中傳播時具有無方向性、頻散及多模式等特點，這就為損傷特征的提取增加了很多困難。同時，目前結構健康監測技術發射的波形固定，檢測環境發生改變時，僅依靠對結果信息的分析難以達到理想的檢測效果。因此激勵信號波形選擇在提升結構的損傷檢測效率方面有著重要的研究意義。

1 Lamb波激勵信號

基于Lamb波的主動結構健康監測方法中，通常采用寬帶激勵或窄帶激勵的方式激發出被測結構中的Lamb波。本文的研究工作使用窄帶Lamb波信號，因此下文中將對窄帶激勵進行闡述。

通常情況下，選擇調制正弦信號作為Lamb波窄帶激勵信號，其表達式[3]如下：

sin2πfct

(1)

式中，n為信號波峰個數；fc為信號中心頻率；A為信號的幅度調制；H為Heaviside階梯函數。

圖1為中心頻率為300 kHz的五周期正弦調制信號的時域波形及其傅里葉變換后的頻域波形。

圖1 五周期激勵信號及其傅里葉變換

根據Lamb波頻散特性，由于部分模態具有截止頻率，為了減少激發的模態數量，避免多模態的混疊，因此需要選擇合適的激勵信號中心頻率，以確保被測結構中的Lamb波只存在S0與A0兩種模式。激勵信號的周期數一般在3.5～13.5之間較為合適[4]。

由于Lamb波的自身特點會造成其傳播過程較為復雜，因此通過數值計算的方法模擬Lamb波激勵信號傳播一段時間后的波形變得十分有意義，本節將對此展開介紹。

設激勵信號為F(t)，激勵信號傳播一段距離后的信號為u(x，t)，其中，x為傳播距離，t為傳播時間。激勵信號在t=0時刻作用于坐標系上x=0的位置，二者之間的關系可由以下初始條件確定：

u(x，0)=0

(2)

u(0，t)=F(t)

(3)

由文獻[5]可推導出結構內任意一點的信號為：

(4)

式中，k(ω)為角頻率為ω時的波數；F(ω)是f(t)的傅里葉變換信號。

當激勵信號遇到類似于邊界或損傷的反射源時，反射后的信號為g(t)，信號接收點接收到的第j個反射源的反射信號為gj(t)，設激勵點到反射源的距離為dj，反射源到接收點的距離為xj，反射系數為Aj(ω)，則gj(t)可表示為：

(5)

在MATLAB里進行編程模擬激勵信號傳播一段時間后的波形，以中心頻率為300 kHz的五周期正弦調制信號為激勵信號，模擬激勵信號傳播一段距離(10.7 cm)，波形如圖2所示，現只討論其中S0直達波的波達時間。

從圖3可以看出，該算法所得數值模擬得到的直達波波形與在ANSYS中進行仿真模擬得到的直達波波形，除個別峰值外，每個周期峰值的波達時間基本吻合，可認為該算法在一定條件下，可準確模擬激勵信號在頻散情況下的傳播。這為后續的波形選擇算法奠定了基礎。

圖2 模擬傳播一段距離后的波形

圖3 數值解與仿真模擬S0模式對比

2 波形選擇

無論是在雷達領域還是在結構健康監測系統，環境噪聲的干擾一直是重要的影響因素，為了排除噪聲干擾，提高對目標的測量精度，近年來雷達領域興起了對認知雷達的研究，而結構健康監測也引入眾多的濾波方法。同時，結構健康監測系統本身與雷達系統具有相似的特點，所以可以在結構健康監測領域中引入雷達系統中的先進方法。

本節旨在給出一種提高靜態目標檢測精度的方法，將新興的認知雷達跟蹤動態目標[6]轉化為結構損傷檢測的靜態目標。首先從信號的回波模型出發，由信號估計理論推導測量噪聲協方差的卡拉美羅下界(Cramer-Rao Low Bound，CRLB)，得出波形參數對其影響與意義，然后在Kalman濾波的基礎上增加波形選擇，選擇合適波形使濾波協方差最小，從而提高對目標狀態的檢測精度。

2.1 信號參量估計

2.1.1 回波信號模型

Lamb波在板結構中傳播形式較為復雜，為了降低計算量，本節提出激勵信號傳播一段距離后波形的簡化計算模型。激發出的Lamb波信號f(t)，所激發的場量用u(x，t)表示，其中t與x分別為激勵信號傳播時間與傳播距離，激勵器處可表示為：

u(x，t)=f(t)

(6)

板中任意一點處：

(7)

根據圖3數值解與仿真模擬S0模式對比，MATLAB仿真得出的傳播一段距離后的S0波形與ANSYS仿真情況下得出的S0波形相似度極高，即MATLAB仿真波形與及ANSYS仿真波形的頻散特性基本相同。這里我們使用MATLAB數值模擬傳播一段距離后的S0波形，并與激勵信號作對比。圖4為傳播0.27 cm的S0波形與激勵信號的對比。

圖4 MATLAB數值模擬的S0波形與激勵信號的對比

同時，使用ANSYS有限元軟件仿真模擬激勵信號傳播一段距離后的S0波形與激勵信號對比如圖5所示。

圖6顯示的是在損傷長度為4 cm裂縫的條件下，激勵信號與ANSYS仿真模擬的損傷散射信號進

圖5 ANSYS仿真模擬S0波形與激勵信號的對比

行的對比。在對幅值進行歸一化處理后，能夠看出損傷散射信號與激勵信號的波形較為相似。

圖6 ANSYS仿真模擬S0波形與激勵信號的對比

通過圖5、圖6對比可以得出，在小頻散的情況下，傳播一段距離后S0模式波形和損傷散射波與激勵信號的波形十分相似。這里提出一個簡化的計算模型，接下來在本文的研究中我們使用式(8)代替式(7)。

(8)

其中，θ1表示損傷或者邊界反射后幅值的衰減程度，θ2為反射信號的時延。

2.1.2 無偏估計量的CRLB

在信號處理中，我們通常要根據觀測信號來估計信號的某些特征參量，這就需要用到信號參量估計的相關知識。在結構的損傷檢測中，本文將著重研究信號的振幅、時延的聯合估計，即多參量估計。本文采用的是最大似然估計，最大似然估計對于解決復雜估計問題具有較好的性能。

在本文中，對接收到的信號進行觀測時，我們只研究信號的幅值及時延，故設觀測矢量x=[x1x2]T，由于噪聲的干擾造成的觀測結果不準確，待估計參量設為θ=[θ1θ2]T，則以θ為參量的概率密度函數為f(x|θ)。在已觀測到的觀測值x=x1的條件下，似然函數f(x=x1|θ1)反映了θ1取各個值的可能性大小。估計問題其本質是根據觀測值確定未知量，換言之，根據已知觀測值x=x1來估計未知量θ1的值。

(9)

若似然函數是可導函數，那么最大似然估計的必要條件為：

(10)

我們假定真實信號為s(t，θ)，s(t，θ)是關于t的連續函數，θ=[θ1，θ2，…，θM]T，包含M個未知參量，θi可以是信號的頻率、時延、幅度等。在一段觀測時間[0，T]內，觀測到的混入噪聲的信號為：

x(t)=s(t，θ)+n(t)，0≤t≤T

(11)

式中，n(t)是功率譜密度為N0/2的高斯白噪聲。

在得到一個估計量后，為了衡量其性能是否達到最佳，且使參量估計為有效估計，我們引入克拉美羅下界[7]來揭示無偏估計量的估計誤差的最小值。

假定概率密度f(x|θ)滿足正則條件：

(12)

(13)

式中，

(14)

(15)

I-1(θ)即為無偏估計量的CRLB，其倒數I(θ)稱之為費希爾信息。

如果有效估計量存在，即存在達到CRLB的估計量，那么這個估計量必定是最大似然估計，此時的參量估計具有很好的性能。這里CRLB會隨發射的波形參數的變化而改變，是下文中波形選擇的關鍵。

2.1.3 測量噪聲方差的CRLB

由信號參量估計理論，我們接下來推導基于簡化模型的測量噪聲方差的卡拉美羅下界。

將(8)式帶入(11)式，信噪比表示為SNR。

x(t)=

(16)

由(10)式：

(17)

(18)

由(14)式：

(19)

則無偏估計量的CRLB為：

(20)

計算可得測量噪聲協方差的CRLB：

(21)

根據測量噪聲協方差的CRLB可以得出以下結論：

1)信噪比越大，信號時域脈沖寬度(這里對應于n/fc)越小，信號的時延估計方差就越小，精度也就越高。

2)對于信號的幅值，信噪比越大，信號的帶寬越小(時域脈沖寬度越大，即n/fc越大)，其估計方差就越小，精度也就越高。

2.2 基于Kalman濾波的波形選擇

Kalman濾波理論是由Wiener波理論的發展而發展起來的。隨著Kalman濾波理論的發展，它在隨機過程參數估計中的應用越來越廣泛，其后來廣泛的應用于解決各種最優濾波問題[8]。在本節中，在Kalman濾波基礎上添加目標狀態的一步預測是算法的關鍵內容。

2.2.1Kalman濾波系統模型

首先引入一個離散控制過程的系統，其狀態方程可表示為：

X(k)=AX(k-1)+W(k)

(22)

離散時間系統的測量方程可表示為：

Z(k)=HX(k)+V(k)

(23)

其中，A為狀態轉移矩陣，由于本文的研究對象為檢測結構的損傷狀況，且結構的損傷狀態不隨時間而改變，故A為單位矩陣。H為測量矩陣，由于測量值反映了結構的損傷狀況，故同樣為單位矩陣。W(k)和V(k)對應的協方差分別為Q和R，分別表示預測和測量過程中的噪聲。

首先利用系統的過程模型預測下一狀態，基于系統的上一狀態預測現在狀態：

X(k|k-1)=AX(k-1|k-1)

(24)

X(k-1|k-1)是上一狀態最優的結果。用P表示X(k|k-1)的預測協方差為：

P(k|k-1)=P(k-1|k-1)+Q

(25)

Kalman增益K為：

K(k)=P(k|k-1)/(P(k|k-1)+R)

(26)

當前時刻目標狀態的最優化估算值X(k|k)：

X(k|k)=

X(k|k-1)+K(k)(Z(k)-X(k|k-1))

(27)

更新X(k|k)的協方差以為保證Kalman濾波器繼續運行下去，則：

P(k|k)=(E-K(k))P(k|k-1)

(28)

而如果要在Kalman濾波的基礎上增加波形選擇部分，我們需要更進一步計算出濾波協方差的一步預測。

預測協方差的一步預測為：

P′(k+1|k)=P(k|k)+Q

(29)

Kalman增益的一步預測為：

K′(k+1)=P′(k+1|k)/(P′(k+1|k)+R)

(30)

濾波協方差的一步預測為：

P′(k+1|k+1)=(E-K′(k+1))P′(k+1|k)

(31)

濾波協方差的一步預測滿足一定波形選擇準則后，所選取的波形參數會對下一時刻的目標狀態的最優化估算值產生影響，進而使得基于Kalman濾波的目標狀態的最優化估計循環計算下去。

2.2.2 波形選擇準則

為使每一時刻狀態估計誤差的均方差最小，引入均方誤差最小準則作為波形選擇的依據，這一準則使狀態空間每個維數上誤差的平方和最小。即：

(32)

Zk為k時刻的測量值集合，濾波協方差定義為：

P(k|k)=

(33)

等號兩邊同時取跡，則：

Tr{P(k|k)}

(34)

則：

(35)

由(34)式與(35)式，可得：

Tr{Pk|k(αk)}=

Tr{Pk|k-1-Pk|k-1(Pk|k-1+R(αk))-1Pk|k-1}

(36)

根據(36)式可以得出，由預測協方差和測量噪聲方差能夠求得k時刻的濾波協方差的跡。由此可知，在均方誤差最小準則這一波形選擇的依據中，測量噪聲方差是與波形的參數相關的變量，也正是前文中所提到的是波形選擇的關鍵所在。

2.2.3 基于Kalman濾波的波形選擇流程圖

如圖7所示，圖中左邊的兩部分為常規Kalman濾波，右邊的部分為所增加的波形選擇模塊。

圖7 基于Kalman濾波的波形選擇流程圖

根據流程圖，可以看出在整個波形選擇過程中，包含有兩個循環計算過程，其一是狀態估計協方差的循環計算過程，這其中包括了協方差的一步預測與波形選擇；其二是目標狀態的估計，包括了目標狀態的預測方程與測量方程，其狀態估計的對象為前文中所提到的信號的時延與幅值，信號的時延與幅值精度是反映本算法是否可行的關鍵指標。

3 基于MATLAB的模擬計算

本節通過在MATLAB中進行模擬計算，驗證基于Kalman濾波的波形選擇這一算法的有效性，并與不帶有波形選擇模塊的Kalman濾波結果對比，驗證上文所述算法的計算精度。

在模擬計算環境時，我們將上一時刻最優值作為該時刻的預測，并且假定測量目標的真實狀態對應的信號時延為2.4×10-5s，對應的幅值變化為0.4。目標初始狀態設為：時延3.5×10-5s、幅值變化為0.5。為了檢驗所研究的波形選擇方法在不同強度的噪聲下的濾波效果，信噪比分別設為15 dB與25 dB，初始發射的波形為五周期漢寧窗調制函數，其波形參數為n=5以及fc=300 kHz。根據Lamb波的傳播特性以及頻散曲線，波形庫的建立周期數的取值區間為3～8，中心頻率取值區間為[250 kHz，350 kHz]，步長為10 kHz，整個波形庫包含有66種波形。

在經過50次迭代計算后，將不含有波形選擇的Kalman濾波與帶有波形選擇模塊的Kalman濾波作對比，比較信號的時延誤差與幅值變化的誤差。如圖8與圖9所示，虛線、點線、實線分別表示帶有波形選擇的Kalman濾波、不帶有波形選擇的Kalman濾波、參量的真實值。其中圖10為信號的時延估計誤差的對比，圖11為信號的幅值變化估計誤差的對比。

圖8 時延估計對比圖(SNR=25dB)

圖9 幅值估計對比圖(SNR=25dB)

圖10 多次平均的時延估計對比圖(SNR=25dB)

圖11 多次平均的幅值估計對比圖(SNR=25dB)

為了檢驗算法在不同強度噪聲下的適用性，當信噪比為15 dB時，信號的時延估計與幅值估計顯示在圖12和圖13中。

圖12 時延估計對比圖(SNR=15dB)

圖13 幅值估計對比圖(SNR=15dB)

根據圖14和圖15，在經過Kalman濾波后，信號的時延估計與幅值估計的精度都較初始值有了較明顯的提高。在經過一定次數的迭代計算之后，在初始值誤差比較大的情況下，Kalman濾波仍能將誤差顯著的減小。并且，無論是信號的時延估計還是信號的幅值估計，虛線(帶有波形選擇的Kalman濾波)比點線(不帶有波形選擇的Kalman濾波)更接近實線(真實值)，這表明了帶有波形選擇的Kalman濾波的計算結果要好于不帶有波形選擇的Kalman濾波。

圖14 多次平均的時延估計對比圖(SNR=15dB)

圖15 多次平均的幅值估計對比圖(SNR=15dB)

根據圖8與圖12、圖9與圖13的對比，在噪聲強度改變時，經過50次的迭代計算之后，帶有波形選擇的Kalman濾波精度均較初值有所提高，但是當噪聲強度比較大時，波形選擇后結果的誤差相對就會大一些。信號的時延估計與幅值估計精度的提升驗證了波形選擇算法在信噪比不同時均具有良好的效果。

4 結論

研究基于Kalman濾波的波形自適應選擇，從建立簡化的回波信號模型出發，經過信號估計理論中的最大似然估計，詳細推導了基于簡化回波模型的測量噪聲協方差的卡拉美羅下界，然后在Kalman濾波的基礎上增加波形選擇模塊，并以均方誤差最小準則作為波形選擇的依據，最后通過模擬計算驗證了經過波形選擇后信號的時延與幅值精度有了較大的提高。波形選擇依據的是均方誤差最小準則，由于是在原有激勵信號波形的基礎上進行波形選擇，且在波形選擇的準則下，所選擇的新波形對應的Kalman濾波的協方差會比原有激勵信號波形對應的Kalman濾波協方差更小，這是帶有波形選擇模塊的Kalman濾波精度更高的根本原因。證明了波形選擇方法是提高噪聲干擾下的結構損傷定位精度的另外一種思路。