MATLAB窗函数法实现FIR的高通,带通和低通滤波器的程序(2)

figure(4);

plot(y);xlabel('时间/s');ylabel('振幅');title('滤波前振幅特性'); figure(5);

subplot(2,2,1),plot(yy1);xlabel('时间/s');ylabel('振幅');title('hanning窗函数滤波振幅特性'); subplot(2,2,2),plot(yy2);xlabel('时间/s');ylabel('振幅');title('hamming窗函数滤波振幅特性'); subplot(2,2,3),plot(yy3);xlabel('时间/s');ylabel('振幅');title('blackman窗函数滤波振幅特性'); subplot(2,2,4),plot(yy4);xlabel('时间/s');ylabel('振幅');title('kaiser窗函数滤波振幅特性'); %滤波前后的信号的相位对比 figure(6);

plot(angle(Y));xlabel('时间/s');ylabel('相位');title('滤波前的相位特性'); figure(7);

subplot(2,2,1),plot(angle(YY1));xlabel('时间/s');ylabel('相位');title('hanning窗函数滤波相位特性');

subplot(2,2,2),plot(angle(YY2));xlabel('时间/s');ylabel('相位');title('hamming窗函数滤波相位特性');

subplot(2,2,3),plot(angle(YY3));xlabel('时间/s');ylabel('相位');title('blackman窗函数滤波相位特性');

subplot(2,2,4),plot(angle(YY4));xlabel('时间/s');ylabel('相位');title('kaiser窗函数滤波相位特性');

2.设计一个hanning hamming blackman kaiser窗函数bandpass_FIR %设计一个hanning hamming blackman kaiser窗函数bandpass_FIR function bandpassfilter Fs=100;%采样频率

fp1=15;%通带下限截止频率 fp2=20;%通带上限截止频率 fs1=10; fs2=25;

wp1=2*pi*fp1/Fs;%将通带下限截止频率转换为数字滤波器频率 wp2=2*pi*fp2/Fs;%将通带上限截止频率转换为数字滤波器频率 ws1=2*pi*fs1/Fs;%将通带下限截止频率转换为数字滤波器频率 ws2=2*pi*fs2/Fs;%将通带上限截止频率转换为数字滤波器频率

Bt=wp1-ws1;

N0=ceil(6.2*pi/Bt); N=N0+mod(N0+1,2);

wn=[(wp1+ws1)/2/pi,(wp2+ws2)/2/pi];

window1=hanning(N);%使用hanning窗函数 window2=hamming(N);%使用hamming窗函数 window3=blackman(N);%使用blackman窗函数 %设过渡带宽度为5Hz

[n,Wn,beta,ftype]=kaiserord([10 15 20 25],[0 1 0],[0.01 0.01 0.01],100);%求阶数n以及参数beta window4=kaiser(n+1,beta);%使用kaiser窗函数 %设计加窗函数fir1

b1=fir1(N-1,wn,window1);

b2=fir1(N-1,wn,window2); b3=fir1(N-1,wn,window3);

b4=fir1(n,Wn,window4,'noscale'); %求取频率响应

[H1,W1]=freqz(b1,1,512,2); [H2,W2]=freqz(b2,1,512,2); [H3,W3]=freqz(b3,1,512,2); [H4,W4]=freqz(b4,1,512,2);

figure(1);

subplot(2,2,1),plot(W1,20*log10(abs(H1)));%绘制频率响应图形 axis([0,1,-100,100]);

title('带通hanning窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

subplot(2,2,2),plot(W2,20*log10(abs(H2)));%绘制频率响应图形 axis([0,1,-100,100]);

title('带通hamming窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

subplot(2,2,3),plot(W3,20*log10(abs(H3)));%绘制频率响应图形 axis([0,1,-100,100]);

title('带通blackman窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

subplot(2,2,4),plot(W4,20*log10(abs(H4)));%绘制频率响应图形 axis([0,1,-100,100]);

title('带通kaiser窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

T=1/Fs;

L=100;%信号长度

t=(0:L-1)*T;%定义时间范围和步长

y=sin(2*pi*5*t)+5*sin(2*pi*15*t)+8*sin(2*pi*40*t);%滤波前的图形 NFFT = 2^nextpow2(L); % Next power of 2 from length of y Y = fft(y,NFFT)/L;%将时域信号变换到频域 f = Fs/2*linspace(0,1,NFFT/2+1);%频域采样 figure(2);

plot(f,2*abs(Y(1:NFFT/2+1)));xlabel('frequency/Hz');ylabel('Amuplitude') ;%滤波前频谱 title('滤波前的频谱');

%滤波后频谱

%采用hanning窗滤波器

yy1=filter(b1,1,y);%调用滤波函数

YY1=fft(yy1,NFFT)/L;%进行傅里叶变换，下同。 f1=Fs/2*linspace(0,1,NFFT/2+1); figure(3);

subplot(2,2,1),plot(f1,2*abs(YY1(1:NFFT/2+1))) ;xlabel('frequency/Hz');ylabel('Amuplitude'); title('hanning窗的滤波效果'); %采用hammning窗滤波器 yy2=filter(b2,1,y);

YY2=fft(yy2,NFFT)/L;

f1=Fs/2*linspace(0,1,NFFT/2+1);

subplot(2,2,2),plot(f1,2*abs(YY2(1:NFFT/2+1))) ;xlabel('frequency/Hz');ylabel('Amuplitude'); title('hamming窗的滤波效果'); %采用blackman窗滤波器 yy3=filter(b3,1,y);

YY3=fft(yy3,NFFT)/L;

f1=Fs/2*linspace(0,1,NFFT/2+1);

subplot(2,2,3), plot(f1,2*abs(YY3(1:NFFT/2+1))) ;xlabel('frequency/Hz');ylabel('Amuplitude'); title('blackman窗的滤波效果'); %采用kaiser窗滤波器 yy4=filter(b4,1,y);

YY4=fft(yy4,NFFT)/L;

f1=Fs/2*linspace(0,1,NFFT/2+1);

subplot(2,2,4),plot(f1,2*abs(YY4(1:NFFT/2+1))) ;xlabel('frequency/Hz');ylabel('Amuplitude'); xlabel('frequency/Hz');ylabel('Amuplitude'); title('kaiser窗函数滤波效果');

%滤波前后的信号的时域对比 figure(4);

plot(y);xlabel('时间/s');ylabel('振幅');title('滤波前振幅特性'); figure(5);

subplot(2,2,1),plot(yy1);xlabel('时间/s');ylabel('振幅');title('hanning窗函数滤波振幅特性'); subplot(2,2,2),plot(yy2);xlabel('时间/s');ylabel('振幅');title('hamming窗函数滤波振幅特性'); subplot(2,2,3),plot(yy3);xlabel('时间/s');ylabel('振幅');title('blackman窗函数滤波振幅特性'); subplot(2,2,4),plot(yy4);xlabel('时间/s');ylabel('振幅');title('kaiser窗函数滤波振幅特性'); %滤波前后的信号的相位对比 figure(6);

plot(angle(Y));xlabel('时间/s');ylabel('相位');title('滤波前的相位特性'); figure(7);

subplot(2,2,1),plot(angle(YY1));xlabel('时间/s');ylabel('相位');title('hanning窗函数滤波相位特性');

subplot(2,2,2),plot(angle(YY2));xlabel('时间/s');ylabel('相位');title('hamming窗函数滤波相位特

性');

subplot(2,2,3),plot(angle(YY3));xlabel('时间/s');ylabel('相位');title('blackman窗函数滤波相位特性');

subplot(2,2,4),plot(angle(YY4));xlabel('时间/s');ylabel('相位');title('kaiser窗函数滤波相位特性');

3.分别设计hanning hamming blackman kaiser窗函数highpass_FIR function highpassfilter clc;

clear all;

Fs=100;%采样频率

fs=35;%高通阻带模拟截止频率 fp=40;%高通通带模拟起始频率 ws=2*pi*fs/Fs; wp=2*pi*fp/Fs; wn=(wp+ws)/2/pi; Bt=wp-ws;

N0=ceil(55*pi/Bt); N=N0+mod(N0+1,2); %调用窗函数

window1=hanning(N); window2=hamming(N); window3=blackman(N);

[n,Wn,beta,ftype]=kaiserord([35,40],[0 1],[0.01 0.01],100); window4=kaiser(n+1,beta); %设计加窗函数fir1

b1=fir1(N-1,wn,'high',window1); b2=fir1(N-1,wn,'high',window2); b3=fir1(N-1,wn,'high',window3);

b4=fir1(n,Wn,'high',window4 ,'noscale'); %求取频率响应

[H1,W1]=freqz(b1,1,512,2); [H2,W2]=freqz(b2,1,512,2); [H3,W3]=freqz(b3,1,512,2); [H4,W4]=freqz(b4,1,512,2);

figure(1);

subplot(2,2,1),plot(W1,20*log10(abs(H1)));%绘制频率响应图形 axis([0,1,-100,100]);

title('高通hanning窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

subplot(2,2,2),plot(W2,20*log10(abs(H2)));%绘制频率响应图形 axis([0,1,-100,100]);

title('高通hamming窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

subplot(2,2,3),plot(W3,20*log10(abs(H3)));%绘制频率响应图形 axis([0,1,-100,100]);

title('高通blackman窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

subplot(2,2,4),plot(W4,20*log10(abs(H4)));%绘制频率响应图形 axis([0,1,-100,100]);

title(' 高通kaiser窗的频率响应图形'); xlabel('频率（Hz)'); ylabel('幅值');

T=1/Fs;

L=100;%信号长度

t=(0:L-1)*T;%定义时间范围和步长

y=sin(2*pi*5*t)+5*sin(2*pi*15*t)+8*sin(2*pi*40*t);%滤波前的图形 NFFT = 2^nextpow2(L); % Next power of 2 from length of y Y = fft(y,NFFT)/L;%将时域信号变换到频域 f = Fs/2*linspace(0,1,NFFT/2+1);%频域采样 figure(2);

plot(f,2*abs(Y(1:NFFT/2+1)));xlabel('frequency/Hz');ylabel('Amuplitude') ;%滤波前频谱 title('滤波前的频谱');

%滤波后频谱

%采用hanning窗滤波器

yy1=filter(b1,1,y);%调用滤波函数

YY1=fft(yy1,NFFT)/L;%进行傅里叶变换，下同。 f1=Fs/2*linspace(0,1,NFFT/2+1); figure(3);

subplot(2,2,1),plot(f1,2*abs(YY1(1:NFFT/2+1))) ;xlabel('frequency/Hz');ylabel('Amuplitude'); title('hanning窗的滤波效果'); %采用hammning窗滤波器 yy2=filter(b2,1,y);

YY2=fft(yy2,NFFT)/L;

f1=Fs/2*linspace(0,1,NFFT/2+1);

subplot(2,2,2),plot(f1,2*abs(YY2(1:NFFT/2+1))) ;xlabel('frequency/Hz');ylabel('Amuplitude'); title('hamming窗的滤波效果');