% Introduction to Fourier series
%
% Samuli Siltanen Nov 2015

% Create evaluation points
Nth = 1024;
th  = 2*pi*[0:(Nth-1)]/Nth;
Dth = th(2)-th(1);

% Construct a function for studying
%f = th.*(2*pi-th);
%f = exp(1i*4*th);
f = zeros(size(th));
f((th>pi/2)&(th<pi)) = 1;

% Choose level of approximation
N = 42;

% Evaluate the inner products, or in other workds,
% Fourier coefficients
nvec  = [-N:N];
nvec  = nvec(:);
coefs = zeros(size(nvec));
for iii = 1:length(nvec)
    n     = nvec(iii);
    fii_n = exp(1i*n*th)/sqrt(2*pi);
    coefs(iii) = Dth*sum(f.*conj(fii_n));
end

% Approximate f by fN using Fourier series
fN = zeros(size(f));
for iii = 1:length(nvec)
    n     = nvec(iii);
    fii_n = exp(1i*n*th)/sqrt(2*pi);
    fN    = fN+coefs(iii)*fii_n;
end

% Show results
figure(1)
clf
plot(th,f,'b')
hold on 
plot(th,fN,'r')
xlim([0 2*pi])
%ylim(1.2*[-max(abs(f)) max(abs(f))])
ylim([-.2 1.2*max(abs(f))])
axis square