CN115050386A - 一种中华白海豚哨叫声信号自动检测和提取方法 - Google Patents

Abstract

本发明公开了一种中华白海豚哨叫声信号自动检测和提取方法，步骤1：对采集到的包含中华白海豚哨叫声信号的数据进行预处理；步骤2：得到的包含中华白海豚哨叫声信号的数据进行截取，切分成若干个时间片段，对全部时间片段逐个执行如下操作，得到含有中华白海豚哨叫声信号的信号起止时间结果；步骤3：根据步骤2中自动检测并保存的含有中华白海豚哨叫声信号的信号起止时间结果，重新对步骤1预处理后的信号进行切分，然后对切分后得到的包含中华白海豚哨叫声信号的信号生成时频图，并自动将每段信号生成的时频图保存。本发明实现针对在海豚研究中进行长期声呐信号采集获得的海量数据的全自动检测和提取，降低对人工干预的依赖，提升效率。

Description

一种中华白海豚哨叫声信号自动检测和提取方法

技术领域

本发明属于水生生物学、水声学和水声信号处理领域，涉及一种中华白海豚哨叫声信号自动检测和提取方法。

背景技术

在遗传或者形态学数据无法获得的情况下，中华白海豚的哨叫声可以揭示不同地区种群之间关系，确定独立种群，特别是在地理上相邻的中华白海豚种群，对于预测其种群状况和调节其保护是至关重要的，因此哨声信号受到了广大研究人员的重视。但是从茫茫大海中采集并快速提取到高质量的中华白海豚哨叫声信号目前还有很大的难度，主要在于，包含中华白海豚哨叫声的信号在所有采集信号中占比非常小，并且对于其拾取目前主要通过人工手段对海量的信号进行筛选，对于人工成本和时间成本都是极大的开销。因此，如果能使用自动手段对走航和原位式方式采集的信号进行分析，从中自动检测并提取出具有高质量中华白海豚哨叫声信号的片段，对于研究效率将具有极大的提升。

发明内容

针对上述现有技术，本发明要解决的技术问题是提供一种中华白海豚哨叫声信号自动检测和提取方法，自动从海量声纳数据中检测并提取中华白海豚哨叫声信号，自动对走航和原位式方式采集的信号进行分析，从中自动检测并提取出具有高质量中华白海豚哨叫声信号的片段，提升工作效率，克服目前仅通过人工手段对海量的声纳信号进行筛选的问题。

为解决上述技术问题，本发明的一种中华白海豚哨叫声信号自动检测和提取方法，包括以下步骤：

步骤1：对采集到的包含中华白海豚哨叫声信号的数据进行预处理；

步骤2：对步骤1得到的包含中华白海豚哨叫声信号的数据进行截取，切分成若干个时间片段，对全部时间片段逐个执行如下操作，得到含有中华白海豚哨叫声信号的信号起止时间结果：

步骤A:对时间片段逐帧进行短时功率谱计算，形成二维时频图；

步骤B：对时间片段的二维时频图，设置相关样本长度CorrNum和时间窗长度SegNum，满足SegNum＞CorrNum，CorrNum表示每一个采样点计算与之相关的后续点个数，SegNum表示每个时间片段被切分为若干个时间窗的长度，将时间片段按时间窗长度SegNum划分为若干个时间窗，对每个时间窗执行以下操作：对时间窗内的所有采样点对应的时频图进行互相关计算，得到互相关矩阵secCoord_sup，对矩阵secCoord_sup沿行方向求和，得到时间窗内的积分列向量；

步骤C：求取时间片段的互相关积分曲线；

步骤D：根据设定的阈值判断时间片段的互相关积分曲线中是否存在中华白海豚哨叫声信号，并得到信号起止时间；

步骤3：根据步骤2中自动检测并保存的含有中华白海豚哨叫声信号的信号起止时间结果，重新对步骤1预处理后的信号进行切分，然后对切分后得到的包含中华白海豚哨叫声信号的信号按照步骤2中方法生成时频图，并自动将每段信号生成的时频图保存。

进一步的，步骤1中预处理包括：降采样处理、高通滤波和时域平滑滤波。

进一步的，步骤B中对时间窗内的所有采样点对应的时频图进行互相关计算，得到互相关矩阵secCoord_sup包括：

对时间窗内的所有采样点对应的功率谱进行互相关计算：

其中，ρ_XY为第X列和第Y列的互相关值，Cov(X,Y)代表时频图第X列和第Y列的协方差，σ_X代表时频图第X列的标准差，σ_Y代表时频图第Y列的标准差；ρ_XY作为互相关矩阵secCoord第X行、第Y列的元素，得到互相关矩阵secCoord，所述互相关矩阵secCoord为SegNum维方阵；取secCoord的上三角矩阵secCoord_utri，每行保留从secCoord_utri主对角线元素开始的前CorrNum个元素取值，其他位置全部置0，令q＝SegNum-CorrNum+1，从第q+i开始，每行补充i个元素ρ_{(q+i)(SegNum+a)}得到矩阵secCoord_sup，其中i＝1,2,3,…,CorrNum-1，a＝1,2，…，i。

进一步的，步骤C中求取时间片段的互相关积分曲线包括：

将时间片段内全部时间窗的积分列向量依次连接组成一个列向量矩阵，所述列向量矩阵中第m行元素表示时间片段内第m个采样点对应的互相关值，然后将列向量矩阵转置，建立以时间为横坐标轴、以互相关值为纵坐标轴的二维坐标系，以采样点的采样时间为横坐标，以采样点的互相关值为纵坐标，得到时间片段各个采样点的互相关积分曲线。

进一步的，步骤D中根据设定的阈值判断时间片段的互相关积分曲线中是否存在中华白海豚哨叫声信号，并计算其起止时间包括：

当采样点的互相关值大于设定的阈值，则判断所述采样点包含中华白海豚哨叫声信号，在所述采样点前的采样点中，距离所述采样点最近且互相关值与阈值重合的采样点对应的采样时间为中华白海豚哨叫声信号的起始时间，记为t_start；

当采样点的互相关值小于设定的阈值，则判断所述采样点中华白海豚哨叫声信号消失，在所述采样点前的采样点中，距离所述采样点最近且互相关值与阈值重合的采样点对应的采样时间为中华白海豚哨叫声信号的终止时间，记为t_end；

当时间片段中包含中华白海豚哨叫声信号起始时间但不包含其终止时间时，记录中华白海豚哨叫声信号的起始时间，并在下个时间片段时记录中华白海豚哨叫声信号的终止时间。

本发明的有益效果：本发明利用中华白海豚的哨叫声信号相比于click信号及其他噪声信号具有显著的周期性和强相关性的特点，采用适当的方式去除不同机理产生的干扰噪声，从而突出哨叫声信号特征，提升了检测精度，并实现了全自动化检测和提取。本发明可以实现针对在海豚研究中进行长期声纳信号采集获得的海量数据的全自动检测和提取，降低对人工干预的依赖，提升效率。

附图说明

图1是本发明的整体流程框图；

图2是自动提取的中华白海豚哨叫声结果；

图3是自动提取的中华白海豚哨叫声保存结果。

具体实施方式

下面结合说明书附图和实施例对本发明做进一步说明。

由于实际采集到鲸目动物哨声观测信号中并不是每一时刻都含有哨声，因此需要一个合适的哨声信号检测方法，从观测信号中检测并确定哨声信号位置，通过检测出的信号位置重新切分并生成时频图。结合图1，本发明包含以下步骤:

步骤1：中华白海豚哨叫声信号预处理。对采集到的包含中华白海豚哨叫声信号的数据进行降采样处理，提高处理速度；针对中华白海豚的哨叫声信号的频域特点，选择合适频带对降采样的信号进行高通滤波，滤除掉显著低于哨叫声信号频带内的杂波，去除海洋背景噪声的干扰，提高信噪比；针对中华白海豚在浅海时会同时发出哨叫声和click信号(回声定位信号串)的特点，使用最小二乘平滑滤波对信号进行时域平滑滤波，剔除掉click信号和同类脉冲信号干扰，降低短时脉冲信号对哨叫声频谱分析的干扰，提高信噪比。

步骤2：中华白海豚哨叫声信号自动检测。针对中华白海豚哨叫声信号持续时间有一定范围，对包含中华白海豚哨叫声信号的数据进行截取，切分成若干个时间片段，每次取1个时间片段逐帧(帧是指每个时间点)进行短时功率谱计算，形成二维时频图；根据中华白海豚哨叫声在一个周期内的各帧功率谱是连续的且彼此具有强相关性的特点，对于各个时间片段的二维时频图，设置相关样本长度和时间窗长度(时间片段长度大于时间窗长度大于相关样本长度)，相关样本长度用以决定对于每一个采样点需要计算与之相关的后续点的个数，时间窗长度用以决定每个时间片段被切分为若干个时间窗的长度，之后逐个对时间窗内的所有采样点对应的时频图进行互相关计算，得到互相关矩阵。为简化计算，取其上三角矩阵，沿对角线方向每一行的保留的长度等于相关样本长度，其他位置全部置0。由于是上三角矩阵，会出现保留长度小于相关样本长度的情况，因此要补充对应缺失的值。对于补充值后的矩阵，沿行方向对其求和，得到时间窗内的积分列向量。对每一个时间窗重复上述步骤，将全部时间窗的互相关积分组成一个新的矩阵，该向量是当前时间片段内所有采样点对应的互相关值，即第一行代表第一个采样点的互相关值，然后将该向量转置，按照每个采样点的顺序依次将互相关值画在二维坐标系上，得到当前时间片段各个时间点的互相关积分曲线。

相关积分阈值分析。对当前时间片段的相关积分按时间点形成的互相关积分曲线，通过设置合理的阈值对其进行判决，自动从互相关积分曲线中判断是否存在海豚的哨叫声信号，并计算其起止时间，某个时间点的互相关积分结果大于阈值，则判断为当前时间点有海豚哨叫声信号，在大于阈值之前的那个重合点即为哨叫声信号的起始时间；某个时间点的互相关积分结果小于阈值，则判断为当前时间点海豚哨叫声信号消失，在小于阈值之前的那个重合点即为哨叫声信号的终止时间。由于每次切片时间为1-2秒，当前时间片段白海豚哨叫声会出现有始无终的情况，此时，仅记录哨叫声起点，根据中华白海豚哨叫声信号持续时间范围是60-1972ms，在下个时间片段记录哨叫声终点，并将最后将起止时间记录并保存。

之后对每个时间片段进行以上的步骤，得到全部信号的起止时间。

步骤3：中华白海豚哨叫声信号自动提取。

根据步骤2中自动检测并保存的含有中华白海豚哨叫声信号的起止时间结果，重新对全部时间信号切分，对这些信号使用步骤2中的参数生成时频图，并自动将每段信号生成的时频图保存。

下面结合具体参数给出实施例。

结合图1，本发明包括以下步骤：

步骤1：中华白海豚哨叫声信号预处理。

波形信号降采样，提高处理速度。由于中华白海豚(中国海南地区)的哨叫声信号频率范围0.17-15.80kHz，因此首先将信号采样率降至36-42kHz，得到数据dataRaw_ds。

波形信号降噪，提高信噪比和检测精度。由于信号采集过程中，在700Hz以下的频段存在大量背景噪声，且中华白海豚哨叫声信号基频范围是0.71-17.30kHz，因此对时序信号dataRaw_ds进行高通滤波，设置通带频率fpass为700Hz，得到滤波后信号data_hipass。针对夹杂在中华白海豚哨叫声信号中的click信号，其频带和强度均会对信号的检测产生显著的干扰，但是click信号属于典型的脉冲信号，可以使用时域平滑进行降噪，在本专利中，对data_hipass进行基于时域的最小二乘平滑滤波，用一阶Savitzky-Golay滤波器平滑，得到平滑后的时序信号data_dn。

步骤2：中华白海豚哨叫声信号自动检测。

针对中华白海豚哨叫声信号持续时间范围是60-1972ms，对降噪后的观测信号进行截取，切分成若干个长度为2s时间片段，信号长度为Ns，设置采样率为Fs，即每个时间片段包含2Fs个采样点；然后每次取出1个时间片段：S_i(n),n＝1,2,3,...,2Fs,

其傅里叶变换结果表示为X_i(k),k＝1,2,3,...,2Fs,

再计算信号的功率谱

k＝1,2,3,...,2Fs,

对该时间片段应用平滑的汉宁窗、680点的FFT、98％窗口重叠率、500Hz频率分辨率计算0-20kHz的功率谱并生成时频图。

计算互相关矩阵。对该时间片段内的时频图,设置相关样本长度CorrNum和时间窗长度SegNum(SegNum＞CorrNum)，CorrNum用以决定对于每一个采样点需要计算与之相关的后续点的个数，SegNum用以判断时间窗的长度。依次取SegNum个采样点集合SecMeta(j),j＝1,...,M,M＝2Fs/SegNum；之后对SecMeta(j)中的所有采样点对应的功率谱进行互相关计算：

其中，ρ_XY代表时频图的第X列和第Y列的互相关值即Pearson相关系数，Cov(X,Y)代表时频图第X列和第Y列的协方差，σ_X代表时频图第X列的标准差。之后得到该SecMeta集合的互相关矩阵secCoord，其中secCoord为SegNum维方阵，secCoord矩阵为：

互相关矩阵简化处理。secCoord为对称矩阵，为了简化计算，取其上三角矩阵，记为secCoord_utri：

由于每一次计算只关心与每一个采样点临近的CorrNum个采样点的相关性，因此，矩阵secCoord_utri沿对角线方向每一行只保留前CorrNum个值，其后的全部位置的点均被置为0，记为secCoord_trap：

但当CorrNum+q-1＝SegNum时，即q＝SegNum-CorrNum+1时，第q行以后(不含q行)的每行不足CorrNum个值(共CorrNum-1行不足CorrNum个值)，故第q+1行补一个值为ρ_{(q+1)(SegNum+1)}，第q+2行补两个值为ρ_{(q+2)(SegNum+1)}和ρ_{(q+2)(SegNum+2)}，依此类推，第SegNum行(最后一行)补CorrNum-1个值。补值后的矩阵记为secCoord_sup:

生成互相关积分曲线。对于矩阵secCoord_sup，沿行方向对其求和，得到积分向量

之后对每个时间窗重复上述步骤，即对每一个SecMeta(j),j＝1,2,3,...,2Fs/SegNum均计算积分向量CCAvg(j)，并将所有CCAvg(j)按顺序连接起来，得到第i个时间片段的总相关积分向量AllCCAvg(i)，AllCCAvg(i)的维度是1列，2Fs行。该向量是当前时间片段内所有采样点对应的互相关值，即第一行代表第一个采样点的互相关值，然后将该向量转置，按照每个采样点的顺序依次将互相关值画在二维坐标系上，得到当前时间片段各个时间点的互相关积分曲线。

相关积分阈值分析。对于由总相关积分向量AllCCAvg(i)按时间点形成的互相关积分曲线，通过设置合理的阈值shd对其进行判决，根据设置好的阈值，自动从互相关积分曲线中判断是否存在海豚的哨叫声信号，并计算其起止时间，一个时间点的互相关积分结果大于阈值，则判断为当前时间点有海豚哨叫声信号，在大于阈值之前的那个重合点即为哨叫声信号的起始时间，记为t_start；一个时间点的互相关积分结果小于阈值，则判断为当前时间点海豚哨叫声信号消失，在小于阈值之前的那个重合点即为哨叫声信号的终止时间，记为t_end，由于每次切片时间为2秒，当前时间片段中华白海豚哨叫声会出现有始无终的情况，此时，仅记录哨叫声起点，根据中华白海豚哨叫声信号持续时间范围是60-1972ms，当处理下个时间片段时记录哨叫声终点。之后对每个时间片段进行以上的步骤，得到全部时间信号的起止时间，并记录保存为

步骤3：中华白海豚哨叫声信号自动提取。

根据步骤2中自动检测并保存的含有中华白海豚哨叫声信号的起止时间结果，重新对全部时间信号切分，对这n段信号使用步骤2中的参数生成时频图，并自动将每段信号生成的时频图保存，保存后的信号如图2和图3所示。

Claims (5)

1.一种中华白海豚哨叫声信号自动检测和提取方法，其特征在于，包括以下步骤：

步骤1：对采集到的包含中华白海豚哨叫声信号的数据进行预处理；

步骤2：对步骤1得到的包含中华白海豚哨叫声信号的数据进行截取，切分成若干个时间片段，对全部时间片段逐个执行如下操作，得到含有中华白海豚哨叫声信号的信号起止时间结果：

步骤A:对时间片段逐帧进行短时功率谱计算，形成二维时频图；

步骤B：对时间片段的二维时频图，设置相关样本长度CorrNum和时间窗长度SegNum，满足SegNum＞CorrNum，CorrNum表示每一个采样点计算与之相关的后续点个数，SegNum表示每个时间片段被切分为若干个时间窗的长度，将时间片段按时间窗长度SegNum划分为若干个时间窗，对每个时间窗执行以下操作：对时间窗内的所有采样点对应的时频图进行互相关计算，得到互相关矩阵secCoord_sup，对矩阵secCoord_sup沿行方向求和，得到时间窗内的积分列向量；

步骤C：求取时间片段的互相关积分曲线；

步骤D：根据设定的阈值判断时间片段的互相关积分曲线中是否存在中华白海豚哨叫声信号，并得到信号起止时间；

步骤3：根据步骤2中自动检测并保存的含有中华白海豚哨叫声信号的信号起止时间结果，重新对步骤1预处理后的信号进行切分，然后对切分后得到的包含中华白海豚哨叫声信号的信号按照步骤2中方法生成时频图，并自动将每段信号生成的时频图保存。

2.根据权利要求1所述的一种中华白海豚哨叫声信号自动检测和提取方法，其特征在于：步骤1所述预处理包括：降采样处理、高通滤波和时域平滑滤波。

3.根据权利要求1所述的一种中华白海豚哨叫声信号自动检测和提取方法，其特征在于：步骤B所述对时间窗内的所有采样点对应的时频图进行互相关计算，得到互相关矩阵secCoord_sup包括：

对时间窗内的所有采样点对应的功率谱进行互相关计算：

其中，ρ_XY为第X列和第Y列的互相关值，Cov(X,Y)代表时频图第X列和第Y列的协方差，σ_X代表时频图第X列的标准差，σ_Y代表时频图第Y列的标准差；ρ_XY作为互相关矩阵secCoord第X行、第Y列的元素，得到互相关矩阵secCoord，所述互相关矩阵secCoord为SegNum维方阵；取secCoord的上三角矩阵secCoord_utri，每行保留从secCoord_utri主对角线元素开始的前CorrNum个元素取值，其他位置全部置0，令q＝SegNum-CorrNum+1，从第q+i开始，每行补充i个元素ρ_{(q+i)(SegNum+a)}得到矩阵secCoord_sup，其中i＝1,2,3,…,CorrNum-1，a＝1,2，…，i。

4.根据权利要求1所述的一种中华白海豚哨叫声信号自动检测和提取方法，其特征在于：步骤C所述求取时间片段的互相关积分曲线包括：

将时间片段内全部时间窗的积分列向量依次连接组成一个列向量矩阵，所述列向量矩阵中第m行元素表示时间片段内第m个采样点对应的互相关值，然后将列向量矩阵转置，建立以时间为横坐标轴、以互相关值为纵坐标轴的二维坐标系，以采样点的采样时间为横坐标，以采样点的互相关值为纵坐标，得到时间片段各个采样点的互相关积分曲线。

5.根据权利要求1所述的一种中华白海豚哨叫声信号自动检测和提取方法，其特征在于：步骤D所述根据设定的阈值判断时间片段的互相关积分曲线中是否存在中华白海豚哨叫声信号，并计算其起止时间包括：

当采样点的互相关值大于设定的阈值，则判断所述采样点包含中华白海豚哨叫声信号，在所述采样点前的采样点中，距离所述采样点最近且互相关值与阈值重合的采样点对应的采样时间为中华白海豚哨叫声信号的起始时间，记为t_start；

当采样点的互相关值小于设定的阈值，则判断所述采样点中华白海豚哨叫声信号消失，在所述采样点前的采样点中，距离所述采样点最近且互相关值与阈值重合的采样点对应的采样时间为中华白海豚哨叫声信号的终止时间，记为t_end；

当时间片段中包含中华白海豚哨叫声信号起始时间但不包含其终止时间时，记录中华白海豚哨叫声信号的起始时间，并在下个时间片段时记录中华白海豚哨叫声信号的终止时间。

Images