CN108387887A - 一种水声信号的混合降噪方法 - Google Patents

Abstract

本发明公开一种水声信号的混合降噪方法，对含噪水声信号进行VMD分解，其中VMD的分解数量由EMD的分解结果确定；计算由VMD分解出的每个变分模态分量和原始信号的相关系数，根据相关系数阈值筛选出噪声变分模态分量；将所有噪声变分模态分量进行求和，然后进行SVD差分谱降噪；重构信号并进行Savitzky‑Golay滤波，得到降噪后的水声信号；计算水声信号降噪前后的关联维度和噪声强度，给出降噪前后的时序波形和二维相图；本发明能够进一步消除噪声，降噪后水声信号的相图具有更光滑的轨道，吸引子具有更好的几何结构规律性。

Description

一种水声信号的混合降噪方法

技术领域

本发明涉及水声信号处理技术领域，具体为一种水声信号的混合降噪方法。

背景技术

文献“Research on ship-radiated noise denoising using secondaryvariational mode decomposition and correlation coefficient,Sensors,2018,18(1):48(1-17)”公开了一种基于变分模态分解和相关系数的水声信号降噪方法，首先对原始水声信号进行变分模态分解(Variational Mode Decomposition，VMD)，然后使用相关系数求出每个变分模态分量u_k和原始水声信号的相关程度，将相关系数小于0.2的变分模态分量标记为噪声成分并予以剔除，最后重构信号，取得较好的降噪效果。但是需要说明的是，该方法存在将噪声变分模态分量直接剔除而导致有用信息被去除的问题。

发明内容

本发明的目的在于提供一种水声信号的混合降噪方法，解决了现有的水声降噪方法，在去除噪声的时候，存在会将有用信息直接被去除的缺陷。

为了达到上述目的，本发明采用的技术方案是：

步骤1，对原始含噪水声信号进行VMD分解得到一系列变分模态分量u_k；

步骤2，计算由VMD分解出的每个变分模态分量u_k和原始含噪水声信号的相关系数ρ，然后将所得的相关系数ρ与噪声阈值进行比较，其中，若所得的相关系数ρ小于噪声阈值，则筛选出的为噪声变分模态分量u_k'，否则，筛选出的为纯洁变分模态分量u_k”；

步骤3，将所有筛选出的噪声变分模态分量u_k'进行求和，然后再进行SVD差分谱降噪，得到降噪后的噪声成分；

步骤4，将降噪后的噪声成分与纯洁变分模态分量u_k”进行重构，并进行Savitzky-Golay滤波，得到降噪后的水声信号。

优选地，步骤1中，首先将原始含噪水声信号进行EMD分解，得到一系列固有模态函数，接着将固有模态函数的数量设置为VMD的分解数目。

优选地，步骤2中，按以下公式计算每个变分模态分量u_k和原始信号的相关系数ρ：

式中，f和u_k分别代表原始信号和变分模态分量，D和Cov分别代表数学期望和协方差。

优选地，步骤3中，对所有噪声变分模态分量u_k'之和进行SVD差分谱降噪的具体方法是：

首先，假设含噪的离散时间序列X＝[x(1),x(2),…,x(N)]；

然后，通过对离散时间序列X＝[x(1),x(2),…,x(N)]进行相空间重构，得到Hankel矩阵A_Hankel：

式中，N＝m+n-1，N为Hankel矩阵A_Hankel中的数据个数；m为Hankel矩阵A_Hankel中的的行数，且m的数值由N的奇偶性来确定：若N是偶数，则m＝N/2；否则m＝(N+1)/2；

接着，对A_Hankel进行奇异值分解，可得

式中，U、V是正交矩阵，Σ是非负对角阵；S＝diag(σ₁,σ₂,…,σ_r)，r是矩阵A_Hankel的秩，σ₁,σ₂,…,σ_r是A_Hankel的奇异值；

接着，将奇异值σ₁,σ₂,…,σ_r的正向差分b(i)按下式定义，进而得到奇异值的差分谱序列B＝(b₁,b₂,…,b_r-1):

b(i)＝|σ_i+1－σ_i|,i＝1,2，…,r-1

接着，通过奇异值的差分谱序列B确定SVD差分谱的最大值所对应的点K，同时将小于等于K的对应的奇异值进行保留，将大于K的对应的奇异值置为0，得到S₁＝diag(σ₁,σ₂,…,σ_k)，

最后，对Σ₁进行SVD反重构，即可实现信号降噪。

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

本发明公开一种水声信号的混合降噪方法，对含噪水声信号进行VMD分解，其中VMD的分解数量由EMD的分解结果确定；计算由VMD分解出的每个变分模态分量和原始信号的相关系数，根据相关系数阈值筛选出噪声变分模态分量；将所有噪声变分模态分量进行求和，然后进行SVD差分谱降噪；重构信号并进行Savitzky-Golay滤波，得到降噪后的水声信号；计算水声信号降噪前后的关联维度和噪声强度，给出降噪前后的时序波形和二维相图。本发明能够进一步消除噪声，降噪后水声信号的相图具有更光滑的轨道，吸引子具有更好的几何结构规律性。

附图说明

图1为本发明提供的一种水声信号的混合降噪方法的流程图；

图2为本发明水声信号的实测数据与本发明降噪方法结果数据的对比图；

图3为本发明水声信号的实测数据与文献降噪方法结果数据的对比图；

图4为本发明降噪方法与文献降噪方法降噪前后相图对比图，其中，图4a是降噪前的相图、图4b是文献降噪方法降噪后的相图、图4c是本发明降噪方法降噪后的相图。

具体实施方式

下面结合具体的实施例对本发明做进一步的详细说明，所述是对本发明的解释而不是限定。

下面结合附图和具体实施方式对本发明做进一步详细描述。

如图1所示，本发明公开一种水声信号的混合降噪方法，包括如下步骤：

步骤1，对原始含噪水声信号进行EMD分解为一系列固有模态函数，统计该一系列固有模态函数的数量，并将其设置为VMD的分解数量，对含噪水声信号进行VMD分解为一系列变分模态分量u_k；

步骤2，计算由VMD分解出的每个变分模态分量u_k和原始含噪水声信号的相关系数ρ，然后将所得的相关系数ρ与噪声阈值进行比较，用于从一系列变分模态分量u_k中筛选出噪声变分模态分量u_k'和纯洁变分模态分量u_k”，其中，若所得的相关系数ρ小于阈值，则筛选出的为噪声变分模态分量u_k'，剩余的为纯洁变分模态分量u_k”，同时，噪声阈值设置为0.2；

步骤3，将所有筛选出的噪声变分模态分量u_k'进行求和，然后再进行SVD差分谱降噪，得到降噪后的噪声成分；

步骤4，将降噪后的噪声成分与纯洁变分模态分量u_k”进行重构，并进行Savitzky-Golay滤波，得到降噪后的水声信号。

具体的，本发明提供的一种水声信号的混合降噪方法，包括如下步骤：

步骤1，对原始含噪水声信号进行EMD分解为一系列固有模态函数，统计固有模态函数的数量，并将其设置为VMD的分解数量，对含噪水声信号进行VMD分解为一系列变分模态分量u_k；

本实施例中，对舰船辐射噪声信号进行采样处理，舰船辐射噪声信号的频率为20kHz，数据样本长2048点，将其作为实测数据，其时域波形如图2所示。

将2048个数据进行EMD分解得到一系列固有模态函数，同时求出分解得到的固有模态函数的数量，并将其设置为VMD的分解数量，将原始信号进行VMD分解,得到一系列变分模态分量u_k；

步骤2，计算由VMD分解出的每个变分模态分量u_k和原始含噪水声信号的相关系数，根据相关系数的阈值筛选出噪声变分模态分量u_k'，一系列变分模态分量u_k筛选后剩余纯洁变分模态分量u_k”；其中，阈值设置为0.2；

根据以下公式求得每个变分模态分量u_k和原始信号的相关系数ρ：

式中，f和u_k分别代表原始信号和变分模态分量，D和Cov分别代表数学期望和协方差，ρ代表了变分模态分量u_k和原始信号的相关系数，设置噪声阈值为0.2，若某个变分模态分量u_k和原始信号的相关系数小于0.2，即可判定变分模态分量u_k为噪声变分模态分量u_k'，剩余的为纯洁变分模态分量u_k”。

步骤3，将所有筛选出的噪声变分模态分量u_k'进行求和，然后再进行SVD差分谱降噪，得到降噪后的噪声成分；

将所有噪声变分模态分量进行求和，奇异值降噪的阶次选择方式有中值、均值、差分谱等，本发明采用奇异值分解差分谱降噪，具体的方法是：

首先，假设含噪的离散时间序列X＝[x(1),x(2),…,x(N)]，通过对离散时间序列X＝[x(1),x(2),…,x(N)]进行相空间重构得到Hankel矩阵A_Hankel：

式中，N＝m+n-1，N为Hankel矩阵A_Hankel中的数据个数；m为Hankel矩阵A_Hankel中的行数，m的数值由N的奇偶性来确定：如果N是偶数，则m＝N/2；否则m＝(N+1)/2。

其次，对A_Hankel进行奇异值分解，可得

式中，U、V是正交矩阵，Σ是非负对角阵，S＝diag(σ₁,σ₂,…,σ_r)，r是矩阵A_Hankel的秩，σ₁,σ₂,…,σ_r是A_Hankel的奇异值；

接着，将奇异值σ₁,σ₂,…,σ_r的正向差分b(i)定义为:

b(i)＝|σ_i+1－σ_i|,i＝1,2，…,r-1

然后，序列B＝(b₁,b₂,…,b_r-1)被称为奇异值的差分谱，它合理地描述了两个相邻奇异值之间的变化状态；

通过B确定SVD差分谱的最大值的对应点K，同时将小于等于K的对应的奇异值进行保留，将大于K的对应的奇异值置为0，得到S₁＝diag(σ₁,σ₂,…,σ_k)，

最后，对Σ₁进行SVD反重构，即可实现信号降噪。

步骤4，将降噪后的噪声成分与纯洁变分模态分量u_k”进行重构，并进行Savitzky-Golay滤波，得到降噪后的水声信号；

将降噪后的噪声成分与纯洁变分模态分量u_k”叠加进行重构，然后将重构后的信号进行Savitzky-Golay滤波，得到降噪后的水声信号。

试验结束后，通过水声信号降噪前后的关联维度和噪声强度，给出降噪前后的时序波形和二维相图。

降噪结果如图2所示，黑色虚线代表原始的水声信号，黑色实线代表降噪后的水声信号。

为了显示算法的优越性，画出了降噪前后的二维相图，并与文献降噪方法作比较，结果如图4所示，可以看出本发明降噪后的相图具有更光滑的轨道，吸引子具有更好的几何结构规律性。同时计算水声信号降噪前后的关联维度和噪声强度，并与文献降噪方法作比较：

对于离散时间序列X＝[x(1),x(2),…,x(N)]，它的噪声强度(Noise Intensity，NI)可以近似表示为：

式中：为时间序列的均值。

关联维度D描述如下：

考虑某一相空间的序列X＝[x(1),x(2),…,x(N)]，把所有可能在距离上比给定的距离ε小的点对数占总的点对数的比重称为关联系数，它的表达式为

当N→∞，并且距离ε比较小时，若存在一个指数规律C(ε)∝ε^D，则定义关联维度D为：

式中：首先对离散时间序列进行相空间重构，然后对于给定的距离ε，计算关联系数C(ε,N)，并求解lnC(ε,N)关于lnε的线性区域的斜率，从而得到时间序列的关联维度D。

水声信号降噪前后的各特征参数比较如表1所示。基于本发明提供的水声信号混合降噪方法的关联维度和噪声强度分别为1.5355和0.1955，均小于文献降噪方法的降噪结果。由图4可知，可以看出本发明降噪后的相图具有更光滑的轨道，吸引子具有更好的几何结构规律性，进一步证明本发明能够更好地消除噪声。所以本发明提供的水声信号混合降噪方法能够很好的对水声信号时间序列进行降噪，是一种比较好的降噪方法。

文献“Research on ship-radiated noise denoising using secondaryvariational mode decomposition and correlation coefficient,Sensors,2018,18(1):48(1-17)”公开了一种水声信号降噪方法，从其降噪后的信号可以看出，降噪效果不明显，画出其相图可以看出自相似结构有些杂乱。由表1和图3可知，采用本发明提出的水声信号降噪方法对水声信号进行降噪处理后，水声信号的关联维和噪声强度同降噪前相比较均有明显的改善，降噪效果优于文献降噪方法。

表1水声信号降噪前后的各特征参数比较

通过对实测水声信号进行降噪处理，得到了较为满意的结果。结果表明，采用本发明提供的水声信号混合降噪方法对水声信号进行降噪处理是一种行之有效的方法，在水声信号的预测、提取与分类等领域具有重要作用。

本发明所举的具体实施例仅是对此发明精神的诠释，本发明技术领域的技术人员可以对描述的具体实施例进行修改或类似的方法替代，并不偏离本发明的精神。

Claims (4)

1.一种水声信号的混合降噪方法，其特征在于，包括以下步骤：

步骤1，对原始含噪水声信号进行VMD分解得到一系列变分模态分量u_k；

步骤2，计算由VMD分解出的每个变分模态分量u_k和原始含噪水声信号的相关系数ρ，然后将所得的相关系数ρ与噪声阈值进行比较，其中，若所得的相关系数ρ小于噪声阈值，则筛选出的为噪声变分模态分量u_k'，否则，筛选出的为纯洁变分模态分量u_k”；

步骤3，将所有筛选出的噪声变分模态分量u_k'进行求和，然后再进行SVD差分谱降噪，得到降噪后的噪声成分；

步骤4，将降噪后的噪声成分与纯洁变分模态分量u_k”进行重构，并进行Savitzky-Golay滤波，得到降噪后的水声信号。

2.根据权利要求1所述的一种水声信号的混合降噪方法，其特征在于，步骤1中，首先将原始含噪水声信号进行EMD分解，得到一系列固有模态函数，接着将固有模态函数的数量设置为VMD的分解数目。

3.根据权利要求1所述的一种水声信号的混合降噪方法，其特征在于，步骤2中，按以下公式计算每个变分模态分量u_k和原始信号的相关系数ρ：

式中，f和u_k分别代表原始信号和变分模态分量，D和Cov分别代表数学期望和协方差。

4.根据权利要求1所述的一种水声信号的混合降噪方法，其特征在于，步骤3中，对所有噪声变分模态分量u_k'之和进行SVD差分谱降噪的具体方法是：

首先，假设含噪的离散时间序列X＝[x(1),x(2),…,x(N)]；

然后，通过对离散时间序列X＝[x(1),x(2),…,x(N)]进行相空间重构，得到Hankel矩阵A_Hankel：

式中，N＝m+n-1，N为Hankel矩阵A_Hankel中的数据个数；m为Hankel矩阵A_Hankel中的的行数，且m的数值由N的奇偶性来确定：若N是偶数，则m＝N/2；否则m＝(N+1)/2；

接着，对A_Hankel进行奇异值分解，可得

式中，U、V是正交矩阵，Σ是非负对角阵；S＝diag(σ₁,σ₂,…,σ_r)，r是矩阵A_Hankel的秩，σ₁,σ₂,…,σ_r是A_Hankel的奇异值；

接着，将奇异值σ₁,σ₂,…,σ_r的正向差分b(i)按下式定义，进而得到奇异值的差分谱序列B＝(b₁,b₂,…,b_r-1):

b(i)＝|σ_i+1－σ_i|,i＝1,2，…,r-1

接着，通过奇异值的差分谱序列B确定SVD差分谱的最大值所对应的点K，同时将小于等于K的对应的奇异值进行保留，将大于K的对应的奇异值置为0，得到S₁＝diag(σ₁,σ₂,…,σ_k)，

最后，对Σ₁进行SVD反重构，即可实现信号降噪。