一种可配置多模式水声信号特征提取方法与流程

本发明涉及水声信号处理，具体涉及一种可配置多模式水声信号特征提取方法。

背景技术：

1、随着海洋的战略地位日益突显，各国都在积极开发利用海洋资源和空间。声波是目前在海洋中唯一能够进行远距离传播的能量形式，水声目标识别对于海洋开发、国防安全有着重大意义，现已成为水声领域的研究热点之一。水下目标自动识别主要包括特征提取与构建分类器两大部分。当前主流的特征提取方法包括时域波形结构分析、频域谱估计以及时频域分析3个方面。频域谱估计可提取信号的频率、功率、包络等特征，以及利用高阶谱分析非高斯信号的特征。这类方法原理简单、易于实现，仅通过采集到的原始水声信号即可获得，但是提取的特征需要一定的先验知识进行信号预处理，在时变的海洋环境下泛化性较弱。时频分析方法提供了时间域与频率域的联合分布信息，可以清楚地描述信号频率随时间变化的关系，是目前应用效果最好、应用最广的特征提取方法。在特征提取之后再训练隐马尔可夫模型(hidden markov model,hmm)、支持向量机、k近邻、神经网络等分类器以实现水下目标的识别。

2、具体在特征提取方面，研究者们多年研究发展出了时域特征，频谱特征，非线性特征和听觉特征等等特征方法。时域分析本身直观方便理解，但往往反映的是信号统计特性随时间变化的整体特征，因此需要对整段音频的局部特征做精准的描述。由此研究者发展出了短时傅里叶变换(short-time fourier transform,stft),小波变换(wavelettransform,wt),魏格纳-威尔分布(wigner-ville distribution,wvd)等联合时域与频域特性的二维时域空间算法。以混沌，分形等数学理论为基础发展的非线性特征在发展中尚未足够成熟，且相关计算较为复杂。频谱特征作为常用的信号处理概念，广泛用于船舶噪声特征提取，且由原来的傅里叶频谱、功率谱不断丰富扩展，发展到lofar谱，听觉谱，mel频率倒谱等等，越来越契合人耳的听觉感知模型。除此之外，近年图像处理技术和机器视觉的迅猛发展也催生了将音频切片，将一维的声音信号转化为二维的频谱图从而提取图像特征并使用机器视觉进行分类的方法，特征提取越来越多样化，精准化。

3、针对一系列的特征提取，大多数是使用软件方法实现，其中最多的是使用librosa库来进行特征提取，librosa是一个python模块，用于分析一般的音频信号，是一个非常强大的python语音信号处理的第三方库。调用librosa库，可以加载音频文件，读取采样率，并且librosa库也提供了stft、mfcc等特征提取工具。但是librosa提供的库并不能看到源代码。

4、目前针对水声信号时频域特征提取的硬件加速方法通常只针对特定的特征提取算法进行优化，导致了硬件固定化。这限制了在不同应用场景下适应不同特征提取需求的能力。其次，现有的硬件加速方法往往缺乏足够的灵活性，无法在运行时进行算法修改和参数调整。这对于应对实时变化的需求或不同信号类型的适应性较差。由于硬件资源有限，例如对于fpga而言，乘法器资源较为有限，难以满足较多的乘法运算，也是目前导致性能瓶颈和限制性能的原因。

5、可见现有的方案因各种限制因素导致特征处理技术频谱少，参数不可配置。

技术实现思路

1、有鉴于此，本发明提供了一种可配置多模式水声信号特征提取方法，能够克服现有特征处理技术频谱少，参数不可配置的缺点，提供一种可配置水声信号特征预处理方法，从而实现多频谱，可配置的提取水声信号特征。

2、为达到上述目的，本发明的技术方案为：所提取的水声信号特征包含提取的特征包括stft谱，lofar谱、mel功率谱以及mfcc系数，具体包括如下步骤：

3、接收水声信号，对所述水声信号进行预处理，将预处理后的水声信号直接进行fft变换，再进行功率谱计算，结果得到stft谱。

4、采用mel尺度的三角函数滤波器组针对stft功率谱进行mel滤波，获得mel功率谱。

5、将得到的mel功率谱进行一个dct变换处理最终得到mfcc系数。

6、将预处理后的水声信号进行归一化和中心化处理，将经过归一化和中心化处理的水声信号进行fft变换得到lofar谱。

7、进一步地，对所述水声信号进行预处理，具体过程为：

8、s101：采用预加重滤波器针对所述水声信号进行预加重处理；

9、s102：将预加重处理后的水声信号进行分帧处理，分好的帧组成帧集；

10、s103：针对所述帧集中的每一个帧乘以一个窗口函数，获得加窗后的新帧，即为预处理结果。

11、进一步地，预加重滤波器具体为：

12、pre_mphasized(t)＝signal(t)-α×signal(t-1)

13、其中signal(t)为t时刻的水声信号，pre_mphasized(t)为针对t时刻的水声信号进行预加重处理后的信号，α为预加重滤波器系数。

14、进一步地，将预加重处理后的水声信号进行分帧处理，分好的帧组成帧集，具体为：

15、s1021：进行分帧参数设置，所设置的分帧参数包括帧长frame_length，帧移frame_step，信号的长度signal_length，帧的数目number_frames，信号补零的长度pad_signal_length，补零后的信号pad_signal，帧集frames；

16、s1022：取信号前frame_length个采样点为第一帧；

17、根据帧移frame_step第二次从第frame_step个采样点开始往后取frame_length个采样点为第二帧；

18、第三次从frame_length+frame_step个采样点开始往后取frame_length个采样点，为第三帧；

19、以此类推，直到取完最后frame_length点作为最后一帧，这些帧组成帧集frames；

20、若信号的长度signal_length减去帧长frame_length之后能够被帧移frame_step整除，则不需要补零，否则对信号进行补零使得信号的长度signal_length减去帧长frame_length之后能够被帧移frame_step整除。

21、进一步地，针对所述帧集中的每一个帧乘以一个窗口函数，其中窗口函数为hamming窗口，采用查表法实现加窗处理，具体地，将hamming窗函数的各个值存在ram里面，依次查找各值对每一帧进行加窗处理。

22、进一步地，采用mel尺度的三角函数滤波器组针对stft功率谱进行mel滤波，获得mel功率谱，具体为：

23、mel滤波的过程包括s201～s206，每帧的stft功率谱经过以上s201～s206的mel滤波过程处理之后，再将输出能量结果取对数得到每一帧的mel功率谱；

24、s201～s206分别为：

25、s201：首先确定最低频率fmax、最高频率fmax和滤波器的个数m；所述最低频率fmax是梅尔刻度的起点，对应于线性频率的0hz，最高频率fmax是采样频率的一半，基于奈奎斯特定理nyquist theorem确定；通过确定最低频率和最高频率得到梅尔刻度的范围；

26、s202：使用从频率转换为梅尔刻度的公式计算最高频率和最低频率对应的梅尔刻度mel(f)，具体为：

27、

28、将最低频率fmax、最高频率fmax分别作为f代入上述从频率转换为梅尔刻度的公式中分别得到梅尔刻度的最小值lowmel和最大值highmel。

29、s203：在梅尔刻度上均匀地分布滤波器，首先根据滤波器个数m来计算相邻两个梅尔滤波器的频率间距dmel，根据滤波器数量和间隔，在起点和终点的梅尔刻度值之间均匀地选择中心频率的梅尔刻度值；其中计算频率间距dmel的公式如下：

30、

31、其中dmel为相邻两个梅尔滤波器的频率间距，highmel为梅尔刻度的最大值，lowmel为梅尔刻度的最小值；由此得到第一个中心频率的梅尔刻度值为dmel，第二个中心频率的梅尔刻度值为2dmel，以此类推，最后一个中心频率的梅尔刻度值为mdmel。

32、s204：对于每个滤波器的中心频率的梅尔刻度值，使用梅尔刻度的逆变换公式将其转换为对应的线性频率f_center，具体为

33、

34、由此得到第一个滤波器的线性中心频率为第二个中心频率为以此类推，最后一个中心频率为

35、s205：得到每个梅尔滤波器的线性中心频率之后，将这些频率映射到fft的频率分量上，以便在频谱上应用滤波器响应；

36、计算出fft频率分辨率δf，具体为：其中fs为采样频率，fft_length为离散傅里叶变换(fft)的点数，δf为fft频率分辨率；

37、然后对于每个滤波器的线性中心频率f_center，将其除以δf并取整得到滤波器在fft频谱分量上的对应下标k：

38、s206：梅尔滤波器的频率响应采用三角形滤波器的形式表示，计算第m个滤波器的频率响应，其中心频率为并且已知相邻滤波器的中心频率为和计算相邻滤波器中心频率之间的频率偏移

39、然后计算左边界频率和右边界频率

40、将得到的中心频率左边界频率右边界频率代入三角函数滤波器的频率响应表达式，由此计算到第m个滤波器在频率轴上的每个频率f对应的频率响应hm(f)；

41、最后对于每个滤波器，将滤波器的频率响应值hm(f)与fft频谱数据在对应的频率分量上的值相乘，得到每个滤波器对应的输出能量：输出能量＝hm(f)×fft频谱数据[k]。

42、进一步地，将预处理后的水声信号进行归一化和中心化处理，将经过归一化和中心化处理的水声信号进行fft变换得到lofar谱，具体分为如下步骤：

43、s301：根据加窗后得到的新帧，找出每一帧的最大值和最小值并且求得每一帧的平均值，归一化公式如下：

44、

45、其中frame为加窗后的帧，mean_frame为加窗后的帧的平均值，max_frame为加窗后的帧的最大值，min_frame为加窗后的帧的最小值，normalization_frame为归一化后的帧；

46、s302：将归一化之后的每一帧求平均值，然后每个归一化后的帧中的数据减去这个平均值，得到中心化之后的帧，公式如下：

47、centralization_frame＝normalization_frame-mean_normalization_frame

48、其中centralization_frame为中心化后的帧，mean_normalization_frame为归一化后的帧的平均值；

49、s303：将经过归一化和中心化处理的新帧centralization_frame进行快速傅里叶变换-fft变换得到lofar谱。

50、进一步地，s303中，采用的快速傅里叶变换是按时间抽取的基2-fft算法，其计算模式为定点数计算，快速傅里叶变换输入的点数为2的n次幂，在硬件实现中将所使用的旋转因子存在ram中，实部和虚部分开计算，其中每一帧的长度是给定的，fft计算的点数跟帧的长度一致。

51、有益效果：

52、1：本发明提供了一种可配置多模式水声信号特征提取方法，所提取的水声信号特征包含提取的特征包括stft谱，lofar谱、mel功率谱以及mfcc系数，对水声信号进行预处理，将预处理后的水声信号直接进行fft变换，再进行功率谱计算，结果得到stft谱；采用mel尺度的三角函数滤波器组针对stft功率谱进行mel滤波，获得mel功率谱；将得到的mel功率谱进行一个dct变换处理最终得到mfcc系数；将预处理后的水声信号进行归一化和中心化处理，将经过归一化和中心化处理的水声信号进行fft变换得到lofar谱。以上步骤均可在fpga硬件中实现，无需进行硬件固定化，因此能够适用不同应用场景下进行不同特征提取的需求，从而实现多频谱，可配置的提取水声信号特征。

53、2：本发明在预加重过程中，将除法运算化简为移位运算，降低了计算复杂度。在加窗操作中，采用查表法实现加窗处理，将窗函数的计算预先存储在查找表中，避免了在运行时进行复杂的数学计算，节省了fpga资源；并且查表法可以在一个时钟周期内完成窗函数值的查找，从而实现快速的加窗操作。