% Filter basics

%convolution

conv([1 2 3],[1 2 1])
conv([2 2 1],[1 1 2])
conv([0.3 0.3 0.3], [5 10 -4 6 -2 9 -8])

% Convlution Theorem
a = [1 2 3 4];
b = [4 5 6 7];
A = fft(a, 8);
B = fft(b, 8);
C = A .* B;
round(real(ifft(C)))
conv(a, b)

% averaging filter
x = randn(10,1)   % A random vector of length 5
h = [1 1 1 1]/4;  % Length 4 averaging filter
y = conv(h,x);
plot(1:10, x, 1:10, y(1:10));
figure;
yy = filter(h, 1, x);
plot(1:10, x, 1:10, yy);


fs = 8192;
t = 0:1/fs:1 - 1/fs;
n = length(t);
x = sin(2*pi*200*t) + sin(2*pi*1000*t) + sin(2*pi*2000*t);
y = fft(x); 
m = abs(y(1:n/2));
p = unwrap(angle(y(1:n/2)));
subplot(2,1,1);
plot((1:n/2)/n*fs, 20 * log10(m))
%semilogy((1:n/2)/n*fs, m)
ylabel('Abs. Magnitude'), grid on
subplot(2,1,2), plot((1:n/2)/n*fs, p * 180/pi)
ylabel('Phase [Degrees]'), grid on
xlabel('Frequency [Hertz]')

% y(n) = (x(n) + x(n-1)) / 2
zplane([0.5 0.5], 1);   % zplane(b, a)  row vectors
zplane([-1]', 0);       % zplane(zeros, poles)  column vectors

% FIR
h = [1 1 1 1]/4;  % Length 4 averaging filter
freqz(h, 1);
freqz(h, 1, 512, 8192);

figure;
xx = filter(h, 1, x);
y = fft(xx);
m = abs(y(1:n/2));
p = unwrap(angle(y(1:n/2)));
subplot(211);
plot((1:n/2)/n*fs, 20*log10(m))
ylabel('Magnitude'), grid on
subplot(212), plot((1:n/2)/n*fs, p * 180/pi)
ylabel('Phase [Degrees]'), grid on
xlabel('Frequency [Hertz]')

% IIR
b = 1;  %no zero
a = [1 -0.9]  % one pole

impz(b, a);
freqz(b, a);
freqz(b, a, 512, 8192);


#FFT of impulse response == freq response
[h t] = impz(b, a);
H = fft(h);
%Mag
figure, subplot(211)
plot(20*log10(abs(H)))
grid on, axis([0 50 -10 20])
%Phase
subplot(212)
plot(180 * (angle(H) ./ (pi)))
grid on, axis([0 50 -80 0])


figure;
xx = filter(b, a, x);
y = fft(xx);
m = abs(y(1:n/2));
p = unwrap(angle(y(1:n/2)));
subplot(211);
plot((1:n/2)/n*fs, 20*log10(m))
ylabel('Magnitude'), grid on
subplot(212), plot((1:n/2)/n*fs, p * 180/pi)
ylabel('Phase [Degrees]'), grid on
xlabel('Frequency [Hertz]')


% Butterworth
[b,a] = butter(1, 0.1);
freqz(b,a,256, 8192);
figure;
[b,a] = butter(3, 0.1);
freqz(b,a,256, 8192);
figure;
[b,a] = butter(5, 0.1);
freqz(b,a,256, 8192);


[b,a] = butter(2, 0.1);
[h1, f, s] = freqz(b,a,256, 8192);
[b,a] = butter(5, 0.1);
h2 = freqz(b,a,256, 8192);
[b,a] = butter(10, 0.1);
h3 = freqz(b,a,256, 8192);

h = [h1 h2 h3];
s.plot = 'mag';
freqzplot(h, f, s);


figure;
[b,a] = butter(5, 0.1, 'high');
freqz(b,a,256, 8192);

figure;
[b,a] = butter(5, [0.4 0.6], 'stop');
freqz(b,a,256, 8192);

figure;
[b,a] = butter(5, [0.4 0.6]);
freqz(b,a,256, 8192);

% FIR filter Design
b = fir1(48,[0.35 0.65]);
freqz(b,1,512)

% IIR filter design
% Butterworth   - no ripple
% Chebycheff I  - ripple in pass band
% Chebycheff II - ripple in stop band
% Elliptic      - ripple in pass and stop band; narrow transition bands
