CN101793955B - 超声波首波检测与定位的方法和装置 - Google Patents

Abstract

本发明公开了一种超声波首波检测与定位的方法和装置，该方法包括：构造超声波信号序列的轨道矩阵，对轨道矩阵进行奇异值分解获取超声波信号序列的奇异谱熵，对奇异谱熵进行中值滤波，去除信号毛刺，对已经过中值滤波的奇异谱熵进行阈值化处理，则已经过阈值化处理的奇异谱熵的下降沿即为超声波信号序列的首波，在轨道矩阵合适的情况下，超声波信号序列与噪声区的奇异谱熵落差大，边缘陡峭，奇异谱熵对超声波信号序列的幅值具有很强的鲁棒性，奇异谱熵不会因为超声波信号序列的中某一片段的信噪比的强弱而有较大变化，因此，本发明的超声波首波检测与定位的方法抗干扰能力强，既使对于淹没在噪声中极弱的首波也能够较为正确的定位。

Description

超声波首波检测与定位的方法和装置

技术领域

本发明涉及一种信号处理技术领域，尤其涉及一种超声波首波检测与定位的方法和装置。

背景技术

首波是也称为头波，是超声和其它类型信号脉冲式检测领域中一个有效信号的开始。首波的可靠检测和准确定位对这类测量的精度有关键作用。首波提取的最大困难是其幅度一般很小，非常易于和系统噪声混淆，噪声的存在对其测量的准确性影响很大。而且一般情况下首波波幅随着信号频率的增大衰减很快，因此在频率较高的测量中，其首波检测更加困难。

目前，常用的首波检测方法主要有阈值法、加窗相关峰值法、长短时窗能量比法等。阈值法是设置一个门限(即阈值)，若信号幅值小于该门限则认为该信号是噪声信号，若信号幅值大于门限值，认为该信号是声波信号，而第一个大于该门限的信号即是首波。为防止个别噪声变异大于门限而导致误判，可采用连续多个信号幅值大于门限值来确认首波。加窗相关峰值法是检测波形与标准波形进行加窗相关对比，计算探索范围内的相关系数，最大相关系数的位置即是首波的波至点。长短时窗能量比法是在首波初至的检测范围内计算短步进时窗与长累积时窗的能量比，长累积时窗是第一个到当前位置处短步进时窗的并集，在首波到达前，长累积时窗的能量与短步进时窗的能量均很小，相差不大；在首波到达时，长累积时窗的能量仍很小，而短步进时窗的能量却突然变大，此时的长短时窗能量比出现大的跃变，需要选择合适的能量阈值，才能可以检测到弱初至首波。

上述几种方法在信号信噪比较强时，尤其是在首波附近的局部信噪比较强的时候，检测效果较好，反之，若首波附近局部区域的信噪比较小其检测效果将恶化，其中恶化最严重的是阈值法，其次是长短时窗能量比法。

发明内容

本发明的目的在于提供一种抗干扰能力、能准确检测首波位置的超声波首波检测与定位的方法和装置。

为达到上述目的，本发明一方面提供了一种超声波首波检测与定位的方法，包括以下步骤：

构造超声波信号序列的轨道矩阵；

对所述轨道矩阵进行奇异值分解获取所述超声波信号序列的奇异谱熵；

对所述奇异谱熵进行中值滤波，去除信号毛刺；

对已经过中值滤波的奇异谱熵进行阈值化处理，则已经过阈值化处理的奇异谱熵的第一个下降沿即为所述超声波信号序列的首波。

本发明的超声波首波检测与定位的方法，所述构造超声波信号序列的轨道矩阵，具体包括如下步骤：

对一个长度为N点的超声波信号序列，截取以n点开始长度为L点的超声波信号子序列，构造出所述超声波信号子序列的轨道矩阵A_n为：

$A_n=\left[\begin{array}{cccc}x\left(n\right)& x(n+1)& \·\·\·& x(n+L-M)\\ x(n+1)& x(n+2)& \·\·\·& x(n+L-M+1)\\ .& .& .& .\\ .& .& .& .\\ .& .& .& .\\ x(n+M-1)& x(n+M)& \·\·\·& x(n+L-1)\end{array}\right]$

其中，n＝1，2，……，N-L+1，L小于N而大于超声波信号序列的周期。

本发明的超声波首波检测与定位的方法，所述对轨道矩阵进行奇异值分解获取所述超声波信号序列的奇异谱熵，具体包括如下步骤：

对所述轨道矩阵A_n进行奇异值分解，获取其r个奇异值为σ_ni(i＝1，2，……r)；

根据公式：得到所述超声波信号序列的归一化奇异谱概率密度函数p_n(k)，其中，σ_nk、σ_nj均为奇异值，k＝1，2，…，r；

根据公式：

得到所述超声波信号序列的奇异谱熵，其中，p_n(i)为归一化奇异谱概率密度函数。

本发明的超声波首波检测与定位的方法，所述对已经过中值滤波的奇异谱熵进行阈值化处理，具体包括：

根据公式：Th＝a×[max(H(n))-min(H(n))]+min(H(n))，得到所述超声波信号序列的奇异谱熵的上升沿和下降沿，其中，max(H(n))、max(H(n)分别为所述超声波信号序列的奇异谱熵的最大值、最小值，a为比例因子。

本发明的超声波首波检测与定位的方法，在所述获取超声波信号序列的奇异谱熵之前，还包括：

根据所述超声波信号序列的短时能量对所述超声波信号序列进行预处理，得到所述超声波信号序列出现的概率区间。

另一方面，本发明还提供了一种超声波首波检测与定位的装置，所述超声波首波检测与定位的装置执行上述超声波首波检测与定位的方法。

本发明的超声波首波检测与定位的方法：首先构造超声波信号序列的轨道矩阵，然后对轨道矩阵进行奇异值分解获取超声波信号序列的奇异谱熵，接着，对奇异谱熵进行中值滤波，去除信号毛刺，最后，对已经过中值滤波的奇异谱熵进行阈值化处理，则已经过阈值化处理的奇异谱熵的下降沿即为超声波信号序列的首波，在轨道矩阵合适的情况下(即超声波信号子序列的长度大于超声波信号序列的周期而小于超声波信号序列的长度时)，超声波信号序列与噪声区的奇异谱熵落差大，边缘陡峭，而且奇异谱熵对超声波信号序列的幅值具有很强的鲁棒性，奇异谱熵不会因为超声波信号序列的中某一片段的信噪比的强弱而有较大变化，因此，本发明的超声波首波检测与定位的方法具有较强的抗干扰能力，既使对于淹没在噪声中极弱的首波也能够较为正确的定位，这对远距离测量时的弱信号波识别和定位至关重要。

附图说明

图1为本发明的超声波首波检测与定位的方法流程图；

图2为本发明的超声波首波检测与定位的装置的结构框图；

图3为信噪比为27.4dB的超声波信号序列的幅值图；

图4为信噪比为27.4dB的超声波信号序列的奇异谱熵图；

图5为信噪比为4.4dB的超声波信号序列的幅值图；

图6为信噪比为4.4dB的超声波信号序列的奇异谱熵图。

具体实施方式

下面结合附图对本发明的具体实施方式进行详细描述：

由于奇异谱熵具有如下特征：

1)信号的奇异谱熵不同于噪声的奇异谱熵。

2)理论上，如果谱的分布保持不变，信号幅值的大小不会影响奇异谱熵。而实际上，奇异谱熵会随信号随机性而变化，但与能量特征相比，奇异谱熵的变化很小。

3)奇异谱熵对噪声具有一定的稳健性，当信噪比(SNR)下降时奇异谱熵的形状保持不变。本发明的核心思想基于奇异值分解理论，将由信号构成的具有一定嵌入维度的轨道矩阵进行奇异值分解，再利用信息熵的统计特性对奇异值集合进行不确定度分析，从而对原始信号的复杂程度给出一个确定的量度。

对于一个L点的信号序列x(i)，根据相空间重构理论，选择一个小于L的正整数M，就可以由信号序列构造出轨道矩阵A，表示为：

$A=\left[\begin{array}{cccc}x\left(1\right)& x\left(2\right)& \·\·\·& x(L-M+1)\\ x\left(2\right)& x\left(3\right)& \·\·\·& x(L-M+2)\\ .& .& .& .\\ .& .& .& .\\ .& .& .& .\\ x\left(M\right)& x(M+1)& \·\·\·& x\left(L\right)\end{array}\right]$

式中的轨道矩阵也称为Hankle矩阵。在矩阵构造中，常把相邻两行/列的重叠元素个数称为嵌入维度，Hankle矩阵嵌入维度为M-1，也可以构造出其他嵌入维度的矩阵。

由矩阵理论知，对于一个m×n的矩阵A，则存在一个m×m的酉矩阵U和n×n的酉矩阵V，使得A可分解为：

A＝UPV^H

其中，上标“H”表示矩阵的共轭转置；当m＝n时，P＝∑；当m＞n时，

当m＜n时，

O_(m-n)×n为(m-n)×n的零矩阵，O_m×(n-m)为m×(n-m)的零矩阵，而∑＝diag(σ₁，σ₂，σ₃，……σr)，r＝min(m，n)。σ_i(i＝1，2，……，r)称为矩阵A的奇异值，为非负数。它们按照从大到小的顺序排列，即：

σ₁＞σ₂＞σ₃＞……＞σ_r-1＞σ_r

U和V分别称为左右奇异矩阵，而U＝[u₁，u₂，…，u_m]^T和V＝[v₁，v₂，…，v_n]^T中的列向量u_i和v_j分别称为A的左右奇异向量。可以证明，u_i是方阵AA^T的特征向量，v_j是方阵A^TA的特征向量；奇异值σ_i(i＝1，2，……，r)是方阵AA^T或A^TA的特征值λ_i(i＝1，2，……，r)的平方根。

基于上述思想，参考图1，本发明的超声波首波检测与定位的方法，包括以下步骤：

步骤S101、对超声波信号序列从某点开始截取一定长度的超声波信号子序列，构造其轨道矩阵。具体的，对一个长度为N点的超声波信号序列，截取以n点开始长度为L点的超声波信号子序列，构造出超声波信号子序列的轨道矩阵A_n为：

其中，n＝1，2，……，N-L+1，L小于N而大于超声波信号的周期，此外，为了能够缩小搜索范围，减少计算量，提高搜索效率在构造轨道矩阵之前可以先对超声波信号序列进行预处理，由超声波信号序列的短时能量初步得到信号可能出现的概率区间，在这些概率区间上再进行基于奇异谱熵的首波定位。

步骤S102、对轨道矩阵A_n进行奇异值分解，计算出其奇异谱熵，按照此方法依次获得整个超声波信号序列的奇异谱熵H(n)。

对轨道矩阵A_n进行奇异值分解，获取其r个奇异值为σ_ni(i＝1，2，……r)；

根据公式：

得到超声波信号序列的归一化奇异谱概率密度函数p_n(k)，其中，σ_nk、σ_ni均为奇异值，k＝1，2，…，r；

根据公式：

得到超声波信号序列的奇异谱熵，其中，p_n(i)为归一化奇异谱概率密度函数。需要说明的是，信号的奇异谱熵一般定义为：

$H\left(n\right)=-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\{p_n\left(i\right)\×\log [p_n\left(i\right)\left]\right\}$

为了避免某点的奇异谱概率密度为零而导致上式运算错误，把上式改为：

$H\left(n\right)=-\underset{i=1}{\overset{r}{\mathrm{\Σ}}}\{p_n\left(i\right)\×\log [p_n\left(i\right)+1\left]\right\}.$

步骤S103、对奇异谱熵H(n)进行中值滤波，去除信号毛刺。具体的，假设是对H(n)进中第i点数值的7点中值滤波，则把H(j)，(j＝i-3，i-2，i-1，i，i+1，i+2，i+3)，7个点进行从从大到小或者从小到大排列，中值滤波median(H(i))等于排列于中间位置上的数值。

步骤S104、对已经过中值滤波的奇异谱熵进行阈值化处理，找出已经过阈值化处理的奇异谱熵的下降沿，即为超声波信号序列的首波位置。具体的，

根据公式：Th＝a×[max(H(n))-min(H(n))]+min(H(n))，得到超声波信号序列的奇异谱熵的上升沿和下降沿，其中，max(H(n))、max(H(n)分别为超声波信号序列的奇异谱熵的最大值、最小值，a为比例因子。通过调整a可以调整阈值的高度，a越大阈值越大。

本领域的技术人员可以清楚地了解到本发明可以以软件产品的形式体现出来，该软件产品可以存储在一个非易失性存储介质(可以是CD-ROM，U盘，移动硬盘等)中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明实施例所述的方法。基于这样的理解，本发明的技术方案还可以通过硬件实现，也可以可借助软件加必要的通用硬件平台的方式来实现，如图2所示，本发明的超声波首波检测与定位的装置，包括轨道矩阵构造模块21，用于构造超声波信号序列的轨道矩阵；奇异谱熵获取模块22，用于对轨道矩阵构造模块21构造出的轨道矩阵进行奇异值分解获取超声波信号序列的奇异谱熵；中值滤波模块23，用于对奇异谱熵获取模块22得到的奇异谱熵进行中值滤波，去除信号毛刺；阈值处理模块24，用于对已经过中值滤波模块23处理后的奇异谱熵进行阈值化处理，则已经过阈值化处理的奇异谱熵的下降沿即为所述超声波信号序列的首波。

本发明的超声波首波检测与定位的方法：首先构造超声波信号序列的轨道矩阵，然后对轨道矩阵进行奇异值分解获取超声波信号序列的奇异谱熵，接着，对奇异谱熵进行中值滤波，去除信号毛刺，最后，对已经过中值滤波的奇异谱熵进行阈值化处理，则已经过阈值化处理的奇异谱熵的下降沿即为所述超声波信号序列的首波。在轨道矩阵合适的情况下(即超声波信号子序列的长度大于超声波信号序列的周期而小于超声波信号序列的长度时)，超声波信号序列与噪声区的奇异谱熵落差大，边缘陡峭，这正符合首波准确定位的需要，而且奇异谱熵对超声波信号序列的幅值具有很强的鲁棒性，奇异谱熵不会因为超声波信号序列的中某一片段的信噪比的强弱而有较大变化，如图3所示，信噪比为27.4dB的超声波信号序列，其奇异谱熵如图4所示，超声波信号序列与噪声区的奇异谱熵落差大，边缘陡峭；而如图5所示，信噪比为4.4dB的超声波信号序列，其频谱谱熵如图6所示，超声波信号序列与噪声区的奇异谱熵依然落差大，边缘陡峭，其中，图3和图5的横坐标表示采样点数，其纵坐标表示信幅值，图4和图6的横坐标表示采样点数，其纵坐标表示奇异谱熵值的相反数，所以从图4或图6上看到的上升沿实际上是该奇异谱熵的下降沿。由此可见，本发明的超声波首波检测与定位的方法具有较强的抗干扰能力，既使对于淹没在噪声中极弱的首波也能够较为正确的定位，这对远距离测量时的弱信号波识别和定位至关重要。

以上的实施例仅仅是对本发明的优选实施方式进行描述，并非对本发明的范围进行限定，在不脱离本发明设计精神的前提下，本领域普通工程技术人员对本发明的技术方案作出的各种变形和改进，均应落入本发明的权利要求书确定的保护范围内。

Claims (4)

1.一种超声波首波检测与定位的方法，其特征在于，包括以下步骤：

构造超声波信号序列的轨道矩阵；

对所述轨道矩阵进行奇异值分解获取所述超声波信号序列的奇异谱熵；

对所述奇异谱熵进行中值滤波，去除信号毛刺；

对已经过中值滤波的奇异谱熵进行阈值化处理，则已经过阈值化处理的奇异谱熵的下降沿即为所述超声波信号序列的首波；其中，所述构造超声波信号序列的轨道矩阵，具体包括如下步骤：

对一个长度为N点的超声波信号序列，截取以n点开始长度为L点的超声波信号子序列，构造出所述超声波信号子序列的轨道矩阵A_n为：

$A_n=\left[\begin{array}{cccc}x\left(n\right)& x(n+1)& \·\·\·& x(n+L-M)\\ x(n+1)& x(n+2)& \·\·\·& x(n+L-M+1)\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ x(n+M-1)& x(n+M)& \·\·\·& x(n+L-1)\end{array}\right]$

其中，n＝1，2，……，N-L+1，L小于N而大于超声波信号序列的周期，M为小于L的正整数；

其中，所述对轨道矩阵进行奇异值分解获取所述超声波信号序列的奇异谱熵，具体包括如下步骤：

对所述轨道矩阵A_n进行奇异值分解，获取其r个奇异值为σ_ni，其中i＝1，2，…，r；

根据公式：

得到所述超声波信号序列的归一化奇异谱概率密度函数p_n(k)，其中，σ_nk、σ_nj均为奇异值，k＝1，2，…，r；

根据公式：

得到所述超声波信号序列的奇异谱熵，其中，p_n(i)为归一化奇异谱概率密度函数。

2.根据权利要求1所述的超声波首波检测与定位的方法，其特征在于，所述对已经过中值滤波的奇异谱熵进行阈值化处理，具体包括：

根据公式：Th＝a×[max(H(n))-min(H(n))]+min(H(n))，得到所述超声波信号序列的奇异谱熵的上升沿和下降沿，其中，max(H(n))、min(H(n)分别为所述超声波信号序列的奇异谱熵的最大值、最小值，a为比例因子。

3.根据权利要求1或2所述的超声波首波检测与定位的方法，其特征在于，在所述获取超声波信号序列的奇异谱熵之前，还包括：

根据所述超声波信号序列的短时能量对所述超声波信号序列进行预处理，得到所述超声波信号序列出现的概率区间。

4.一种超声波首波检测与定位的装置，其特征在于，所述超声波首波检测与定位的装置执行包括权利要求3所述的超声波首波检测与定位的方法，其包括：

轨道矩阵构造模块，用于构造超声波信号序列的轨道矩阵；

奇异谱熵获取模块，用于对所述轨道矩阵构造模块构造出的轨道矩阵进行奇异值分解获取超声波信号序列的奇异谱熵；

中值滤波模块，用于对所述奇异谱熵获取模块得到的奇异谱熵进行中值滤波，去除信号毛刺；

阈值处理模块，用于对已经过所述中值滤波模块处理后的奇异谱熵进行阈值化处理，则已经过阈值化处理的奇异谱熵的下降沿即为所述超声波信号序列的首波。

Images