% nichols_loci_v2.m
% Adds the loci of constant gain of the complementary sensitivity 
% function to a Nyquist plot
%   kja 071212
aminit;
%
% Create patch vectors
M=1.1;cx1=-M^2/(M^2-1);r1=M/(M^2-1);
dpi=pi/200;xl=200; yl=r1*dpi; yh=xl;
X1=[-xl -xl];Y1=[-yh -yl];
for ff=-pi+dpi:dpi:pi-dpi;    
    X1=[X1 cx1+r1*cos(ff)];Y1=[Y1 r1*sin(ff)];
end
X1=[X1 -xl -xl];Y1=[Y1 yl yh];
M=1/M;cx2=-M^2/(M^2-1);r2=-M/(M^2-1);
X2=[xl xl];Y2=[yh yl];
for ff=dpi:dpi:2*pi-dpi,
    X2=[X2 cx2+r2*cos(ff)];Y2=[Y2 r2*sin(ff)];
end
X2=[X2 xl xl];Y2=[Y2 -yl -yh];
X=[X1 X2 -xl];Y=[Y1 Y2 -yh];
Z=X+i*Y;LZ=log(Z);
patch(unwrap(imag(LZ)),real(LZ),[0.7,0.7,0.7]);
hold on;
%
% Range of M-values
MMt=[1.02 1.05 1.1 1.2 1.5 2];
s='b-';
%
% Constant gain loci
for ii=1:6,
M=MMt(ii);
cx=M^2/(M^2-1);cy=0;r=M/(M^2-1);
nicplot(-cx,cy,r,s);hold on;
M=1/M;
cx=-M^2/(M^2-1);cy=0;r=M/(M^2-1);
nicplot(cx,cy,r,s);hold on;
end
%
% Constant phase loci
PP=pi*[10 20 40 60 80]/180;s=':';
for jj=1:5,
  fi=PP(jj);
  cx=-0.5;cy=0.5/tan(fi);r=0.5/sin(fi);
  nicplot(cx,cy,r,s);
  nicplot(cx,-cy,r,s);
end  
%
%   Labels and axes
axis([-2*pi 0 -pi/2 3*pi/2])
axis square
line([-2*pi 0],[0 0]);
line([-pi -pi],[-pi/2 3*pi/2]);
xxx=logspace(-2,2);
z=-0.5+i*xxx;
pl=plot(-unwrap(imag(log(z))),real(log(z)),'b--');
set(pl,'Linewidth',1.5);
xlabel('arg L(iw)');ylabel('log L(iw)');

amprint('nichols-loci.eps');
