fsamp=10000;
Tsamp=1/fsamp;
fp1=1250;   %band pass edges
fp2=3255;
fs1=1750;   %band stop edges
fs2=3000;

%Digital Frequencies
wp1=2*pi*fp1*Tsamp;
wp2=2*pi*fp2*Tsamp;
ws1=2*pi*fs1*Tsamp;
ws2=2*pi*fs2*Tsamp;

%Attenuation stop and pass

as=46;          %20log(0.005)=-46
ap=0.91515;     %20log(1-.1)=-0.91515
T=1;
%Prewarping
Wp1=(2/T)*tan(wp1/2);
Wp2=(2/T)*tan(wp2/2);
Wp=[Wp1,Wp2];
Ws1=(2/T)*tan(ws1/2);
Ws2=(2/T)*tan(ws2/2);
Ws=[Ws1,Ws2];

[N,Wn]=ellipord(Wp,Ws,ap,as,'s');
[B,A]=ellip(N,ap,as,Wn,'Stop','s');
[b,a]=bilinear(B,A,1/T);

w=[0:pi/1024:pi];
H=freqz(b,a,w);
magH=abs(H);
phaseH=angle(H);
figure(1);
subplot(2,1,1);
plot(w,magH);
xlabel('Omega');
ylabel('Amplitude');
title('Magnitude Response -PartD');
grid on;


subplot(2,1,2);
plot(w,phaseH);
xlabel('Omega');
ylabel('Angle in Radians');
title('Phase Resonse');
grid on;

figure(2);
grd=grpdelay(b,a,w);
plot(w,grd);
xlabel('Omeag');
ylabel('Samples');
title('Groud dealy of Filter');
grid on

