f1 = 880e6;
f2 = 13.56e6;
r = 50;
l1 = 1e-8;
l2 = 1e-6;
c1 = (1 / (2 * pi * f1))^2 / l1;
c2 = (1 / (2 * pi * f2))^2 / l2;

% Impedanzen der Spulen und Kondensatoren
f = [1:1e3:1000e6];
xl1 = i * (2 * pi * f * l1);
xl2 = i * (2 * pi * f * l2);
xc1 = -i*(1 ./ (2 * pi * f * c1));
xc2 = -i*(1 ./ (2 * pi * f * c2));

% Zwei parallele Schwingkreise an 50 Ohm
zs1 = xl1 + xc1 + r;
zs2 = xl2 + xc2 + r;
zp = 1 ./ (1 ./ zs1 + 1 ./ zs2);
z = zp ./ (r + zp);

subplot(2,1,1);
g = mag2db(abs(z));
plot(f, mag2db(abs(z)), 'b', f1, g(f1 / 1e3), 'ro', f2, g(f2 / 1e3), 'ro');
title(['Bandselektion von f1=' num2str(f1, 3) ' Hz und f2=' num2str(f2, 3) 'Hz']);
ylabel('Dämpfung [dB]');
xlabel('Frequenz [Hz]');
subplot(2,1,2);
a = angle(z) / pi * 180;
plot(f, a, 'b', f1, a(f1 / 1e3), 'ro', f2, a(f2 / 1e3), 'ro');
ylabel('Phase [grad]');
xlabel('Frequen [Hz]');

b1 = f1 / (2 * pi * f1 * l1 / r);
b2 = f2 / (2 * pi * f2 * l2 / r);
disp(['l1=' num2str(l1) ' c1=' num2str(c1) ' b1=' num2str(b1, 3)]);
disp(['l2=' num2str(l2) ' c2=' num2str(c2) ' b1=' num2str(b2, 3)]);

%z1 = zs1 ./ (r + zs1);
%z2 = zs2 ./ (r + zs2);
