杆状一维构件在工业建筑中应用很广。在这些杆状构件的生产过程中，这类构件主要起承重作用，面临的压力较大，如果存在瑕疵很容易破损，酿成事故。检测过程不能对结构体造成损坏，因此对杆状构件进行无损检测就具有重要意义^[1]。

时间反转法从光学中的相位共轭演变而来^[2]，能实现先到的信号后发，后到的信号先发。由于各信号能够同时同相到达目标位置，实现时间与空间上的聚焦，此项技术可以在不损坏被测物体的基础上实现无损检测。我国关于时反的研究取得了极具价值的研究成果。中科院声学研究所的汪承灏院士等^[3]对各种介质下的时反进行了研究，证明其聚焦能力。汪承灏、魏炜等^[4]对固体时反进行了科研，验证了时反抑噪能力与固体介质中的聚焦能力。邓菲等^[5]率先进行了先进性的管道聚焦实验。

本文提出了基于小波包分析的时间反转聚焦算法的改进算法，在噪声的情况下，有效地分离出奇异点，实现损伤信号的聚焦和重构，并能提高监测系统的信噪比。

1 小波去噪

含噪信号分解后^[6]，信号和噪声的小波系数的幅值将会不同。通过阈值的设定可以实现对噪声和信号的分离。小于阈值的系数去除，大于阈值的系数保留，所保留的即为去噪后小波系数。将保留下来的的有用信号的小波系数进行重构，即可获得去噪后的信号。

假设有如下观测信号：

$f\left( t \right) = s\left( t \right) + n\left( t \right)$

式中： $f\left( t \right)$ 表示原始信号， $n\left( t \right)$ 是方差为 ${\sigma ^2}$ 的高斯白噪声，服从正态分布 $N\left( {0,{\sigma ^2}} \right)$ 。

现有一维信号 $f\left( t \right)$ 。首先对信号进行离散采样，可以得N点离散信号 $f\left( n \right)$ ， $n = 0,1, \cdots ,N - 1$ ，则 $f\left( t \right)$ 小波变换系数为

${w_f}\left( {j,k} \right) = {2^{\frac{j}{2}}}\sum\limits_{n = 0}^{N - 1} {f\left( n \right)} \psi \left( {{2^j}n - k} \right)$

在实际应用中，通过双尺度变换得到递归方程。

$\begin{array}{*{20}{c}} {{s_f}\left( {j + 1,k} \right) = {s_f}\left( {j,k} \right) {\cdot} h\left( {j,k} \right)} \\ {{w_f}\left( {j + 1,k} \right) = {s_f}\left( {j,k} \right) {\cdot} g\left( {j,k} \right)} \end{array}$

式中： ${S_f}\left( {0,k} \right)$ 代表原始信号 $f\left( k \right)$ 、 ${S_f}\left( {j,k} \right)$ 代表j尺度上的逼近系数， ${w_f}\left( {j,k} \right)$ 则代表小波系数。那么，小波变换重构公式可以写成

$\begin{gathered} {S_f}\left( {j - 1,k} \right) = {S_f}\left( {j,k} \right) {\cdot} {\overset{\frown} h} \left( {j,k} \right) + {w_f}\left( {j,k} \right){\cdot}{\overset{\frown} g} \left( {j,k} \right) \\ \end{gathered} $

当混有白噪声时：

${\omega _{j,k}} = {u_{j,k}} + {v_{j,k}}$

阈值的选择^[7]：

$T = \sigma \sqrt {2\ln n} /\exp (\left( {l - 1} \right)/2)$

阈值函数的选择：

${{\overset{\frown} \omega} _{j,k}} = \left\{ {\begin{array}{*{20}{c}} {{\omega _{j,k}} + T - \displaystyle\frac{{2T}}{{\exp \left( {\left| {\displaystyle\frac{{{\omega _{j,k}}}}{T}} \right| - 1} \right) + 1 + n}}},\;\;{{\omega _{j,k}} \leqslant - T} \\ {\displaystyle\frac{{2\operatorname{sgn} \left( {{\omega _{j,k}}} \right){{\left| {{\omega _{j,k}}} \right|}^{n + 1}}}}{{\left[ {\exp \left( {\left| {\displaystyle\frac{{{\omega _{j,k}}}}{T}} \right| - 1} \right) + 1 + n} \right]{T^n}}}},\;\;{\left| {{\omega _{j,k}}} \right| < T} \quad\quad\quad\\ {{\omega _{j,k}} - T + \displaystyle\frac{{2T}}{{\exp \left( {\left| {\displaystyle\frac{{{\omega _{j,k}}}}{T}} \right| - 1} \right) + 1 + n}}},\;\;{{\omega _{j,k}} \geqslant T} \end{array}} \right.$

小波阈值去噪的基本步骤如下：

1）利用小波变换分解含噪信号 $f\left( k \right)$ ，可得到一组小波系数 ${\omega _{j,k}}$ ；

2）用阈值处理小波系数 ${\omega _{j,k}}$ ，确定小波系数的估计值，使 $\left\| {{{{\overset{\frown} \omega} }_{j,k}} - {\omega _{j,k}}} \right\|$ 最小；

3）用小波逆变换对 ${{\overset{\frown} \omega} _{j,k}}$ 进行重构，得到估计的损伤信号 ${\overset{\frown} f} \left( k \right)$ ，即为去噪后的信号。

2 小波包分解

对小波重构信号 ${\overset{\frown} f} \left( k \right)$ 进行小波包分析，小波包分析是小波分析的改进^[8]，它使信号的低频分量和高频分量都得以保留，处理更细致，分辨率也更高。

${L^2}\left( R \right) = \mathop \oplus \limits_{j \in Z} {W_j}$ 是 ${L^2}\left( R \right)$ 按照不同空间尺度j下分解而成的子空间 ${W_j}\left( {j \in Z} \right)$ 的正交和。

按二进制细分频域中的 ${W_j}$ ，用新的子空间 $M_j^n$ 统一表征尺度子空间 ${V_j}$ 和小波子空间 ${W_j}$ ，若令

$\left\{ {\begin{array}{*{20}{c}} \begin{gathered} M_j^0 = {V_j} \\ M_j^1 = {W_j} \\ \end{gathered} &({j \in {\bf Z}} )\end{array}} \right.$

则Hilbert空间的正交分解 ${V_{j + 1}} = {V_j} \oplus {W_j}$ 可用 $M_j^n$ 分解为

$\begin{array}{*{20}{c}} {M_{j + 1}^0 = M_j^0 \oplus M_j^1}&({j \in {\bf Z}}) \end{array}$

若 $M_j^n$ 为 ${m_n}\left( t \right)$ 的闭包空间，则 ${m_n}\left( t \right)$ 满足双尺度方程：

$\left\{ {\begin{array}{*{20}{c}} {{m_{2n}}\left( t \right) = \sqrt 2 \displaystyle\sum\limits_{k \in Z} {h\left( k \right)} {m_n}\left( {2t - k} \right)} \\ {{m_{2n + 1}}\left( t \right) = \sqrt 2 \displaystyle\sum\limits_{k \in Z} {g\left( k \right)} {m_n}\left( {2t - k} \right)} \end{array}} \right.$

(1)

式中 $\left\{ {{m_n}\left( t \right)} \right\}\left( {n \in {Z_ + }} \right)$ 是由基函数在 ${m_0}\left( t \right) = \varphi \left( t \right)$ 的条件下的正交小波包，并在 $g\left( k \right) = {\left( { - 1} \right)^k}h\left( {1 - k} \right)$ 条件下，2个系数成正交关系。 $n = 0$ 时，式（1）可作进一步的简化：

$\left\{ {\begin{array}{*{20}{c}} {{m_0}\left( t \right) = \displaystyle\sum\limits_{k \in Z} {{h_k}} {m_0}\left( {2t - k} \right)} \\ {{m_1}\left( t \right) = \sqrt 2 \displaystyle\sum\limits_{k \in Z} {{g_k}} {m_0}\left( {2t - k} \right)} \end{array}} \right.$

(2)

在多分辨率分析中 $\varphi \left( t \right)$ 和 $\psi \left( t \right)$ 的双尺度方程为

$\left\{ {\begin{array}{*{20}{c}} {\varphi \left( t \right) = \displaystyle\sum\limits_{k \in Z} {{h_k}\varphi \left( {2t - k} \right),{{\left\{ {{h_k}} \right\}}_{k \in Z}}} \in {l^2}} \\ {\psi \left( t \right) = \displaystyle\sum\limits_{k \in Z} {{g_k}\varphi \left( {2t - k} \right),{{\left\{ {{g_k}} \right\}}_{k \in Z}}} \in {l^2}} \end{array}} \right.$

(3)

二者联立即可将 ${m_0}\left( t \right)$ 和 ${m_1}\left( t \right)$ 变为尺度函数 $\varphi \left( t \right)$ 和小波函数 $\psi \left( t \right)$ 。式（2）为式（3）的等价形式。把这种等价表示推广到（非负整数），可得到式（1）的等价形式：

$M_{j + 1}^n = M_j^n \oplus M_j^{2n + 1}(j \in {\bf Z},n \in {{\bf Z}_ + })$

小波包重构：

$g_j^n\left( t \right) \in M_j^n,g_j^n\left( t \right) = \sum\limits_l {d_l^{j,n}} {m_n}\left( {{2^j}t - l} \right)。$

3 时反聚焦

随机噪声设为 $n\left( t \right)$ 且系统噪声与信号不相干，若波源S到损伤散射点D之间的系统传递响应函数为

$h\left( t \right) = \sum\limits_{i = 1}^\infty {{A_i}} \left( {t - {\tau _i}} \right)$

则根据声互易性原理^[9]，因为模型中激励元与接收元在同端，用信号处理技术模拟信号的接收（即被动时间反转技术）^[10]。首波 $h\left( t \right)$ 作用信号2次，反射系数为 ${r_c}$ ，最终接收阵元接收信号：

$y\left( t \right) = {r_c}s\left( t \right) \text{·} h\left( t \right) + n\left( t \right)$

(4)

将式（4）时反后再次发射，第2次接收的信号可表示为

$y\left( t \right) = {r_c}s\left( t \right)\cdot h\left( t \right) \cdot h\left( t \right) + n\left( t \right)$

(5)

将式（5）时反后再次发射，第二次接收的信号可表示为

$y\left( t \right) = {r_c}y\left( { - t} \right) \cdot h\left( t \right) \cdot h\left( t \right) + n\left( t \right)$

(6)

对式（6）作傅里叶变换得到

$\begin{aligned} z\left( \omega \right) =& {r_c}{Y^ * }\left( \omega \right)H\left( \omega \right)H\left( \omega \right) + N\left( \omega \right) = \\ & r_c^2 \cdot X \cdot \left( \omega \right){\left| {H\left( \omega \right)} \right|^4} + r_c^2{N^ * }\left( \omega \right){\left| {H\left( \omega \right)} \right|^2} + N\left( \omega \right) \\ \end{aligned} $

(7)

式（7）可推导出 $H\left( \omega \right) > 1$ 时，信号将被增强，出现聚焦相关主峰，噪声得到了抑制和削弱。

4 改进算法

利用小波阈值去噪算法和小波包分解算法增强与时间反转算法相结合，增强了时间反转算法对桩基损伤的定位效果以及聚焦抑噪效果，具体的过程如下：

1）在symN小波基上对检测信号 $y\left( t \right)$ 进行小波分解；对分解后的小波系数进行阈值去噪处理；

2）小波重构，获得去噪信号 $s\left( t \right)$ ；

3）对去噪后的信号 $s\left( t \right)$ 进行小波包分解，从分解后的信号中得到M个奇异点，截取时宽为 ${T_{\rm{i}}} = {t_{i2}} - {t_{i1}}$ ， $i = 1, 2, \cdots, M$ ， $i$ 为第 $i$ 个奇异点；

4）桩基应力波检测信号数据和激励信号的函数表达都是已经确定的，分别为 $y\left( t \right)$ 和 $x\left( t \right)$ 。将 $y\left( t \right)$ 和 $x\left( t \right)$ 变换到频域，变换后的信号分别为 $Y\left( \omega \right)$ 和 $X\left( \omega \right)$ ，则桩基检测信号的信道冲击响应函数为

$H\left( \omega \right) = Y\left( \omega \right)/X\left( \omega \right)$

5）用宽度为 ${T_{{i}}}$ 的时间窗对信号进行截取，窗内乘一，窗外置零，得到M个截取信号，时域和频域分别为为 ${s_{ic}}\left( t \right)$ 和 ${S_{ic}}\left( \omega \right)$ 的信号；

6）将每个截取信号 ${s_{ic}}\left( t \right)$ 在整个信号时间T内进行时间反转，时反后得到的M个时域为 ${t_i}\left( t \right) = {s_{cr}}\left( {T - t} \right)$ 和 ${r_i}\left( t \right) = {S_{cr}}\left( \omega \right) \cdot \exp \left( { - {\rm i}{ \omega t}} \right)$ 的时反信号。

含噪信号经过小波包分解后，奇异点难以分辨；经过小波去噪后，使得奇异点容易分辨，有利于时间反转镜对损伤信号的聚焦，新的时反算法不但可以更好地识别损伤信号，而且使信号检测系统的信噪比有了进一步的提高。

5 实验与仿真

单音频信号为幅值 $A = 1\;000$ 、中心频率ω=0.5 kHz脉宽 $\tau = \omega /{\text π} $ 的窄带脉冲信号，其时域波形与频域波形如图1、2所示，利用MATLAB软件进行仿真：

$s\left( t \right) = A \cdot \sin \left( {\omega t} \right)$

在图2中 $[2,4]$ ms区间内也有一个峰值相对较大的波包，这个回波信号比较明显，但其他波形部分无法识别是否包含损伤信息，这是在桩基类构件中信号的衰减和多径效应共同作用下，检测回波信号发生时间拓展引发失真。

对检测信号进行3层sym6小波包分解，得到如图3所示的子波，图中完整低信号区间为 $[0,10.23]$ ms，纵坐标表示接受阵元处质点的运动幅值，可以看出在 $\left[ {1.13,1.34} \right]$ 、 $\left[ {1.29,1.48} \right]$ 、 $\left[ {1.59,1.83} \right]$ 、 $\left[ {2.00,2.26} \right]$ 、 $\left[ {3.72,3.88} \right]$ 、 $\left[ {4.62,4.76} \right]$ 、 $\left[ {4.99,5.14} \right]$ ms，奇异点在以上时区内出现，用矩形窗进行截取，使能量更为集中。

截取的损伤信号如图4所示。将yc1W~yc7W分别在时反窗内进行时间反转，为了更好地研究检测信号，利用时间反转算法进一步增强损伤回波信号。激励信号与首次回波信号已知，分别为 $X\left( \omega \right)$ 和 $Y\left( \omega \right)$ ，则估计出的虚拟信道 $H\left( \omega \right) = Y\left( \omega \right)/X\left( \omega \right)$ ，将时反后的信道虚拟重发到信道中去，得到的二次回波聚焦信号如图5所示，得到的平均聚焦增益18.496 6 dB，聚焦增益以及聚焦后的最大幅值如表1所示。

由表1可知，表中7处无噪损伤信号聚焦增益的平均值为18.496 6 dB。

当信号加入−0.67 dB高斯白噪声是时，首次回波信号为 ${y_s}\left( t \right) = s\left( t \right) \cdot h\left( t \right) + n\left( t \right)$ ，含噪回波信号如图6所示，小波包分解后的信号如图7所示，截取的信号以及聚焦信号分别如图8、9所示。

对 ${y_s}\left( t \right)$ 进行3层小波阈值去噪，之后进行小波包分解，分解后的回波信号如图10所示。

从去噪后小波分解的回波信号中，可以更有效地分辨出奇异点，对奇异点进行截取，信道传递函数仍为 $H\left( \omega \right)$ ，对截取的时反信号进行时反聚焦实验，所得结果如图10、11所示。损伤点不同处理过程的信噪如表2所示。

经计算， ${y_s}\left( t \right)$ 通过时间反转得到的增益为5.1 dB，经本文新处理方法处理后信噪比提高至9.34 dB。

6 结论

从以上实验结果可以得出如下结论：

1）本文基于小波分析的时间反转改进算法，可以在含有噪声的情况下，识别出损伤点；

2）可以对损伤信号在损伤点处进行聚焦，并具有更强的抑制噪声能力，对损伤信号的信噪比平均提升了9.4 dB；

3）对每个损伤点都具有聚焦、抑制噪声和提高信噪比的能力，改进算法可以对一维构件实测中损伤的模式识别与定位以及信号的处理有一定的益处。

今后可以利用迭代算法进一步提升系统的信噪比以及聚焦增益。