pkg load signal

function [b,a] = bandstop(f)

        fstop = [f-5, f+5];
        fpass = [0, 10000];
        Rpass = 1;
        Rstop = 3;
        Wpass = 2 * pi * fpass;
        Wstop = 2 * pi * fstop;
        [n, Wn_p, Wn_s] = cheb2ord (Wpass, Wstop, Rpass, Rstop, "s")
        [b, a] = cheby2 (n, Rstop, Wn_s, "stop", "s");

endfunction

function H = double_bandstop(f,b1,a1,b2,a2)

        W = 2 * pi * f;

        H1 = freqs (b1, a1, W);
        H2 = freqs (b2, a2, W);

        H = H1.*H2

endfunction

f1 = 2765
f2 = 5529
[ b1,a1 ] = bandstop(f1)
[ b2,a2 ] = bandstop(f2)


f_range = 2400:6200;
H = double_bandstop(f_range,b1,a1,b2,a2)

subplot(3,1,1)
plot (f_range, 20 * log10 (abs (H)));
title("Dubbelt chebyshev typ 2 bandsparr filter")
xlabel("Frekvens (Hz)");
ylabel ("Dampning (dB)");
ylim ([-80, 0]);

f_range1 = 2700:2820
H_z1 = double_bandstop(f_range1,b1,a1,b2,a2)

subplot(3,1,2)
plot (f_range1, 20 * log10 (abs (H_z1)));
title("inzoomning vid 2765 Hz")
xlabel("Frekvens (Hz)");
ylabel ("Dampning (dB)");
ylim ([-80, 0]);

f_range2 = 5470:5590;
H_z2 = double_bandstop(f_range2,b1,a1,b2,a2)

subplot(3,1,3)
plot (f_range2, 20 * log10 (abs (H_z2)));
title("inzoomning vid 5529 Hz")
xlabel("Frekvens (Hz)");
ylabel ("Dampning (dB)");
ylim ([-80, 0]);
print -dpdf bandstop.pdf