clear;

% Especificaciones
ws = [1000 8000];
wp = [1500 5000];
Ap = 1;As = 60;

% Paso a especificaciones de pasobajo
% Tablas 1 y 2
Bw1 = wp(2)-wp(1);
% Fijamos las frecuencias de pasabanda
wx2 = wp(1)*wp(2);wx=sqrt(wx2);
if ws(1)*ws(2) < wx2,
	ws(2) = wx2/ws(1);
else
	ws(1) = wx2/ws(2);
end
Bw2 = ws(2)-ws(1);

% Especificiones de pasobajo
Wp = Bw1;
Ws = Bw2;

k1 = 10^(0.1*Ap)-1;
k2 = 10^(0.1*As)-1;
% Filtro de Butterworth
e2b = k1;
nb = ceil((log10(sqrt(k2/k1)))/log10(Ws/Wp));
v3 = (1/k1)^(1/(2*nb));

% Filtro de Chebyshev I
e2c1 = k1;
nc1 = ceil((acosh(sqrt(k2/k1)))/acosh(Ws/Wp));

% Filtro de Chebyshev II
e2c2 = 1/k2;
nc2 = ceil((acosh(sqrt(k2/k1)))/(acosh(Ws/Wp)));

% Filtro elíptico
[ne, Wn]=ellipord(Wp, Ws, Ap, As,'s'); 


% Diseńo del prototipo de filtro filtro pasobajo
% Filtro de Butterworth
[Z,P,K] = buttap(nb);
[Bbp,Abp] = zp2tf(Z,P,K);

% Filtro de Chebyshev I
[Z,P,K] = cheb1ap(nc1,Ap);
[Bc1p,Ac1p] = zp2tf(Z,P,K);

% Filtro de Chebyshev I
[Z,P,K] = cheb2ap(nc2,As);
[Bc2p,Ac2p] = zp2tf(Z,P,K);

%Filtro elíptico
[Z,P,K] = ellipap(ne,Ap,As);
[Be,Ae] = zp2tf(Z,P,K);


% Convertir a Pasabanda
% Filtro de Butterworth
[Bb,Ab] = lp2bp(Bbp,Abp,wx,Bw1*v3);

% Filtro de Chebyshev I
[Bc1,Ac1] = lp2bp(Bc1p,Ac1p,wx,Bw1);

% Filtro de Chebyshev II
[Bc2,Ac2] = lp2bp(Bc2p,Ac2p,wx,Bw2);

%Filtro elíptico
[Be,Ae] = lp2bp(Be,Ae,wx,Bw1);

fprintf(1,'Orden del filtro de Butterworth = %d\n',nb);
fprintf(1,'Orden del filtro de Chebyshev I = %d\n',nc1);
fprintf(1,'Orden del filtro de Chebyshev II = %d\n',nc2);
fprintf(1,'Orden del filtro elíptico = %d\n',ne);


% Respuestas frecuenciales de los tres filtros obtenidos
W = 1:10:10000;
[Hb]=freqs(Bb,Ab,W);
Hb = 20*log10(abs(Hb));
[Hc1]=freqs(Bc1,Ac1,W);
Hc1 = 20*log10(abs(Hc1));
[Hc2]=freqs(Bc2,Ac2,W);
Hc2 = 20*log10(abs(Hc2));
[He]=freqs(Be,Ae,W);
He = 20*log10(abs(He));



plot(W,Hb,'r',W,Hc1,'g',W,Hc2,'m',W,He,'b');hold;
plot([wp ws],-[Ap Ap As As],'*y');grid;zoom;