实验三采用快速卷积法对先前的生成的混杂噪声的音频信号进行滤波

Forest Frame

已于 2022-11-21 16:27:38 修改

阅读量1.9k

点赞数 4

分类专栏： MATLAB数字信号处理文章标签：音视频 matlab

于 2022-11-17 21:46:07 首次发布

本文链接：https://blog.csdn.net/SnowyForest___/article/details/127912610

版权

MATLAB数字信号处理专栏收录该内容

5 篇文章 132 订阅

订阅专栏

本文介绍了如何使用快速卷积方法对含有噪声的音频信号进行滤波。通过MATLAB实现，确定单位取样响应，利用FFT进行卷积运算，比较了滤波前后的信号波形和效果，展示了快速卷积与线性卷积的差异。

摘要由CSDN通过智能技术生成

实验三采用快速卷积法对先前的生成的混杂噪声的音频信号进行滤波

要求：选择子作业１中的音频信号，自行给定滤波器的单位取样响应，采用快速卷积实现对音频信号的滤波，比较滤波前后信号的波形和回放的效果。
最终整体效果如下图：

alt

一、快速卷积的原理

快速卷积主要有FFT卷积和分块卷积两种模式，以下讨论FFT卷积形式。
设语音信号为 x(n) ，单位采样相应为 h(n) 。时域上做卷积即为频域上做乘积，因此根据这个原理，将 x(n) 和 h(n) 分别通过快速傅里叶变换FFT得到频域信号 $H(\omega)和X(\omega)$ 然后相乘得到 $Y(\omega)$ ，最后将其通过IFFT从频域变换为时域即可。其原理框图如下：

alt

其中 h(n) 与 x(n) 长度通常是不同的，但在频域相乘中 $H(\omega)和X(\omega)$ 必须是点数相同的，因此需要通过尾部补0的方式，将 x(n) 和 h(n) 的长度填充到L，L大于两者点数之和减一。

二、卷积滤波的实现

1、单位取样响应的确定（以下方法均为自行摸索得出，CSDN原创）

为了获得单位取样相应函数，这次我选用了无限长滤波器（因为有限长滤波器点数太多，且点幅度的主要部分分布在中间段，不方便取，无限长滤波器可以采取选择单位采样信号的部分点，且其点幅值主要部分在采样开始部分），在MATLAB滤波器设计工具中，我看到了所设计滤波器的Impulse Response即单位冲激响应波形，如下图所示：

alt

接着我放大图像，从起始部分开始选择了部分点。并将这些点的值添加进了MATLAB的工作区间（光标选中点后右键选择Export Cursor Data to Workspace）。数据点的添加如下图所示：

alt

将数据导入工作区间后会有个变量（我命名为hn），其为一个结构体变量，包括一个绘图对象Target和一个表示数据点位置和幅值信息的Position元胞数组。选中Position右键选择根据所选内容新建工作区变量，再选择新建元胞数组可得到一个新的关于单位取样相应的变量。操作如下图所示：

alt

最后将此变量导出为hn.mat数据文件, 供后续调用时使用。

2、卷积滤波的代码实现

①单位取样信号的获取

②快速卷积

③线性卷积

附录：MATLAB代码

clear;
clc;
format long;
close all;

% ==========原始信号========== %
[x, fs] = audioread('./实验三语音信号/小提琴.wav');
x = x(:, 1);
x = x';

N = length(x);  % 整个图由N1个样点构成
dt = 1 / fs;
tscale = dt * N;  % X轴显示的时间长度，单位为秒
t = 0 : dt : tscale - tscale / N;

subplot(1, 2, 1);
% subplot(2, 4, 1);
plot(t .* 1000, x);
title('小提琴信号时域图');
xlabel('t/ms', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;

y = fft(x);
realy = 2 * abs(y(1 : length(x))) / length(x);
realf = (0 : length(x) - 1) * (fs / length(x)); 
subplot(1, 2, 2);
% subplot(2, 4, 5);
stem(realf, realy, '.');
title('小提琴信号频谱图');
axis([0, 4000, 0, 0.04]);
xlabel('f/Hz', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;

figure;
% ==========小提琴混杂噪声信号========== %
[Xn, fs1] = audioread('./实验三语音信号/小提琴混杂声音_缩混.wav');
Xn = Xn(:, 1);
Xn = Xn';

N1 = length(Xn);  % 整个图由N1个样点构成
dt1 = 1 / fs1;
tscale1 = dt1 * N1;  % X轴显示的时间长度，单位为秒
t1 = 0 : dt1 : tscale1 - tscale1 / N1;

subplot(1, 2, 1);
% subplot(2, 4, 2);
plot(t1 .* 1000, Xn);
title('小提琴混杂噪声信号时域图');
xlabel('t/ms', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;

Y1 = fft(Xn);
realy = 2 * abs(Y1(1 : length(Xn))) / length(Xn);
realf = (0 : length(Xn) - 1) * (fs1 / length(Xn)); 
subplot(1, 2, 2);
% subplot(2, 4, 6);
stem(realf, realy, '.');
title('小提琴混杂噪声信号频谱图');
axis([0, 4000, 0, 0.04]);
xlabel('f/Hz', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;

% ==========卷积滤波========== %
hn = load('hn.mat');  %此文件数据是前面通过滤波器工具箱生成的滤波器中的得到的单位取样相应的部分数据点
for i = 1 : 122  %从保存数据的元胞数组中提取单位取样相应的幅值
    Hn(i) = hn.Position{i, 1}(1, 2);
end

L1 = pow2(nextpow2(length(Xn) + length(Hn) - 1));  %确定FFT快速卷积的点数
Xk = fft(Xn, L1);  %计算Xn的L点FFT,结果为Xn
Hk = fft(Hn, L1);  %计算Hn的L点FFT,结果为Hk
Yk = Xk .* Hk;  %计算YK,频域相乘即为时域相卷
y1n = ifft(Yk, L1);  %对YK调用IFFT，求得y1(n)
y2n = conv(Xn, Hn);  %计算y2(n)的卷积，此为线性卷积

figure;
% =====快速卷积===== %
dt2 = 1 / fs1;
tscale2 = dt2 * L1;  % X轴显示的时间长度，单位为秒
t2 = 0 : dt2 : tscale2 - tscale2 / L1;
subplot(1, 2, 1);
% subplot(2, 4, 3);
plot(t2 .* 1000, y1n);
title('快速卷积滤波后信号时域图');
xlabel('t/ms', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;
audiowrite('./实验三语音信号/快速卷积滤波后信号.wav', y1n, 8000);

Y2 = fft(y1n);
realy = 2 * abs(Y2(1 : length(y1n))) / length(y1n);
realf = (0 : length(y1n) - 1) * (fs1 / length(y1n)); 
subplot(1, 2, 2);
% subplot(2, 4, 7);
stem(realf, realy, '.');
title('快速卷积滤波后信号频谱图');
axis([0, 4000, 0, 0.04]);
xlabel('f/Hz', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;

figure;
% =====线性卷积===== %
L2 = length(y2n);
dt3 = 1 / fs1;
tscale3 = dt3 * L2;  % X轴显示的时间长度，单位为秒
t3 = 0 : dt3 : tscale3 - tscale3 / L2;
subplot(1, 2, 1);
% subplot(2, 4, 4);
plot(t3 .* 1000, y2n);
title('线性卷积滤波后信号时域图');
xlabel('t/ms', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;
audiowrite('./实验三语音信号/线性卷积滤波后信号.wav', y2n, 8000);

Y3 = fft(y2n);
realy = 2 * abs(Y3(1 : length(y2n))) / length(y2n);
realf = (0 : length(y2n) - 1) * (fs1 / length(y2n)); 
subplot(1, 2, 2);
% subplot(2, 4, 8);
stem(realf, realy, '.');
title('线性卷积滤波后信号频谱图');
axis([0, 4000, 0, 0.04]);
xlabel('f/Hz', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
ylabel('电压/V', 'FontName', '宋体', 'FontWeight', 'normal', 'FontSize', 14);
grid on;