%   Script ZFBANAL.M zur Herstellung von Bildern
%
%   Verhalten von zeitdiskreten Systemen im Frequenzbereich
%
% verwendet EPSBILD.M
%
% J. Lunze
% 6.1.1997                Daten wie 27.11.1996
% für 2. Auflage: 15.12.2001
% für 3. Auflage: 15.3.2004
% für 5. Auflage: 18.10.2007
% für 6. Auflage: 5.1.2010
% für 7. Auflage: 25.5.2012
% für 8. Auflage: 4.4.2014
% für 9. Auflage: 9.5.2016
% für 10. Auflage: 9.9.2019
%
%
clear
close all

A = [-0.8 2;
     -2 -0.8];
b = [2;
     2];
c = [1 1];
d = 0.6;
x0=[1, 1];
Ta=0.4;
[Ad, bd]=c2d(A, b, Ta);

%%  Modelleigenschaften (Beispiel 13.1)
%
Strecke=ss(A, b, c, d);
%  Regler  u(k) = -K x(k) + kw w(k)
K=1/(1+d)*c;
kw=1/(1+d);
%  Systemmatrix des geschlossenen Kreises:
Aquer=Ad-bd*K;
bquer=bd*kw;
%  offene Kette = Strecke (da Einheitsr\"uckf\"uhrung)
dStrecke=ss(Ad, bd, c, d);
[Z, N]=tfdata(dStrecke,'v');


%%  Ortskurve eine zeitdiskreten Systems
%
%
figure(1);
subplot(2,3,2)
[Real, Imag, W]=dnyquist(Ad, bd, c, d, Ta);
plot(Real, Imag, 'b', Real, -Imag, 'b--');
hold on;
plot([0 0], [-4 4], 'k:', [-4 4], [0 0], 'k:');
plot([-1 -1], [0.2 -0.2], 'b');
plot(ddcgain(Ad, bd, c, d), 0, 'bo');
text(1.5, -2.8, '\leftarrow');
hold off;
axis([-3 4 -3.5 3.5]);
xticks([-3 0 3]);
yticks([-3 0 3]);
axis('square');
latexxlabel('Re');
latexylabel('Im');
rotateY;
latextitle('Ortskurve eines zeitdiskreten Systems');
subplot(2,3,1);
System=ss(A, b, c, d);
[R, I, Omega]=nyquist(System);
plot(R(:), I(:), 'b', R(:), -I(:), 'b--');hold on;
plot([0 0], [-4 4], 'k:', [-4 4], [0 0], 'k:');
plot([-1 -1], [0.2 -0.2], 'k');
plot(dcgain(System), 0, 'bo');
text(1.5, -2, '\leftarrow');
hold off;
latexxlabel('Re');
latexylabel('Im');
rotateY;
axis([-3 4 -3.5 3.5]);
xticks([-3 0 3]);
yticks([-3 0 3]);
axis('square');

%%  Ortskurve bei Verdopplung der Abtastzeit
%
%
figure(1);
subplot(2, 3, 3);
Ta=0.8;
[Ad, bd]=c2d(A, b, Ta);
[Real, Imag]=dnyquist(Ad, bd, c, d, Ta);
plot(Real, Imag, 'b', Real, -Imag, 'b--');
hold on;
plot([0 0], [-5 5], 'k:', [-5 5], [0 0], 'k:');
plot([-1 -1], [0.1 -0.1], 'k');
plot(ddcgain(Ad, bd, c, d), 0,'bo');
text(0.5, -2.2, '\leftarrow');
hold off;
axis([-2.5 3.5 -3 3]);
xticks([-3 0 3]);
yticks([-3 0 3]);
axis('square');
latexxlabel('Re');
latexylabel('Im');
rotateY;
latextitle('Ortskurve bei doppelter Abtastzeit');

%%
%%
%    Näherungsbetrachtung: kontinuierliches System mit Totzeit T/2
%
%
figure(3)
subplot(3,3,[1, 2, 4, 5]);
W=logspace(-1.5, 1.5, 300);
System=ss(A, b, c, d);
[Betrag, Phase, W]=bode(System, W);
for k=1:max(size(W))
    omega=W(k);
    PhaseT(k)=-(omega*Ta/2)*360/2/pi;   % Umrechnung in Grad
    Phase(1,1,k)=Phase(1,1,k)+PhaseT(k);
end;
%Phase(:)=Phase(:)+PhaseT(:);
Phase=Phase*2*pi/360;        % Umrechnung in rad
G = Betrag.*(cos(Phase)+i*sin(Phase));
plot(real(G(:)), imag(G(:)), 'b', [-3 3], [0,0], 'k:', [0,0], [-4 4], 'k:');
hold on;
plot(dcgain(System), 0, 'bo');
plot(-1, 0, 'b+');
text(0.5, -2.8, '\leftarrow');
latextitle('Kontinuierliches System mit Totzeit $$T_t=T/2$$');
latexxlabel('Re');
latexylabel('Im');
rotateY;
axis('square');
axis([-1.5, 3, -3.5, 1]);
xticks([-1 0 1 2 3]);
yticks([-3 -2 -1 0 1]);
hold off;


%%
%   Bodediagramm des zeitdiskreten Systems
%
%
figure(4);
%   kontinuierliches System aus Abschnitt RTI-6.5.3
zk=[1 -2];
nk=[1 2 5];
[A, B, C, D]=tf2ss(zk, nk);
Ta=0.4;
[Ad, Bd]=c2d(A, B, Ta);
[zs, ns]=ss2tf(Ad, Bd, C, D);
[Betrag, Phase, W]=dbode(Ad, Bd, C, 0, Ta);
subplot(4,1,[1,2]);;
semilogx(W, 20*log10(Betrag(:)), 'b');
latextitle('Bodediagramm des zeitdiskreten Systems');
xleer;
latexylabel('$$|G({\rm e}^{{\rm j}\omega T})|$$');
rotateY;
yticks([-10 -5]);
grid('on');
%
subplot(4,1,[3,4]);;
semilogx(W, Phase(:), 'b');
latexxlabel('$$\omega$$ in $${\rm rad/s}$$'); %  in $$\frac{\rm rad}{\rm s}$$
latexylabel('$${\rm arg}(G)$$');
rotateY;
grid('on');
xticks([0.01 0.1 1 10 100]);
yticks([-200 -100 0 100 200]);

%%
%%
%   Vergrößerung des Frequenzbereiches \"uber wT hinaus
%
%
figure(5);
W=logspace(-1, 2, 700);
[Betrag, Phase, W]=dbode(Ad, Bd, C, 0, Ta, 1, W);
subplot(4,1,[1,2]);;
semilogx(W, 20*log10(Betrag(:)), 'b');
latextitle('Bodediagramm des zeitdiskreten Systems');
xticks([0.01 0.1 1 10 100]);
yticks([-10 -5]);
latexylabel('$$|G({\rm e}^{{\rm j}\omega T})|$$');
rotateY;
grid('on');
%subplot(4,1,[3,4]);;
%semilogx(W, Phase(:));
%latexxlabel('$$\omega$$ in $$\frac{\rm rad}{\rm s}$$');
%ylabel('arg G in Grad');
%rotateY;
%axis([0.1, 100, -1000, 270]);
%grid('on');


%%
%     Ausgabe der Bilder als .eps-Dateien
%
Epsbildc15('zfbanal');
