%   Script ROBSTAB.M
%
%   Toleranzband für die Ortskurve eines Systems
%
% J. Lunze
% 17.11.1995
% für 2. Auflage: 2.11.1998
% für 3. Auflage: 3.4.2001
% für 4. Auflage: 18.9.2003
% für 5. Auflage: 18.4.2005
% für 6. Auflage: 4.12.2006
% für 7. Auflage: 28.4.2008
% für 9. Auflage: 10.2.2012
% für 10. Auflage: 27.1.2014
% für 12. Auflage: 15.07.2019
%
%   
clear
close all
z=[1];

%%  Ortskurve des Systems 1/(1 + 2dTs + T^2s^2)(1+sT3)
%   für T=1, d=0.5, T3=2
%
%   mit Fehlerabschätzung |dG| < |0.5/(1+2 dr Tr s + Tr^2s^2)|
%   mit Tr = 1,77, dr=0,333
%   
figure(1);
d=0.5;
n=[1 2*d 1];
r=roots(n);
r(3) = -0.5;
n=poly(r);
System=tf(z, n);
[Real, Imag]=nyquist(System);
plot(Real(:), Imag(:), 'b', Real(:), -Imag(:), 'b--');
grid('on');
latextitle('Toleranzband f\"ur die Ortskurve');
latexxlabel('Re');
latexylabel('Im');
rotateY;
axis([-2 3 -3 2]);
hold on;
plot([-2 3],[0 0],'k:',[0 0],[-3 2],'k:');
%                                            Robustheitsabschätzung
zr=[0.5];
nr=[3.49 1 1];
Regler=tf(zr, nr);
O=[0:0.1:2];
kO=length(O);
[Betrag, Phase] = bode(Regler, O);
[BetragG, PhaseG] = bode(System, O);
for k=1:kO
    k1 = 0;
    RealG = BetragG(k)*cos(PhaseG(k)*pi/180);
    ImagG = BetragG(k)*sin(PhaseG(k)*pi/180);
    for x1=-Betrag(k):Betrag(k)/100:Betrag(k)
        k1=k1+1;
        K1(k1, 1)= x1+RealG;
        K1(k1, 2)= ImagG + sqrt(Betrag(k)^2 - x1^2);
        K1(402-k1, 1) = x1+RealG;
        K1(402-k1, 2) = ImagG - sqrt(Betrag(k)^2 - x1^2);
    end;
    plot(K1(:, 1), K1(:, 2), 'b');
end;
axis([-1.3 2.7 -2.7 1.3]);
xticks([-1 0 1 2 3]);
yticks([-3 -2 -1 0 1 2 3]);
axis('square');

%%     Bilder auf .EPS- und .PS-Dateien
%
Epsbildc15('robstab');
