- 用窗函数法完成低通、高通滤波器的设计;
(1)低通滤波器指标:通带截止频率wp=π/4rad,阻带截止频率ws=π/2 rad,通带最大衰减ap=1 dB,阻带最小衰减a s=40 dB。
解:设计代码如下所示:
wp=0.25*pi;wr=0.5*pi;
tr_width=wr-wp; %过渡带宽度
N=ceil(6.2*pi/tr_width)+1 %滤波器的长度, N=奇数为 1 型; N=偶数为 2 型
n=0:1:N-1;
wc=(wr+wp)/2; %理想低通的截止频率
hd=ideal_lp(wc,N); %理想低通的单位脉冲响应
w_ham=(hanning(N))'; %汉宁窗
h=hd.*w_ham; %截取得到实际单位脉冲响应
[db,mag,pha,w]=freqz_m(h,1); %计算实际滤波器的幅度响应
delta_w=2*pi/1000;
Ap=-(min(db(1:1:wp/delta_w+1))); %实际通带波动
Ar=-round(max(db(wr/delta_w+1:1:501))); %最小阻带衰减
subplot(221);stem(n,hd);title('理想单位脉冲响应 hd(n)')
subplot(222);stem(n,w_ham);title('汉宁窗)w(n)')
subplot(223);stem(n,h);title('实际单脉冲响应 h(n)')
subplot(224);plot(w/pi,db);title('幅度响应(db)')
axis([0,1,-100,10]);
运行结果如下图所示:
(2)高通滤波器指标:通带截止频率wp=π/2 rad,阻带截止频率ws=π/4 rad,通带最大衰减a p=1 dB,阻带最小衰减a s=40 dB。
解:设计代码如下所示:
wp=0.5*pi;wr=0.25*pi;
tr_width=wp-wr; %过渡带宽度
N=ceil(6.2*pi/tr_width)+1 %滤波器的长度, N=奇数为 1 型; N=偶数为 2 型
n=0:1:N-1;
wc=(wr+wp)/2/pi; %理想高通的截止频率
hd=ideal_hp(wc,N); %理想高通的单位脉冲响应
w_ham=(hanning(N))'; %汉宁窗
h=hd.*w_ham; %截取得到实际单位脉冲响应
[db,mag,pha,w]=freqz_m(h,1); %计算实际滤波器的幅度响应
delta_w=2*pi/1000;
Ap=-(min(db(1:1:wp/delta_w+1))); %实际通带波动
Ar=-round(max(db(wr/delta_w+1:1:501))); %最小阻带衰减
subplot(221);stem(n,hd);title('理想单位脉冲响应 hd(n)')
subplot(222);stem(n,w_ham);title('汉宁窗)w(n)')
subplot(223);stem(n,h);title('实际单脉冲响应 h(n)')
subplot(224);plot(w/pi,db);title('幅度响应(db)')
axis([0,1,-100,10]);
运行结果如下图所示:
(3)带通滤波器的指标为wp1 = 0.35p , wp 2 = 0.65p ,ws1 = 0.2p , ws 2 = 0.8p ,as = 60dB至少用三种窗函数完成 FIR 数字滤波器的设计,并对各类滤波器的单位脉冲响应 h(n)和频率响应 H(ejω)进行比较和分析。
解:设计代码如下所示:
% 数字滤波器指标
ws1=0.2*pi; wp1=0.35*pi;
ws2=0.8*pi; wp2=0.65*pi;
As=60;
tr_width=abs(min((wp1-ws1),(wp2-ws2)));
M=ceil(11*pi/tr_width)+1;
n=0:1:M-1;
wc1=(ws1+wp1)/2;
wc2=(ws2+wp2)/2;
hd=ideal_lp(wc2,M)-ideal_lp(wc1,M);
% 生成 blackman 窗
w_bla=(blackman(M))';
h=hd.*w_bla
% 频域图像的绘制
freqz(h,1)
figure(2);
subplot(2,2,1),stem(n,hd);title('idael impulse response')
axis([0 M-1 -0.4 0.5]);xlabel('n');ylabel('hd(n)')
xa=0.*n;
hold on
plot(n,xa,'k');
hold off
subplot(2,2,2),stem(n,w_bla);title('blackman window')
axis([0 M-1 0 1.1]);xlabel('n');ylabel('w(n)')
subplot(2,2,3),stem(n,h);title('actual impulse response')
axis([0 M-1 -0.4 0.5]);xlabel('n');ylabel('h(n)')
hold on
plot(n,xa,'k');
hold off
运行结果如下图所示:
比较和分析:以上低通、高通、带通三种滤波器
- 利用 MATLAB 软件读入一段语音信号,对语音信号进行低通、高通、带通、带阻滤波,并对滤波前后的频谱进行分析,播放处理前后的语音,并对其实际效果从算法上加以分析。
(1)对信号进行低通滤波
[audio,fs]=audioread('C:\Users\DELL\Desktop\数字信号实验\caiq\ly.wav');%声音读取
audio = audio(:,1); %双通道变单通道
n=length(audio);
s=audio;
T = 1/fs;%采样间隔
t = (0:n-1)*T;%时间轴
f = (0:n-1)/n*fs;%频率轴样
%快速傅里叶变换
audio_fft=fft(audio,n)*T;
%设计低通滤波器
wp=0.25*pi;ws=0.5*pi;
tr_width=ws-wp; %过渡带宽度
N=ceil(6.2*pi/tr_width)+1;
wn=(ws+wp)/2; %理想低通的截止频率
[b,a]=butter(N,wn,'low','s'); %S域频率响应的参数即:滤波器的传输函数
[bz,az]=bilinear(b,a,0.5); %利用双线性变换实现频率响应S域到Z域的变换
figure(2);%低通滤波器特性
[h,w]=freqz(bz,az);
title('低通滤波器');
plot(w*fs/(2*pi),abs(h));
grid;
%滤波
z=filter(bz,az,s);
z_fft=fft(z); %滤波后的信号频谱
figure(1);
%绘出原始音频时域波
subplot(2,2,1);
plot(t,audio);
xlabel('时间/s');
ylabel('幅度');
title('初始信号波形');
grid;
%绘出原始音频频域频谱
subplot(2,2,2);
audiof = abs(audio_fft);
plot(f(1:(n-1)/2),audiof(1:(n-1)/2));
title('初始信号频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
%绘出滤波音频时域波
subplot(2,2,3);
plot(t,z);
title('低通滤波后的信号波形');
xlabel('时间/s');
ylabel('幅度');
grid;
%绘出滤波音频频域波
subplot(2,2,4);
zf = abs(z_fft);
plot(f(1:(n-1)/2),zf(1:(n-1)/2));
title('低通滤波后信号的频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
audio_final = [audio;s;z]; %原始语音,加噪语音,滤波语音的合成音频矩阵
sound(audio_final,fs); %播放语音
(2)对信号进行高通滤波
[audio,fs]=audioread('C:\Users\DELL\Desktop\数字信号实验\caiq\ly.wav');%声音读取
audio = audio(:,1); %双通道变单通道
n=length(audio);
s=audio;
T = 1/fs;%采样间隔
t = (0:n-1)*T;%时间轴
f = (0:n-1)/n*fs;%频率轴样
%快速傅里叶变换
audio_fft=fft(audio,n)*T;
%设计高通滤波器
wp=0.5*pi;ws=0.25*pi;
tr_width=wp-ws; %过渡带宽度
N=ceil(6.2*pi/tr_width)+1;
wn=(ws+wp)/2/pi; %理想高通的截止频率
[b,a]=butter(N,wn,'high','s'); %S域频率响应的参数即:滤波器的传输函数
[bz,az]=bilinear(b,a,0.5); %利用双线性变换实现频率响应S域到Z域的变换
figure(2);%低通滤波器特性
[h,w]=freqz(bz,az);
title('高通滤波器');
plot(w*fs/(2*pi),abs(h));
grid;
%滤波
z=filter(bz,az,s);
z_fft=fft(z); %滤波后的信号频谱
figure(1);
%绘出原始音频时域波
subplot(2,2,1);
plot(t,audio);
xlabel('时间/s');
ylabel('幅度');
title('初始信号波形');
grid;
%绘出原始音频频域频谱
subplot(2,2,2);
audiof = abs(audio_fft);
plot(f(1:(n-1)/2),audiof(1:(n-1)/2));
title('初始信号频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
%绘出滤波音频时域波
subplot(2,2,3);
plot(t,z);
title('高通滤波后的信号波形');
xlabel('时间/s');
ylabel('幅度');
grid;
%绘出滤波音频频域波
subplot(2,2,4);
zf = abs(z_fft);
plot(f(1:(n-1)/2),zf(1:(n-1)/2));
title('高通滤波后信号的频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
audio_final = [audio;s;z]; %原始语音,加噪语音,滤波语音的合成音频矩阵
sound(audio_final,fs); %播放语音
(3)对信号进行带阻滤波
[audio,fs]=audioread('C:\Users\DELL\Desktop\数字信号实验\caiq\ly.wav');%声音读取
audio = audio(:,1); %双通道变单通道
n=length(audio);
s=audio;
T = 1/fs;%采样间隔
t = (0:n-1)*T;%时间轴
f = (0:n-1)/n*fs;%频率轴样
%快速傅里叶变换
audio_fft=fft(audio,n)*T;
rp = 1;rs = 60;
wp1 = 0.35*pi;
wp2 = 0.65*pi;
ws1 = 0.2*pi; % 带阻截止频率1
ws2 = 0.8*pi; % 带阻截止频率2
[N,wn] = ellipord([wp1, wp2], [ws1, ws2], rp, rs, 's');
[b, a] = ellip(N, rp, rs, wn, 'stop', 's');
[bz,az] = bilinear(b, a, 0.5);
% 绘制带阻滤波器的频率响应
figure;
[h,w] = freqz(bz, az);
title('带阻滤波器');
plot(w * fs / (2 * pi), abs(h));
grid;
%滤波
z=filter(bz,az,s);
z_fft=fft(z); %滤波后的信号频谱
figure(1);
%绘出原始音频时域波
subplot(2,2,1);
plot(t,audio);
xlabel('时间/s');
ylabel('幅度');
title('初始信号波形');
grid;
%绘出原始音频频域频谱
subplot(2,2,2);
audiof = abs(audio_fft);
plot(f(1:(n-1)/2),audiof(1:(n-1)/2));
title('初始信号频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
%绘出滤波音频时域波
subplot(2,2,3);
plot(t,z);
title('带阻滤波后的信号波形');
xlabel('时间/s');
ylabel('幅度');
grid;
%绘出滤波音频频域波
subplot(2,2,4);
zf = abs(z_fft);
plot(f(1:(n-1)/2),zf(1:(n-1)/2));
title('带阻滤波后信号的频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
audio_final = [audio;s;z]; %原始语音,加噪语音,滤波语音的合成音频矩阵
sound(audio_final,fs); %播放语音
(4)对信号进行带通滤波
[audio,fs]=audioread('C:\Users\DELL\Desktop\数字信号实验\caiq\ly.wav');%声音读取
audio = audio(:,1); %双通道变单通道
n=length(audio);
s=audio;
T = 1/fs;%采样间隔
t = (0:n-1)*T;%时间轴
f = (0:n-1)/n*fs;%频率轴样
%快速傅里叶变换
audio_fft=fft(audio,n)*T;
rp = 1;
rs = 60;
wp1 = 0.35*pi;
wp2 = 0.65*pi;
ws1 = 0.2*pi; % 带通截止频率1
ws2 = 0.8*pi; % 带通截止频率2
[N,wn] = ellipord([wp1, wp2], [ws1, ws2], rp, rs, 's');
[b, a] = ellip(N, rp, rs, wn, 'bandpass', 's');
[bz,az] = bilinear(b, a, 0.5);
% 绘制带通滤波器的频率响应
figure;
[h,w] = freqz(bz, az);
title('带通滤波器');
plot(w * fs / (2 * pi), abs(h));
grid;
%滤波
z=filter(bz,az,s);
z_fft=fft(z); %滤波后的信号频谱
figure(1);
%绘出原始音频时域波
subplot(2,2,1);
plot(t,audio);
xlabel('时间/s');
ylabel('幅度');
title('初始信号波形');
grid;
%绘出原始音频频域频谱
subplot(2,2,2);
audiof = abs(audio_fft);
plot(f(1:(n-1)/2),audiof(1:(n-1)/2));
title('初始信号频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
%绘出滤波音频时域波
subplot(2,2,3);
plot(t,z);
title('带通滤波后的信号波形');
xlabel('时间/s');
ylabel('幅度');
grid;
%绘出滤波音频频域波
subplot(2,2,4);
zf = abs(z_fft);
plot(f(1:(n-1)/2),zf(1:(n-1)/2));
title('带通滤波后信号的频谱');
xlabel('频率/Hz');
ylabel('幅度');
grid;
audio_final = [audio;s;z]; %原始语音,加噪语音,滤波语音的合成音频矩阵
sound(audio_final,fs); %播放语音
分析:在算法上,低通滤波器去除高频部分,可能导致失真,适用于保留基本音调。高通滤波器去除低频部分,可能使语音尖锐,适用于强调清晰度。带通滤波器突出特定频率范围,适用于强调特定音素或特征。带阻滤波器去除特定频带内噪声,但可能损失语音信息。