%   Script DMOTOR.M zur Herstellung von Bildern
%
%   Diskrete Realisierung eines Reglers für den Gleichstrommotor
%   Reglerparameter aus FKLMOTR.M (Band 1)
%
% J. Lunze
% für 3. Auflage: 19.4.2004
% für 5. Auflage: 19.10.2007
% für 7. Auflage: 25.5.2012
% für 10. Auflage: 9.9.2019
%
%
clear
close all
%%    Beispiel 2: Entwurf einer Geschwindigkeitsregelung
%                für einen Gleichstrommotor
Bspmotor;
%
%  Drehzahlgeber
Td = 0.008;       % s
Drehzahlgeber=tf([1], [Td 1]);
MotorMitDrehzahl=series(MotorUF, Drehzahlgeber);
O=logspace(-1, 3, 300);

%%   Verwendung eines PI-Reglers
%    mit kp=4, TN=0.0777
%    (also Nullstelle bei 7.08; |Gpi(s)| = 40dB für s->oo)
%
figure(1);
subplot(4,1,[1,2]);
kP=4;
TI=0.0777;
zr=[kP*TI kP];
nr=[TI 0];
PIRegler=tf(zr, nr);
offeneKette=series(MotorMitDrehzahl, PIRegler);
%
%   Übergangsfunktion des PI-geregelten Motors
%
Regelkreis=feedback(offeneKette, 1);
[Y, T]=step(Regelkreis);
plot(T, Y, 'b');
latextitle('F\"uhrung\"ubergangsfunktion des Motors');
latexxlabel('$$t$$ in s');
latexylabel('$$h_{\rm w}$$ in $$\frac{1}{\rm s}$$');
rotateY;
grid('on');
axis([0 1.5 0 1.1]);
xticks([0 0.5 1 1.5]);
yticks([0 1]);


%%  zeitdiskrete Realisierung des PI-Reglers
%
%
figure(2);
%  Grenzfrequenz der offenen Kette  2 rad/s
%  also nach Abtasttheorem:  2 pi/T > 4 rad/s, also T < 1.6
Ta=0.02;
dRegelstrecke=c2d(MotorMitDrehzahl, Ta);
%  Reglerparameter nach Abschn. 14.1
%  K_PI(z)=(c2 z + c1)/(z-1)
c2=kP + kP*Ta/2/TI;
c1=-kP + kP*Ta/2/TI;
dPIRegler=tf([c2 c1],[1 -1],Ta);
% Vergleich der Bodediagramme der Regler
subplot(4,1,[1,2]);
W=logspace(0, 3);
[Betrag, Phase]=bode(PIRegler, W);
[dBetrag, dPhase, dW]=dbode([c2 c1], [1 -1], Ta);
semilogx(W, 20*log10(Betrag(:)), 'b--', dW, 20*log10(dBetrag(:)), 'b');
latexylabel('$$|K|$$ in dB');
rotateY;
xleer;
latextitle('Kontinuierliche und zeitdiskrete Drehzahlregelung ($$T=0,02 s$$)');
axis([1 1000 0 50]);
yticks([0 40]);
grid('on');
% Vergleich der Bodediagramme der offenen Ketten
subplot(4,1,[3,4]);
doffeneKette=series(dRegelstrecke, dPIRegler);
dRegelkreis=feedback(doffeneKette, 1);
[Betrag, Phase, W]=bode(offeneKette,W);
[zo, no]=tfdata(doffeneKette,'v');
[dBetrag, dPhase, dW]=dbode(zo, no, Ta);
semilogx(W, 20*log10(Betrag(:)), 'b--', dW, 20*log10(dBetrag(:)), 'b');
latexxlabel('$$\omega$$ in $${\rm rad/s}$$'); %  in $$\frac{\rm rad}{\rm s}$$
yticks([-100 0]);
latexylabel('$$|G|$$ in dB');
rotateY;
grid('on');
axis([1 1000 -120 20]);
xticks([1 10 100 1000]);
yticks([-100 0]);


%%   Vergleich des Führungsverhaltens
%
%
figure(3)
subplot(4,1,[1,2]);
[dz, dn]=tfdata(dRegelkreis, 'v');
dY=dstep(dz, dn);
T=Ta*[0:length(dY)-1];
stairs(T, dY, 'Color','b','LineWidth',1.5);
hold on
%plot(T, dY, 'b*', [0 10],[0 0],'k:');
plot([0 10],[0 0],'k:');
Y=step(Regelkreis, T);
plot(T, Y, 'b--');
axis([0 1.5 -0.1 1.1]);
xticks([0 0.5 1 1.5]);
yticks([0 1]);
latexxlabel('$$t$$ in s');
latexylabel('$$h_{\rm w}$$ in $$\frac{1}{\rm s}$$');
rotateY;
latextitle('Vergleich der F\"uhrungs\"ubergangsfunktionen (T=0,02 s)');
hold off;


%%
%
%   zeitdiskrete Realisierung des PI-Reglers bei vergrößerter Abtastzeit
%
%
figure(4);
% neue Abtastzeit
Ta=0.2;
dRegelstrecke=c2d(MotorMitDrehzahl, Ta);
%  Reglerparameter nach Abschn. 14.1
%  K_PI(z)=(c2 z + c1)/(z-1)
c2=kP + kP*Ta/2/TI;
c1=-kP + kP*Ta/2/TI;
dPIRegler=tf([c2 c1],[1 -1],Ta);
% Vergleich der Bodediagramme der Regler
subplot(4,1,[1,2]);
[Betrag, Phase, W]=bode(PIRegler);
[dBetrag, dPhase, dW]=dbode([c2 c1], [1 -1], Ta);
semilogx(W, 20*log10(Betrag(:)), 'b--', dW, 20*log10(dBetrag(:)), 'b');
latexylabel('$$|K|$$ in dB');
rotateY;
xleer;
latextitle('Kontinuierliche und zeitdiskrete Drehzahlregelung (T=0,2 s)');
axis([1 1000 0 50]);
yticks([0 40]);
grid('on');
% Vergleich der Bodediagramme der offenen Ketten
subplot(4,1,[3,4]);
doffeneKette=series(dRegelstrecke, dPIRegler);
dRegelkreis=feedback(doffeneKette, 1);
[Betrag, Phase, W]=bode(offeneKette);
[zo, no]=tfdata(doffeneKette,'v');
[dBetrag, dPhase, dW]=dbode(zo, no, Ta);
semilogx(W, 20*log10(Betrag(:)), 'b--', dW, 20*log10(dBetrag(:)), 'b');
latexxlabel('$$\omega$$ in $${\rm rad/s}$$'); %  in $$\frac{\rm rad}{\rm s}$$
latexylabel('$$|G|$$ in dB');
rotateY;
axis([1 1000 -120 20]);
xticks([1 10 100 1000]);
yticks([-100 0]);
grid('on');


%%   Vergleich des F\"uhrungsverhaltens
%
%
figure(5)
subplot(4,1,[1,2]);
[dz, dn]=tfdata(dRegelkreis, 'v');
dY=dstep(dz, dn);
T=Ta*[0:length(dY)-1];
stairs(T, dY, 'Color','b','LineWidth',1.5);
hold on
plot(T, dY, 'bo', [0 10],[0 0],'k:');
%plot([0 10],[0 0],'k:');
Y=step(Regelkreis, T);
plot(T, Y, 'b--');
axis([0 1.5 -0.1 1.1]);
xticks([0 0.5 1 1.5]);
yticks([0 1]);
latexxlabel('$$t$$ in s');
latexylabel('$$h_{\rm w}$$ in $$\frac{1}{\rm s}$$');
rotateY;
latextitle('Vergleich der F\"uhrungs\"ubergangsfunktionen (T=0,2 s)');
hold off;


%%   Abgetastete Signale mit und ohne Messrauschen
%
%
figure(6);
% neue Abtastzeit
Ta=0.05;
Tend=1.5;
omegaR=1.1*2*pi/Ta;
kend=Tend/Ta;
T=Ta*[0:kend-1];
dt=0.001;
T2=0:dt:Tend;
%   Darstellung der Störung vor und hinter Abtaster
subplot(4,1,[1,2]);;
kend2=length(T2);
for k=1:kend2
   R2(k)=sin(T2(k)*omegaR+0.1);
end;
R3=[];
for k=1:kend
   R3(k)=sin((k-1)*Ta*omegaR+0.2);
end;
plot(T2, R2, 'b', T, R3, 'bo', [0 Tend],[0 0],'k:');
axis([0 Tend -1.1 1.1]);
yticks([-1 0 1]);
latexylabel('$$r$$ in $$\frac{1}{\rm s}$$');
rotateY;
xleer;
latextitle('St\"orung des drehzahlgeregelten Gleichstrommotors (T=0,05 s)');
subplot(4,1,[3,4]);
stairs(T, R3, 'Color','b','LineWidth',1.5);
hold on
plot(T, R3, 'bo', [0 Tend],[0 0],'k:');
latexylabel('$$r$$ in $$\frac{1}{\rm s}$$');
rotateY;
latexxlabel('$$t$$ in s');
axis([0 Tend -1.1 1.1]);
xticks([0 0.5 1 1.5]);
yticks([-1 0 1]);
hold off;


%%  gestörter Regelkreis
%
%
figure(7)
subplot(4,1,[1,2]);
Y=lsim(Regelkreis, -R2, T2);
plot(T2, Y, 'b', [0 Tend],[0 0],'k:');
latextitle('Gest\"orter kontinuierlicher und diskreter Regelkreis (T=0,05 s)');
axis([0 Tend -0.03 0.01]);
xleer;
latexylabel('$$y$$ in $$\frac{1}{\rm s}$$');
rotateY;
yticks([-0.03 0]);
subplot(4,1,[3,4]);
% Vergleich mit diskretem Regelkreis
plot(T2, Y, 'b', [0 Tend],[0 0],'k:');
hold on
dRegelstrecke=c2d(MotorMitDrehzahl, Ta);
%  Reglerparameter nach Abschn. 14.1
%  K_PI(z)=(c2 z + c1)/(z-1)
c2=kP + kP*Ta/2/TI;
c1=-kP + kP*Ta/2/TI;
dPIRegler=tf([c2 c1],[1 -1],Ta);
doffeneKette=series(dRegelstrecke, dPIRegler);
dRegelkreis=feedback(doffeneKette, 1);
[dz, dn]=tfdata(dRegelkreis,'v');
Yd=dlsim(dz, dn, -R3);
plot(T, Yd,'bo');
stairs(T, Yd, 'Color','b','LineWidth',1.5);
latexylabel('$$y$$ in $$\frac{1}{\rm s}$$');
rotateY;
latexxlabel('$$t$$ in s');
axis([0 Tend -0.4 0.4]);
xticks([0 0.5 1 1.5]);
yticks([-0.4 0 0.4]);



figure(8)
% Einbau eines Antialiasing-Filters
% Grenzfrequenz entsprechend Abtastzeit Ta bei omega_Gr=pi/Ta
subplot(4,1,[1,2]);
omegaGr=pi/Ta/3;
Filter=tf([1],[1/omegaGr 1]);   % PT1-Glied
%  Näherung: Messstörung R wird erst gefiltert, dann an den
%  zeitdiskreten Regelkreis angelegt
R2F=lsim(Filter, R2, T2);
% Abtasten des gefilterten Messrauschens
dk=ceil(Ta/dt);
for k=1:kend
    R2d(k,1)=R2F((k-1)*dk+1,1);
end
YdF=dlsim(dz, dn, -R2d);
plot(T, YdF,'bo', [0 Tend],[0 0],'k:');
hold on
stairs(T, YdF, 'Color','b','LineWidth',1.5);
latexylabel('$$y$$ in $$\frac{1}{\rm s}$$');
rotateY;
latexxlabel('$$t$$ in s');
axis([0 Tend -0.4 0.4]);
xticks([0 0.5 1 1.5]);
yticks([-0.4 0 0.4]);


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