% E04_2.m               05apr17  12apr21
%   Transfer functions and freq. response of ANALOG LTI systems
clear all; clf;

 zz = [];   pp = [-1];      q0=3;
% zz = [];   pp = [-1+j*1.3   -1-j*1.3];    q0=3;
% zz = [0];   pp = [-1+j*1.3  -1-j*1.3];    q0=3;  % note: [] is different from [0] !!
% zz = [-1];   pp = [-1+j*1.3  -1-j*1.3];    q0=3;
% zz = [0];   pp = [-1+j*.7  -1-j*.7  -1+j*3  -1-j*3];   q0=5;
% zz = [j*1.8  -j*1.8];   pp = [-1+j*1.8   -1-j*1.8];     q0=3;

qq = [-q0:.077:q0]; % avoid "divide by 0" warning
[sig,ome] = meshgrid(qq,qq);
s = sig + j*ome;
NUM = 1;    for k=1:max(size(zz));    NUM = NUM .* (s - zz(k));    end;
DEN = 1;    for k=1:max(size(pp));    DEN = DEN .* (s - pp(k));   end;
W = NUM ./ DEN;
Wm = abs(W);
Wm(Wm>20) = 20;  % clip the peaks

subplot(221);  mesh(sig,ome,Wm);   xlabel('sigma'); ylabel('Omega'); zlabel('|W|');
axis([-q0,q0,-q0,q0,0,5]);
subplot(222); contour(sig,ome,Wm,50);  xlabel('sigma'); ylabel('Omega'); title('|W|');
hold on; plot(linspace(0,0,2),linspace(0,q0,2),'r'); hold off;

Ome = [0:.1:q0];
NUM = 1;    for k=1:max(size(zz));    NUM = NUM .* (j*Ome - zz(k));    end;
DEN = 1;    for k=1:max(size(pp));    DEN = DEN .* (j*Ome - pp(k));   end;
W = NUM ./ DEN;
Wm = abs(W);
subplot(223);  plot(Ome,Wm);   xlabel('Omega   (lin)'); ylabel('|W|   (lin)');
subplot(224);  loglog(Ome,Wm);   xlabel('Omega   (log)'); ylabel('|W|   (log)');

% set(gcf,'PaperPositionMode','auto');
% print -deps   fdtexample.eps
return;
