%%%%%%%% MATLAB CODE for evaluation of `multiperiodic' function F %%%%%%%%%
%(in notation of http://www.ippt.gov.pl/~edanicki/danickibook.pdf Chapt. 3)
g=[1 2;3 6;7 8;11 13]/16;% 0<|g|<1 - left&right strip edges in 2-col. matrx
N=512;%number of harmonics [0,N-1]
l=size(g);ne=l(1);               %ATT: ne : required EVEN number of strips!
g=g-sum(sum(g')'/2)/ne;h=1i*pi*sum(g')';y=[1i;zeros(N-1,1)];
for m=1:ne;t=cos(pi*(g(m,2)-g(m,1)));a=1i;b=1i*t;z=[zeros(1,N-2) b*exp(h(m)) a];
for l=2:N-1;c=(t*(2*l-1)*b-(l-1)*a)/l;a=b;b=c;z(N-l)=c*exp(l*h(m));end;
for i=N:-1:1;y(i)=2*z(1+N-i:N)*y(1:i);end;end;y=[zeros(ne/2,1);y];F=fft(y(1:N));
plot([real(F) imag(F)]);%this can be verified by applying Fourier transformation
 
