clearvars 
close all
clc;

%% exercise 8.1
Kp = 1;
Ti = 1.5;
Td = 1;
N = 10;
h = 0.007;
tau = 0.7;

Kd = Kp*Td;
Ki= Kp/Ti;

% Define the transfer function of the given process 
G = tf(1,[1 0.8 0.5], 'InputDelay', tau);
H = c2d(G, h, 'zoh');

% Define the continuous-time PID-controller
Gpid = tf([Kp*Td Kp/Ti Kp],[0 1 0]);

% Define the dicrete-time PID-controller (backward euler)
% H_pid_bw = tf([K*(1+Td/h + h/Ti) K*(-1-2*Td/h) (K*Td/h)],[1 -1 0],h)
H_pid_bw  = pid(Kp, Ki, Kd,'Ts', h,'IFormula','BackwardEuler','DFormula', 'BackwardEuler');

% Define the dicrete-time PID-controller (Tustin)
opt2 = c2dOptions('Method','tustin');
H_pid_tustin = c2d(Gpid,h,opt2)

%{ 
Using the approach in line 36 does not allow you to use 'Trapezoidal' for
'DFormula', due to controller becoming unstable. 
%}

% H_pid_tustin = pid(Kp, Ki, Kd,'Ts', h,'IFormula','Trapezoidal')

% Define the dicrete-time PID-controller (Practical)
H_pid_practical = pidstd(Kp, Ti, Td, N, h,'IFormula','ForwardEuler','DFormula', 'BackwardEuler')

% Closed-loop systems
Gcl = feedback(G*Gpid, 1);
Hcl_bw = feedback(H*H_pid_bw, 1);
Hcl_tustin = feedback(H*H_pid_tustin, 1);
Hcl_practical = feedback(H*H_pid_practical, 1);

figure(1); clf;
t=15;

step(t, Gcl, 'b--')
hold on;
step(t, Hcl_bw,'g-.')
step(t, Hcl_tustin,'r:')
step(t, Hcl_practical,'k-')
legend('Continuous-time','Backward Euler', 'Tustin','Practical','Location','SouthEast')

%% exercise 8.2

Ts = 1; % sampling time

% Transfer function of the process
G_process = tf([1],[10 1 0]);
H = c2d(G_process, Ts, 'zoh');

% Transfer function of the controller
Hc = tf([13 -0.88*13],[1 0.5], 1);

% closed-loop system
Hcl = feedback(H*Hc, 1);

%{
Returns the values of Mp and ts calculated in the exercise among others. 
More information can be found in MATLAB's documentation.
%} 

stepinfo(Hcl) 

% figure(2);clf;
% step(Hcl)
% hold on;
% step(Hc/10)