/* 2010.11.25 A.Date 
  Fourier expantion; f(t) = |t|  see text p.70 
*/

#include	<stdio.h>
#include	<stdlib.h>
#include	<math.h>

int N = 20;
int SHOW_EVERY_N_STEP = 5;

double start_t = -3.0*M_PI;
double end_t =  3.0*M_PI;
double step = 0.01;


int main(int argc, char *argv[])
{
	int i;
	double pi = M_PI;
	double t, y;
	char filename[256] = "2010am2.dat";
	FILE *fp;

	if((fp = fopen(filename, "wt")) == NULL){
		printf("NOT open file (%s)\n", filename);
		exit(1);
	}       

	for ( t = start_t; t < end_t; t=t+step ){
	  y = pi/2.0;
	  fprintf(fp, "%.2lf \t",t);
	  for ( i = 1; i <= N; i++ ){
	    y = y - 4.0/pi/(2.0*i-1.0)/(2.0*i-1.0)*cos( (2.0*i-1.0)*t );
		if (i % SHOW_EVERY_N_STEP == 0){
	  		fprintf(fp, "%.5lf ",y);
		}
	  }
	  fprintf(fp, "\n");
	}

	fclose(fp);
}