Thiết kế bộ lọc FIR phương pháp tần số lấy mẫu (MATLAB CODE - design FIR filter)

Bởi
Thiết kế bộ lọc thông cao theo phương pháp tần số lấy mẫu với các thông số sau :
            ws= 0.6 pi,     As= 50 dB
            wp= 0.8 pi,     Rp= 1 dB 
M = 33, T1 = 0.1095; T2 = 0.598.
% Freq. Samp. Tech.: Highpass, Optimum method T1
% ws=0.6pi, wp=0.8pi, Rp=1dB, As=50dB
% M=33, T1 = 0.1095; T2 = 0.598;
M = 33; alpha = (M-1)/2; l = 0:M-1; wl = (2*pi/M)*l;
T1 = 0.1095; T2 = 0.598;
Hrs = [zeros(1,11),T1,T2,ones(1,8),T2,T1,zeros(1,10)];
Hdr = [0,0,1,1]; wdl = [0,0.6,0.8,1];
k1 = 0:floor((M-1)/2); k2 = floor((M-1)/2)+1:M-1;
angH = [-alpha*(2*pi)/M*k1, alpha*(2*pi)/M*(M-k2)];
H = Hrs.*exp(j*angH);
h = real(ifft(H,M));
[db,mag,pha,grd,w] = freqz_m(h,1);
[Hr,ww,a,L] = Hr_Type1(h);
subplot(1,1,1)
subplot(2,2,1);plot(wl(1:17)/pi,Hrs(1:17),'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Highpass: M=33,T1=0.1095,T2=0.598')
xlabel('frequency in pi units'); ylabel('Hr(k)')
set(gca,'XTickMode','manual','XTick',[0;.6;.8;1])
set(gca,'XTickLabelMode','manual','XTickLabels',[' 0';'.6';'.8';' 1'])
set(gca,'YTickMode','manual','YTick',[0,0.109,0.59,1]); grid
subplot(2,2,2); stem(l,h); axis([-1,M,-0.4,0.4])
title('Impulse Response'); xlabel('n'); ylabel('h(n)');
subplot(2,2,3); plot(ww/pi,Hr,wl(1:17)/pi,Hrs(1:17),'o');
axis([0,1,-0.1,1.1]); title('Amplitude Response')
xlabel('frequency in pi units'); ylabel('Hr(w)')
set(gca,'XTickMode','manual','XTick',[0;.6;.8;1])
set(gca,'XTickLabelMode','manual','XTickLabels',[' 0';'.6';'.8';' 1'])
set(gca,'YTickMode','manual','YTick',[0,0.109,0.59,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-100,10]); grid
title('Magnitude Response'); xlabel('frequency in pi units');
ylabel('Decibels');
set(gca,'XTickMode','manual','XTick',[0;.6;.8;1])
set(gca,'XTickLabelMode','manual','XTickLabels',[' 0';'.6';'.8';' 1'])
set(gca,'YTickMode','Manual','YTick',[-50;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['50';' 0'])
 Design the bandpass filter for the following specifications using the frequency sampling technique.
        lower stopband edge:   w1s = 0.2 pi,   As = 60 dB
        lower passband edge:   w1p = 0.35 pi,  Rp = 1 dB
        upper passband edge:  w2p = 0.65 pi,  Rp = 1 dB
        upper stopband edge:   w2s = 0.8 pi,   As = 60 dB
% ws1=0.2pi, wp1=0.35pi, wp2=0.65pi, ws2=0.8pi, Rp=1dB, As=60dB
% T2 = 0.59417456, T1=0.109021
M = 40; alpha = (M-1)/2; l = 0:M-1; wl = (2*pi/M)*l;
T1 = 0.109021; T2 = 0.59417456;
Hrs = [zeros(1,5),T1,T2,ones(1,7),T2,T1,zeros(1,9),T1,T2,ones(1,7),T2,T1,zeros(1,4)];
Hdr = [0,0,1,1,0,0]; wdl = [0,0.2,0.35,0.65,0.8,1];
k1 = 0:floor((M-1)/2); k2 = floor((M-1)/2)+1:M-1;
angH = [-alpha*(2*pi)/M*k1, alpha*(2*pi)/M*(M-k2)];
H = Hrs.*exp(j*angH);
h = real(ifft(H,M));
[db,mag,pha,grd,w] = freqz_m(h,1);
[Hr,ww,a,L] = Hr_Type2(h);
subplot(1,1,1)
subplot(2,2,1);plot(wl(1:21)/pi,Hrs(1:21),'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Bandpass: M=40,T1=0.5941, T2=0.109')
xlabel('frequency in pi units'); ylabel('Hr(k)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1])
set(gca,'YTickMode','manual','YTick',[0,0.59,0.109,1]); grid
subplot(2,2,2); stem(l,h); axis([-1,M,-0.4,0.4])
title('Impulse Response'); xlabel('n'); ylabel('h(n)');
subplot(2,2,3); plot(ww/pi,Hr,wl(1:21)/pi,Hrs(1:21),'o');
axis([0,1,-0.1,1.1]); title('Amplitude Response')
xlabel('frequency in pi units'); ylabel('Hr(w)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1])
set(gca,'YTickMode','manual','YTick',[0,0.59,0.109,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-100,10]); grid
title('Magnitude Response'); xlabel('frequency in pi units');
ylabel('Decibels');
set(gca,'XTickMode','Manual','XTick',[0,0.2,0.35,0.65,0.8,1]);
set(gca,'YTickMode','Manual','YTick',[-60;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['60';' 0'])
 Thiết kế bộ lọc thông thấp theo phương pháp tần số lấy mẫu.
    wp= 0.2 pi,     Rp= 0.25 dB
    ws= 0.3 pi,     As= 50 dB 
%  (a) naive approach :  20 samples to have value 1 at 0.2*pi and value 0 at 0.3*pi
%
M = 20; alpha = (M-1)/2; l = 0:M-1; wl = (2*pi/M)*l;
Hrs = [1,1,1,zeros(1,15),1,1];
Hdr = [1,1,0,0]; wdl = [0,0.25,0.25,1];
k1 = 0:floor((M-1)/2); k2 = floor((M-1)/2)+1:M-1;
angH = [-alpha*(2*pi)/M*k1, alpha*(2*pi)/M*(M-k2)];
H = Hrs.*exp(j*angH);
h = real(ifft(H,M));
[db,mag,pha,grd,w] = freqz_m(h,1);
[Hr,ww,a,L] = Hr_Type2(h);
subplot(1,1,1)
subplot(2,2,1);plot(wl(1:11)/pi,Hrs(1:11),'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Frequency Samples: M=20')
xlabel('frequency in pi units'); ylabel('Hr(k)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,1]); grid
subplot(2,2,2); stem(l,h); axis([-1,M,-0.1,0.3])
title('Impulse Response'); xlabel('n'); ylabel('h(n)');
subplot(2,2,3); plot(ww/pi,Hr,wl(1:11)/pi,Hrs(1:11),'o');
axis([0,1,-0.2,1.2]); title('Amplitude Response')
xlabel('frequency in pi units'); ylabel('Hr(w)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-60,10]); grid
title('Magnitude Response'); xlabel('frequency in pi units');
ylabel('Decibels');
set(gca,'XTickMode','Manual','XTick',[0;0.2;0.3;1]);
set(gca,'YTickMode','Manual','YTick',[-16;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['16';' 0']) 
 % (b) Increase No. of samples to 40 and transition sample T1 = 0.5
%
M = 40; alpha = (M-1)/2; l = 0:M-1; wl = (2*pi/M)*l;
Hrs = [ones(1,5),0.5,zeros(1,29),0.5,ones(1,4)];
Hdr = [1,1,0,0]; wdl = [0,0.25,0.25,1];
k1 = 0:floor((M-1)/2); k2 = floor((M-1)/2)+1:M-1;
angH = [-alpha*(2*pi)/M*k1, alpha*(2*pi)/M*(M-k2)];
H = Hrs.*exp(j*angH);
h = real(ifft(H,M));
[db,mag,pha,grd,w] = freqz_m(h,1);
[Hr,ww,a,L] = Hr_Type2(h);
subplot(1,1,1)
subplot(2,2,1);plot(wl(1:21)/pi,Hrs(1:21),'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Frequency Samples: M=40,T1=0.5')
xlabel('frequency in pi units'); ylabel('Hr(k)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,0.5,1]); grid
subplot(2,2,2); stem(l,h); axis([-1,M,-0.1,0.3])
title('Impulse Response'); xlabel('n'); ylabel('h(n)');
subplot(2,2,3); plot(ww/pi,Hr,wl(1:21)/pi,Hrs(1:21),'o');
axis([0,1,-0.1,1.1]); title('Amplitude Response')
xlabel('frequency in pi units'); ylabel('Hr(w)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,0.5,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-100,10]); grid
title('Magnitude Response'); xlabel('frequency in pi units');
ylabel('Decibels');
set(gca,'XTickMode','Manual','XTick',[0;0.2;0.3;1]);
set(gca,'YTickMode','Manual','YTick',[-30;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['30';' 0'])
% (c) Using optimum value for the transion sample as T1 = 0.3904(from linear programing or table)
%
M = 40; alpha = (M-1)/2; l = 0:M-1; wl = (2*pi/M)*l;
Hrs = [ones(1,5),0.3904,zeros(1,29),0.3904,ones(1,4)];
Hdr = [1,1,0,0]; wdl = [0,0.25,0.25,1];
k1 = 0:floor((M-1)/2); k2 = floor((M-1)/2)+1:M-1;
angH = [-alpha*(2*pi)/M*k1, alpha*(2*pi)/M*(M-k2)];
H = Hrs.*exp(j*angH);
h = real(ifft(H,M));
[db,mag,pha,grd,w] = freqz_m(h,1);
[Hr,ww,a,L] = Hr_Type2(h);
subplot(1,1,1)
subplot(2,2,1);plot(wl(1:21)/pi,Hrs(1:21),'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Frequency Samples: M=40,T1=0.39')
xlabel('frequency in pi units'); ylabel('Hr(k)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,0.39,1]); grid
subplot(2,2,2); stem(l,h); axis([-1,M,-0.1,0.3])
title('Impulse Response'); xlabel('n'); ylabel('h(n)');
subplot(2,2,3); plot(ww/pi,Hr,wl(1:21)/pi,Hrs(1:21),'o');
axis([0,1,-0.1,1.1]); title('Amplitude Response')
xlabel('frequency in pi units'); ylabel('Hr(w)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,0.39,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-100,10]); grid
title('Magnitude Response'); xlabel('frequency in pi units');
ylabel('Decibels');
set(gca,'XTickMode','Manual','XTick',[0;0.2;0.3;1]);
set(gca,'YTickMode','Manual','YTick',[-43;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['43';' 0']) 
 % Freq. Samp. Tech. Optimum method T1 & T2  % wp = 0.2pi, ws=0.3pi, Rp=0.25dB, As=50dB
% T1 = 0.5925, T2=0.1099
M = 60; alpha = (M-1)/2; l = 0:M-1; wl = (2*pi/M)*l;
Hrs = [ones(1,7),0.5925,0.11,zeros(1,43),0.11,0.5925,ones(1,6)];
Hdr = [1,1,0,0]; wdl = [0,0.2,0.3,1];
k1 = 0:floor((M-1)/2); k2 = floor((M-1)/2)+1:M-1;
angH = [-alpha*(2*pi)/M*k1, alpha*(2*pi)/M*(M-k2)];
H = Hrs.*exp(j*angH);
h = real(ifft(H,M));
[db,mag,pha,grd,w] = freqz_m(h,1);
[Hr,ww,a,L] = Hr_Type2(h);
subplot(1,1,1)
subplot(2,2,1);plot(wl(1:31)/pi,Hrs(1:31),'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Lowpass: M=60,T1=0.59, T2=0.109')
xlabel('frequency in pi units'); ylabel('Hr(k)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,0.59,0.109,1]); grid
subplot(2,2,2); stem(l,h); axis([-1,M,-0.1,0.3])
title('Impulse Response'); xlabel('n'); ylabel('h(n)');
subplot(2,2,3); plot(ww/pi,Hr,wl(1:31)/pi,Hrs(1:31),'o');
axis([0,1,-0.1,1.1]); title('Amplitude Response')
xlabel('frequency in pi units'); ylabel('Hr(w)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])
set(gca,'YTickMode','manual','YTick',[0,0.59,0.109,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-100,10]); grid
title('Magnitude Response'); xlabel('frequency in pi units');
ylabel('Decibels');
set(gca,'XTickMode','Manual','XTick',[0;0.2;0.3;1]);
set(gca,'YTickMode','Manual','YTick',[-63;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['63';' 0'])

Nhận xét

Đăng nhận xét