clear
function Z=dFt(Y)
    N=length(Y);
    Z=zeros(Y);
    for k=1:N
        s=0;
        for l=1:N
            s=s+Y(l)*exp(-2*%i*(k-1)*(l-1)*%pi/N);
        end
        Z(k)=s/N;
    end
endfunction
dFt([1 2 3 4])

function q=pN(x,Y)
    N=length(Y);
    M=N/2;
    Z=dFt(Y);
   // Z=fft(Y)/length(Y);// si par hasard dFt ne marche pas
    q=Z(M+1)*cos(M*x);
    for k=1:M
        q=q+Z(k)*exp(%i*(k-1)*x);
    end
    for k=1:M-1
        q=q+Z(2*M-k+1)*exp(-%i*k*x);
    end
endfunction
function y=f(x)
    //y=exp(cos(x));
    y=(1+0.9*cos(x))^-1;
endfunction
N=6;
Xint=(2*%pi/N)*[0:N-1];
Yint=f(Xint);
x=[0:0.01:2*%pi];
y=real(pN(x,Yint));
testimag=norm(imag(y),'inf')
clf
plot2d(x,f(x),style=1)
plot2d(x,real(pN(x,Yint)),style=5)
plot2d(Xint,Yint,style=-2)
xtitle('f et pN pour N=6','x','y')
legends(['f';'pN'],[1,2],opt='ul')
test=max(abs(f(x)-pN(x,Yint)))

halt()
    TM=[];
    Terreur=[];
    for M=2:40
        N=2*M;
        Xint=(2*%pi/N)*[0:N-1];
        Yint=f(Xint);
        erreur=max(abs(f(x)-pN(x,Yint)));
        TM=[TM M];
        Terreur=[Terreur erreur];
    end
clf
plot2d(TM,Terreur,logflag='nl')
xtitle('Convergence','M','erreur')
// une droite en échelle 'normal-logarithme' signifie  une convergence lineaire (ou géométrique) en C/q^M avec q>1
