function intDFT = intFFT(a,b,n,f)
%Integral of a function using truncated fourier series and fft
%Code by: James Hateley
%
%INPUT
%   a = left endpoint
%   b = right endpoint
%   n = number of points for partition must be 2^k
%   f = function wanting to integrate
%
%OUTPUT
%   intDFT = value of integral at each point of partition
%
%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%   Initialize
%
x=linspace(a,b,n+1);
x(n+1)=[]; N=n/2;
Jvec = 1:n;
%
%   Load coefficients
%
beta=fft(f(x),n);
%
%   Calculate coeffients for integral
%
Bbeta=zeros(size(beta));
Bbeta(1)=Bbeta(1) + 1i*sum((beta(2:N)-beta(n:-1:N+2))./(Jvec(1:N-1)));
Bbeta(2:N) = -1i*(beta(2:N)./Jvec(1:N-1));
Bbeta(n:-1:N+2) = 1i*(beta(n:-1:N+2)./Jvec(1:N-1));
%
%   Compute integral
%
intDFT = ifft(Bbeta,n)+ beta(1)*x+beta(1,N+1)/N*sin(N*x);
clear Bbeta
clear Beta
clear x
clear Jvec
end
