CN104568113A - 一种基于模型的海洋声传播调查爆炸波自动截取方法 - Google Patents

Abstract

本发明公开了一种基于模型的海洋声传播调查爆炸波自动截取方法，针对发射船按测线走航投爆炸声源，垂直接收阵在固定点接收信号的声传播测量实验配置，根据噪声和爆炸波的频谱特性，处理前对接收数据高通滤波，然后，利用模型预报不同发射-接收阵距离的声传播损失值，结合近距离处实际接收的爆炸波样本能量统计值，估计实验中不同发射-接收阵距离处接收爆炸波信号的能量值，进而设计随发射-接收阵距离变化的检测门限，使得近发射-接收阵距离处门限值高，远距离处检测门限低，保证检测概率的条件下尽量降低虚警概率。该方法可以有效地自动截取爆炸波，对声传播测量实验数据自动化处理、降低人工干预度有重要意义。

Description

一种基于模型的海洋声传播调查爆炸波自动截取方法

技术领域

本发明属于水下瞬态信号检测领域，具体涉及一种基于模型的海洋声传播调查爆炸波自动截取方法。

背景技术

海洋声传播调查研究相关海域水声信道的能量衰减特性，既是水声信道物理规律性研究的基础，又为声纳设备的设计和研制提供依据。随着海洋声学调查的发展，从中可以发现其发展趋势是：海洋水声调查的要素越来越多、越来越精细、数据量越来越大，从而对数据处理技术的要求也越来越高。

爆炸声信号的检测和截取是实现高效率数据处理的关键因素。传统上对爆炸波的信号手动截取：将数据分成小段，观察每小段的波形或计算每小段的频谱，当观察到有类似爆炸波的信号，对信号进行截取。这样做效率低，虚警高，而且对数据处理人员专业性要求比较高。所以，开发自动信号截取的方法是做好声学数据处理的关键技术。数据处理过程中，对爆炸波信号自动检测与截取要求检测器有高检测概率和低虚警概率，分别对应信号漏报低和信号虚警低。这直接影响到后期将截取到的爆炸声信号和爆炸声源记录对应的相关工作量的大小。

海洋声传播调查实验中，一般采用双船配置、单船和锚底潜标配置，他们的共同点是，发射船沿着固定的测线走航，同时在不同的位置投放声源，声源可能采用爆炸波、等离子体源或有规发射信号等，接收垂直阵在固定的位置接收声源信号，这里考虑爆炸声源。

从大量的实验数据分析可知，噪声由海洋背景噪声和脉冲状干扰组成，前者可以近似认为是高斯噪声，后者是具有重拖尾的非高斯噪声。噪声的存在使爆炸波信号检测时引起虚警。由于爆炸波信号的波形并不确知，需要看作为一类随机信号，导出的检验统计量为能量检测器。当检验统计量大于设定的门限值时，我们认为检测到信号。然而由于一条测线通常跨越数十到数百公里，信号的能量存在3～5个量级的变化，导致几十分贝信噪比的变化，因而若采用单个门限值，在发射-接收阵距离很大时，爆炸波能量衰减大信噪比低，导致很高的漏报。若降低门限值，在发射-接收距离小时，又由于脉冲状干扰的存在，导致高的虚警。因此单一门限，很难满足实现海洋环境中近距离处的低虚警概率和远距离处的高检测概率的需求，有必要采用随着距离变化的自适应能量检测门限，近发射-接收阵距离处提高门限，以有效降低脉冲状噪声引起的虚警，远发射-接收阵距离处降低门限，以有效提高爆炸波的检测概率。基于此，提出了一种基于模型的海洋声传播调查爆炸波自动截取方法。

发明内容

针对现有技术存在的技术问题，本发明提供了一种基于模型的海洋声传播调查爆炸波自动截取方法，使得对同一条测线上声传播数据爆炸波检测有较高的检测概率和较低的虚警概率。

一种基于模型的海洋声传播调查爆炸波自动截取方法，所述测量方法包括如下步骤：

(1)对接收阵接收的数据，在声源-接收阵距离r₀附近人工截取爆炸波数据样本，计算每个样本所有通道滤波后爆炸波能量和平均值

(2)根据实际测量的海深、声速剖面、声源和接收阵位置、海底声学特性等信息，利用声传播模型预报水声调查环境下声源-接收阵距离r₀处的声传播损失值

(3)读入一帧接收阵接收的数据；

(4)对每个通道数据进行滤波处理；

(5)计算各个通道能量，和所有通道能量累加和，找到能量最大位置；

(6)以最大能量位置为基准，加窗w_n，利用窗外数据估计该帧数据滤波后单通道噪声成分的平均方差σ²；

(7)以最大能量位置为基准，加窗w_s，利用窗内数据计算检验统计量T；

(8)计算窗w_s的窗位置对应的时间，计算该时刻声源-接收阵之间的距离r；

(9)根据实际测量的海深、声速剖面、声源和接收阵位置、海底声学特性等信息，利用声传播模型预报水声调查环境下声源-接收阵距离r处的声传播损失值TL_r；

(10)根据(1)(2)(6)(9)的结果预估的爆炸波截取窗w_s内的总能量并设定检测门限γ'_max(r)；

(11)将噪声分解为独立的高斯海洋背景噪声和脉冲状干扰的线性叠加，根据给定的虚警概率P'_FA，计算只有高斯海洋背景噪声存在时的检测门限γ'(r)；

(12)根据以下公式计算门限：γ(r)＝max(γ'_max(r),γ'(r))

(13)如果检验统计量T大于门限值，则截取窗w_s内的信号作为爆炸波信号，否则，重复(3)～(13)。

所述的步骤(1)或(2)中发射-接收阵距离r₀在采用爆炸声源的前提下参考值为5～10倍接收阵布阵位置海深的距离。

所述的步骤(2)或(9)中利用声场传播模型预报不同声源-接收阵距离处的声传播损失值的具体方法为：首先，从接收数据中，人工截取n个爆炸波信号样本，分别做离散傅里叶变换，频谱归一化后求平均得到声源信号频谱；利用声传播模型计算并合成该声源在不同接收位置的宽带信号，计算宽带声源信号在相应声源-接收位置处的声传播损失，并对不同接收水听器深度的声传播损失取平均。

所述的步骤(1)或步骤(4)中滤波器为高通滤波器，滤波器的截止频率根据噪声或干扰占据的主要频带范围确定，截止频率取为有效噪声或干扰频带的上限，例如取100Hz或200Hz。

所述的步骤(5)中所有通道能量累加和计算方法为：首先，对每个通道信号按同样的时间窗分段，计算每小段内能量，分段长度根据数据特征设置，然后，将不同通道同一时间窗对应的段内能量累加，得到所有通道能量累加和。

所述的步骤(6)中平均方差的计算方法为：首先，以最大能量位置为基准加窗w_n，窗长应使得当存在爆炸波时将爆炸波全部包含在窗内；然后，将窗外完全不存在爆炸波的数据作为噪声，计算每个通道噪声的方差，再取平均得到单通道噪声平均方差。

所述的步骤(7)或(8)或(10)或(13)中窗w_s的长度应考虑信道传播的时延扩展，设置为大于近场(例如100～200米)测量到的爆炸波的长度，长度与声源-接收阵距离有关，距离越远，长度越大。

所述的步骤(6)中的窗w_n比所述的步骤(7)中的窗w_s的长，以保证爆炸波信号不会进入噪声方差的统计。

所述的步骤(7)中能量统计量的计算方法为下式：

$T=\underset{c=1}{\overset{M}{\mathrm{\Σ}}}\underset{i=i_0}{\overset{i_0+l_s-1}{\mathrm{\Σ}}}{x}_{\mathrm{ci}}^2$

其中，M表示接收水听器通道数，i₀表示加窗的起始位置坐标，x_ci表示第c个通道接收数据中坐标为i的值，l_s表示窗w_s的长度。

所述的步骤(8)中根据时间计算发射-接收阵距离的具体方法为：找出该时刻对应的发射点GPS(全球定位系统)和接收点GPS，计算出声源-接收阵距离。

所述的步骤(10)中同时存在高斯海洋背景噪声和脉冲状干扰时检测门限的计算方法为下式：

${\mathrm{\γ}}_{\max }^{\′}\left(r\right)=\mathrm{\α}\·{\stackrel{\‾}{E}}_{r_0}\·{10}^{(({\mathrm{TL}}_{r_0}-{\mathrm{TL}}_r)/10)}+l_s\·M\·{\mathrm{\σ}}^2$

其中，为发射-接收阵距离r₀处爆炸波样本滤波后所有通道能量累加和平均值，TL_r与分别为模型计算的声源-接收阵距离r时与声源-接收阵距离r₀时声传播损失，α∈(0,1)为调节系数，σ²为单通道噪声平均方差，M表示接收水听器通道数，l_s表示窗w_s的长度。

所述的步骤(11)中检测门限的计算方法为下式：

${\mathrm{\γ}}^{\′}\left(r\right)=l_s\·(Q^{-1}\left(P_{\mathrm{FA}}^{\′}\right)\·\sqrt{\frac{2{\left({\mathrm{\σ}}^2\right)}^2}{l_s}}+{\mathrm{\σ}}^2)$

其中，P'_FA为设定的虚警概率，Q^-1(P)为标准正态分布随机变量右尾函数的逆函数，γ'(r)为只有海洋背景噪声时由虚警概率确定的门限值，σ²为单通道噪声平均方差。

与现有技术相比，本发明方法具有的优点为：利用声场传播模型预报测线上不同发射-接收阵位置的声传播损失，进而预报不同发射-接收阵位置的爆炸波能量值。利用预报的能量值设置不同发射-接收阵距离时接收数据的能量检测器门限值。近发射-接收阵距离时，信噪比高，门限也设置较高，以有效降低脉冲状干扰引起的虚警概率；远发射-接收阵距离时，信噪比低，门限也相应降低，以有效提高检测概率。当然，门限值降低，随着检测概率的提高，虚警也会随着增大。但综合整条测线上爆炸波的检测和截取，可以同时获得较高的检测概率和较低的虚警概率。在后期的与投弹表对应过程中，可有效减少人工处理工作量，提高数据处理效率。

附图说明

图1为基于模型的海洋声传播调查爆炸波自动截取流程图。

图2为针对某次海洋声传播实验中一条测线数据的能量检测门限设置结果。

具体实施方式

本实施例中，声场传播模型采用Kraken简正波模型，模型输入参数为声速剖面、海深、海底声学特性(如沉积层厚度、密度、衰减等)、声源频率和发射与接收阵位置等信息，模型输出为单频信号从发射到接收阵的信道频率响应。计算频段(通常为三分之一倍频程宽)内所有频率点的信道频率响应，然后根据爆炸波信号的频谱特性，合成宽带信号。声场传播模型也可采用射线模型、抛物方程模型和快场模型等。

图1为基于模型的海洋声传播调查爆炸波自动截取流程图，主要步骤如下：

(1)读入垂直阵接收数据文件，显示所有通道波形，人工截取5～10倍接收阵布阵位置海深的发射-接收阵距离r₀附近的所有接收通道爆炸波样本。对每一接收通道数据做高通滤波(截止频率f_c根据噪声和干扰能量所占据的主要频带确定，取其频带的上限值，例如采用100Hz、200Hz)，然后，计算每个样本所有通道接收爆炸波数据的能量和，在样本集上求平均，记为

(2)利用声场传播模型预报水声调查环境下发射-接收阵距离r₀时的声传播损失值具体方法为：

(a)对发射-接收阵距离r₀时爆炸波样本单通道接收数据做截止频率为f_c的高通滤波；

(b)令表示高通滤波后的第q个样本，对其做离散傅里叶变换，有：

$X_{r_{0}q}\left(k\right)=\underset{j=0}{\overset{N-1}{\mathrm{\Σ}}}{x}_{r_{0}q}\left(j\right){\mathrm{\ω}}_N^{\mathrm{jk}},k=0,...,N-1.$

其中ω_N＝e^(-2πi)/N，N为离散傅里叶变换的点数，频率采样间隔为第k个频点的频率为

(c)计算Q个爆炸波样本归一化平均频谱：

${\stackrel{\‾}{X}}_{r_0}\left(k\right)=\frac{1}{Q}\underset{q=1}{\overset{Q}{\mathrm{\Σ}}}\left(\frac{X_{r_{0}q}\left(k\right)}{\underset{k\∈[1,N]}{\max }\left|\right|X_{r_{0}q}\left(k\right)\left|\right|}\right)$

(d)平均频谱采样值为其中，j＝k₁,k₂,k₃,...,k_J，对应需要计算的频段。利用Kraken声场模型计算声源-接收距离r₀时，不同接收水听器深度z_m处各频率点的信道传递函数，并合成宽带信号，计算各水听器位置处声传播损失，进而计算不同水听器深度的平均声传播损失

(3)读入一帧数据，记为X^M×P，M表示通道数，P表示每个通道的数据长度。

(4)对X^M×P的每行数据做截止频率为f_c的高通滤波，滤波后记为

(5)计算所有通道能量累加和，并找到能量最大位置。c通道第i段能量为

$E_{X_f}(c,i)=\underset{p=(i-1)\·L+1}{\overset{i\·L}{\mathrm{\Σ}}}{X}_f^2(c,p)$

其中，L表示单通道分段求能量时每小段的长度，c∈[1,M]。所有通道的能量累加和由下式求得：

$E_{X_{\mathrm{fa}}}\left(i\right)=\underset{c=1}{\overset{M}{\mathrm{\Σ}}}{E}_{X_f}(c,i)$

(6)找到最大值处位置，即以此为基准加窗w_n，那么窗对应的坐标范围为：

I_n＝[max(0,(II-1)·L-ε·l_n+1),min(P,(II-1)·L+(1-ε)·l_n)]

其中，ε∈(0,1)为窗位置相对于能量最大位置的调节系数，利用窗外的信号估计单通道噪声数据平均方差σ²；

(7)在最大值位置处，加窗w_s，该窗对应的坐标范围为I_w＝[max(0,(II-1)·L-β·l_s+1),min(P,(II-1)·L+(1-β)·l_s)]，其中，β∈(0,1)为窗位置相对于能量最大位置的调节系数，计算检验统计量T；

(8)计算窗w_s的窗位置对应的时间，找到该时刻的发射GPS和接收GPS，计算该时刻声源-接收阵之间的距离r。

(9)根据实际测量的海深、声速剖面、声源和接收阵位置、海底声学特性等信息，利用声传播模型预报水声调查环境下声源-接收阵距离r处的声传播损失值TL_r，具体处理方法如(2)(d)中所阐述。

(10)估计窗w_s内的信号能量为：

$\hat{E}=M\·l_s\·{\mathrm{\σ}}^2+{\stackrel{\‾}{E}}_{r_0}\·{10}^{(({\mathrm{TL}}_{r_0}-{\mathrm{TL}}_r)/10)}$

(11)根据预估的窗w_s内的总能量设定检测门限γ'_max(r)：

${\mathrm{\γ}}_{\max }^{\′}\left(r\right)=\mathrm{\α}\·{\stackrel{\‾}{E}}_{r_0}\·{10}^{(({\mathrm{TL}}_{r_0}-{\mathrm{TL}}_r)/10)}+l_s\·M\·{\mathrm{\σ}}^2$

(12)根据只有海洋背景噪声时的虚警概率P'_FA，计算检测门限：

${\mathrm{\γ}}^{\′}\left(r\right)=l_s\·(Q^{-1}\left(P_{\mathrm{FA}}^{\′}\right)\·\sqrt{\frac{2{\left({\mathrm{\σ}}^2\right)}^2}{l_s}}+{\mathrm{\σ}}^2)$

(13)计算门限γ(r)＝max(γ'_max(r),γ'(r))，将检验统计量T与门限γ(r)比较，如果T≥γ(r)，则将窗w_s内的数据截取为爆炸波，否则读取下一帧数据，继续进行(3)～(13)判决过程。

本方法的核心在于能量检测器门限的自适应设置，结合图1流程图，进行了随距离变化的判决门限仿真。图2给出了某次海洋声传播调查时，某一条测线的能量检测门限随距离变化曲线。图中只给出发射-接收阵距离20km以外的门限，20km以内由于信噪比很高，不会出现漏报和虚警情况。在仿真中，爆炸波采样点的发射-接收阵距离为20km，调节系数取为α＝0.01，窗长l_s＝10000，窗长l_n＝50000，采样率f_s＝10000Hz，当噪声成分只含有海洋背景噪声时，虚警概率设置为P′_FA＝0.005。图2表明：(1)当只有海洋背景噪声的时候，由虚警概率P′_FA决定的检测门限值低于实际采用的检测门限，说明实际数据处理过程中，虚警概率主要是由瞬态干扰引起的，这与实际情况是吻合的。(2)如果采用传统的单一门限，当选择门限1时，60km以外的爆炸波信号将会出现漏报，当选择门限2时，由于整条测线上门限值都很低，那么容易将瞬态干扰信号误判为爆炸波，造成很高的虚警，而采用随距离变化的门限，则可以克服门限1与门限2的缺点，兼顾了高检测概率和低虚警概率的要求。

Claims (10)

1.一种基于模型的海洋声传播调查爆炸波自动截取方法，其特征在于，包括如下步骤：

(1)对接收阵接收的数据，在声源-接收阵距离r₀附近人工截取爆炸波数据样本，计算每个样本所有通道滤波后爆炸波能量和平均值

(2)利用声传播模型预报水声调查环境下声源-接收阵距离r₀处的声传播损失值

(3)读入一帧接收阵接收的数据，并对每个通道数据进行滤波处理；

(4)计算各个通道能量，以及所有通道能量累加和，找到能量最大位置；

(5)以最大能量位置为基准，加窗w_n，利用窗外数据估计该帧数据滤波后单通道噪声成分的平均方差σ²；

(6)以最大能量位置为基准，加窗w_s，利用窗内数据计算检验统计量T；

(7)计算窗w_s的窗位置对应的时间，计算该时刻声源-接收阵之间的距离r；

(8)利用声传播模型预报水声调查环境下声源-接收阵距离r处的声传播损失值TL_r；

(9)根据σ²、TL_r预估爆炸波截取窗w_s内的总能量并计算检测门限γ'_max(r)；

(10)将噪声分解为独立的高斯海洋背景噪声和脉冲状干扰的线性叠加，根据给定的虚警概率P'_FA，计算只有高斯海洋背景噪声存在时的检测门限γ'(r)；

(11)根据公式γ(r)＝max(γ'_max(r),γ'(r))计算门限：

(12)如果检验统计量T大于门限值，则截取窗w_s内的信号作为爆炸波信号，否则，重复步骤(3)～步骤(12)。

2.如权利要求1所述的爆炸波自动截取方法，其特征在于，所述的步骤(1)和/或步骤(2)中发射-接收阵距离r₀在采用爆炸声源的前提下参考值为5～10倍接收阵布阵位置海深的距离。

3.如权利要求1所述的爆炸波自动截取方法，其特征在于，所述的步骤(2)和/或步骤(8)中，利用声场传播模型预报不同声源-接收阵距离处的声传播损失值的具体方法为：

首先从接收数据中截取n个爆炸波信号样本，分别做离散傅里叶变换，频谱归一化后求平均得到声源信号频谱；再利用声传播模型计算并合成该声源在不同接收位置的宽带信号，计算宽带声源信号在相应声源-接收位置处的声传播损失，并对不同接收水听器深度的声传播损失取平均。

4.如权利要求1所述的爆炸波自动截取方法，其特征在于，所述的步骤(4)中所有通道能量累加和计算方法为：首先，对每个通道信号按同样的时间窗分段，计算每小段内能量，分段长度根据数据特征设置；然后，将不同通道同一时间窗对应的段内能量累加，得到所有通道能量累加和。

5.如权利要求1所述的爆炸波自动截取方法，其特征在于，所述的步骤(5)中平均方差的计算方法为：以最大能量位置为基准加窗w_n，窗长应使得当该帧数据中存在爆炸波时爆炸波全部包含在窗内；再将窗外完全不存在爆炸波的数据作为噪声，计算每个通道噪声的方差，再取平均得到单通道噪声平均方差。

6.如权利要求1和5所述的爆炸波自动截取方法，其特征在于，所述窗w_s的长度应大于近场测量到的爆炸波的长度，且窗w_n的窗长大于窗w_s。

7.如权利要求1所述的爆炸波自动截取方法，其特征在于，所述的步骤(6)中检验统计量的计算方法为下式：

$T=\underset{c=1}{\overset{M}{\mathrm{\Σ}}}\underset{i=i_0}{\overset{i_0+l_s-1}{\mathrm{\Σ}}}{x}_{\mathrm{ci}}^2$

其中，M表示接收水听器通道数，i₀表示加窗的起始位置坐标，x_ci表示第c个通道接收数据中坐标为i的值，l_s表示窗w_s的长度。

8.如权利要求1所述的爆炸波自动截取方法，其特征在于，所述的步骤(7)中根据时间计算发射-接收阵距离的具体方法为：找出该时刻对应的发射点GPS和接收点GPS，计算出声源-接收阵距离。

9.如权利要求1所述的爆炸波自动截取方法，其特征在于，所述的步骤(9)中，同时存在高斯海洋背景噪声和脉冲状干扰时检测门限的计算方法为下式：

${\mathrm{\γ}}_{\max }^{\′}\left(r\right)=\mathrm{\α}\·{\stackrel{\‾}{E}}_{r_0}\·{10}^{(({\mathrm{TL}}_{r_0}-{\mathrm{TL}}_r)/10)}+l_s\·M\·{\mathrm{\σ}}^2$

其中，为发射-接收阵距离r₀处爆炸波样本滤波后所有通道能量累加和平均值，TL_r与分别为声源-接收阵距离r时与声源-接收阵距离r₀时声传播损失，α∈(0,1)为调节系数，σ²为单通道噪声平均方差，M表示接收水听器通道数，l_s表示窗w_s的长度。

10.如权利要求1所述的爆炸波自动截取方法，其特征在于：所述的步骤(10)中采用下式计算检测门限γ'(r)：

${\mathrm{\γ}}^{\′}\left(r\right)=l_s\·(Q^{-1}\left(P_{\mathrm{FA}}^{\′}\right)\·\sqrt{\frac{2{\left({\mathrm{\σ}}^2\right)}^2}{l_s}}+{\mathrm{\σ}}^2)$

其中，P'_FA为设定的虚警概率，Q^-1(P)为标准正态分布随机变量右尾概率的逆函数，γ'(r)为只有海洋背景噪声时由虚警概率确定的门限值，σ²为单通道噪声平均方差。