% response.m - MATLAB script to generate step and frequency response
%
% This section uses a single linear system to generate the step and
% frequency response plots used in the chapter.  The system used is a
% damped, second order oscillator with a first order lag.
%
% Note that some MATLAB commands can take the A, B, C, D matrices directly,
% both others require us to pass a "system" as the argument.  You can create
% a system by using the "ss" function, as shown below.
%

clear fun;				% remove any extraneous functions
aminit;					% initialize MATLAB

% Generate the system matrices
A = [-0.2 2 0; -0.5 -0.2 4; 0 0 -10];
B = [0; 0; 0.5]; C = [0.25 0 0]; D = [0];

% Set up as a system object and normalize the response
sys = ss(A,B,C,D);
sys = sys / freqresp(sys, 0);

%
% Step response
%
[y,time] = step(sys, 30);
figure(1); clf; subplot(211)
plot(time, y, 'Linewidth', AM_data_linewidth);
amaxis([0 30 0 2]);
xlabel('Time [sec]'); ylabel('Output');

% Add labels to the plot showing the features
% 95% and 105% lines
yfn = y(length(y)); tfn = max(time);
line([0 tfn], [0.95*yfn 0.95*yfn], 'Color', 'black', ...
  'LineWidth', AM_ref_linewidth);
line([0 tfn], [1.05*yfn 1.05*yfn], 'Color', 'black', ...
  'LineWidth', AM_ref_linewidth);

% Rise time
line([0.55 0.55], [0.05 1.15*yfn], 'Color', 'black', ...
  'LineWidth', AM_ref_linewidth);
line([1.7 1.7], [0.05 1.15*yfn], 'Color', 'black', ...
  'LineWidth', AM_ref_linewidth);
arrow([-0.3, 1.8*yfn/3], [0.55, 1.8*yfn/3], AM_small_arrowsize);
arrow([1.7+0.55+0.3, 1.8*yfn/3], [1.7, 1.8*yfn/3], AM_small_arrowsize);
amtext(2.8, 1.8*yfn/3, 'Rise time Tr');

% Settling time
line([13.6 13.6], [0.05 1.15*yfn], 'Color', 'black', ...
  'LineWidth', AM_ref_linewidth);
arrow([0 yfn/3], [13.6 yfn/3], AM_small_arrowsize);
text(14, yfn/3, 'Settling time Ts');

% Overshoot
[ypk, pki] = max(y); tpk = time(pki);
line([tpk-1 tpk+1], [ypk+0.001 ypk+0.001], 'Color', 'black', ...
  'LineWidth', AM_ref_linewidth);
line([0 tfn], [yfn yfn], 'Color', 'black', 'LineStyle', '--', ...
  'LineWidth', AM_ref_linewidth);
arrow([tpk+1 1], [tpk+1, ypk], AM_small_arrowsize);	% up
arrow([tpk+1 ypk], [tpk+1, 1], AM_small_arrowsize);	% down
amtext(tpk+1.5, (ypk+yfn)/2 + 0.04, 'Overshoot Mp');

% Steady state
arrow([25 0.6], [29.5, 0.99], AM_small_arrowsize);
amtext(20, yfn-0.5, 'Steady state value yss');

% Redraw the plot so that it's on top
hold on; plot(time, y);

amprint('response-step.eps');

%
% Frequency response
%
tvec=0:0.1:60;
u1 = sin(0.75*tvec);
[y1, x1] = lsim(sys*0.75, u1, tvec);

% Shift to steady state
ssi = ceil(length(tvec)/2):length(tvec);
tvec = tvec - tvec(ssi(1));

figure(2); clf; subplot(211);
plot(tvec(ssi), u1(ssi), '--', tvec(ssi), y1(ssi), '-', ...
  'Linewidth', AM_data_linewidth);
amaxis([0 20 -2 2]);
set(gca, 'DataAspectRatio', [5 1.33 1]);
xlabel('Time [sec]'); ylabel('Input, output');
lgh = legend('Input', 'Output', 'Orientation', 'horizontal');
legend(lgh, 'boxoff');

% Add annotation
line([0 30], [0 0], 'Color', 'black', 'LineWidth', AM_ref_linewidth);

% Input amplitude
amtext(3.7, 0.5, 'Au', 'BackgroundColor', 'white');
arrow([5.55 0.1], [5.55, 0.95], AM_small_arrowsize);
arrow([5.55 0.95], [5.55, 0], AM_small_arrowsize);

% Output amplitude
amtext(6.6, 1, 'Ay', 'BackgroundColor', 'white', 'Margin', 1);
arrow([6.4 0], [6.4, 1.35], AM_small_arrowsize);
arrow([6.4 1.35], [6.4, 0], AM_small_arrowsize);

% Phase
line([11.8 11.8], [-1.3 0.2]);
line([12.72 12.72], [-0.8 0.2]);
arrow([10.8 -0.3], [11.8 -0.3], AM_small_arrowsize);
arrow([13.72 -0.3], [12.72 -0.3], AM_small_arrowsize);
amtext(14, -0.3, 'DT');

% Period
line([11.8-8*pi/3 11.8-8*pi/3], [-1.3 0.2]);
arrow([11.8-8*pi/3 -1], [11.8 -1], AM_small_arrowsize);
arrow([11.8 -1], [11.8-8*pi/3 -1], AM_small_arrowsize);
amtext(7.5, -1.3, 'T');

amprint('response-sinusoids.eps');

%
% Bode plot
% 
% In addition to the normal step response, we also generate a Bode plot
% so that we can illustrate some of the properties.
%

figure (3); clf; subplot(221);

% Get the frequency response
[magh, phaseh] = ambode(0.75*sys, {0.05, 10});

% Add some annotation
subplot(magh); hold on;
loglog([0.05 10], [1 1], 'k-', 'Linewidth', AM_ref_linewidth);

subplot(phaseh); hold on;
semilogx([0.05 10], [-180 -180], 'k-', 'Linewidth', AM_ref_linewidth);

amprint('response-magfreq.eps');

% Local Variables:
% mode: Matlab
% End:
