clc; clear
close all

% En dominio continuo
s=tf('s');
tc=0:0.001:6;
Gp=1/(s*(s+2));
Gplc=feedback(Gp,1);
Gc=9*(s+2)/(s+3);
Gla=Gc*Gp;
Glc=feedback(Gla,1);

step(Glc,tc);
legend('Sistema Compensado en Tiempo Continuo')

figure; margin(Gla)
figure; bode(Glc)

% Selección del periodo de muestreo
wn=3; xita=0.5; wd=wn*sqrt(1-xita^2);
Nd=5; % muestras x periodo de oscilación
Td=wd/2/pi;
T=Td/Nd;

% En el dominio de tiempo discreto
% T=0.1;
td=0:T:5-T;
ref1=ones(length(td),1);
z=tf('z',T);
Gpd=c2d(Gp,T);

Gpdsc=feedback(Gpd,1);
figure; step(Gpdsc)

% Gcd=c2d(Gc,T);
Gcda=(9*z-6*exp(-3*T)-3)/(z-exp(-3*T));
Glad=Gcda*Gpd;
Glcd=feedback(Glad,1);

S1=stepinfo(Glcd);
Mp1=S1.Overshoot/100;
xita1=-log(Mp1)/sqrt(log(Mp1)^2 + pi^2);

figure; step(Glc,Glcd,td); 
legend('Sistema Compensado en Tiempo Continuo','Sistema Compensado en Tiempo Discreto,T=0.0827')
figure; axis square; pzmap(Glcd)

[y1,td]=step(Glcd,td);
error_01=ref1 - y1; 
[yc,td]=step(Glc,td);
error_c=ref1 - yc; 
figure; lsim(Gcda,error_01,td); hold on
lsim(Gc,error_c,td);

figure; margin(Glad)
figure; bode(Glcd)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Nd=5; % muestras x periodo de oscilación
Td=wd/2/pi;
T=Td/Nd;

td=0:T:5-T;
ref2=ones(length(td),1);
z=tf('z',T);
Gpd=c2d(Gp,T);
Gcd=(9*z-6*exp(-3*T)-3)/(z-exp(-3*T));
Glad=Gcd*Gpd;
Glcd=feedback(Glad,1);

S2=stepinfo(Glcd);
Mp2=S2.Overshoot/100;
xita2=-log(Mp2)/sqrt(log(Mp2)^2 + pi^2);

figure; step(Glc,Glcd,td); 
legend('Sistema Compensado en Tiempo Continuo','Sistema Compensado en Tiempo Discreto,T=0.0827')
figure; axis square; pzmap(Glcd)

[y2,td]=step(Glcd,td);
error_02=ref2 - y2; 
figure; lsim(Gcd,error_02,td);

figure; margin(Glad)
figure; bode(Glcd)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Nd=10; % muestras x periodo de oscilación
Td=wd/2/pi;
T=Td/Nd;

td=0:T:5-T;
ref3=ones(length(td),1);
z=tf('z',T);
Gpd=c2d(Gp,T);
Gcd=(9*z - 6*exp(-3*T) -3)/(z-exp(-3*T));
Glad=Gcd*Gpd;
Glcd=feedback(Glad,1);

S3=stepinfo(Glcd);
Mp3=S3.Overshoot/100;
xita3=-log(Mp3)/sqrt(log(Mp3)^2 + pi^2);

figure; step(Glc,Glcd,td); 
legend('Sistema Compensado en Tiempo Continuo','Sistema Compensado en Tiempo Discreto,T=0.0413')
figure; axis square; pzmap(Glcd)

[y3,td]=step(Glcd,td);
error_03=ref3 - y3; 
figure; lsim(Gcd,error_03,td);

figure; margin(Glad)
figure; bode(Glcd)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% close all

T=0.001;
td=0:T:5-T;
ref4=ones(length(td),1);
z=tf('z',T);
Gpd=c2d(Gp,T);
Gcd=(9*z-6*exp(-3*T)-3)/(z-exp(-3*T));
Gcdr =3.3821*(z-0.9988)/(z-0.9986);

% Gcd=Gcdr
Glad=Gcd*Gpd;
Glcd=feedback(Glad,1);
figure; axis square; pzmap(Glcd)

S4=stepinfo(Glcd);
Mp4=S4.Overshoot/100;
xita4=-log(Mp4)/sqrt(log(Mp4)^2 + pi^2);

figure; step(Glc,Glcd,td); 
legend('Sistema Compensado en Tiempo Continuo','Sistema Compensado en Tiempo Discreto,T=0.001')
figure; axis square; pzmap(Glcd)

[y4,td]=step(Glcd,td);
error_04=ref4 - y4; 
[yc,td]=step(Glc,td);
error_c=ref4 - yc; 

figure; lsim(Gcd,error_04,td); hold on
lsim(Gc,error_c,td);

figure; margin(Glad)
figure; bode(Glcd)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%