clear all;
RD=0.060; L= 33e-6;             % main inductor + parasitics
RC=0.080; C=100e-6; Lp=500e-12; % main capacitor + parasitics
fg = 1/(2*pi*sqrt(L*C)); % approximative cutoff frequency
fnC = 1/(2*pi*RC*C);     % approximative 1. zero
fnLp= RC/(2*pi*Lp);      % approximative 2. zero
fqn3= RD/(2*pi*L);       % additional zero of QTF

% computation of PTF(s)
ap0=1; ap1=RC*C;      ap2=Lp*C;
bp0=1; bp1=(RC+RD)*C; bp2=(L+Lp)*C;
% computation of QTF(s)
aq0=RD;  aq1=RC*RD*C+L; aq2=RD*Lp*C+RC*L*C; aq3=Lp*L*C;
bq0=bp0; bq1=bp1;       bq2=bp2;            bq3=1e-14;

% computing accurate 1. and 2. zeors
spn1=-RC/(2*Lp)+sqrt((RC/(2*Lp))^2-1/(Lp*C)); 
spn2=-RC/(2*Lp)-sqrt((RC/(2*Lp))^2-1/(Lp*C)); 
fpn1=-spn1/(2*pi); % accurate 1. zero (compare to fnC)
fpn2=-spn2/(2*pi); % accurate 2. zero (compare to fnLp)
% computing accurate cutoff frequency fpp
spp1=-(RC+RD)/(2*(L+Lp)) + sqrt(((RC+RD)/(2*(L+Lp)))^2 - 1/((L+Lp)*C));
spp2=-(RC+RD)/(2*(L+Lp)) - sqrt(((RC+RD)/(2*(L+Lp)))^2 - 1/((L+Lp)*C));
fpp = sqrt(spp1*spp2)/(2*pi);

% implementing time-continuous transfer functions
PTF_s = tf([ap2 ap1 ap0],[bp2 bp1 bp0]);
QTF_s = tf([aq3 aq2 aq1 aq0],[bq3 bq2 bq1 bq0]);
figure(1);
o=bodeoptions; o.FreqUnits='Hz';
bodeplot(PTF_s,QTF_s,o); grid on;
figure(2); 
subplot(211); step(PTF_s,QTF_s); grid on;
subplot(212); impulse(PTF_s,QTF_s); grid on;




