% Calculating margins of process, hp=k*exp(-tau*s)/s
% with SIMC PI controller tuning with Kp=1/(2*k*tau)
% and Ti=8*tau

k=1;
tau=1;

Kp=1/(2*k*tau); Ti=8*tau;

s=tf('s');

hp=k*exp(-tau*s)/s;
hc=Kp*(1+Ti*s)/(Ti*s);

disp('The loop transfer function, h_0(s)=hc*hp')
h0=hc*hp

% Solution using the margin function
[Gm,Pm,W180,Wc] = margin(h0)

% Solution using the definitions in Ch 9 for: w180, GM, wc and Pm

% phase=atan(Ti*w)-tau*w-pi
x=fsolve('-x+atan(8*x)',1);  % solve for 
w180=x

mag_h0=Kp*k*sqrt(1+(Ti*x)^2)/(Ti*x^2);  % mag_h0=| h_0(j*w180) |
gm=1/mag_h0                            % Gain Margin

x=fsolve('0.5*sqrt(1+(8*x)^2)/(8*x^2)-1',0.2); % Solve for w_c=x
wc=x
pm=atan(Ti*wc)-tau*wc;
pm=pm*180/pi

