%   Script AIRBAG.M zur Herstellung von Bildern
%
%   Analyse und Entwurf der Kompensationsrückführung
%   des Airbag-Beschleunigungssensors
%
% J. Lunze
% 5.11.1998
% für 3. Auflage: 2.4.2001
% für 4. Auflage: 18.9.2003
% für 5. Auflage: 11.6.2005
% für 6. Auflage: 4.12.2006
% für 7. Auflage: 28.4.2008
% für 9. Auflage: 16.2.2012
% für 10. Auflage: 27.1.2014
% für 12. Auflage: 5.8.2019
%
%
clear
close all
%
%   Maßeinheiten
%
% Zeit in Sekunden
% Spannung in Volt
%
%   Problembeschreibung
%
c=3;
a=c*9.81;       % die zu messende Beschleunigung ist
            % der c-fache Wert der Erdbeschleunigung
            % c=3...50
%
%   Daten
%
mm=0.25*1e-3;       % Masse (in Kilogramm) der Membran
dm=0.5;         % Dämpfung der Membran
omegam=50000;       % Grenzfrequenz (in Hz) der Membran
Fa=50*9.81*mm;      % Kraft (in Newton) bei Beschleunigung 50g
kc=10/Fa;       % "Verstärkungsfaktor" der Membran
            % einschließlich Demodulator bei Beschleunigung
            % 50g, so daß Signalwert 10 (Volt) entsteht
U=10;           % Maximalspannung (in Volt) für Kompensation
alpha=Fa/U/U;       % Proportionalitätsfaktor in quadratischer
            % Abhängigkeit von Kraft und Spannung
s1=-omegam*dm+omegam*sqrt(1-dm*dm)*i;
s2=-omegam*dm-omegam*sqrt(1-dm*dm)*i;
n=poly([s1, s2]);
z=s1*s2*kc;


%%    Linearisierung der quadratischen Abhängigkeit F(U) bei 3g
%
%
figure(1);
subplot(2,1,1);
Fa3=3*9.81*mm;
U3g=sqrt(Fa3/alpha);
kprop=sqrt(3*9.81*mm*alpha);
k=0;
for u=0:0.01:10
    k=k+1;
    Un(k)=u;
    Fn(k)=alpha*u*u;
end
plot(Un, Fn, 'b');
latextitle('Nichtlineares \"Ubertragungsglied');
hold on;
plot([0 U U], [Fa Fa 0], 'k:');
plot([0 U3g U3g], [Fa3 Fa3 0], 'k:');
plot([0 2*U3g], [0 2*Fa3], 'b--');
Fa20=20*9.81*mm;
U20g=sqrt(Fa20/alpha);
kprop20=sqrt(20*9.81*mm*alpha);
plot([0 1.2*U20g], [0 1.2*Fa20], 'b--');
latexxlabel('$$u_{\rm c}$$ in V');
latexylabel('$$F_{\rm e}$$  in N');
rotateY;
axis([0 11 0 0.13]);
xticks([0 2 4 6 8 10]);
yticks([0 0.05 0.1]);
hold off

%%  I-erweiterte Strecke
%
%
figure(2)
ni=conv(n, [1 0]);
W=logspace(3, 6, 200);
Strecke=tf(z, ni);
[Betrag0, Phase0]=bode(kprop*Strecke, W);
%  Entwurf der Rückführung
%
kI=75000;
[Betrag1, Phase1]=bode(kI*kprop*Strecke, W);
subplot(2,1,1);
semilogx(W, 20*log10(Betrag0(:)), 'b', W, 20*log10(Betrag1(:)), 'b--');
latextitle('Offene Kette der Kompensation');
axis([1e+4 1e+6 -200 50]);
yticks([-200 -100 0]);
xleer;
latexylabel('$$|G|$$ in dB');
rotateY;
grid on;
subplot(2,1,2);
semilogx(W, Phase1(:), 'b');
latexylabel('$$\phi$$ in Grad');
rotateY;
axis([1e+4 1e+6 -300 0]);
yticks([-360 -270 -180 -90 0]);
latexxlabel('$$\omega$$ in rad/s');
xticks([10000 100000 1000000]);
grid on;



%%  Verhalten des geschlossenen Kreises
%
%
figure(3)
subplot(3,1,[1,2]);
Regler=tf(kI*kprop, [1 0]);
Strecke1=tf(z, n);
Regelkreis=feedback(Strecke1, Regler);
T=linspace(0, 0.0004, 300);
Y=step(Regelkreis, T);
plot(T, Y, 'b', [-10 10],[0 0],'k:');
axis([0 0.0004 -20 80]);
yticks([-20 0 20 40 60 80]);
xticks([0 0.0002 0.0004]);
latexxlabel('$$t$$ in s');
latexylabel('\begin{minipage}{2cm}Demodulierte Spannung\end{minipage}');
rotateY;

%%
figure(4)
subplot(2,1,1)
offeneKette=series(Strecke1, Regler);
Regelkreis1=feedback(offeneKette, 1);
Y=step(Regelkreis1, T);
plot(T, 3*Y, 'b', [0 1],[3 3], 'k:');   % für Beschleunigung 3g
latextitle('Messergebnis bei Beschleunigung 3g');
latexxlabel('$$t$$ in s');
latexylabel('$$a$$ in $$g$$');
rotateY;
axis([0 0.0004 0 3.5]);
xticks([0 0.0002 0.0004]);
yticks([0 1 2 3 4]);


%%  Analyse des Meßsystems für Beschleunigung 20g
%
%
figure(5)
U20g=sqrt(20*9.81*mm/alpha);
kprop20=20*9.81*mm/U20g;
[Betrag3, Phase3]=bode(kI*kprop*Strecke, W);
[Betrag20, Phase20]=bode(kI*kprop20*Strecke, W);
subplot(2,1,1);
semilogx(W, 20*log10(Betrag3(:)), 'b--', W, 20*log10(Betrag20(:)), 'b');
axis([1e+4 1e+6 -80 20]);
grid on;
latextitle('Kompensationskreis bei Beschleunigung 20g');
latexxlabel('$$\omega$$ in rad/s');
latexylabel('$$|G|$$ in dB');
rotateY;
yticks([-80 -40 0]);
xticks([10000 100000 1000000]);


%%
figure(6)
subplot(2,1,1);
Regler2=tf(kI*kprop20, [1 0]);
offeneKette2=series(Strecke1, Regler2);
Regelkreis2=feedback(offeneKette2, 1);
Y2=step(Regelkreis2, T);
plot(T, 20*Y, 'b--', T, 20*Y2, 'b', [0 0.0004], [20 20], 'k:');
latextitle('Messergebnis bei 20g');
latexxlabel('$$t$$ in s');
latexylabel('$$a$$ in $$g$$');
rotateY;
xticks([0 0.0002 0.0004]);
yticks([0 20 40]);




%%   Herstellung von .ps und .eps-Bildern
%
Epsbildc15('airbag');
