/* gcc fft_fftb.c -L/home/tama/FFT/lib -lfft -lm	*/

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

#include <f2c.h>

int ezffti_(long int *n, float *wsave);
int ezfftf_(long int *n, float *r__, float *azero, float *a, float *b, float *wsave);
int ezfftb_(long int *n, float *r__, float *azero, float *a, float *b, float *wsave);

void fft_fftb( float *x1, int num );

#define MAXDATA		100000
#define MAXLINE		1024



int main( int argc, char *argv[] )
{
/* read data from   stdin      x0 y0
                               x1 y1
                               x2 y2
			       ....
   output 		       avr   0
                               a1    b1
			       .....
		a[k] = 2./n*r(i)*cos(k*(i-1)*2*pi/n)
		b[k] = 2./n*r(i)*sin(k*(i-1)*2*pi/n)

*/


	float x[MAXDATA], y1[MAXDATA], y2[MAXDATA];
	float cspc[MAXDATA], co[MAXDATA];
	float freq[MAXDATA], p1[MAXDATA], p2[MAXDATA];
	int kmax;

	long int num = 0 ;


	char buf[MAXLINE];
	int i;
	float dF ;


	while(   fgets(buf, MAXLINE, stdin) != NULL ) {
		sscanf(buf, "%f", &y1[num]);
		y2[num] = y1[num];
		num++;

		if(num >= MAXDATA-1) {
			fprintf( stderr, "Too much data!\n");
			exit(EXIT_FAILURE);
		}
		
	}

	if(num < 2) exit(EXIT_FAILURE);

	if( num%2 == 0 ) {
		kmax = num/2 -1 ;
	} else {
		kmax = (num-1)/2;
	}

	fft_fftb( y1, num);
	for( i=0; i<num; i++){
		fprintf(stdout, "%e %e\n", y1[i], y2[i]);
	}

}

void fft_fftb( float *x1, int num )
{
	/*  Cospectral density btwn x1 & x2 ==>cosp, normalized in cosp_n */

	float *wsave, *a1, *b1, avr1;
	int i, i1, j,  kmax, jmax;
	long int lnum;

	if(num%2 == 0) {
		kmax = num/2-1;
	} else {
		kmax =(num-1)/2;
	}

	wsave = (float *)malloc( sizeof(float)*(3*num+15) );
	a1 = (float *)malloc( sizeof(float)*kmax );
	b1 = (float *)malloc( sizeof(float)*kmax );
	if( wsave == NULL || a1 == NULL || b1 == NULL) {
		fprintf( stderr, "Can nor malloc in fft_fftb\n");
		exit( EXIT_FAILURE);
	}


	lnum = num;
	ezffti_( &lnum, wsave);
	ezfftf_( &lnum, x1, &avr1, a1, b1, wsave);
	ezfftb_( &lnum, x1, &avr1, a1, b1, wsave);

	


}