Lamb波波包混疊分離方法的結構損傷定位?

郭 禹，張 超，季宏麗，吳義鵬，裘進浩，王 勇

（1.南京航空航天大學機械結構力學及控制國家重點實驗室 南京，210016）

（2.上海宇航系統工程研究所結構系統研究室 上海，201109）

引言

結構健康監測（structural health monitoring，簡稱SHM）是一種實時、在線的監測技術，通過獲取與結構損傷相關的多種參數，識別結構中存在的損傷，進而預測結構的使用壽命，防止損傷的擴展，減小安全事故發生的概率［1-3］。在眾多SHM技術中，Lamb波能夠在大尺度板狀結構中遠距離傳輸，信號的衰減量小，并且對于結構中的損傷具有很高的靈敏度，因此被廣泛用于損傷監測［4-6］。

在基于Lamb波的損傷監測技術中，主要有損傷因子監測技術和利用損傷因子進行損傷成像的監測技術。損傷因子旨在提取結構在時域、頻域、時頻域及波數域的特征變化來衡量結構的健康狀態［7-9］。而損傷成像技術通過成像的方法表征結構中的損傷，快速地確定損傷在結構中的位置。因為損傷成像直觀、損傷信息多等優點，成為了許多學者研究的熱點。Wang等［10］提出了時間反轉成像法（time-reversal method，簡稱TR），通過時間反轉并重新激發傳感器陣列采集到的信號完成損傷的聚焦成像。這種方法不需要提前知道結構的基線信號就可以完成損傷定位。Zhao等［11-12］研究了損傷概率重構（reconstruction algorithm for probabilistic inspection of defects，簡稱RAPID）的方法，計算不同信號對于損傷位置的概率的貢獻度，重構出損傷的位置。Michaels等［13］使用時間延時-累加定位法精準地定位了鋁板上的缺口和腐蝕損傷，延時-累加（delay-andsum，簡稱DAS）定位法基于殘余信號的理論波達時間完成對損傷的定位成像。

對于復雜結構和傳感器距離邊界較近的情況，損傷散射的Lamb波遇到邊界會產生反射，實際得到的殘余信號會出現多個波包疊加的情況。當損傷靠近邊界時甚至會發生損傷信號和邊界反射波包混疊的情況，從而導致DAS定位法出現偏離和偽損傷。Shan等［14］針對復雜結構的DAS定位提出了一種自適應有效數據提取的方法，自適應地截取信號的波包來提高DAS定位的精度。然而當損傷散射波包和反射波在時域上重疊時，直接截掉反射波的同時也會截斷大部分的損傷波包時域信號，導致定位結果變形。在成像損傷誤判的問題上，Sharifkhodaei等［15］提出了一種改進的延時-累加方法（windowed energy arrival method，簡稱WEAM），對殘余信號能量的包絡加上對數正態分布窗函數來捕捉首個波峰作為損傷的實際散射信號，從而提升DAS成像的精度。James等［16］提出了最小方差無失真響應的方法，通過給傳統的DAS定位增加加權系數來減少噪聲的投影。Lu等［17］使用邊界反射系數來減少邊界反射信號造成的偽損傷投影，反射系數由換能器的分布和反射信號的強度來決定。但是這些方法都沒能完全提取信號中的反射波波包。另一方面，Lamb波由于頻散效應會使得波包變寬，時域信號會產生變形，如果直接使用殘余信號進行延時-累加定位，將使得損傷區域面積變大，降低定位的精度。傳統的信號處理方法，如希爾伯特變換、黃氏變換、短時傅里葉變換和小波變換都不能對Lamb波進行頻散補償。

針對以上反射波混疊和Lamb波頻散補償的問題，筆者研究了Lamb波波包混疊分離的結構損傷定位方法，提升DAS定位算法的成像精度。首先，對復雜邊界條件下反射波和損傷散射波混疊的問題，建立了含頻散效應的Lamb波波包混疊和函數模型；其次，提出了基于隱變量參數求解的波包分離方法，重構了每個波包的分布情況，消除了反射波投影引起的傳統DAS定位中的偽損傷；然后，對重構出的每個波包進行了頻散補償，抑制了波包隨傳播距離產生的變形，提高了定位成像的分辨率；最后，在飛機復合材料加筋壁板上實驗驗證了改進DAS方法的可靠性。從結果來看，波包分離法能夠解決反射波混疊和頻散補償的問題，提升損傷定位的精度。

1 含頻散效應的Lamb波波包混疊模型

Lamb波信號vac（t）會因為傳播時間變化發生頻散，導致波包參數發生改變。激勵信號vac（t）為高斯窗調制的窄帶信號，vac（t）的頻域表示為Vac（f）

其中：F為傅里葉變換。

不同波達時間τk的Lamb波波包信號的頻域和時域表達式為

其中：IF為逆傅里葉變換；rk為波達時間τk下Lamb波的傳播距離；Cp（f）為A0模式的相速度。

通過復Morlet小波ψ（t）提取頻散后信號vk（t）的包絡Φk（t）為

其中：ω0為小波的中心角頻率；γ為高斯寬度；abs為取模；a為尺度因子；CWT為小波變換。

圖1 殘余信號在無反射、有反射下的傳播路徑Fig.1 The residual signal in the non-reflective and reflective path

在不考慮邊界的板結構中（如圖1（a）所示），驅動器i到傳感器j之間的損傷散射信號包絡為Φij（t）。然而，在典型的航空加筋結構中（如圖1（b）所示），加筋引起殘余信號的反射，導致獲取的殘余信號Φij（t）中包含多個不同時延的損傷反射波包。實際的損傷散射信號Φij（t）是由直達波Φ0ij（t）和后續的反射波AkΦkij（t）混疊而成，如式（6）所示

其中：Ak為第k個反射波信號幅值系數；Φkij（t）為波達時間τk的反射波信號；實際的包絡信號Φij（t）為含頻散效應的Lamb波波包和函數。

在延時累加算法中，每個點（x，y）的像素值I（x，y）由散射信號包絡Φij（t）在延時tij（x，y）之后得到，tij（x，y）是任意位置散射點（x，y）在驅動器（xi，yi）到傳感器（xj，yj）路徑上的理論傳播時間，表達式為

其中：N為監測區域內的傳感器數量；toff為激勵偏置時間；cg為群速度。

混疊的時域信號投影到定位圖上會對定位結果產生影響。針對該問題，文獻［14］中采用了如下的自適應提取的方法對信號Φ進行截取，信號長度表達式top為

其中：tfa為首波峰的時間；tth為信號截取閾值。

但是直接截取不能有效地分離混疊的波包，需要對混疊波包函數Φij（t）進行解耦分離。

2 波包分離方法及數值驗證

2.1 含頻散效應波包的高斯擬合

以厚度為1 mm的碳纖維復合材料板為例，考慮A0模式的相速度曲線（如圖2所示），根據式（2）和（3）可以計算得到不同波達時間τk的Lamb波波包信號，如圖3所示。其中：圖3（a）為50 kHz的窄帶激勵下的響應信號，波達時間最大為850μs；圖3（b）為對應信號的包絡Φk（t），τk越大，信號包絡越寬。激勵的窄帶信號vac（t）的表達式為

其中：fc為激勵信號的中心頻率。

由式（10）可知，窄帶信號vac（t）是由高斯窗調制的，采用式（4）中的小波變換，可以推得激勵信號的包絡Φac（t）是一個嚴格的高斯函數

圖2 碳纖維復合材料板A0模式相速度曲線Fig.2 A0 mode phase velocity curve of CFRP panel

因此采用高斯函數對頻散后的波包信號進行擬合重構。單個波包的波達時間τk已知，用方差σk2來表示Φk的寬度，擬合的表達式為

其中：tr為信號Φk（t）的樣本；R為樣本數；τk為高斯函數的均值。

重構得到的頻散波包的包絡（如圖3（c）所示）和原始信號的包絡一致，并且前后包絡信號的相對平均誤差不超過0.1%，完成對單個頻散波包參數的估計。

波包的寬度隨波達時間的增大而增大，對波達時間τk=［50∶1∶1 000］μs的序列進行方差σk2的擬合，得到如圖4所示的差值函數σk2=σ（τk）2。代入式（12），單個頻散波包的函數式Φk=Φk（t|τk）只由參數波達時間τk決定。混疊波包模型Φ的表達式為

圖3 不同波達時間Lamb波Fig.3 Lamb wave at different arrival times

圖4 方差關于波達時間的插值函數σk2=σ(τk)2Fig.4 Interpolation functionσk2=σ(τk)2 of variance on arrival time

其中：θ=（A，τ）；A=（A0，…，A k，…，AK）；τ=（τ0，…，τk，…，τK）。

2.2 波包函數的分離補償方法

為了識別出波包混疊函數Φ對應最優的參數項θ=（A，τ），將模型的對數似然函數L（T|θ）作為優化的性能參數目標。對數似然函數L（T|θ）的計算公式為

其中：T={t1，…，tr，…，tR}為信號Φ（t）對應的樣本集。

采用隱變量迭代估計的方法求解極大似然估計L，輸入參數θ（0）=（A（0），τ（0））通過K-mean算法對樣本集T初始化得到，θ（0）作為第0次迭代的參數［18-19］。其中波包信號樣本T表示給定觀測變量的數據，此時反映觀測數據tr來自分波包函數Φk的數據是未知的，用隱變量zrk表示，其定義為

其 中：zrk組 成 隱 隨 機 變 量 的 數 據Z={z10，…，zrk，…，zRK}。

T和Z連在一起稱為完全數據，得到完全數據下的對數似然函數L*為

已知初始參數θ（0）后，每一次迭代的求解都分為E步和M步［20］。E步需要確定Q函數Q（θ，θ（l）），即完全數據下的對數似然函數L*關于觀測變量T和第l次 迭 代 的 參 數θ（l）下 對 未 觀 測 數 據Z的 期 望，由式（17）得到

其中：E為對Z的數學期望。

迭代的M步是求解在θ（l）下Q（θ，θ（l））對θ的極大值，得到新一輪的迭代的模型參數為

重復以上迭代，直到式（19）中對數似然函數L的相對變化量小于一個極小值β。L（l）為第l次迭代的對數似然函數的值，L（l+1）為第l+1次迭代的值。收斂條件為

通過最終收斂的參數θ重構出的每個波包函數Φk，隨著τk的增加，波包的寬度也在增加。此時根據圖4中得到的方差的插值函數σk2=σ（τk）2，定義頻散補償系數為

其中：σ（τk）為τk對應的標準差；σ（toff）為激勵信號的標準差；p0為固定頻散補償系數。

補償后的波包函數Φk的表達式為

其中：Akod和σkod為補償前的參數。

從Φk中篩選重組出原本信號Φ的直達波包Φ*，計算公式為

其中：Ath用來控制幅值的閾值，為信號包絡最大值Amax的0.5，篩選出能量較大的波包；top為控制波達時間的閾值，通過文獻［14］中自適應截取信號長度的方法計算得到（如式（9））。

篩選出波達時間較早的波包，得到分離反射波后的直達波波包Φ*。而式（9）中的閾值tth的計算公式為

其中：rmr為以陣列對角線為長軸且經過長邊中點橢圓的監測距離；cg為波的群速度。

2.3 Lamb波波包分離方法的數值驗證

為了驗證波包分離法對于混疊波包信號的重構效果，對數值模擬的混疊Lamb波信號進行波包分離重構。激勵的信號選取50 kHz的窄帶信號。隨機生成圖5（a）中6個不同波達時間τk和幅值系數Ak的Lamb波信號。由式（4）小波變換求解得到各個波包的原始信號包絡分布如圖5（b）所示。將圖5（c）中信號包絡對應的觀測樣本T作為算法的輸入，計算樣本的初始最大似然估計L（0）和新一輪的模型參 數θ（1），然 后 計 算 參 數θ（1）對 應 的 最 大 似 然 估 計L（1），直到重構參數θ（l）對應的L（l）收斂。

由最終重構得到的參數θ，得到分離后的波包分布Φk如圖5（d）所示。重構后的波包分布收斂于正確的波包初始參數，成功對混疊波包信號Φ實現了分離。

對比多組不同情況下的隨機仿真信號，對比輸入的混疊波包參數θ和分離重構得到的波包參數θ。重構前后參數θ=（τ，A）的多組誤差的平均值如表1所示。算法前后重構得到參數的誤差都在5%以內，波包分離法的參數估計的可靠性較好。

表1 重構前后參數θ的誤差Tab.1 Error of parameterθafter reconstruction %

圖5 波包分離的數值驗證Fig.5 Numerical validation of wave packet separation

3 實驗驗證

3.1 實驗建立

實驗的對象是飛機復合材料加筋壁板（如圖6所示）。結構尺寸為950 mm×1 000 mm×1 mm，8層纖維的鋪層方向是按照［0°/90°/45°/-45°］S的順序排布的。圖6展示了實驗的監測系統。其中NI-PXI-5412發波卡用于產生激勵信號，通過Trek-2100HF功率放大器給到通道切換電路，NI-PXI-5105采集傳感器的信號。

實驗中選擇的監測區域為圖6中ABCD四傳感器陣列圍成的虛線區域。監測區域的尺寸為220 mm×130 mm，傳感器距離加強筋邊界的距離小于40 mm，響應信號包含驗證的Lamb波反射。激勵信號選用的是窄帶5波峰信號。激勵中心頻率fc為50 kHz。在該頻率下，Lamb波的A0模式在信號中占據主導地位，其他模式的信號可以忽略不計。激勵信號通過功率放大器放大到50 V作用在復材加筋板上，每個通道的采樣頻率為1 ms/s。單損傷和雙損傷的半徑為10 mm，采用吸波介質來模擬孔狀損傷。

圖6 實驗系統和結構Fig.6 Experimental system and structure

3.2 單損傷實驗結果

在單個損傷的情況下，對監測區域ABCD進行掃描監測。傳感器B和C之間殘余信號的串擾部分已經置零（如圖7（a）所示）。利用傳統DAS定位方法（將波包信號帶入式（7）），殘余信號的包絡Φ（t）（如圖7（b）所示）投影得到的損傷定位結果如圖8（a）所示。其中：白色符號“×”表示損傷實際的位置；黑色符號“+”表示偽損傷的位置。損傷的實際位置是坐標（75，45），對應時域信號中330μs處的直達波。偽損傷的位置發生在左下角（40，0）上，偽損傷對應圖7（b）時域信號中400μs左右的強反射波。反射波的幅值甚至超過了直達波的幅值，導致定位結果在（40，0）處出現了偽損傷。

采用文獻［14］中自適應提取信號的方法，計算得到信號長度為380μs，截取時域信號波包（如圖7（b））。由于直達波（330μs）和反射波混疊（400μs）的原因，截取的殘余信號中仍包含一部分反射波（400μs），并且截取了混疊信號中的直達波信號，在帶入DAS算法（式（7））投影后，破壞原有的定位結果（如圖8（b））。

將圖7（b）中的殘余信號波包Φ（t）的樣本T作為算法的輸入，對波包進行隱變量概率模型重構，得到波包函數Φ的分布，并對分離后波包進行頻散補償，得到分離補償之后的各個波包分布（如圖7（c））。根據邊界距離設置反射波達到的時間閾值，最終分離反射波后得到直達波波包Φ*（圖7（d））。在對時域信號波包進行波包分離和頻散補償之后，得到定位結果如圖8（c）所示。定位結果在（75，45）左右與實際損傷位置符合，像素圖只有直達波包（330μs）對應的投影，這是因為分離了反射波包（400μs），所以消除了波包在定位圖上（40，0）處的投影。另一方面，由于頻散效應得到了一定的補償，損傷定位圖的分辨率得到了提高。

圖7 單損傷傳感器路徑B-C間波包分離Fig.7 Wave packet separation of single damage on sensor BC path

Lamb波波包混疊分離方法的目的是為了消除反射波的成像投影和頻散補償，提升定位成像的信噪比和空間分辨率。為了定量評估損傷定位成像的效果，定義了損傷圖像的信噪比SNR（dBs）為

其中：變量mean（Id）和mean（Ia）分別為損傷定位圖在實際損傷區域和健康區域內的圖像強度均值。

SNR的效果和信號處理、圖像處理中的一樣，SNR的值越高代表損傷圖像的分辨率越高，定位噪聲越小。

在單損傷的情況下，原始圖像對應的mean（Id）為0.953，mean（Ih）為0.586，得到波包分離前的信噪比SNR為2.11。改 進的DAS方 法 中，mean（Id）為0.920，mean（Ih）為0.413，得到SNR值為3.48。選取不同位置的5組損傷進行成像，得到的損傷圖像信噪比如表2所示，損傷定位成像的結果得到了提升。

圖8 單損傷信號定位結果Fig.8 Localization result of single damage signal

表2 損傷成像評估Tab.2 Damage imaging assessment dB

3.3 雙損傷實驗結果

對監測區域ABCD進行雙損傷定位實驗，損傷的實際位置在（150，90）和（140，30）。圖9（a）中的殘余信號在傳感器B和C之間，兩個損傷的散射波包的波達時間分別在200和220μs左右，由于波達時間距離太近，在時域上出現了大面積重疊。對雙損傷包絡信號Φ（t）定位的結果如圖10（a）所示，損傷的實際位置在（150，90）和（140，30）（對應200和220μs）。然而，相比于損傷1，損傷2的幅值由于分辨率太低，難以從健康區域的幅值中突顯出來。從時域信號上看，圖9（b）中2個損傷在時域上的波包距離非常近，受到頻散影響，殘余信號的波包變寬，時域中波包重疊量大，時域信號投影到像素圖中只有1個峰，導致了損傷的丟失。

圖9 雙損傷傳感器路徑B-C間波包分離Fig.9 Wave packet separation of double damage on sensor BC path

對信號的包絡Φ（t）進行波包分離（如圖9（c）所示），時域上200和220μs的波包被分離出來，篩選直達波后映射到圖10（b）中的定位圖上識別出了2個損傷。在此基礎上，對分離后的波包信號Φk進行頻散補償（如圖9（d）所示），得到圖9（e）中的直達波。因為波包分離和頻散補償的原因，使得時域上的200和220μs的波包重疊度變小，2個損傷對應的峰明顯分離，映射到圖10（c）的定位圖上，計算得到改進后的DAS成像的SNR為3.16。而原始圖像的SNR為2.68，雙損傷定位的分辨率和精度得到了提高。

圖10 雙損傷信號定位結果Fig.10 Localization result of double damage signal

4 結束語

為了解決飛機結構中復雜邊界條件對Lamb波產生反射而引起的損傷定位偏差，筆者提出了基于Lamb波波包混疊分離的損傷定位識別方法。根據Lamb波頻散特性，建立了含有頻散效應的Lamb波波包混疊模型，采用隱變量的概率估計方法實現了波包和函數的分離，并對重構出的波包和函數進行頻散補償和重組。以碳纖維復合材料板為例，通過數值仿真驗證了波包分離方法的有效性。在飛機復合材料加筋壁板中進行了單損傷和多損傷的定位實驗。結果表明，該方法可以分離直達波和反射波波包，并對各個直達波包進行頻散補償，結合自適應波包截取算法，實現了單損傷分辨率2.11～3.48，雙損傷2.68～3.16的提升，相比于傳統的延時累加定位方法具有更好的定位精度和抗混疊干擾的能力。