1 简介
窗函数法是设计数字滤波器的一种方法,本文采用各种矩形窗、汉宁窗、哈明窗、布莱克曼窗等不同窗形式对滤波器性能的影响.最后运用Matlab进行仿真比较.
2 完整代码
function y = plot_freq(Fs,L,x,s,i)
% 绘制频域幅度波形图
% Fs:采样频率
% L:序列长度
% x:幅值
% s:标题名称
% i:标识序号
figure(i);
y_fft = fft(x);
P2 = abs(y_fft/L);
P1 = P2(1:floor(L/2)+1);
P1(2:end-1) = 2*P1(2:end-1);
f = Fs*(0:(L/2))/L;
plot(f,P1);
title(s);
xlabel('f (Hz)');
ylabel('|X(k)|');
axis([0 3000 0 0.005]);
%{
y_fft = (fft(x)/L*2); % 幅度图 /(N/2) 恢复真实幅值
y2 = abs(fftshift(y_fft)); % 平移
f = (Fs * (0:L-1)/L)';
f_s = (-L/2:L/2-1)*(Fs/L);
plot(f_s,y2);
%plot(f(1:floor(length(f)/2)),y_fft(1:floor(length(f)/2))); % 取正主值区间
title(s);
axis([0 2000 0 max(y2)]); % 取主频部分
xlabel('f (Hz)');
ylabel('|X(k)|');
%}
function y = plot_time(t,x,s,i)
% 绘制时域波形图
% t:时间轴
% x:幅值
% s:标题名称
% i:标识序号
figure(i);
plot(t,x);
title(s);
xlabel('t (seconds)');
ylabel('x(k)');
%% 声音的采集
clc;
clear;
close all;
filename = '三狮.wav';
[y,Fs] = audioread(filename);
y = y(:,1); % 取单声道
L = length(y); % 采样点
T = 1/Fs; % 周期
t = (0:L-1)*T;
info = audioinfo(filename); % wav详细信息
disp(info); % 输出音频信号 采样频率,采样位数 等信息
sound(y,Fs); %以采样率 Fs 向扬声器发送音频信号 y
i = 0;
%% 绘图
i = i+1;
s = '原音频信号时域波形图';
plot_time(t,y,s,i); % 绘制时域波形图
i = i+1;
s = '原音频信号频域波形图';
plot_freq(Fs,L,y,s,i); % 绘制频域波形图
%% 窗函数对频谱分析的影响
%% 矩形窗
win_r = boxcar(length(y));
y_win_r = y.*win_r; % 时域乘积
sound(y_win_r,Fs);
i = i+1;
s = '矩形窗时域波形图';
plot_time(t,y_win_r,s,i);
i = i+1;
s = '矩形窗频率域波形图';
plot_freq(Fs,L,y_win_r,s,i);
% 矩形窗 1/2
win_r_2 = boxcar(floor(length(y)/2)); % 取语音长度1/2为窗口长度
win_r_2 = [win_r_2;zeros(ceil(length(y)/2),1)]; % 对齐补0
y_win_r_w = y.*win_r_2;
sound(y_win_r_w,Fs);
i = i+1;
s = '矩形窗(N = L/2)时域波形图';
plot_time(t,y_win_r_w,s,i);
i = i+1;
s = '矩形窗(N = L/2)频率域波形图';
plot_freq(Fs,L,y_win_r_w,s,i);
%% 汉宁窗
win_h = hann(length(y));
y_win_h = y.*win_h;
sound(y_win_h,Fs);
i = i+1;
s = '汉宁窗时域波形图';
plot_time(t,y_win_h,s,i);
i = i+1;
s = '汉宁窗频率域波形图';
plot_freq(Fs,L,y_win_h,s,i);
% 汉宁窗 1/2
win_h_2 = hann(floor(length(y)/2));
win_h_2 = [win_h_2;zeros(ceil(length(y)/2),1)]; % 对齐补0
y_win_h_2 = y.*win_h_2;
sound(y_win_h_2,Fs);
i = i+1;
s = '汉宁窗(N = L/2)时域波形图';
plot_time(t,y_win_h_2,s,i);
i = i+1;
s = '汉宁窗(N = L/2)频率域波形图';
plot_freq(Fs,L,y_win_h_2,s,i);
%% 哈明窗
win_hm = hamming(length(y));
y_win_hm = y.*win_hm;
sound(y_win_hm,Fs);
i = i+1;
s = '哈明窗时域波形图';
plot_time(t,y_win_hm,s,i);
i = i+1;
s = '哈明窗频率域波形图';
plot_freq(Fs,L,y_win_hm,s,i);
% 哈明窗 1/2
win_hm_2 = hamming(floor(length(y)/2));
win_hm_2 = [win_hm_2;zeros(ceil(length(y)/2),1)]; % 对齐补0
y_win_hm_2 = y.*win_hm_2;
sound(y_win_hm_2,Fs);
i = i+1;
s = '哈明窗(N = L/2)时域波形图';
plot_time(t,y_win_hm_2,s,i);
i = i+1;
s = '哈明窗(N = L/2)频率域波形图';
plot_freq(Fs,L,y_win_hm_2,s,i);
%% 布莱克曼窗
win_blk = blackman(length(y));
y_win_blk = y.*win_blk;
sound(y_win_blk,Fs);
i = i+1;
s = '布莱克曼窗时域波形图';
plot_time(t,y_win_blk,s,i);
i = i+1;
s = '布莱克曼窗频率域波形图';
plot_freq(Fs,L,y_win_blk,s,i);
% 布莱克曼窗 1/2
win_blk_2 = blackman(floor(length(y)/2));
win_blk_2 = [win_blk_2;zeros(ceil(length(y)/2),1)]; % 对齐补0
y_win_blk_2 = y.*win_blk_2;
sound(y_win_blk_2,Fs);
i = i+1;
s = '布莱克曼窗(N = L/2)时域波形图';
plot_time(t,y_win_blk_2,s,i);
i = i+1;
s = '布莱克曼窗(N = L/2)频率域波形图';
plot_freq(Fs,L,y_win_blk_2,s,i);
%% 截取2s
time = 2;
samples = [1,floor(L/2)];
[y_cut,Fs_cut] = audioread(filename,samples);
y_cut = y_cut(:,1); % 取单声道
L_cut = length(y_cut); % 采样点
T_cut = 1/Fs_cut; % 周期
t_cut = (0:L_cut-1)*T_cut;
sound(y_cut,Fs_cut);
i = i+1;
s = '音频信号截取2s时域波形图';
plot_time(t_cut,y_cut,s,i); % 绘制时域波形图
i = i+1;
s = '音频信号截取2s频域波形图';
plot_freq(Fs_cut,L_cut,y_cut,s,i); % 绘制频域波形图
%% 降采样对频谱分析的影响
Fs_res = Fs/2;
y_res=resample(y,Fs_res,Fs);%信号降采样处理
% y_res = downsample(y,2);
L_res = length(y_res);
sound(y_res,Fs_res);
t_res = 1/(Fs_res) * (0:L_res-1);
i = i+1;
s = '降采样时域波形图';
plot_time(t_res,y_res,s,i);
i = i+1;
s = '降采样频率域波形图';
plot_freq(Fs_res,L_res,y_res,s,i);
3 仿真结果
4 参考文献
[1]汪依帆, 王敏, 李海龙. 利用FFT进行信号频谱分析及MATLAB实现[J]. 数码设计, 2019, 8(7):1.
[2]余会娟. 基于窗函数法的数字滤波器语音信号处理[J]. 电脑知识与技术:学术版, 2018(1X):4.
博主简介:擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真,相关matlab代码问题可私信交流。
部分理论引用网络文献,若有侵权联系博主删除。