% Problem 9.19
clear
% (a) System model for T=0
num=[0 1 2];
den=conv([1 1],[1 10]);

diary sol9_19.dat
w=logspace(-1,2,100);
[mag,phase]=bode(num,den,w);
figure(1)
margin(mag,phase,w)
subplot(211)
grid
subplot(212)
grid
print -deps sol9_19a
pause

% (b) T=0.1
T=0.1;
w1=w*180/pi;
dph=-T*w1';
phase01=phase+dph;
figure(2)
[Gm,Pm,Wcp,Wcg]=margin(mag,phase01,w);
fprintf('\n\nStability margins for system with desired ess:\n\n')
fprintf('Gain Margin              = %f dB\n',20*log10(Gm))
fprintf('Phase Margin             = %f degrees\n',Pm)
fprintf('Gain Crossover Frequency = %f rad/sec\n',Wcg)
fprintf('Phase Crossover Frequency= %f rad/sec\n',Wcp)
%margin(mag,phase01,w)

magdB=20*log10(mag);
subplot(211)
semilogx(w,magdB)
grid
xlabel('Frequency (rad/sec)')
ylabel('Gain dB')
subplot(212)
semilogx(w,phase01)
grid
xlabel('Frequency (rad/sec)')
ylabel('Phase deg')
print -deps sol9_19b
pause

% (c) T=0.2
T=0.2;
dph=-T*w1';
phase02=phase+dph;
figure(3)
[Gm,Pm,Wcp,Wcg]=margin(mag,phase02,w);
fprintf('\n\nStability margins for system with desired ess:\n\n')
fprintf('Gain Margin              = %f dB\n',20*log10(Gm))
fprintf('Phase Margin             = %f degrees\n',Pm)
fprintf('Gain Crossover Frequency = %f rad/sec\n',Wcg)
fprintf('Phase Crossover Frequency= %f rad/sec\n',Wcp)
%margin(mag,phase02,w)
magdB=20*log10(mag);
subplot(211)
semilogx(w,magdB)
grid
xlabel('Frequency (rad/sec)')
ylabel('Gain dB')
subplot(212)
semilogx(w,phase02)
grid
xlabel('Frequency (rad/sec)')
ylabel('Phase deg')
print -deps sol9_19c
pause

% T=0.5
T=0.5;
dph=-T*w1';
phase05=phase+dph;
figure(4)
[Gm,Pm,Wcp,Wcg]=margin(mag,phase05,w);
fprintf('\n\nStability margins for system with desired ess:\n\n')
fprintf('Gain Margin              = %f dB\n',20*log10(Gm))
fprintf('Phase Margin             = %f degrees\n',Pm)
fprintf('Gain Crossover Frequency = %f rad/sec\n',Wcg)
fprintf('Phase Crossover Frequency= %f rad/sec\n',Wcp)
magdB=20*log10(mag);
subplot(211)
semilogx(w,magdB)
grid
xlabel('Frequency (rad/sec)')
ylabel('Gain dB')
subplot(212)
semilogx(w,phase05)
grid
xlabel('Frequency (rad/sec)')
ylabel('Phase deg')
print -deps sol9_19d
diary off