%truncated FIR low pass filter demonstration
%
%
clear
clf
%set alpha = 4 (frequency cutoff equal to r/4,
%where r is the sampling rate).
alpha = 4.;
%
%  Use 128 point filters.
%
N = 128;
%assign normalized frequency range, f = k/N
%where f=1 is sampling rate
df=1/1024;
f=0:df:1/2;

%
%  First, use a rectangular window.
%
figure(1)

for ii=1:N
j=ii-N/2;
if j==0
y(ii)=2/alpha;
else
y(ii) = (2/alpha)*sin(pi*2*j/alpha)/(pi*2*j/alpha);
end
end

for k=1:length(f)
s(k)=0;
for l=1:N
s(k) = s(k) + (y(l)*exp(i*2*pi*(l-1)*f(k)));
end
end

spec=abs(s);
%
%  Plot the windowed filter.
%
subplot(3,1,1)
plot(1:N,y)
title('Rectangular Window; N=128')
axis([0 N-1 -2/alpha 2/alpha]);
%
% Plot the frequency response.
%
subplot(3,1,2)
semilogy(f,spec);
xlabel('f/r')
ylabel('|W_k|');
axis([0 1/2 1.0e-5 10]);
%
% Plot the phase response.
% 
phase=unwrap(angle(s),3)*180/pi;
subplot(3,1,3);
plot(f,phase);
xlabel('f/r');
ylabel('\theta (degrees)');

%
%  Now, the Bartlett window.
%

ystore = y;

y=y.*(bartlett(length(y)))';

for k=1:length(f)
s(k)=0;
for l=1:N
s(k) = s(k) + (y(l)*exp(i*2*pi*(l-1)*f(k)));
end
end

spec=abs(s);

figure(2)
subplot(3,1,1)
plot(1:N,y)
title('Bartlett Window; N=128')
axis([0 N-1 -2/alpha 2/alpha]);
subplot(3,1,2)
semilogy(f,spec);
xlabel('f/r')
ylabel('|W_k|');
axis([0 1/2 1.0e-5 10]);
%
% Plot the phase response.
% 
phase=unwrap(angle(s),3)*180/pi;
subplot(3,1,3);
plot(f,phase);
xlabel('f/r');
ylabel('\theta (degrees)');

y=y.*(hamming(length(y)))';

for k=1:length(f)
s(k)=0;
for l=1:N
s(k) = s(k) + (y(l)*exp(i*2*pi*(l-1)*f(k)));
end
end

spec=abs(s);

%
% Now, the Hamming window.
%
figure(3)
subplot(3,1,1)
plot(1:N,y)
title('Hamming Window; N=128')
axis([0 N-1 -2/alpha 2/alpha]);
subplot(3,1,2)
semilogy(f,spec);
xlabel('f/r');
ylabel('|W_k|');
axis([0 1/2 1.0e-5 10]);
%
% Plot the phase response.
% 
phase=unwrap(angle(s),3)*180/pi;
subplot(3,1,3);
plot(f,phase);
xlabel('f/r');
ylabel('\theta (degrees)');