#include <stdio.h>
#include <sys/time.h>		/* struct timeval */
#include <sys/resource.h>	/* struct rusage  */

/* time declarations */
struct rusage start_struct, stop_struct;
float print_times(struct rusage, struct rusage, int);
void print_vector(double *, int, int);
void finish(struct rusage, struct rusage, int, int, int);
void test_original (int, int);
void test_ct_spec(void (*)(), int, int);
void test_rt_spec(int, int);

#define ORIGINAL 0
#define CTSPEC   1
#define RTGEN    2
#define RTSPEC   3

extern void fft(int, double *, double *);
extern void fft_16();
extern void fft_32();
extern void fft_64();
extern void fft_128();

/*
** FFT tester
**      Loads random numbers into complex array, ffts it, and shows result.
**	Beats me what the output should be; but any two different ffts
**      should return the same result.
**
**      % tester output 10000
*/

#define MAX 128
#define REAL double

REAL init_real[140001], real[140001];
REAL init_imag[140001], imag[140001];

void (*spec)(REAL *, REAL *);

main(argc,argv)
int argc;
char **argv;
{
 int i,j, repetitions;
 FILE *o;
 char *output;

 if (argc < 2) repetitions = 1; else repetitions = atoi(argv[1]);

 printf ("fft_duhamel, %d repetitions\n\n", repetitions);

 for (i=0;i<MAX;i++)  
   {
     init_real[i]=random()%10;
     init_imag[i]=random()%10;
   }

 /* 16 points */ 
 test_original (16, repetitions);
 test_ct_spec(&fft_16, 16, repetitions);
 test_rt_spec(16, repetitions);

 /* 32 points */
 test_original (32, repetitions);
 test_ct_spec(&fft_32, 32, repetitions);
 test_rt_spec(32, repetitions);

 /* 64 points */
 test_original (64, repetitions);
 test_ct_spec(&fft_64, 64, repetitions);
 test_rt_spec(64, repetitions);

 /* 128 points */
 test_original (128, repetitions);
 test_ct_spec(&fft_128, 128, repetitions);
 test_rt_spec(128, repetitions);
}

void test_original (int elements, int repetitions)
{
  int i, j;

 for (j=0; j < 1000; j++) /* for cache */
   {
     for (i=0;i<elements;i++)
       {
	 real[i]=init_real[i];
	 imag[i]=init_imag[i];
       }
     fft(elements, real, imag);
   }
 getrusage(RUSAGE_SELF, &start_struct);
 for (j=0; j < repetitions; j++)
   {
     for (i=0;i<elements;i++)
       {
	 real[i]=init_real[i];
	 imag[i]=init_imag[i];
       }
     fft(elements, real, imag);
   }
 getrusage(RUSAGE_SELF, &stop_struct);
 finish(start_struct, stop_struct, repetitions, elements, ORIGINAL);
 for (i=0;i<elements;i++) {
     real[i]=0.0;
     imag[i]=0.0;
   }
}

void test_ct_spec(void (*f)(), int elements, int repetitions)
{
  int i, j;

 for (j=0; j < 1000; j++) /* for cache */
   {
     for (i=0;i<elements;i++)
       {
	 real[i]=init_real[i];
	 imag[i]=init_imag[i];
       }
     f(real, imag);
   }
 getrusage(RUSAGE_SELF, &start_struct);
 for (j=0; j < repetitions; j++)
   {
     for (i=0;i<elements;i++)
       {
	 real[i]=init_real[i];
	 imag[i]=init_imag[i];
       }
     f(real, imag);
   }
 getrusage(RUSAGE_SELF, &stop_struct);
 finish(start_struct, stop_struct, repetitions, elements, CTSPEC);
 for (i=0;i<elements;i++) {
     real[i]=0.0;
     imag[i]=0.0;
   }
}

void test_rt_spec(int elements, int repetitions)
{
  int i, j;

 for (j=0; j < 1000; j++) /* for cache */
   {
     spec = (void (*)(REAL *, REAL *))(rts_fft_1(elements));
   }
 getrusage(RUSAGE_SELF, &start_struct);
 for (j=0; j < repetitions; j++)
   {
     spec = (void (*)(REAL *, REAL *))(rts_fft_1(elements));
   }
 getrusage(RUSAGE_SELF, &stop_struct);
 finish(start_struct, stop_struct, repetitions, elements, RTGEN);

 for (j=0; j < 1000; j++) /* for cache */
   {
     for (i=0;i<elements;i++)
       {
	 real[i]=init_real[i];
	 imag[i]=init_imag[i];
       }
     spec(real, imag);
   }
 getrusage(RUSAGE_SELF, &start_struct);
 for (j=0; j < repetitions; j++)
   {
     for (i=0;i<elements;i++)
       {
	 real[i]=init_real[i];
	 imag[i]=init_imag[i];
       }
     spec(real, imag);
   }
 getrusage(RUSAGE_SELF, &stop_struct);
 finish(start_struct, stop_struct, repetitions, elements, RTSPEC);
}


void finish(struct rusage start_struct, 
       struct rusage stop_struct,
       int repetitions, int elements, int which)
{
  static float or_time, ct_time, gen_time, rt_time;

  switch (which)
    {
    case ORIGINAL:
      printf ("Original(%d)\n", elements);
      or_time = print_times(start_struct, stop_struct, repetitions);
      break;
    case CTSPEC:
      printf ("CT Spec(%d)\n", elements);
      ct_time = print_times(start_struct, stop_struct, repetitions);
      break;
    case RTGEN:
      printf ("RT Gen(%d)\n", elements);
      gen_time = print_times(start_struct, stop_struct, repetitions);
      break;
    case RTSPEC:
      printf ("RT Spec(%d)\n", elements);
      rt_time = print_times(start_struct, stop_struct, repetitions);
      break;
    }
  
  print_vector(real, 0, 10);
  print_vector(imag, 0, 10);
  printf("\n");

  if (which == RTSPEC) {
    printf ("CT speedup(%d) \t%f\n",  elements, (or_time/ct_time));
    printf ("RT speedup(%d) \t%f\n",  elements, (or_time/rt_time));
    printf ("RT brkeven(%d) \t%f\n\n",  elements, 
	     (gen_time/(or_time - rt_time)));
  }
}