%% Ejemplo filtro pasabanda de chebyshev-II analógico
clear all;

% Especif.
fp1=45e3;fp2=90e3;
Ap=0.5;
fs1=20e3;fs2=150e3;
As=50;

%% Convertir filtro a pasobajo

wp1 = 2*pi*fp1;
wp2 = 2*pi*fp2;
ws1 = 2*pi*fs1;
ws2 = 2*pi*fs2;

%% Fijamos frecuencias de parabanda
wx2 = ws1*ws2;

if (wp1*wp2 > wx2) 
  wp1 = wx2/wp2;
else
  wp2 = wx2/wp1;
end

wlpp = wp2-wp1;
wlps = ws2-ws1;


%% Prototipo Pasobajo

% Calculo de e2
e2 = 1/(10^(0.1*As)-1);

% Orden del filtro
n = acosh(sqrt(1/(10^(0.1*Ap)-1)/e2))/acosh(wlps/wlpp);
n = ceil(n);

% Frecuencia de corte
wc = wlps;

theta = (1/n)*asinh(1/sqrt(e2));
k=1:n;
alpha = (2*k-1)*pi/(2*n);


% Determinacion de los polos y ceros en el prototipo
pk = 1./(-sinh(theta)*sin(alpha) + j*cosh(theta)*cos(alpha));
zk = -j./cos(alpha);


[B,A]=zp2tf(zk',pk',1);

%% Determinacion de la ganancia K

K=A(length(A))/B(length(B));



[B1,A1] = lp2bp(K*B,A,sqrt(wx2),wc); 

N=500;
wf=2*pi*160e3;

W=(5:N-1)*wf/N;
H1=freqs(B1,A1,W);

subplot(3,1,1); plot(W/(2*pi),20*log10(abs(H1)));grid;
axis([0 wf/(2*pi) -(As+20) inf]);
ylabel(['|H(w)|']);
xlabel(['Frecuencia (Hz)']);
subplot(3,1,2); plot(W/(2*pi),unwrap(angle(H1)));grid;
axis([0 wf/(2*pi) -inf inf]);
ylabel(['Fase H(w) (rad)']);
xlabel(['Frecuencia (Hz)']);

% Retraso de grupo
dW=W(2)-W(1);
phase=unwrap(angle(H1));
dphase = -diff(phase)/dW;
subplot(3,1,3);plot(W(1:length(W)-1)/(2*pi),dphase);grid;
axis([0 wf/(2*pi) 0 inf]);
ylabel(['Retraso de Grupo (s)']);
xlabel(['Frecuencia (Hz)']);
print -depsc2 bp_cheb2_ej