简正波频散曲线提取与声源定位方法和系统

技术领域

本发明涉及水声信号处理技术领域，具体地，涉及一种简正波频散曲线提取与声源定位方法和系统。尤其地，涉及一种基于二次重采样的简正波频散曲线精准提取与声源定位方法。

背景技术

浅海波导瞬态声源被动定位技术一直是水声领域的热点和难点问题。在浅海中，瞬态声源接收信号可以描述为一系列具有复杂频散特征的简正波模态的叠加，各阶简正波频散曲线包含声源传播距离信息，是现有研究进行声源被动定位、声场参数反演的理论基础。

浅海波导海面海底情况复杂，水文条件随时空显著变化，水声传播存在显著的频散效应和多途效应，使得接收信号呈现出复杂的多模态耦合、强非线性频散、强噪声干扰特点。以上复杂信号特征对常规信号处理方法造成了严重挑战，传统的时频分析或模态分解方法很难从复杂的接收信号中同时分离和高精度表征各阶频散曲线，且分离的频散曲线时频分辨率低、频散曲线估计精度不高，难以与拷贝声场频散曲线准确匹配，从而极大地限制了瞬态声源被动定位的精度与可靠性。例如，Warping变换类方法是当前单水听器接收信号进行模态分解的主要方法，将接收信号在时域(或频域)进行一次重采样，可以将各阶频散模态同时变换到截止频率附近，从而通过带通滤波实现模态分解。然而，Warping变换分解后的各阶简正波仍然需要进行时频表征，在时频图上各阶频散曲线时频分辨率较差、时频能量分布不集中，从中提取频散曲线较为杂乱，需要人工干预选择提取效果较好的片段，导致频散曲线提取精度差、效率低，最后影响与拷贝声场的匹配程度，使得声源被动定位精度低且稳定性差。

专利文献CN115166817A公开了一种基于冰层模态群慢度差特征的冰声定位方法，利用i个加速度计分别采集冰层中声源发出的A0模态冲击信号；对每个加速度计采集的A0模态冲击信号采用图像去模糊时频分析提取A0模态频散曲线；对每个加速度计对应A0模态频散曲线均执行如下操作，得到声源到每个加速度计的距离：在A0模态频散曲线的频率点内选取两个不同频率点构成频点组合，所述频点组合为全部可能的组合，计算每个频点组合的到达时间差；计算每个频点组合的群慢度差；根据群慢度差和到达时间差求解得到声源到加速度计的距离；最后根据声源到每个加速度计的距离求解得到声源坐标。然而该专利无法完全解决目前存在的技术问题，也无法满足本发明的需求。

发明内容

针对现有技术中的缺陷，本发明的目的是提供一种简正波频散曲线提取与声源定位方法和系统。

根据本发明提供的简正波频散曲线提取与声源定位方法，包括：

步骤1：利用单个水听器接收宽带脉冲声源在浅海中的传播信号，并对接收信号进行时间截取与带通滤波处理：选取第一阶简正波到达时刻作为接收信号的时间起点，选择待分析的带宽对接收信号进行带通滤波；

步骤2：设定时窗大小，对接收信号在各个时窗内进行第一次重采样变换，将各阶简正波紧邻密集频散曲线的强非线性频散特征同时消除；

步骤3：对各时窗的第一次重采样信号进行傅里叶变换，提取频谱峰值频率作为中心频率构建带通滤波器组，然后对所有时窗内的重采样信号进行带通滤波，从而实现各阶模态分量同步分离；

步骤4：对各时窗的各阶模态分量进行第二次重采样变换，将第一次重采样消除的强非线性频散曲线恢复至原始频散曲线，同时实现各时窗内紧邻密集强非线性频散曲线的时频能量聚集，提高时频分辨率，获得时频精细表征的频散曲线；

步骤5：对二次重采样后各阶频散模态的强非线性频散曲线进行精细提取；

步骤6：基于海洋环境参数，利用Kraken模型或Pekeris理论公式计算理论频散曲线值，并与上述精细提取的各阶频散模态的强非线性频散曲线值进行匹配，最后获取声源传播距离。

优选的，所述步骤2包括：

步骤2.1：设定时窗大小为tw，时窗中心点为tk，则依次在时窗[t

其中，Γ

以t

其中，τ

步骤2.2：将式(2)带入式(1)，得到N个时窗第一次重采样变换结果为：

经过上述第一次重采样后，重采样信号的频散曲线变为：

其中，φ[·]表示信号相位，

优选的，所述步骤3包括：

对式(3)的N个时窗的重采样信号进行傅里叶变换，获得N个傅里叶频谱，进行幅值叠加后获得频谱如下：

其中，

然后采用峰值选取方法，提取频谱中的若干峰值点所在频率f

优选的，所述步骤4包括：

步骤4.1：对第k个时窗的第m阶模态分量进行第二次重采样变换，表达式为：

其中，Γ

其中，

步骤4.2：将式(7)带入式(6)，第二次重采样变换进一步表达为：

在此过程中，第一次重采样时刻点τ

经过第二次重采样后，第k个时窗第m阶模态分量的频散曲线变为：

因此，经过第二次重采样后，重采样信号的第m阶频散曲线在第k个时窗的中心t

优选的，所述步骤5包括：

将式(8)中第m阶频散模态沿时间变量θ进行相同点数的傅里叶变换，获得第m阶频散模态的时频分布：

然后，按照能量最大原则从第m阶频散模态的时频分布中，提取时频脊线作为频散曲线提取值：

遍历各阶频散模态，精准提取各阶频散模态的频散曲线值；

所述步骤6包括：

频散曲线匹配准则构建如下：

其中，

根据本发明提供的简正波频散曲线提取与声源定位系统，包括：

模块M1：利用单个水听器接收宽带脉冲声源在浅海中的传播信号，并对接收信号进行时间截取与带通滤波处理：选取第一阶简正波到达时刻作为接收信号的时间起点，选择待分析的带宽对接收信号进行带通滤波；

模块M2：设定时窗大小，对接收信号在各个时窗内进行第一次重采样变换，将各阶简正波紧邻密集频散曲线的强非线性频散特征同时消除；

模块M3：对各时窗的第一次重采样信号进行傅里叶变换，提取频谱峰值频率作为中心频率构建带通滤波器组，然后对所有时窗内的重采样信号进行带通滤波，从而实现各阶模态分量同步分离；

模块M4：对各时窗的各阶模态分量进行第二次重采样变换，将第一次重采样消除的强非线性频散曲线恢复至原始频散曲线，同时实现各时窗内紧邻密集强非线性频散曲线的时频能量聚集，提高时频分辨率，获得时频精细表征的频散曲线；

模块M5：对二次重采样后各阶频散模态的强非线性频散曲线进行精细提取；

模块M6：基于海洋环境参数，利用Kraken模型或Pekeris理论公式计算理论频散曲线值，并与上述精细提取的各阶频散模态的强非线性频散曲线值进行匹配，最后获取声源传播距离。

优选的，所述模块M2包括：

模块M2.1：设定时窗大小为t

其中，Γ

以t

其中，τ

模块M2.2：将式(2)带入式(1)，得到N个时窗第一次重采样变换结果为：

经过上述第一次重采样后，重采样信号的频散曲线变为：

其中，φ[·]表示信号相位，

优选的，所述模块M3包括：

对式(3)的N个时窗的重采样信号进行傅里叶变换，获得N个傅里叶频谱，进行幅值叠加后获得频谱如下：

其中，

然后采用峰值选取方法，提取频谱中的若干峰值点所在频率f

优选的，所述模块M4包括：

模块M4.1：对第k个时窗的第m阶模态分量进行第二次重采样变换，表达式为：

其中，Γ

其中，

模块M4.2：将式(7)带入式(6)，第二次重采样变换进一步表达为：

在此过程中，第一次重采样时刻点τ

经过第二次重采样后，第k个时窗第m阶模态分量的频散曲线变为：

因此，经过第二次重采样后，重采样信号的第m阶频散曲线在第k个时窗的中心t

优选的，所述模块M5包括：

将式(8)中第m阶频散模态沿时间变量θ进行相同点数的傅里叶变换，获得第m阶频散模态的时频分布：

然后，按照能量最大原则从第m阶频散模态的时频分布中，提取时频脊线作为频散曲线提取值：

遍历各阶频散模态，精准提取各阶频散模态的频散曲线值；

所述模块M6包括：

频散曲线匹配准则构建如下：

其中，

与现有技术相比，本发明具有如下的有益效果：

(1)本发明无需任何环境敏感参数，可同时分离和精细表征各阶强非线性频散曲线，且分离的各阶频散曲线具有很高的时频分辨率和能量聚集性，能够大幅提升各阶频散曲线的提取精度；

(2)本发明提取的频散曲线精度高，能够提升与拷贝声场的匹配程度，从而显著提高声源被动定位的精度和可靠性；

(3)算法无需人工手动选择合适的频散曲线片段，自动化和工程实用性强。

附图说明

通过阅读参照以下附图对非限制性实施例所作的详细描述，本发明的其它特征、目的和优点将会变得更明显：

图1为本发明基于二次重采样的简正波频散曲线精准提取与声源定位方法流程图；

图2为仿真的宽带脉冲声源单水听器接收信号及其第一次重采样结果图：图2中的(a)为仿真声源接收信号时域波形；图2中的(b)为仿真声源接收信号短时傅里叶变换；图2中的(c)为第一次重采样后消除频散的时频图；图2中的(d)为各时窗频谱叠加与峰值频率选取；

图3为本发明针对仿真信号进行第二次重采样后频散曲线精准提取与精准定位结果图：图3中的(a)为本发明频散曲线时频分解与精细表征结果；图3中的(b)为本发明频散曲线精准提取结果；图3中的(c)为本发明声源精准定位结果；

图4为采用传统Warping变换方法进行频散曲线提取与声源定位的对比结果图：图4中的(a)为Warping变换频散曲线时频分解与表征结果；图4中的(b)为Warping变换频散曲线提取结果；图4中的(c)为Warping变换声源定位结果；

图5为海试数据CSS声源单水听器接收信号及其第一次重采样结果图：图5中的(a)为CSS声源接收信号时域波形；图5中的(b)为CSS声源接收信号短时傅里叶变换；图5中的(c)为第一次重采样后消除频散的时频图；图5中的(d)为各时窗频谱叠加与峰值频率选取；

图6为本发明针对海试数据CSS声源接收信号进行第二次重采样后频散曲线精准提取与精准定位结果图：图6中的(a)为本发明海试数据CSS声源频散曲线时频分解与精细表征结果；图6中的(b)为本发明海试数据CSS声源频散曲线精准提取结果；图6中的(c)为本发明CSS声源精准定位结果。

具体实施方式

下面结合具体实施例对本发明进行详细说明。以下实施例将有助于本领域的技术人员进一步理解本发明，但不以任何形式限制本发明。应当指出的是，对本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变化和改进。这些都属于本发明的保护范围。

实施例1

本发明提供了一种基于二次重采样的简正波频散曲线精准提取与声源定位方法，其流程图如图1所示，具体步骤包括：

步骤1：仿真利用单个水听器接收宽带脉冲声源在浅海中传播的信号，并对接收信号进行时间截取与带通滤波处理。

仿真海洋环境参数如下：声源和水听器分别位于水深10m和水深80m处，水听器距离声源10km，海水声速和密度分别为1500m/s和1000kg/m

步骤2：在各个时窗内对接收信号进行第一次重采样变换，将各阶简正波紧邻密集频散曲线的强非线性频散特征同时消除。

步骤2a：设定时窗大小t

其中，τ

步骤2b：以各采样时刻点为中心，依次对各个时窗内的信号进行第一次重采样变换，得到N个时窗第一次重采样变换结果为：

经过上述第一次重采样后，能够同时消除各阶简正波频散曲线的强非线性特征，获得频域上相对独立的频谱峰值(或时频域上近似平行于时间轴的水平线，如图2(c)所示)，从而为下一步带通滤波与模态分解奠定基础。

步骤3：对各时窗的第一次重采样信号进行傅里叶变换，获得N个时窗的傅里叶频谱，对所有频谱进行幅值叠加后，提取频谱峰值频率作为中心频率构建带通滤波器组，然后对所有时窗内的重采样信号进行带通滤波，从而实现各阶模态分量同步分离。

如图2(d)所示，提取叠加频谱中的8个峰值点所在频率f

步骤4：对各时窗的各阶模态分量进行第二次重采样变换，将第一次重采样消除的强非线性频散曲线恢复至原始频散曲线，同时实现各时窗内紧邻密集强非线性频散曲线的时频能量聚集，提高时频分辨率，获得时频精细表征的频散曲线。

具体过程如下：

步骤4.1：设定第k个时窗的第二次重采样函数为：

其中，Γ

步骤4.2：依次对第k个时窗的各阶模态分量进行第二次重采样变换，表达式为：

在此过程中，第一次重采样时刻点τ

经过第二次重采样后，第k个时窗第m阶模态分量的频散曲线变为：

因此，经过第二次重采样后，重采样信号的第m阶频散曲线在第k个时窗的中心t

步骤5：对二次重采样后各阶频散模态的强非线性频散曲线进行精细提取。具体过程如下：

将第m阶频散模态沿时间变量θ进行相同点数的傅里叶变换，获得第m阶频散模态的时频分布：

然后，按照能量最大原则从第m阶频散模态的时频分布中，提取时频脊线作为频散曲线提取值：

遍历各阶频散模态，可以精准提取各阶频散模态的频散曲线值。本发明频散曲线提取结果如图3(b)所示，可以看出本发明提取的前8阶频散曲线与理论频散曲线保持高度一致；在局部放大图中可以进一步看出，即使是紧邻密集的前三阶频散曲线，本发明提取结果仍然具有较好的提取精度。

步骤6：基于海洋环境参数，利用Pekeris理论公式计算理论频散曲线值，并与上述精准提取的各阶模态频散曲线值进行匹配，最后精准获取声源传播距离。声源传播距离估计如下：

其中，

声源传播距离估计结果如图3(c)所示，距离估计值为10km，与仿真设定的水听器与声源距离10km一致，说明了本发明提取的频散曲线具有较高精度，并且由此带来的声源定位精度较高。

为了说明本发明在提升频散曲线提取与声源定位精度方面的优越性，采用时域Warping变换方法作为对比，利用Warping变换进行频散曲线提取与声源定位，其结果如图4所示。从图4(a)中可以发现，Warping变换获得的前8阶频散曲线时频分辨率较差，时频能量分散、集中性差，从中提取频散曲线精度较差。如图4(b)所示，Warping变换提取的频散曲线与理论频散曲线有较大误差，特别是前三阶紧邻密集的频散曲线，Warping变换方法提取的频散曲线大幅偏离了理论频散曲线，需要人工手动选择一定范围的频散曲线片段再进行匹配，导致后续声源定位精度低、稳定性差，难以准确匹配声源传播距离，声源定位结果如图4(c)所示，Warping变换方法距离估计值为9.6km，与真实值10km存在一定差距。因此，上述对比结果进一步说明了本发明在简正波紧邻密集、强非线性频散曲线提取精度和声源定位精度方面的优越性。

为了进一步说明本发明在真实海洋环境中简正波频散曲线提取与声源被动定位的有效性和准确性，将本发明应用于公开海试数据燃烧声源(Combustive Sound Source，CSS声源)在浅海中传播的频散曲线提取与声源定位。海试数据CSS声源单水听器接收信号及其第一次重采样结果如图5所示，其中包含声源接收信号时域波形和短时傅里叶变换，选取待分析频带为1-850Hz，设定时窗大小和时间变量与上述仿真步骤相同，进行第一次重采样变换后，能够同时消除各阶频散曲线的强非线性频散特征，并在各时窗叠加频谱中提取前8阶频散曲线。进而，如图6所示，进行第二次重采样变换和频散曲线精细提取，并利用提取的频散曲线精准匹配基于海洋环境参数计算的理论频散曲线，获得声源精准定位结果。可以发现，本发明在海试数据处理中具有较好的有效性，各阶紧邻密集的频散曲线得到了精细表征和分解，其时频分辨率和能量聚集性大幅提高，频散曲线提取准确，声源定位精度4.83km与实际GPS测量值4.8km十分接近，验证了本发明在频散曲线提取和声源定位精度方面的高度准确性。

实施例2

本发明还提供一种简正波频散曲线提取与声源定位系统，所述简正波频散曲线提取与声源定位系统可以通过执行所述简正波频散曲线提取与声源定位方法的流程步骤予以实现，即本领域技术人员可以将所述简正波频散曲线提取与声源定位方法理解为所述简正波频散曲线提取与声源定位系统的优选实施方式。

根据本发明提供的简正波频散曲线提取与声源定位系统，包括：模块M1：利用单个水听器接收宽带脉冲声源在浅海中的传播信号，并对接收信号进行时间截取与带通滤波处理：选取第一阶简正波到达时刻作为接收信号的时间起点，选择待分析的带宽对接收信号进行带通滤波；模块M2：设定时窗大小，对接收信号在各个时窗内进行第一次重采样变换，将各阶简正波紧邻密集频散曲线的强非线性频散特征同时消除；模块M3：对各时窗的第一次重采样信号进行傅里叶变换，提取频谱峰值频率作为中心频率构建带通滤波器组，然后对所有时窗内的重采样信号进行带通滤波，从而实现各阶模态分量同步分离；模块M4：对各时窗的各阶模态分量进行第二次重采样变换，将第一次重采样消除的强非线性频散曲线恢复至原始频散曲线，同时实现各时窗内紧邻密集强非线性频散曲线的时频能量聚集，提高时频分辨率，获得时频精细表征的频散曲线；模块M5：对二次重采样后各阶频散模态的强非线性频散曲线进行精细提取；模块M6：基于海洋环境参数，利用Kraken模型或Pekeris理论公式计算理论频散曲线值，并与上述精细提取的各阶频散模态的强非线性频散曲线值进行匹配，最后获取声源传播距离。

所述模块M2包括：

模块M2.1：设定时窗大小为t

其中，Γ

以t

其中，τ

模块M2.2：将式(2)带入式(1)，得到N个时窗第一次重采样变换结果为：

经过上述第一次重采样后，重采样信号的频散曲线变为：

其中，φ[·]表示信号相位，

所述模块M3包括：

对式(3)的N个时窗的重采样信号进行傅里叶变换，获得N个傅里叶频谱，进行幅值叠加后获得频谱如下：

其中，

然后采用峰值选取方法，提取频谱中的若干峰值点所在频率f

所述模块M4包括：

模块M4.1：对第k个时窗的第m阶模态分量进行第二次重采样变换，表达式为：

其中，Γ

其中，

模块M4.2：将式(7)带入式(6)，第二次重采样变换进一步表达为：

在此过程中，第一次重采样时刻点τ

经过第二次重采样后，第k个时窗第m阶模态分量的频散曲线变为：

因此，经过第二次重采样后，重采样信号的第m阶频散曲线在第k个时窗的中心t

所述模块M5包括：

将式(8)中第m阶频散模态沿时间变量θ进行相同点数的傅里叶变换，获得第m阶频散模态的时频分布：

然后，按照能量最大原则从第m阶频散模态的时频分布中，提取时频脊线作为频散曲线提取值：

/>

遍历各阶频散模态，精准提取各阶频散模态的频散曲线值；

所述模块M6包括：

频散曲线匹配准则构建如下：

其中，

本领域技术人员知道，除了以纯计算机可读程序代码方式实现本发明提供的系统、装置及其各个模块以外，完全可以通过将方法步骤进行逻辑编程来使得本发明提供的系统、装置及其各个模块以逻辑门、开关、专用集成电路、可编程逻辑控制器以及嵌入式微控制器等的形式来实现相同程序。所以，本发明提供的系统、装置及其各个模块可以被认为是一种硬件部件，而对其内包括的用于实现各种程序的模块也可以视为硬件部件内的结构；也可以将用于实现各种功能的模块视为既可以是实现方法的软件程序又可以是硬件部件内的结构。

以上对本发明的具体实施例进行了描述。需要理解的是，本发明并不局限于上述特定实施方式，本领域技术人员可以在权利要求的范围内做出各种变化或修改，这并不影响本发明的实质内容。在不冲突的情况下，本申请的实施例和实施例中的特征可以任意相互组合。