% linearity.m - generate plots to illustrate linearity
% RMM, 12 Sep 06

aminit;

% Use the coupled spring massi system, with normalized steady state gain
sys = springmass();
sys = sys / freqresp(sys, 0);
x0 = [2, -1];

% Create input vectors as a function of time
tvec=0:0.1:60;
u1 = 0 * tvec;
u2 = [1-tripuls(0:0.1:29.9, 60, 0) tripuls(0:0.1:30, 60, 0)];

[y1, t1, x1] = lsim(sys, u1, tvec, x0);
[y2, t2, x2] = lsim(sys, u2, tvec);
[y3, t3, x3] = lsim(sys, u1+u2, tvec, x0);

figure(1); clf; 

subplot(431); pl = plot(t1, u1); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  ylabel('Homogeneous');
  title('Input {\itu}');
subplot(432); pl = plot(t1, x1); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  title('State (x1, x2)');
  hold on; plot([0 0], [-1 2], 'k.');
subplot(433); pl = plot(t1, y1); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  title('Output {\ity}');

subplot(434); pl = plot(t2, u2); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  ylabel('Particular');
subplot(435); pl = plot(t2, x2); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  hold on; plot([0], [0], 'k.');
subplot(436); pl = plot(t2, y2); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);

subplot(437); pl = plot(t3, u1+u2); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  ylabel('Complete');
  xlabel('Time {\itt} [sec]', 'FontName', 'times');
subplot(438); pl = plot(t3, x3); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  hold on; plot([0 0], [-1 2], 'k.');
  xlabel('Time {\itt} [sec]');
subplot(439); pl = plot(t3, y3); amaxis([0 60 -2 2]);
  set(pl, 'LineWidth', AM_data_linewidth);
  xlabel('Time {\itt} [sec]');

amprint('linearity.eps');
