N = 1000;
x = [0:N-1] * 2 * pi / N;
fx = (x > pi/2) & (x <= 3*pi/2);

c1 = (1/pi) * sum( fx .* exp(-i*x) * (2*pi/N) );
fx1 = real(c1) * cos(x) + imag(c1) * sin(x);

n = 10;
c = (1/pi) * sum( ((ones(n+1,1)*fx) .* exp([0:n]' * (-i*x)) * (2*pi/N))' );

y0 = 0.5*c(1)*cos(0*x);
y1 = real(c(2))*cos(1*x)  + imag(c(2))*sin(1*x);
y2 = real(c(3))*cos(2*x)  + imag(c(3))*sin(2*x);
y3 = real(c(4))*cos(3*x)  + imag(c(4))*sin(3*x);
y4 = real(c(5))*cos(4*x)  + imag(c(5))*sin(4*x);
y5 = real(c(6))*cos(5*x)  + imag(c(6))*sin(5*x);
y6 = real(c(7))*cos(6*x)  + imag(c(7))*sin(6*x);
y7 = real(c(8))*cos(7*x)  + imag(c(8))*sin(7*x);
y8 = real(c(9))*cos(8*x)  + imag(c(9))*sin(8*x);
y9 = real(c(10))*cos(9*x) + imag(c(10))*sin(9*x);

z0 = y0;
z1 = z0 + y1;
z2 = z1 + y2;
z3 = z2 + y3;
z4 = z3 + y4;
z5 = z4 + y5;
z6 = z5 + y6;
z7 = z6 + y7;
z8 = z7 + y8;
z9 = z8 + y9;
