% Problem 3.26

disp('Problem 3.26')
disp(' ')

% System model
A=[-16.11 -0.39 27.2 0 0;
   0.01 -16.99 0 0 12.47;
  15.11 0 -53.6 -16.57 71.78;
  -15.36 0 0 -107.2 232.11;
  2.27 60.1 0 2.273 -102.99];

B=[11.12 -2.61 -21.91 -53.5 69.1;
   -12.6 3.36 0 0 0]';

C=[0 0 0 0 1;
   0 1 1 0 0];

D=[0 0;
   0 0];

% Transfer function
[num1,den]=ss2tf(A,B,C,D,1);

disp(' ')
disp('3.26(a) - System transfer function')
disp(' ')
disp('g11(s)=y1(s)/u1(s):')
printsys(num1(1,:),den,'s')
disp(' ')
disp('g21(s)=y2(s)/u1(s):')
printsys(num1(2,:),den,'s')

[num2,den]=ss2tf(A,B,C,D,2);
disp(' ')
disp('g12(s)=y1(s)/u2(s):')
printsys(num2(1,:),den,'s')
disp(' ')
disp('g22(s)=y2(s)/u2(s):')
printsys(num2(2,:),den,'s')

% 3.26(b) - Impulse response
t=0:0.1:5;
[y,x]=impulse(A,B,C,D,1,t);
subplot(211)
plot(t,y(:,1),'-',t,y(:,2),'--')
[xc,yc]=ginput(1)
text(xc,yc,'u2=0')
[xc,yc]=ginput(1)
text(xc,yc,'y1')
[xc,yc]=ginput(1)
text(xc,yc,'y2')
xlabel('t [secs]')
ylabel('outputs')
grid
[y,x]=impulse(A,B,C,D,2,t);
subplot(212)
plot(t,y(:,1),'-',t,y(:,2),'--')
[xc,yc]=ginput(1)
text(xc,yc,'u1=0')
[xc,yc]=ginput(1)
text(xc,yc,'y1')
[xc,yc]=ginput(1)
text(xc,yc,'y2')
xlabel('t [secs]')
ylabel('outputs')
grid

% 3.26(c)
t=0:0.1:20;
u=[exp(-t')+sin(t'),0*t'];
[y,x]=lsim(A,B,C,D,u,t);
figure(2)
subplot(111)
plot(t,y)
xlabel('t [secs]')
ylabel('outputs')
[xc,yc]=ginput(1)
text(xc,yc,'y1')
[xc,yc]=ginput(1)
text(xc,yc,'y2')
grid