clear; close all; clc;
[y, Fs] = audioread('badDolly.mp4');

NFFT = length(y);
f_axis = linspace(-Fs/2,Fs/2,NFFT);
ck = fftshift(fft(y));

[h1,w1] = freqz(y(:,1),1,NFFT,Fs);
[h2,w2] = freqz(y(:,2),1,NFFT,Fs);

figure(2);
%plot(w1/1e3,db(abs(h1)),w2/1e3,db(abs(h2)))
plot(f_axis/1e3,abs(ck))
xlabel('Frekvens [kHz]');
ylabel('|H(e^{j\omega})|');

% Ut i fra figur(2) med FFT, finner dominante frekvenser i:
% 0,3kHz, 0,5kHz, 1kHz
% 300Hz, 500Hz, 1000Hz
% F = f / samplingsfrekvens, cos(2*pi*f/Fs)  = cos(2*pi*F)
r = 0.999;
f_1 = 300/Fs; f_2 = 500/Fs; f_3 = 1000/Fs;
b_1 = [1 (-2*cos(2*pi*f_1)) 1];
b_2 = [1 (-2*cos(2*pi*f_2)) 1];
b_3 = [1 (-2*cos(2*pi*f_3)) 1];

a_1 = [1 (-2*r*cos(2*pi*f_1)) r^2];
a_2 = [1 (-2*r*cos(2*pi*f_2)) r^2];
a_3 = [1 (-2*r*cos(2*pi*f_3)) r^2];

y_filt = filter(b_1,a_1,y);
y_filt = filter(b_2,a_2,y_filt);
y_filt = filter(b_3,a_3,y_filt);
y_filt = y_filt./max(y_filt(:));

% Forskjellen med 0.95,0.99 og 0.999, er at båndbredden
% på filteret minker når radiusen går mot resinatoren.
% e.g, mindre av originalsignalet forsvinner når r øker
f_axis = linspace(-2000,2000,length(y_filt));
plot(f_axis,(abs(fftshift(fft(y_filt)))))
xlabel('Frekvens [Hz]');
ylabel('|H(e^{j\omega})| r = 0.999');

sound(y_filt,Fs)