/* fft_test.f -- translated by f2c (version of 23 April 1993  18:34:30).
   You must link the resulting object file with the libraries:
	-lf2c -lm   (in that order)
*/

#include "f2c.h"

/* Table of constant values */

static integer c__9 = 9;
static integer c__1 = 1;
static integer c__3 = 3;
static integer c_n2 = -2;
static integer c__2 = 2;
static integer c__4 = 4;

/* Main program */ MAIN__()
{
    /* System generated locals */
    integer i__1, i__2, i__3, i__4;
    doublereal d__1;
    complex q__1;

    /* Builtin functions */
    double cos(), sin();
    integer s_wsle(), do_lio(), e_wsle();
    /* Subroutine */ int s_stop();

    /* Local variables */
    extern /* Subroutine */ int harm_();
    static integer isig, iord, ierr;
    extern /* Subroutine */ int c2fft_();
    static complex a[65536];
    static integer j, k, m[3];
    static complex s[1536];
    extern doublereal dtime_();
    static integer iferr, ifset;
    static real timer[2], aa[131072];
    static integer id, nm, is[512], nn, wm[1044], wn[1044], ntimes;
    extern /* Subroutine */ int tprint_(), fft_();
    static real err;

    /* Fortran I/O blocks */
    static cilist io___19 = { 0, 6, 0, 0, 0 };
    static cilist io___20 = { 0, 6, 0, 0, 0 };
    static cilist io___22 = { 0, 6, 0, 0, 0 };
    static cilist io___24 = { 0, 6, 0, 0, 0 };
    static cilist io___25 = { 0, 6, 0, 0, 0 };
    static cilist io___26 = { 0, 6, 0, 0, 0 };


/* 	test the execution speed of FFT */
/* 	Note: dtime is a bsd call used to get the timings */
    ntimes = 10;
/* 	256x256 sizes for harm */
    m[0] = 8;
    m[1] = 8;
    m[2] = 0;
/* 	setup tables */
    ifset = 0;
    for (j = 1; j <= 65536; ++j) {
	i__1 = j - 1;
	d__1 = cos(j * (float).01) * .1 * sin(j * (float).1) * sin(j * (float)
		.0333);
	a[i__1].r = d__1, a[i__1].i = (float)0.;
    }
    harm_(a, m, is, s, &ifset, &iferr);
/* 	initialize timer */
    err = dtime_(timer);
    i__1 = ntimes;
    for (j = 1; j <= i__1; ++j) {
	ifset = -2;
	harm_(a, m, is, s, &ifset, &iferr);
	ifset = 2;
	harm_(a, m, is, s, &ifset, &iferr);
    }
    err = dtime_(timer);
    tprint_(timer, "harm, 256x256", &ntimes, 13L);
/* 	Now do the same for C2FFT */
    for (j = 1; j <= 65536; ++j) {
	i__1 = j - 1;
	d__1 = cos(j * (float).01) * .1 * sin(j * (float).1) * sin(j * (float)
		.0333);
	a[i__1].r = d__1, a[i__1].i = (float)0.;
    }
    id = 256;
    nm = 256;
    nn = 256;
    isig = 0;
    iord = 1;
    c2fft_(a, &id, &nm, &nn, wm, wn, &isig, &iord, is, &ierr);
    if (ierr != 0) {
	s_wsle(&io___19);
	do_lio(&c__9, &c__1, "Error Code ", 11L);
	do_lio(&c__3, &c__1, (char *)&ierr, (ftnlen)sizeof(integer));
	do_lio(&c__9, &c__1, " in C2FFT", 9L);
	e_wsle();
	s_stop("", 0L);
    }
/* 	initialize timer */
    err = dtime_(timer);
    i__1 = ntimes;
    for (j = 1; j <= i__1; ++j) {
	isig = -1;
	c2fft_(a, &id, &nm, &nn, wm, wn, &isig, &iord, is, &ierr);
	if (ierr != 0) {
	    s_wsle(&io___20);
	    do_lio(&c__9, &c__1, "Error Code ", 11L);
	    do_lio(&c__3, &c__1, (char *)&ierr, (ftnlen)sizeof(integer));
	    do_lio(&c__9, &c__1, " in C2FFT", 9L);
	    e_wsle();
	    s_stop("", 0L);
	}
	isig = 1;
/* 	normalize */
	i__2 = nm * nn;
	for (k = 1; k <= i__2; ++k) {
	    i__3 = k - 1;
	    i__4 = k - 1;
	    d__1 = (doublereal) (nm * nn);
	    q__1.r = a[i__4].r / d__1, q__1.i = a[i__4].i / d__1;
	    a[i__3].r = q__1.r, a[i__3].i = q__1.i;
	}
    }
    c2fft_(a, &id, &nm, &nn, wm, wn, &isig, &iord, is, &ierr);
    if (ierr != 0) {
	s_wsle(&io___22);
	do_lio(&c__9, &c__1, "Error Code ", 11L);
	do_lio(&c__3, &c__1, (char *)&ierr, (ftnlen)sizeof(integer));
	do_lio(&c__9, &c__1, " in C2FFT", 9L);
	e_wsle();
	s_stop("", 0L);
    }
    err = dtime_(timer);
    tprint_(timer, "C2FFT 256x256", &ntimes, 13L);
/* 	setup tables */
    ifset = 0;
    for (j = 1; j <= 65536; j += 2) {
	aa[j - 1] = cos(j * (float).01) * .1 * sin(j * (float).1) * sin(j * (
		float).0333);
	aa[j] = (float)0.;
    }
    id = 256;
    nm = 256;
    nn = 256;
/* 	initialize timer */
    err = dtime_(timer);
    i__1 = ntimes;
    for (j = 1; j <= i__1; ++j) {
/* 	back */
	i__2 = nn * nm;
	fft_(aa, &aa[1], &i__2, &nn, &nn, &c_n2);
	i__2 = nn * nm;
	i__3 = nn * nm;
	fft_(aa, &aa[1], &i__2, &nn, &i__3, &c_n2);
/* 	forward */
	i__2 = nm * nn << 1;
	for (k = 1; k <= i__2; ++k) {
	    aa[k - 1] /= nm * nn;
	}
    }
    i__1 = nn * nm;
    fft_(aa, &aa[1], &i__1, &nn, &nn, &c__2);
    i__1 = nn * nm;
    i__2 = nn * nm;
    fft_(aa, &aa[1], &i__1, &nn, &i__2, &c__2);
    err = dtime_(timer);
    tprint_(timer, "fft 256x256", &ntimes, 11L);
/* 	Now repeat with a size of 250x250 */
    for (j = 1; j <= 65536; ++j) {
	i__1 = j - 1;
	d__1 = cos(j * (float).01) * .1 * sin(j * (float).1) * sin(j * (float)
		.0333);
	a[i__1].r = d__1, a[i__1].i = (float)0.;
    }
    id = 250;
    nm = 250;
    nn = 250;
    isig = 0;
    iord = 1;
    c2fft_(a, &id, &nm, &nn, wm, wn, &isig, &iord, is, &ierr);
    if (ierr != 0) {
	s_wsle(&io___24);
	do_lio(&c__9, &c__1, "Error Code ", 11L);
	do_lio(&c__3, &c__1, (char *)&ierr, (ftnlen)sizeof(integer));
	do_lio(&c__9, &c__1, " in C2FFT", 9L);
	e_wsle();
	s_stop("", 0L);
    }
/* 	initialize timer */
    err = dtime_(timer);
    i__1 = ntimes;
    for (j = 1; j <= i__1; ++j) {
	isig = -1;
	c2fft_(a, &id, &nm, &nn, wm, wn, &isig, &iord, is, &ierr);
	if (ierr != 0) {
	    s_wsle(&io___25);
	    do_lio(&c__9, &c__1, "Error Code ", 11L);
	    do_lio(&c__3, &c__1, (char *)&ierr, (ftnlen)sizeof(integer));
	    do_lio(&c__9, &c__1, " in C2FFT", 9L);
	    e_wsle();
	    s_stop("", 0L);
	}
	isig = 1;
/* 	normalize */
	i__2 = nm * nn;
	for (k = 1; k <= i__2; ++k) {
	    i__3 = k - 1;
	    i__4 = k - 1;
	    d__1 = (doublereal) (nm * nn);
	    q__1.r = a[i__4].r / d__1, q__1.i = a[i__4].i / d__1;
	    a[i__3].r = q__1.r, a[i__3].i = q__1.i;
	}
    }
    c2fft_(a, &id, &nm, &nn, wm, wn, &isig, &iord, is, &ierr);
    if (ierr != 0) {
	s_wsle(&io___26);
	do_lio(&c__9, &c__1, "Error Code ", 11L);
	do_lio(&c__3, &c__1, (char *)&ierr, (ftnlen)sizeof(integer));
	do_lio(&c__9, &c__1, " in C2FFT", 9L);
	e_wsle();
	s_stop("", 0L);
    }
    err = dtime_(timer);
    tprint_(timer, "C2FFT 250x250", &ntimes, 13L);
/* 	setup tables */
    ifset = 0;
    for (j = 1; j <= 65536; j += 2) {
	aa[j - 1] = cos(j * (float).01) * .1 * sin(j * (float).1) * sin(j * (
		float).0333);
	aa[j] = (float)0.;
    }
    nm = 250;
    nn = 250;
/* 	initialize timer */
    err = dtime_(timer);
    i__1 = ntimes;
    for (j = 1; j <= i__1; ++j) {
/* 	back */
	i__2 = nn * nm;
	fft_(aa, &aa[1], &i__2, &nn, &nn, &c_n2);
	i__2 = nn * nm;
	i__3 = nn * nm;
	fft_(aa, &aa[1], &i__2, &nn, &i__3, &c_n2);
/* 	normalize */
	i__2 = nm * nn << 1;
	for (k = 1; k <= i__2; ++k) {
	    aa[k - 1] /= nm * nn;
	}
/* 	forward */
	i__2 = nn * nm;
	fft_(aa, &aa[1], &i__2, &nn, &nn, &c__2);
	i__2 = nn * nm;
	i__3 = nn * nm;
	fft_(aa, &aa[1], &i__2, &nn, &i__3, &c__2);
    }
    err = dtime_(timer);
    tprint_(timer, "fft 250x250", &ntimes, 11L);
} /* MAIN__ */

/* Subroutine */ int tprint_(timer, label, ntimes, label_len)
real *timer;
char *label;
integer *ntimes;
ftnlen label_len;
{
    /* Builtin functions */
    integer s_wsle(), e_wsle(), do_lio();

    /* Fortran I/O blocks */
    static cilist io___27 = { 0, 6, 0, 0, 0 };
    static cilist io___28 = { 0, 6, 0, 0, 0 };
    static cilist io___29 = { 0, 6, 0, 0, 0 };
    static cilist io___30 = { 0, 6, 0, 0, 0 };
    static cilist io___31 = { 0, 6, 0, 0, 0 };
    static cilist io___32 = { 0, 6, 0, 0, 0 };


    /* Parameter adjustments */
    --timer;

    /* Function Body */
    s_wsle(&io___27);
    e_wsle();
    s_wsle(&io___28);
    do_lio(&c__9, &c__1, "Timing for routine ", 19L);
    do_lio(&c__9, &c__1, label, label_len);
    e_wsle();
    s_wsle(&io___29);
    do_lio(&c__9, &c__1, "Total User time:", 16L);
    do_lio(&c__4, &c__1, (char *)&timer[1], (ftnlen)sizeof(real));
    e_wsle();
    s_wsle(&io___30);
    do_lio(&c__9, &c__1, "Total System time:", 18L);
    do_lio(&c__4, &c__1, (char *)&timer[2], (ftnlen)sizeof(real));
    e_wsle();
    timer[1] /= (doublereal) (*ntimes << 1);
    timer[2] /= (doublereal) (*ntimes << 1);
    s_wsle(&io___31);
    do_lio(&c__9, &c__1, "User time per FFT:", 18L);
    do_lio(&c__4, &c__1, (char *)&timer[1], (ftnlen)sizeof(real));
    e_wsle();
    s_wsle(&io___32);
    do_lio(&c__9, &c__1, "System time per FFT:", 20L);
    do_lio(&c__4, &c__1, (char *)&timer[2], (ftnlen)sizeof(real));
    e_wsle();
    return 0;
} /* tprint_ */

/* Main program alias */ int fft_test__ () { MAIN__ (); }