% Plant transfer function values from assignment spec
a=3;
b=1;
c=4;

% Controller parameters from PID tuning block in Simulink
Kp=177.582632673963;
Kd=85.3149480379147;
N=3263.86804812559;
Ki=39.0499305396423;

% Numerator coefficients of transfer functions
N1=Kp+(N*Kd);
N2=(Ki+(N*Kp));
N3=N*Ki;

% Denominator coefficients of closed loop transfer function
D1=a;
D2=(a*N)+b;
D3=Kp+(N*Kd)+(b*N)+c;
D4=(N*Kp)+Ki+(c*N);
D5=N*Ki;

% Denominator coefficients of forward transfer function
d1=a;
d2=(a*N)+b;
d3=(b*N)+c;
d4=c*N;
d5=0;

% create forward transfer function
fwd_tf=tf([N1 N2 N3],[d1 d2 d3 d4 d5]);

% Create Nyquist plot of forward transfer function
% This HAS to exclude the feedback as it will be a measure of 
% System stability - how much gain / phase shift needed to make instable
figure;
subplot(2,2,1);
nyquist(fwd_tf);
title('Nyquist Plot of forward TF');

% Create Bode plot of foward transfer function 
% - provides similar info to Nyquist plot
subplot(2,2,2);
margin(fwd_tf);
%title('Bode plot of forward TF')

% Create closed loop transfer function of whole system (inc. feedback)
system_tf=tf([N1 N2 N3],[D1 D2 D3 D4 D5]);

% Plot Pole-Zero diagram of closed loop system
subplot(2,2,3);
pzplot(system_tf);
title('Pole-Zero plot of closed loop TF');

% Create transfer function to represent only the plant
plant_tf = tf(1,[a b c]);

% Plot pole-zero diagram of plant transfer function 
subplot(2,2,4);
pzplot(plant_tf);
title('Pole-Zero plot of plant');

% create transfer function of controller
Controller_tf = tf([(Kp+(N*Kd)) (N*Kp)+Ki N*Ki], [1 N]);

% create new figure
figure;

% with 1 row of 3 graphs, and with the first graph of these....
subplot(1,3,1);
% plot the pole-zero diagram of the plant on its own
pzplot(plant_tf);
title('Plant');

% in the next graph
subplot(1,3,2);
% plot the pole-zero diagram of the controller on its own
pzplot(Controller_tf);
title('Controller');

% and in the last graph
subplot(1,3,3);
% plot the pole-zero diagram of entire system with feedback
pzplot(system_tf);
title('System');