CN114757240B - 一种水声目标运动参数联合估计方法 - Google Patents

一种水声目标运动参数联合估计方法 Download PDF

Info

Publication number
CN114757240B
CN114757240B CN202210677320.5A CN202210677320A CN114757240B CN 114757240 B CN114757240 B CN 114757240B CN 202210677320 A CN202210677320 A CN 202210677320A CN 114757240 B CN114757240 B CN 114757240B
Authority
CN
China
Prior art keywords
grid
cost function
doppler
cpa
search
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202210677320.5A
Other languages
English (en)
Other versions
CN114757240A (zh
Inventor
高大治
高伟
宋文华
李小雷
孙凯
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Ocean University of China
Original Assignee
Ocean University of China
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Ocean University of China filed Critical Ocean University of China
Priority to CN202210677320.5A priority Critical patent/CN114757240B/zh
Publication of CN114757240A publication Critical patent/CN114757240A/zh
Application granted granted Critical
Publication of CN114757240B publication Critical patent/CN114757240B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/02—Preprocessing
    • G06F2218/04—Denoising
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S11/00—Systems for determining distance or velocity not using reflection or reradiation
    • G01S11/14—Systems for determining distance or velocity not using reflection or reradiation using ultrasonic, sonic, or infrasonic waves
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
    • G06F2218/08—Feature extraction

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Signal Processing (AREA)
  • General Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

本申请属于海洋声学测量及水声定位技术领域，提供一种水声目标运动参数联合估计方法，包括以下步骤：构造运动参数
Figure 991851DEST_PATH_IMAGE001
的搜索网格；基于每个网格点对应的Doppler‑warping变换信号计算该网格点处的第一代价函数，遍历所述搜索网格，得到第一代价函数网格；基于第一代价函数网格的多个局部极大值生成多普勒参数耦合曲线；基于每个网格点对应的波导干涉条纹曲线计算该网格点处的第二代价函数，遍历所述搜索网格，得到第二代价函数网格；基于第二代价函数网格生成
Figure 781952DEST_PATH_IMAGE002
耦合直线；基于所述多普勒参数耦合线和所述
Figure 536281DEST_PATH_IMAGE002
耦合直线确定V与
Figure 555184DEST_PATH_IMAGE003
的估计值。本申请的方法能够精确地估计水声目标运动参数。

Description

一种水声目标运动参数联合估计方法
技术领域
本申请属于海洋声学测量及水声定位技术领域，具体地，涉及一种水声目标运动参数联合估计方法。
背景技术
水声无源探测技术一直以来倍受关注，国内外学者根据水中目标所辐射出噪声特性获取目标位置及运动信息。常见的目标无源探测技术有三子阵测距、目标运动分析(TMA)、匹配场处理(MFP)和基于声场干涉结构的方法等。传统三子阵测距和目标运动分析未考虑浅海波导环境的影响，在浅海性能较差；匹配场处理对于环境适配性要求较高，且拷贝场计算量较大因此实际实现难度较大；浅海声场干涉结构方法对环境宽容性较好，不需要水听器之间保持较高的同步性，运算量较小，因此获得大量研究，在使用单水听器或短阵信号，利用干涉条纹进行目标运动参数估计时，只能得到目标最近距离与速度比值，无法将两个参数分别解出。
为解决上述问题，专利201610879384.8提出了一种联合波导不变量和线谱的单水听器被动定位方法，其先利用干涉条纹估计参数耦合值，然后利用线谱互谱分析得到目标径向运动速度，通过多时间点拟合得到速度信息，再结合耦合值最终解出距离信息，但这种方法要获得较好的速度结果需要较长时间的信号积累，在最近点时间处有较大误差；并且要求信号稳定，即信号频率不能有变化，因此对于被动声源以及信号时长较短时方法并不适用。专利201010234991.1公开了一种基于单水听器的被动测距方法，结合引导源信息估计目标径向运动速度以及距离信息，但无法给出绝对运动速度。专利201910790936.1提出了一种浅海单水听器运动目标距离估计方法，使用对宽带信号频段进行fft估计径向速度，然后采用最小二乘法拟合得到距离估计值，该方法同样需要较长时间积累，且在最近点时刻处径向速度估计误差较大从而导致距离估计误差较大。
已知利用线谱多普勒频移可以实现对目标运动速度，正横距离等参数估计。Ferguson和Quinn等(B.G.Ferguson,B.G.Quinn,Application of the Short-timeFourier Transform and the Wigner-Ville Distribution to the AcousticLocalization of Aircraft,Journal of the Acoustical Society of America,96(1994),821-827.)提出使用时频分析提取多普勒信号瞬时频率，然后使用最小二乘法将提取频率与模型进行拟合，从而实现目标运动参数估计。高德洋等(高德洋,高大治,迟静,王良,宋文华.Doppler-warping变换及其应用在声学目标运动速度估计[J].物理学报,2021,70(12):253-260.)提出Doppler-warping变换并应用于水声目标运动参数估计。但水中声信号的多普勒现象相对于空气较弱，单线谱进行多参数估计时会存在较强的参数耦合现象，无法分别解出最近距离与速度值。
发明内容
本申请的目的在于解决上述现有各种水声目标探测技术中存在的缺陷，提供一种能够精确地估计水声目标运动参数的方法。
本申请的实施例可以通过以下技术方案实现：
一种水声目标运动参数联合估计方法，基于水声目标的辐射噪声的多普勒频移线谱及波导干涉条纹联合估计水声目标的运动参数，所述运动参数包括水声目标的速度V及最近距离rCPA，包括以下步骤：
S1：构造运动参数V、rCPA的搜索网格{(vn,rm)}，其中n∈[1,N],m∈[1,M]；
S2：基于每个网格点对应的Doppler-warping变换信号计算该网格点处的第一代价函数，遍历所述搜索网格，得到第一代价函数网格S1；
S3：基于第一代价函数网格S1的多个局部极大值生成多普勒参数耦合曲线；
S4：基于每个网格点对应的波导干涉条纹曲线计算该网格点处的第二代价函数，遍历所述搜索网格，得到第二代价函数网格S2；
S5：基于第二代价函数网格S2生成rCPA/V耦合直线；
S6：基于所述多普勒参数耦合线和所述rCPA/V耦合直线确定V与rCPA的估计值。
进一步地，步骤S2包括以下步骤：
S21：基于所述辐射噪声的原始信号提取多普勒频移线谱并确定水声目标的最近点时间tCPA；
S22：选取所述搜索网格中的任意一个网格点(vi,rj)，其中，i∈[1,N],j∈[1,M]；
S23：构造(vi,rj)对应的Doppler-warping算子
Figure GDA0003765983370000021
其中，c为水中声速；
S24：使用所述h(i,j,t)对辐射噪声的原始信号x(t)进行重采样，得到该网格点(vi,rj)对应的Doppler-warping变换信号xp(i,j,t)；
S25：计算(vi,rj)对应的第一代价函数
Figure GDA0003765983370000022
其中，Xp(i,j,f)为所述xp(i,j,t)的频谱，[f1,f2]为积分上下限；
S26：重复执行步骤S22至S25，直到遍历所述搜索网格；
S27：基于全部网格点处的第一代价函数生成第一代价函数网格S1。
优选地，所述积分上下限[f1,f2]基于所述线谱的多普勒频率展宽情况确定。
进一步地，步骤S4包括以下步骤：
S41：确定所述辐射噪声的时频谱I(f,t)；
S42：选取所述搜索网格中的任意一个网格点(vk,rl)，其中，k∈[1,N],l∈[1,M]；
S43：提取(vk,rl)对应的多条波导干涉条纹曲线
Figure GDA0003765983370000031
其中，β为波导不变量，f0为各条干涉条纹曲线最低点对应的频率值；
S44：计算(vk,rl)对应的第二代价函数
Figure GDA0003765983370000032
其中，ds的积分路径为所述f(k,l,f0,t)；
S45：重复执行步骤S42至S45，直到遍历所述搜索网格；
S46：基于全部网格点处的第二代价函数生成第二代价函数网格S2。
优选地，步骤S41之前还包括步骤S40：选取第一代价函数S1中最大值对应的运动参数v、rCPA组合构造Doppler-warping算子并对辐射噪声的原始信号进行重采样。
进一步地，步骤S5包括以下步骤：
S51：选取所述第二代价函数网格S2的每一列S2(:,m),m∈[1,M]中的最大值；
S52：基于所述多个最大值对应的网格点拟合生成rCPA/V耦合直线。
本申请的实施例提供的一种水声目标运动参数联合估计方法至少具有以下有益效果：
(1)本申请提供的水声目标运动参数联合估计方法，能够结合水声目标的辐射噪声中包含的线谱的多普勒频移特征和时频谱中由于波导效应产生的干涉结构特征，分别生成多普勒参数耦合曲线和基于波导不变量的运动参数耦合直线，从而综合利用了单水听器接收的辐射噪声中包含的多种效应造成的不同信号特征，最大程度地提高了单水听器运动参数估计的精度；
(2)在生成基于波导不变量的运动参数耦合直线之前增加了使用Doppler-warping算子对辐射噪声的原始信号进行重采样的步骤，通过选取第一代价函数S1中最大值对应的运动参数v、rCPA组合构造Doppler-warping算子并对辐射噪声的原始信号进行重采样，能够消除原始信号中的多普勒频移效应，从而有效地提高了拟合运动参数耦合直线的精度。
附图说明
图1为水声目标运动及水听器布设的示意图；
图2为根据本申请实施例的水声目标运动参数联合估计方法的流程图；
图3为一个具体的单频声源辐射噪声的原始信号的频谱；
图4为一个具体的水声目标所产生的辐射噪声的LOFAR图；
图5为根据本申请实施例1的线谱仿真时频图；
图6为根据本申请实施例1的第一代价函数网格的热度图及生成的多普勒参数耦合曲线；
图7为根据本申请实施例1的宽带干涉信号仿真时频图；
图8为根据本申请实施例1的第二代价函数的热度图及生成的rCPA/V耦合直线；
图9为根据本申请实施例1的运动参数估计结果；
图10为根据本申请实施例2的实测线谱信号的时频图；
图11为根据本申请实施例2的第一代价函数网格的热度图及生成的多普勒参数耦合曲线；
图12为根据本申请实施例2的实测宽带干涉信号的时频图；
图13为根据本申请实施例2的第二代价函数的热度图及生成的rCPA/V耦合直线；
图14为根据本申请实施例2的运动参数估计结果。
具体实施方式
以下，基于优选的实施方式并参照附图对本申请进行进一步说明。
此外，为了方便理解，放大或者缩小了图纸上的各种构件，但这种做法不是为了限制本申请的保护范围。
单数形式的词汇也包括复数含义，反之亦然。
在本申请实施例中的描述中，需要说明的是，若出现术语“上”、“下”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，或者是本申请实施例的产品使用时惯常摆放的方位或位置关系，仅是为了便于描述本申请和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本申请的限制。此外，本申请的描述中，为了区分不同的单元，本说明书上用了第一、第二等词汇，但这些不会受到制造的顺序限制，也不能理解为指示或暗示相对重要性，其在本申请的详细说明与权利要求书上，其名称可能会不同。
本说明书中词汇是为了说明本申请的实施例而使用的，但不是试图要限制本申请。还需要说明的是，除非另有明确的规定和限定，若出现术语“设置”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一体地连接；可以是机械连接，可以是直接相连，也可以通过中间媒介间接连接，可以是两个元件内部的连通。对于本领域的技术人员而言，可以具体理解上述术语在本申请中的具体含义。
为更好地说明本申请的技术方案，首先对基于水声目标的辐射噪声进行运动参数估计的原理进行说明。
图1为水声目标运动及水听器布设的示意图，如图1所示，水声目标为船舶或其拖曳的发射换能器等声源1，声源1在一定时间内以等速V进行直线运动，并产生在水声环境中传播的辐射噪声，声源1与布设于水中的水听器2的实时距离为r，在声源1的行进过程中与水听器2的最近距离及对应的最近点时间分别为rCPA与tCPA。在水听器2的位置、海水中声速c以及tCPA已知的情况下，对水声目标的运动参数估计是要是对其速度V和最近距离rCPA进行估计。
常规的使用多水听器进行运动参数估计的方法一般通过在海水中布设不在一条直线上的至少三个水听器，分别基于每个水听器接收到的辐射噪声信号得到各自的运动参数估计值，最后基于多个联立的方程组唯一确定一组运动参数。然而上述水声目标定位方法存在水听器阵元布设复杂、需要对各个水听器精确校时、需要处理的数据量大等问题。
由于浅海海洋环境的复杂性以及水声目标所包含的多种辐射声源的不同频谱特征，使得水听器接收到的辐射噪声信号中叠加了多普勒频移、浅海波导等效应所导致的多种时频谱结构特征，基于浅海中运动的水声目标辐射噪声的上述特性，本申请提出了一种水声目标运动参数联合估计方法，能够分别利用辐射噪声信号所携带的多重有效特征，实现了单水听器联合估计水声目标的运动参数。图2示出了根据本申请实施例提供的一种水声目标运动参数联合估计方法的流程图，如图2所示，包括以下步骤：
S1：构造运动参数V、rCPA的搜索网格{(vn,rm)}，其中n∈[1,N],m∈[1,M]；
S2：基于每个网格点对应的Doppler-warping变换信号计算该网格点处的第一代价函数，遍历所述搜索网格，得到第一代价函数网格S1；
S3：基于第一代价函数网格S1的多个局部极大值生成多普勒参数耦合曲线；
S4：基于每个网格点对应的波导干涉条纹曲线计算该网格点处的第二代价函数，遍历所述搜索网格，得到第二代价函数网格S2；
S5：基于第二代价函数网格S2生成rCPA/V耦合直线；
S6：基于所述多普勒参数耦合线和所述rCPA/V耦合直线确定V与rCPA的估计值。
以下结合附图及具体实施例对步骤S1至S6进行详细的说明。
步骤S1用于构造运动参数V、rCPA的搜索网格，在本申请的一些具体的实施例中，搜索网格的形式为：
Figure GDA0003765983370000051
其中，每个搜索网格的网格点代表了运动参数V、rCPA的一组待搜索的参数组合，[v1,vN]和[r1,rM]分别代表了对V和rCPA进行搜索的上下限。
在完成搜索网格的构造后，分别遍历搜索网格的每一组运动参数组合，基于辐射噪声线谱的多普勒频移特征生成多普勒参数耦合曲线(步骤S2、S3)以及基于水声目标辐射噪声时频谱的波导不变性特征生成rCPA/V耦合直线(S4、S5)，最后通过上述参数耦合曲线和参数耦合直线的交点确定V与rCPA的估计值(步骤S6)。
为更清楚地说明本申请的技术方案，以下分别对基于辐射噪声线谱特征的多普勒参数耦合曲线和基于辐射噪声时频谱的波导不变特征的rCPA/V耦合直线的物理意义及确定方式进行说明。
以下对基于辐射噪声线谱确定多普勒参数耦合曲线的实施方式进行说明。
如图1所示，由于声源1相对于水听器2进行直线运动，因此水听器2接收到的声源1所产生的辐射噪声的原始信号将出现多普勒频移效应，图3示出了一个具体的单频声源辐射噪声的原始信号的频谱，从图3可以看出其单频信号因多普勒频移具有明显的频带展宽现象。
对水声目标所产生的辐射噪声的原始信号进行短时傅里叶变换处理，能够得到其时频谱，一般采用LOFAR(Low Frequency Analysis Recording)图进行表示。图4示出了一个具体的水声目标所产生的辐射噪声的LOFAR图，从图4中能够观察到多条抛物线状的线谱，这些线谱分别由水声目标包含的各种单频噪声源(如船舶上的螺旋桨、机械设备以及发射单频信号的拖曳声源等)产生，其类似抛物线的形态则代表了水声目标与水听器逐渐靠近到再次远离的过程中所出现的多普勒频移现象。
对于待估计的运动参数V、rCPA，如果其取值为水声目标的真实运动参数值(v0,r0),则利用该组真实值构造Doppler-warping算子
Figure GDA0003765983370000061
并利用其对辐射噪声的原始信号x(t)进行重采样后，所得到的重采样信号xp(v0,r0,t)将消除多普勒频移效应，即xp(v0,r0,t)的LOFAR图中的线谱将重新恢复为单频，显然，该组(v0,r0)对应的重采样信号在该单频频率附近的频谱熵的倒数将具有最大值。
利用上述辐射噪声的Doppler-warping变换信号的线谱特征，可以使用不同的运动参数组合构造对应的Doppler-warping算子对原始信号进行重采样，并将重采样后信号的频谱熵的倒数作为代价函数，从中搜索最大值对应的运动参数组合作为运动参数的估计值。然而，由于水声目标运动所产生的多普勒频移远远小于空气中高速运动目标所产生的多普勒频移，运动参数(V,rCPA)之间存在明显的参数耦合现象，导致其代价函数计算结果中存在多个局部极值，因此基于多普勒频移线谱特征进行运动参数估计的结果为一条多普勒参数耦合曲线，水声目标的真实运动参数(v0,r0)即处于该多普勒参数耦合曲线上。利用该条多普勒参数耦合曲线结合联合其他运动参数估计结果，就能够唯一地确定一组运动参数V,rCPA的估计值。
具体地，在本申请的一些实施例中，在步骤S2中遍历搜索网格，通过每个网格点对应的Doppler-warping变换信号计算该网格点处的第一代价函数，最终得到第一代价函数网格S1；在步骤S3中对第一代价函数网格S1中的多个局部极大值进行曲线拟合，从而得到多普勒参数耦合曲线。其中，步骤S2进一步包括以下步骤：
S21：基于所述辐射噪声的原始信号提取多普勒频移线谱并确定水声目标的最近点时间tCPA(在一些具体的实施方式中，如图4所示，可以通过短时傅里叶变换计算水声目标辐射噪声的时频谱并生成LOFAR图，从图中选择线谱并基于线谱的形态确定水声目标的最近点时间tCPA)；
S22：选取所述搜索网格中的任意一个网格点(vi,rj)，其中，i∈[1,N],j∈[1,M]；
S23：构造(vi,rj)对应的Doppler-warping算子
Figure GDA0003765983370000062
其中，c为水中声速；
S24：使用所述h(i,j,t)对辐射噪声的原始信号x(t)进行重采样，得到该网格点(vi,rj)对应的Doppler-warping变换信号xp(i,j,t)；
S25：计算(vi,rj)对应的第一代价函数
Figure GDA0003765983370000071
其中，Xp(i,j,f)为所述xp(i,j,t)的频谱，[f1,f2]为积分上下限，在本申请的一些优选的实施例中，[f1,f2]基于多普勒频移线谱的频率展宽情况确定；
S26：重复执行步骤S22至S25，直到遍历所述搜索网格；
S27：基于全部网格点处的第一代价函数生成第一代价函数网格S1。
在得到第一代价函数网格S1后，在步骤S3中利用第一代价函数网格S1中的多个局部极大值进行曲线拟合，即可得到多普勒参数耦合曲线。利用平面上的多个离散点进行曲线拟合的方法已为本领域技术人员所知晓，在此不再赘述。
以下对基于辐射噪声时频谱的波导不变特征确定rCPA/V耦合直线的实施方式进行说明。
水声目标的宽频辐射噪声信号在浅海水声环境中传播时，由于波导效应的存在，其时频谱具有显著的干涉结构特征。Chuprov等首先提出了波导不变量β的概念用以描述频率、距离以及干涉条纹斜率之间的关系，其表达式为：
Figure GDA0003765983370000072
(1)式中r为声源到水听器距离，ω为声源角频率。将(1)式中的角频率ω换为频率f，并移项整理得到：
Figure GDA0003765983370000073
对(2)式两边积分，可得：
Figure GDA0003765983370000074
如图1所示，当水声目标1在水听器2一侧匀速直线运动时，水声目标1相对于水听器2的距离r可以表示为：
Figure GDA0003765983370000075
将(4)式带入(3)式,可得到如(5)式的干涉条纹方程：
Figure GDA0003765983370000081
(5)式为干涉条纹方程，根据参数f0取值的不同，分别表示了不同的干涉条纹频率f随时间t变化的关系。其中，参数f0为各条干涉条纹曲线最低点对应的频率值，β为波导不变量，一般理想情况下可近似认为β＝1。
Hough变换最早用于直线检测，后来逐渐发展到用于曲线检测问题。其实质上是建立一种曲线到参数空间的映射关系，参数空间中每一点的强度为沿曲线对应参数积分结果。当最近点时间tCPA、波导不变量β已知时，通过遍历运动参数V、rCPA的搜索网格，将各组运动参数带入(5)式进行Hough变换将得到不同参数f0对应的多条干涉条纹曲线并沿各条曲线进行积分。当运动参数的取值正确时，所得到的积分值最大。因干涉条纹方程中两参数为比值形式，因此Hough变换搜索得到两个运动参数的估计值应呈线性关系。
具体地，在本申请的一些实施例中，在步骤S4中遍历搜索网格，通过每个网格点对应的基于波导不变量进行Hough变换得到的干涉条纹方程计算该网格点处的第二代价函数，最终得到第二代价函数网格S2；在步骤S5中对第二代价函数网格S2中每个最近距离搜索值处的最大值进行拟合，从而得到rCPA/V耦合直线。其中，步骤S4进一步包括以下步骤：
S41：确定所述辐射噪声的时频谱I(f,t)；
S42：选取所述搜索网格中的任意一个网格点(vk,rl)，其中，k∈[1,N],l∈[1,M]；
S43：提取(vk,rl)对应的多条波导干涉条纹曲线
Figure GDA0003765983370000082
其中，β为波导不变量，f0为各条干涉条纹曲线最低点对应的频率值；
S44：计算(vk,rl)对应的第二代价函数
Figure GDA0003765983370000083
其中，ds的积分路径为所述f(k,l,f0,t)；
S45：重复执行步骤S42至S45，直到遍历所述搜索网格；
S46：基于全部网格点处的第二代价函数生成第二代价函数网格S2。
在本申请的一些优选的实施例中，步骤S41之前还包括步骤S40：选取第一代价函数S1中最大值对应的运动参数v、rCPA组合构造Doppler-warping算子并对辐射噪声的原始信号进行重采样。通过选取第一代价函数S1中最大值对应的运动参数v、rCPA组合构造Doppler-warping算子并对辐射噪声的原始信号进行重采样，能够消除原始信号中的多普勒频移效应，从而有效地提高了拟合运动参数耦合直线的精度。
在本申请的一些实施例中，步骤S5包括以下步骤：
S51：选取所述第二代价函数网格S2的每一列S2(:,m),m∈[1,M]中的最大值；
S52：基于所述多个最大值对应的网格点拟合生成rCPA/V耦合直线。
实施例1
本实施例为仿真例。其中浅海声场仿真参数设置：Pekeris波导，水体声速1480m/s,水深40m,声源深度2.5m，海底声速1650m/s,密度1.76g/cm3，海底衰减系数0.2dB/λ。声源运动速度V＝9m/s,最近距离rCPA＝900m，到达最近点时间tCPA＝100s，声压场由kraken程序计算生成。
声源线谱中心频率选取400Hz，采样率10kHz，线谱仿真时频图如图5所示。设置运动参数V，rCPA搜索网格：V＝[1：1：20]，rCPA＝[500:10:1500]。对仿真得到的时域信号进行Doppler-warping变换并遍历搜索网格得到第一代价函数网格，图6示出了第一代价函数网格的热度图及生成的多普勒参数耦合曲线(图中黑色曲线)，采用三次多项式拟合得到。
浅海宽带干涉谱仿真参数设置：浅海声场仿真参数与线谱一致，宽带干涉信号带宽为1000Hz-3000Hz，时频图频率分辨率为1Hz，时间步长0.1s，宽带干涉信号仿真时频图如图7所示。根据时频图的干涉结构采用Hough变换，其中β为0.95，遍历搜索网格得到第二代价函数网格，图8示出了第二代价函数的热度图及生成的rCPA/V耦合直线(图中黑色直线)，采用一次多项式拟合得到。
图9示出了多普勒参数耦合曲线与rCPA/V耦合直线的相交的示意图，图中交点即为参数估计结果，可以看到速度估计值为8.99m/s，最近距离估计值为899m，两参数估计误差均为0.11％。
实施例2
本实施例为海上实验数据处理实例。实验使用自容式水听器记录水声信号，水听器自噪声为27dB reμPa2/Hz，灵敏度为-169dB re 1V/μPa，动态范围121dB。接收船抛锚停主辅机，水听器放于水面下5m处,采样率128kHz。
实验目标船为快艇。快艇全速(16m/s)沿直线航行，快艇到达最近点时间为39s，距离水听器最近距离467m，水听器位于其航线一侧，期间航向和速度保持不变,航行中发动机和螺旋桨等辐射噪声为实验声源，总信号时长70s。接收船和快艇均配有差分GPS，实时记录位置、航向、速度信息。实验海域水深约为28m，水体平均声速为1510m/s。
图10为实验接收到的线谱信号时频图，图中可以看到存在较为明显的多普勒频移。设置运动参数V，rCPA搜索网格：V＝[5：0.5：15]，rCPA＝[200:10:800]。对原始接收线谱信号进行Doppler-warping变换，并遍历搜索网格得到第一代价函数网格，图11示出了第一代价函数网格的热度图及生成的多普勒参数耦合曲线(图中黑色曲线)，采用三次多项式拟合得到。
对接收信号进行短时傅里叶变换，时间窗函数为Hanning窗,窗长1s,步长0.1s。选取1000Hz-2000Hz宽带干涉信号得到如图12所示的时频图，根据时频图的干涉结构采用Hough变换，其中β为0.95，遍历搜索网格得到第二代价函数网格，图13示出了第二代价函数的热度图及生成的rCP/V耦合直线(图中黑色直线)，采用一次多项式拟合得到。
图14示出了实验得到的多普勒参数耦合曲线与rCPA/V耦合直线相交的示意图，图中交点即为参数估计结果，速度估计值为16.26m/s,最近距离估计值为477.8m，两参数估计误差分别为1.6％、2.3％。
以上对本申请的具体实施方式作了详细介绍，对于本技术领域的技术人员来说，在不脱离本申请原理的前提下，还可以对本申请进行若干改进和修饰，这些改进和修饰也属于本申请权利要求的保护范围。

Claims (4)

1.一种水声目标运动参数联合估计方法，基于水声目标的辐射噪声的多普勒频移线谱及波导干涉条纹联合估计水声目标的运动参数，所述运动参数包括水声目标的速度V及最近距离rCPA，其特征在于，包括以下步骤：
S1：构造运动参数V、rCPA的搜索网格{(vn,rm)}，其中n∈[1,N],m∈[1,M]；
S2：基于每个网格点对应的Doppler-warping变换信号计算该网格点处的第一代价函数，遍历所述搜索网格，得到第一代价函数网格S1；
S3：基于第一代价函数网格S1的多个局部极大值生成多普勒参数耦合曲线；
S4：基于每个网格点对应的波导干涉条纹曲线计算该网格点处的第二代价函数，遍历所述搜索网格，得到第二代价函数网格S2；
S5：基于第二代价函数网格S2生成rCPA/V耦合直线；
S6：基于所述多普勒参数耦合线和所述rCPA/V耦合直线确定V与rCPA的估计值；
步骤S2进一步包括以下步骤：
S21：基于所述辐射噪声的原始信号提取多普勒频移线谱并确定水声目标的最近点时间tCPA；
S22：选取所述搜索网格中的任意一个网格点(vi,rj)，其中，i∈[1,N],j∈[1,M]；
S23：构造(vi,rj)对应的Doppler-warping算子
Figure FDA0003782876940000011
其中，c为水中声速，t为时间；
S24：使用所述h(i,j,t)对辐射噪声的原始信号x(t)进行重采样，得到该网格点(vi,rj)对应的Doppler-warping变换信号xp(i,j,t)；
S25：计算(vi,rj)对应的第一代价函数
Figure FDA0003782876940000012
其中，Xp(i,j,f)为所述xp(i,j,t)的频谱，[f1,f2]为积分上下限；
S26：重复执行步骤S22至S25，直到遍历所述搜索网格；
S27：基于全部网格点处的第一代价函数生成第一代价函数网格S1；
步骤S4进一步包括以下步骤：
S41：确定所述辐射噪声的时频谱I(f,t)；
S42：选取所述搜索网格中的任意一个网格点(vk,rl)，其中，k∈[1,N],l∈[1,M]；
S43：提取(vk,rl)对应的多条波导干涉条纹曲线
Figure FDA0003782876940000021
其中，β为波导不变量，f0为各条干涉条纹曲线最低点对应的频率值；
S44：计算(vk,rl)对应的第二代价函数
Figure FDA0003782876940000022
其中，ds的积分路径为所述f(k,l,f0,t)；
S45：重复执行步骤S42至S45，直到遍历所述搜索网格；
S46：基于全部网格点处的第二代价函数生成第二代价函数网格S2。
2.根据权利要求1所述的一种水声目标运动参数联合估计方法，其特征在于：
所述积分上下限[f1,f2]基于多普勒频移线谱的频率展宽情况确定。
3.根据权利要求1所述的一种水声目标运动参数联合估计方法，其特征在于，步骤S41之前还包括步骤S40：
选取第一代价函数网格S1中最大值对应的运动参数v、rCPA组合构造Doppler-warping算子并对辐射噪声的原始信号进行重采样。
4.如根据权利要求1所述的一种水声目标运动参数联合估计方法，其特征在于，步骤S5进一步包括以下步骤：
S51：选取所述第二代价函数网格S2的每一列s2(:,m),m∈[1,M]中的最大值；
S52：基于所述多个最大值对应的网格点拟合生成rCPA/V耦合直线。
CN202210677320.5A 2022-06-16 2022-06-16 一种水声目标运动参数联合估计方法 Active CN114757240B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210677320.5A CN114757240B (zh) 2022-06-16 2022-06-16 一种水声目标运动参数联合估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210677320.5A CN114757240B (zh) 2022-06-16 2022-06-16 一种水声目标运动参数联合估计方法

Publications (2)

Publication Number Publication Date
CN114757240A CN114757240A (zh) 2022-07-15
CN114757240B true CN114757240B (zh) 2022-09-16

Family

ID=82336099

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210677320.5A Active CN114757240B (zh) 2022-06-16 2022-06-16 一种水声目标运动参数联合估计方法

Country Status (1)

Country Link
CN (1) CN114757240B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115015839B (zh) * 2022-08-10 2022-11-22 中国海洋大学 一种浅海水下目标被动定位系统

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108226933A (zh) * 2017-12-28 2018-06-29 西北工业大学 一种基于条纹干涉结构的深海宽带目标深度估计方法
CN112415474A (zh) * 2020-10-21 2021-02-26 中国海洋大学 线性差频匹配场处理水中目标声学被动定位方法、介质
WO2022088564A1 (zh) * 2020-10-28 2022-05-05 鹏城实验室 基于运动平台的多载波水声通信多普勒估计与补偿方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108226933A (zh) * 2017-12-28 2018-06-29 西北工业大学 一种基于条纹干涉结构的深海宽带目标深度估计方法
CN112415474A (zh) * 2020-10-21 2021-02-26 中国海洋大学 线性差频匹配场处理水中目标声学被动定位方法、介质
WO2022088564A1 (zh) * 2020-10-28 2022-05-05 鹏城实验室 基于运动平台的多载波水声通信多普勒估计与补偿方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Doppler-warping 变换及其应用在声学目标运动速度估计;高德洋 等;《物理学报》;20210620;第70卷(第12期);摘要、第2-3节 *
孙凯 等.浅海高速目标运动参数获取研究.《2021年浙黑苏鲁沪渝四省二市声学技术学术会议论文集》.2021, *
浅海高速目标运动参数获取研究;孙凯 等;《2021年浙黑苏鲁沪渝四省二市声学技术学术会议论文集》;20211031;引言、第3节 *

Also Published As

Publication number Publication date
CN114757240A (zh) 2022-07-15

Similar Documents

Publication Publication Date Title
CN110146895B (zh) 基于倒置式多波束回声仪的声速剖面反演方法
CN108089155B (zh) 一种深海环境下单水听器声源被动定位方法
CN103048642B (zh) 基于频域最小二乘法的水声脉冲信号匹配场定位方法
CN102749622B (zh) 基于多波束测深的声速剖面及海底地形的联合反演方法
CN104678384B (zh) 一种波束域的声压差互相关谱分析水下目标速度估计方法
CN112083404B (zh) 一种基于多途特征匹配的单矢量水听器声源深度估计方法
CN107179535A (zh) 一种基于畸变拖曳阵的保真增强波束形成的方法
CN103076594B (zh) 一种基于互相关的水声脉冲信号双阵元定位的方法
CN108562891B (zh) 一种深海低信噪比条件下声源深度自主实时跟踪方法
CN113281727B (zh) 一种基于水平线列阵的输出增强的波束形成方法及其系统
CN101915922A (zh) 拖曳线列阵被动测距方法
CN103529441A (zh) 一种被动合成孔径目标信号检测和分辨方法及系统
CN114757241B (zh) 一种多普勒参数耦合线提取方法
CN111273297A (zh) 一种基于ar波数谱的水平阵动目标深度估计方法
CN103176163A (zh) 基于相位模型的舰船线谱噪声源位置识别方法
CN114757240B (zh) 一种水声目标运动参数联合估计方法
Ren et al. Acoustic interferometry for geoacoustic characterization in a soft-layered sediment environment
CN104714235A (zh) 一种低频双矢量水听器阵列的测距方法及系统
CN103344961B (zh) 舰船速度和距离联合测量的被动声多普勒相位方法
CN103513238B (zh) 一种规整化最小二乘子空间相交的目标方位测向方法
CN110488219A (zh) 一种浅海单水听器运动目标距离估计方法
CN103809170A (zh) 一种多普勒频偏修正方法及系统
Koch Proof of principle for inversion of vector sensor array data
CN115236592B (zh) 一种基于单阵元时频曲线匹配的冰声定位方法
Gebbie et al. Aspect-dependent radiated noise analysis of an underway autonomous underwater vehicle

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant