CN115901049A - 基于波形相似度和多频率综合的结构冲击成像定位方法 - Google Patents

Abstract

本发明公开了一种基于波形相似度和多频率综合的结构冲击成像定位方法，属于结构健康监测技术领域；该方法利用信号分解技术将冲击应力波信号分解成多频率信号分量，计算各频率信号分量的包络信号，对监测区域离散成若干网格，引入基于波形的相似度指数，将网格点的相似度指数作为网格点的像素值实现冲击成像，根据峰值信噪比判定冲击成像结果是否有效，将有效图像进行图像融合作为最终的冲击成像结果，将最终融合后的图像像素值最大值的位置作为预测的冲击位置。本发明冲击成像定位方法定位精度高，而且成像聚焦性好，具有较好的工程应用前景。

Description

基于波形相似度和多频率综合的结构冲击成像定位方法

技术领域

本发明属于结构健康监测领域，具体涉及一种基于波形相似度和多频率综合的结构冲击成像定位方法。

背景技术

壁板结构在航空航天、海洋工程、轨道交通和能源动力等领域都有着大量的应用。然而壁板结构在服役期间无法避免地会受到外界冲击载荷产生损伤，导致结构强度严重降低，危害结构健康安全。传统无损检测方法(超声、红外和涡流技术等)维护成本高，且存在漏检的风险。对壁板结构的冲击事件进行实时在线的监测可以识别冲击位置和冲击载荷历程，大大降低维护成本并为线下的结构安全维护提供了指导，对结构的健康监测具有重要意义。

冲击定位是冲击监测技术中的首要任务，冲击在壁板结构中激发的冲击应力波信号可以被压电、应变和加速度等传感器接收。通过对传感器接收到的冲击应力波信号进行特征提取，再结合特定算法可以实现冲击位置的识别。冲击成像定位方法既避免了非线性方程组求解的困难，且可以对结构上各个位置发生冲击事件的概率进行可视化，具有较好的工程应用价值，现有的冲击成像定位算法包括基于时差的冲击成像定位方法、基于相位合成的虚拟时间反转聚焦成像、多重信号分类和基于参考数据库的方法等。其中，基于时差的冲击成像定位方法的关键是要获取精确的波达时刻，受到噪声影响，很难获得精确的波达时刻，定位精度和图像分辨率均受到影响；基于相位合成的虚拟时间反转基于声学互易性原理，通过叠加窄带波的信号幅值来进行定位，受到噪声的干扰冲击成像中容易存在假象，无法获得较好的定位结果；多重信号分类基于信号方向向量与噪声子空间的正交性进行冲击定位，受噪声影响较小，但依赖于传感器线性密集阵列而消耗大量的传感器，且存在扫描盲区和计算度复杂等问题；基于参考数据库的冲击成像定位方法对于一个含有众多子结构的大型结构来说需要存储一个庞大的数据库，目前在实际运用中存在较大困难。此外，目前已有的方法大多基于单一特定频率信号进行冲击定位，在信号的信息利用率上存在不足。

发明内容

为了解决上述问题，本发明提出了一种基于波形相似度和多频率综合的结构冲击成像定位方法，无需校准传感器和建立庞大的数据库，信号的信息利用率较高，冲击位置预测精准。

本发明的技术方案如下：

一种基于波形相似度和多频率综合的结构冲击成像定位方法，包括如下步骤：

步骤1、标记传感器的位置坐标，记录冲击应力波信号；

步骤2、对冲击应力波信号进行分解，获取多频率信号分量；

步骤3、利用希尔伯特变换计算各传感器不同频率信号分量的包络信号；

步骤4、计算不同传感器各频率信号分量的相对时延；

步骤5、计算波形相似度指数，冲击成像；

步骤6、根据有效冲击成像判定结果，进行冲击位置预测。

进一步地，步骤1的具体过程如下：

步骤1.1、标记壁板结构上监测区域内的所有接收冲击应力波信号的传感器的位置坐标，传感器S_m的位置坐标记为(x_m,y_m)，m＝1,2,…,M；M为传感器总数；

步骤1.2、记录所有传感器接收到的冲击应力波信号，将传感器接收到的冲击应力波信号记为h_m(t)，m＝1,2,…,M。

进一步地，步骤2中采用小波变换信号分解技术将冲击应力波信号分解，获得多频率信号分量，具体过程如下：

步骤2.1、对于实数空间下冲击应力波信号h_m(t)的小波变换WT定义如下：

其中，a为尺度因子，b为平移因子；“< >”表示内积，ψ_a,b表示小波函数；R表示实数空间，“ˉ”表示共轭；对于复Morlet母小波函数ψ(·)，定义如下：

其中，i表示虚数单位，t表示时间，f_b表示复Morlet母小波的带宽，f_c为复Morlet母小波的中心频率；

步骤2.2、利用复Morlet小波变换对冲击应力波信号进行分解得到信号分量，该信号分量的中心频率与尺度因子之间的关系如下：

其中，f_n为小波变换分解得到的信号分量的中心频率，f_c为复Morlet母小波的中心频率，f_s为冲击应力波信号的采样频率；

通过调节尺度因子a、复Morlet母小波的中心频率f_c和带宽f_b，根据式(1)-(3)分解得到各传感器多频率信号分量

K为不同频率的总数。

进一步地，步骤3中，求各传感器的不同频率信号分量的包络信号，如下：

其中，

为

的希尔伯特变换。

进一步地，步骤4的具体过程如下：

步骤4.1、将监测区域离散成若干网格；假设网格点(x,y)为冲击源位置，它到第m个传感器的距离为：

步骤4.2、计算不同传感器到冲击源位置的距离差d_mn(x,y)：

d_mn(x,y)＝R_n(x,y)-R_m(x,y)m,n＝1,2,..,M   (6)

其中，n表示第n个传感器；

步骤4.3、计算不同传感器的各频率信号分量的相对时延

其中，C_k为冲击应力波第k个频率信号分量的波速。

进一步地，步骤5的具体过程如下：

步骤5.1、对于网格点(x,y)，根据不同传感器之间各频率信号分量的时延，以包含冲击应力波中直达波信号的时长τ的时间窗口截取包络信号计算余弦相似度

其中：

其中，L为时长为τ的包络信号的总采样点数，l表示第l个采样点；

步骤5.2、对于第k个频率信号分量，定义网格点(x,y)的相似度指数I_k(x,y)如下：

步骤5.3、搜索所有网格点的相似度指数，将相似度指数作为像素值进行冲击成像。

进一步地，步骤6的具体过程如下：

步骤6.1、利用第k个频率信号分量冲击成像结果的峰值信噪比PSNR_k判断该频率信号分量冲击成像结果是否有效：

其中，δ为图像中的某一个像素，max[I_k(δ)]表示第k个频率信号分量冲击成像中的最大像素值，D_k是I_k的子集，由幅值小于某一阈值的像素值组成；Q_k(δ)＝I_k(δ|δ∈D_k)，N为子域D_k中的样本总数；

当某个频率信号分量冲击成像结果的PSNR_k大于给定的阈值，则该冲击成像结果视为有效的冲击成像结果，反之为无效的冲击成像结果；

步骤6.2、图像融合，冲击位置预测；

丢弃无效的冲击成像结果，选择有效的冲击成像结果的乘积作为最终的冲击成像结果I：

其中，v为有效冲击成像的编号；

步骤6.3、最终成像结果I中最大像素的位置代表预测的冲击位置。

本发明所带来的有益技术效果：

(1)本发明利用信号分解技术结合波形相似度进行冲击定位，与已有的基于参考数据库结合波形相似度的方法相比，无需对传感器校准，无需建立庞大的数据库，实施过程更加简单，有利于冲击监测系统地集成。

(2)本发明利用信号分解技术将冲击应力波信号分解，获取多频率信号分量进行特征信息提取，噪声信号幅度在信号分量中降低，有效信号幅度增加，既降低了噪声对冲击定位带来的干扰，又弥补了基于单一特定频率信号进行特征提取的特征信息利用率低的不足，进而提高冲击定位的精度。

(3)本发明利用多频率信号分量分别进行冲击成像，根据峰值信噪比判定有效和无效的冲击成像结果，丢弃无效的冲击成像结果，把有效的冲击成像的图像进行融合能够压制噪声带来的假象，提高了冲击成像的聚焦性及分辨率。

附图说明

图1为本发明基于波形相似度和多频率综合的结构冲击成像定位方法的流程图；

图2为本发明冲击定位原理示意图；

图3为本发明实施例中复合材料加筋壁板及传感器布置示意图；

图4为本发明实施例中复合材料加筋板的典型冲击应力波信号及频谱示意图；

图5为本发明实施例中小波变换将典型冲击应力波信号分解获得的多频率信号分量示意图；

图6为本发明实施例中不同传感器同一频率信号分量对比示意图；

图7为本发明实施例中不添加白噪声的各频率信号分量冲击成像结果图；

图8为本发明实施例中各频率信号分量添加20dB(10％)白噪声的冲击成像结果图；

图9为本发明实施例中添加白噪声的各频率信号分量冲击成像结果的峰值信噪比示意图；

图10为本发明实施例中多频率综合图像融合的冲击定位结果示意图。

具体实施方式

下面结合附图以及具体实施方式对本发明作进一步详细说明：

如图1所示，本发明提出了一种基于波形相似度和多频率综合的结构冲击成像定位方法，具体包括如下步骤：

步骤1、标记传感器的位置坐标，记录冲击应力波信号；

步骤1.1、标记壁板结构上监测区域内的所有接收冲击应力波信号的传感器(包括压电传感器，应变传感器和加速度传感器等)的位置坐标。

以压电传感器为例，如图2所示，对于布置了M个压电传感器的壁板结构监测区域，建立直角坐标系，图中T为冲击源，S₁、S_m、S_n、S_M分别表示不同位置的传感器，例如传感器S_m的位置坐标记为(x_m,y_m)。

步骤1.2、记录所有传感器接收到的冲击应力波信号，即将传感器接收到的冲击应力波信号记为h_m(t)，m＝1,2,…,M。

步骤2、对冲击应力波信号进行分解，获取多频率信号分量；

利用信号分解技术(包括小波变换、S变换和经验模态分解等)将冲击应力波信号进行分解，获得多频率信号分量。

以小波变换为例，具体过程如下：

步骤2.1、对于实数空间下冲击应力波信号h_m(t)的小波变换WT定义如下：

其中，a为尺度因子，b为平移因子，“<>”表示内积，ψ_a,b表示小波函数，R表示实数空间，“ˉ”表示共轭，ψ(·)为母小波函数，本发明选取复Morlet母小波函数：

其中，i表示虚数单位，t表示时间，f_b表示复Morlet母小波的带宽，f_c为复Morlet母小波的中心频率；

步骤2.2、利用复Morlet小波变换对冲击应力波信号进行分解得到信号分量，该信号分量的中心频率与尺度因子之间的关系如下：

其中，f_n为小波变换分解得到的信号分量的中心频率，f_c为复Morlet母小波的中心频率，f_s为冲击应力波信号的采样频率；

所以，通过调节尺度因子a、复Morlet母小波的中心频率f_c和带宽f_b，根据式(1)-(3)即可分解得到各传感器多频率信号分量

K为不同频率的总数。

步骤3、利用希尔伯特变换计算各传感器不同频率信号分量的包络信号。

求各传感器的不同频率信号分量的包络信号，如下：

其中，

为

的希尔伯特变换。

步骤4、计算不同传感器各频率信号分量的相对时延；

将监测区域离散成若干网格，假设网格点(x,y)为冲击源位置，根据传感器到网格点(x,y)的距离和已获得的各频率波速，计算不同传感器间各频率信号分量的相对时延。具体过程如下：

步骤4.1、将监测区域离散成若干网格，假设网格点(x,y)为冲击源位置，它到第m个传感器的距离为：

步骤4.2、计算不同传感器到冲击源位置的距离差d_mn(x,y)：

d_mn(x,y)＝R_n(x,y)-R_m(x,y)m,n＝1,2,..,M   (6)

其中，n表示第n个传感器。

步骤4.3、计算不同传感器各频率信号分量的相对时延

其中，C_k为冲击应力波第k个频率的波速。

步骤5、计算波形相似度指数，冲击成像；

由于不同传感器接收到的信号都来自同一冲击源，各传感器同一频率成分存在相似性，根据相对时延以合适时长的时间窗口截取包络信号的离散时间序列来计算余弦相似度，引入波形相似度指数，利用多频率信号分量对壁板结构监测区域分别进行冲击成像。

由于每个传感器接收到的冲击应力波信号都来自于同一冲击源，不同传感器冲击应力波同一信号分量直达波信号存在相似性，只不过因为不同传感器到冲击源的距离不同导致时延不同，对于两个不同传感器的同一频率信号分量的包络信号，将其相对时延补偿后的离散时间序列看作两个向量，利用向量的相似度计算方法(包括余弦相似度、欧几里得距离和皮尔逊相关系数等)计算相似度会得到最大的相似度值。所以，具体过程如下：

步骤5.1、对于网格点(x,y)，根据不同传感器之间各频率信号分量的时延，本发明以包含直达波信号的时长τ的时间窗口截取包络信号计算余弦相似度

其中：

其中，L为时长为τ的包络信号的总采样点数，l表示第l个采样点。

步骤5.2、对于第k个频率信号分量，定义网格点(x,y)相似度指数I_k(x,y)如下：

步骤5.3、搜索所有网格点的相似度指数，将相似度指数作为像素值即可实现冲击成像。

步骤6、根据有效冲击成像判定结果，进行冲击位置预测；

根据峰值信噪比确定无效和有效的冲击成像结果，判断标准为：当利用某个频率信号分量冲击成像的峰值信噪比大于给定阈值时，该频率信号分量冲击成像的结果视为有效的冲击成像结果；反之，为无效的冲击成像结果。丢弃无效的冲击成像结果，将所有有效的冲击成像结果的乘积作为最终的冲击成像结果，图像像素最大的位置为预测的冲击位置。具体过程如下：

步骤6.1、利用第k个频率信号分量冲击成像结果的峰值信噪比PSNR_k判断该频率信号分量冲击成像结果是否有效：

其中，δ为图像中的某一个像素，max[I_k(δ)]表示第k个频率信号分量冲击成像中的最大像素值，D_k是I_k的子集，由幅值小于某一阈值的像素值组成。Q_k(δ)＝I_k(δ|δ∈D_k)，N为子域D_k中的样本总数。

当某个频率信号分量冲击成像结果的PSNR_k大于给定的阈值，则该冲击成像结果视为有效的冲击成像结果，反之为无效的冲击成像结果。

步骤6.2、图像融合，冲击位置预测；

丢弃无效的冲击成像结果，选择有效的冲击成像结果的乘积作为最终的冲击成像结果I，以此来提高冲击定位精度，并加强冲击成像在冲击源位置的聚焦性，压制噪声所引起的假象：

其中，v为有效冲击成像的编号；

步骤6.3、最终成像结果I中最大像素的位置代表预测的冲击位置。

实施例

为了验证本发明方法的有效性，在如图3所示的总体尺寸为2360mm×1260mm，包含了5根横向加强筋和4根纵向加强筋的碳纤维复合材料加筋壁板结构上进行落锤冲击监测实验。蒙皮铺层顺序为[45/-45/-45/90/45/0/90/0/90/0/45/90/-45/-45/45]，筋条铺层顺序为[45/0/0/-45/90/0/90/-45/0/0/45]。在加筋壁板的蒙皮背面布置了48个压电传感器(S1-S48分别对应1号-48号传感器)，传感器间距为300mm×210mm。以蒙皮正面左下角为原点建立直角坐标系，标记传感器坐标。两端进行固支。利用30通道冲击监测系统实时采集1至30号传感器(S1-S30)的信号，采样率为200kHz。如图4所示，为在加筋壁板上坐标为(649mm，735mm)的位置5焦耳冲击能量落锤冲击，9号压电传感器获取的典型冲击应力波信号(左侧)及傅里叶频谱(右侧)，由频谱可知冲击应力波信号的主要频率成分在0-40kHz。如图5所示，为利用复Morlet小波变换将典型冲击应力波信号分解得到的4kHz，7kHz，10kHz，20kHz，30kHz和40kHz频率从低到高的六个频率信号分量。如图6所示，为证明方法原理的可行性，在16号传感器上方坐标为(649mm，525mm)的位置冲击，利用复Morlet小波变换分解得到15号和17号传感器冲击应力波信号的4kHz，7kHz，10kHz，20kHz，30kHz和40kHz频率的信号分量，对比可以发现两个传感器不同频率信号分量直达波信号高度相似，所以说明了本发明的原理是可行的。

为了进一步验证本发明定位算法，如图7所示，在不施加白噪声的情况下，为利用4kHz，7kHz，10kHz，20kHz，30kHz和40kHz六个频率信号分量分别对坐标为(649mm，525mm)的位置按照步骤2-步骤5进行冲击成像的结果，图中“○”为实际冲击源位置，“×”为本发明方法估计的冲击源位置，对比显示各频率信号分量均能获得较为准确的冲击定位结果。如图8所示，为对信号施加20dB噪声后利用六个频率信号分量分别冲击成像结果，4kHz，7kHz，10kHz的信号分量冲击成像结果仍能准确的估计冲击位置，而20kHz，30kHz，40kHz的信号分量冲击成像结果已经不能准确的估计冲击位置，以冲击成像像素大于等于最大值的95％的元素为子区域计算6号传感器各频率信号分量冲击成像结果的峰值信噪比，结果如图9所示，以-95dB为峰值信噪比阈值，4kHz，7kHz，10kHz信号分量冲击成像的峰值信噪比大于阈值，为有效的冲击成像，20kHz，30kHz，40kHz信号分量冲击成像的峰值信噪比小于阈值为无效的冲击成像。将有效的冲击成像结果进行融合，结果如图10所示，聚焦性及分辨率好，定位误差仅有1.2cm。所以，本发明算法在噪声环境下仍能准确定位，具有较好的工程应用价值。

当然，上述说明并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的技术人员在本发明的实质范围内所做出的变化、改型、添加或替换，也应属于本发明的保护范围。

Claims (7)

1.一种基于波形相似度和多频率综合的结构冲击成像定位方法，其特征在于，包括如下步骤：

步骤1、标记传感器的位置坐标，记录冲击应力波信号；

步骤2、对冲击应力波信号进行分解，获取多频率信号分量；

步骤3、利用希尔伯特变换计算各传感器不同频率信号分量的包络信号；

步骤4、计算不同传感器各频率信号分量的相对时延；

步骤5、计算波形相似度指数，冲击成像；

步骤6、根据有效冲击成像判定结果，进行冲击位置预测。

2.根据权利要求1所述基于波形相似度和多频率综合的结构冲击成像定位方法，其特征在于，所述步骤1的具体过程如下：

步骤1.1、标记壁板结构上监测区域内的所有接收冲击应力波信号的传感器的位置坐标，传感器S_m的位置坐标记为(x_m,y_m)，m＝1,2,…,M；M为传感器总数；

步骤1.2、记录所有传感器接收到的冲击应力波信号，将传感器接收到的冲击应力波信号记为h_m(t)，m＝1,2,…,M。

3.根据权利要求1所述基于波形相似度和多频率综合的结构冲击成像定位方法，其特征在于，所述步骤2中采用小波变换信号分解技术将冲击应力波信号分解，获得多频率信号分量，具体过程如下：

步骤2.1、对于实数空间下冲击应力波信号h_m(t)的小波变换WT定义如下：

其中，a为尺度因子，b为平移因子；“< >”表示内积，ψ_a,b表示小波函数；R表示实数空间，“ˉ”表示共轭；对于复Morlet母小波函数ψ(·)，定义如下：

其中，i表示虚数单位，t表示时间，f_b表示复Morlet母小波的带宽，f_c为复Morlet母小波的中心频率；

步骤2.2、利用复Morlet小波变换对冲击应力波信号进行分解得到信号分量，该信号分量的中心频率与尺度因子之间的关系如下：

其中，f_n为小波变换分解得到的信号分量的中心频率，f_c为复Morlet母小波的中心频率，f_s为冲击应力波信号的采样频率；

通过调节尺度因子a、复Morlet母小波的中心频率f_c和带宽f_b，根据式(1)-(3)分解得到各传感器多频率信号分量

k＝1,2,…,K,K为不同频率的总数。

4.根据权利要求1所述基于波形相似度和多频率综合的结构冲击成像定位方法，其特征在于，所述步骤3中，求各传感器的不同频率信号分量的包络信号

如下：

其中，

为

的希尔伯特变换。

5.根据权利要求1所述基于波形相似度和多频率综合的结构冲击成像定位方法，其特征在于，所述步骤4的具体过程如下：

步骤4.1、将监测区域离散成若干网格；假设网格点(x,y)为冲击源位置，它到第m个传感器的距离为：

步骤4.2、计算不同传感器到冲击源位置的距离差d_mn(x,y)：

d_mn(x,y)＝R_n(x,y)-R_m(x,y) m,n＝1,2,..,M    (6)

其中，n表示第n个传感器；

步骤4.3、计算不同传感器的各频率信号分量的相对时延

其中，C_k为冲击应力波第k个频率信号分量的波速。

6.根据权利要求1所述基于波形相似度和多频率综合的结构冲击成像定位方法，其特征在于，所述步骤5的具体过程如下：

步骤5.1、对于网格点(x,y)，根据不同传感器之间各频率信号分量的时延，以包含冲击应力波中直达波信号的时长τ的时间窗口截取包络信号计算余弦相似度

其中：

其中，L为时长为τ的包络信号的总采样点数，l表示第l个采样点；

步骤5.2、对于第k个频率信号分量，定义网格点(x,y)的相似度指数I_k(x,y)如下：

步骤5.3、搜索所有网格点的相似度指数，将相似度指数作为像素值进行冲击成像。

7.根据权利要求1所述基于波形相似度和多频率综合的结构冲击成像定位方法，其特征在于，所述步骤6的具体过程如下：

步骤6.1、利用第k个频率信号分量冲击成像结果的峰值信噪比PSNR_k判断该频率信号分量冲击成像结果是否有效：

其中，δ为图像中的某一个像素，max[I_k(δ)]表示第k个频率信号分量冲击成像中的最大像素值，D_k是I_k的子集，由幅值小于某一阈值的像素值组成；Q_k(δ)＝I_k(δ|δ∈D_k)，N为子域D_k中的样本总数；

当某个频率信号分量冲击成像结果的PSNR_k大于给定的阈值，则该冲击成像结果视为有效的冲击成像结果，反之为无效的冲击成像结果；

步骤6.2、图像融合，冲击位置预测；

丢弃无效的冲击成像结果，选择有效的冲击成像结果的乘积作为最终的冲击成像结果I：

其中，v为有效冲击成像的编号；

步骤6.3、最终成像结果I中最大像素的位置代表预测的冲击位置。

Images