%====================================================
%====================================================
%====================================================
% Ej de filtro
%====================================================
%====================================================
%====================================================

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ESPECIFICAR LOS DATOS DEL ANALISIS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Vector de frecuencias
vf=linspace(0,1.5,1000);

% Numero de frecuencias de analisis
Nfrec=numel(vf);

% Impedancias de refrencia para los parametros S
Z01=1;
Z02=1;

% Datos filtro
%vg=[0.7654    1.8478    1.8478    0.7654];
vg=[0.9732    1.3723    1.8032    1.3723    0.9732];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CALCULAR LA RED PARA CADA FRECUENCIA
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

for ifrec=1:Nfrec

    w=2*pi*vf(ifrec);

    Yp=j*w*vg(1);
    T1=[1 0; Yp 1];

    Zs=j*w*vg(2);
    T2=[1 Zs; 0 1];

    Yp=j*w*vg(3);
    T3=[1 0; Yp 1];

    Zs=j*w*vg(4);
    T4=[1 Zs; 0 1];

    Yp=j*w*vg(5);
    T5=[1 0; Yp 1];

    T=T1*T2*T3*T4*T5;

    A=T(1,1); B=T(1,2); C=T(2,1); D=T(2,2);
    denom1=A+B/Z02+C*Z01+D*Z01/Z02;
    S11(ifrec)=(A+B/Z02-C*Z01-D*Z01/Z02)/denom1;
    S21(ifrec)=2*sqrt(Z01/Z02)/denom1;
    denom2=D+B/Z01+C*Z02+A*Z02/Z01;
    S22(ifrec)=(D+B/Z01-C*Z02-A*Z02/Z01)/denom2;
    S12(ifrec)=2*sqrt(Z02/Z01)*(A*D-B*C)/denom2;

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PINTAR RESULTADOS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Pintar resultados en unidades naturales
hfig_nat=figure;
hold on
subplot(2,1,1)
plot(vf,abs(S11),'c',vf,abs(S12),'r',vf,abs(S21),'m',vf,abs(S22),'b','Linewidth',1.5);
xlabel('Frecuencia (Hz)'); ylabel('|S_{ij}|'); legend('S_{11}','S_{12}','S_{21}','S_{22}')
subplot(2,1,2)
plot(vf,angle(S11)*180/pi,'c',vf,angle(S12)*180/pi,'r',vf,angle(S21)*180/pi,'m',vf,angle(S22)*180/pi,'b','Linewidth',1.5);
xlabel('Frecuencia (Hz)'); ylabel('<S_{ij} (º)'); legend('S_{11}','S_{12}','S_{21}','S_{22}')

% Pintar resultados en dB´s y retardo de grupo
hfig_nor=figure;
hold on
subplot(2,1,1)
plot(vf,20*log10(abs(S11)),'c',vf,20*log10(abs(S12)),'r',vf,20*log10(abs(S21)),'m',vf,20*log10(abs(S22)),'b','Linewidth',1.5);
xlabel('Frecuencia (Hz)'); ylabel('|S_{ij}| (dB)'); legend('S_{11}','S_{12}','S_{21}','S_{22}')
subplot(2,1,2)
retardo_grupo_12=-diff(unwrap(angle(S12)))./(diff(vf));
retardo_grupo_21=-diff(unwrap(angle(S21)))./(diff(vf));
plot(vf(1:end-1),retardo_grupo_12,'r',vf(1:end-1),retardo_grupo_21,'m','Linewidth',1.5);
xlabel('Frecuencia (Hz)'); ylabel('Retado de grupo (seg)'); legend('\tau_{12}','\tau_{21}')
