Index: /trunk/Ohana/src/opihi/cmd.data/Makefile
===================================================================
--- /trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 16098)
+++ /trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 16099)
@@ -42,4 +42,5 @@
 $(SRC)/erase.$(ARCH).o		\
 $(SRC)/extract.$(ARCH).o	\
+$(SRC)/fft1d.old.$(ARCH).o		\
 $(SRC)/fft1d.$(ARCH).o		\
 $(SRC)/fft2d.$(ARCH).o		\
Index: /trunk/Ohana/src/opihi/cmd.data/fft1d.c
===================================================================
--- /trunk/Ohana/src/opihi/cmd.data/fft1d.c	(revision 16098)
+++ /trunk/Ohana/src/opihi/cmd.data/fft1d.c	(revision 16099)
@@ -50,5 +50,5 @@
   }    
 
-  fft (Ore[0].elements, Oim[0].elements, Npix, Nbit); 
+  fft (Ore[0].elements, Oim[0].elements, Npix, Nbit, TRUE); 
   
   return (TRUE);
Index: /trunk/Ohana/src/opihi/cmd.data/fft1d.old.c
===================================================================
--- /trunk/Ohana/src/opihi/cmd.data/fft1d.old.c	(revision 16099)
+++ /trunk/Ohana/src/opihi/cmd.data/fft1d.old.c	(revision 16099)
@@ -0,0 +1,70 @@
+# include "data.h"
+
+int fft1dold (int argc, char **argv) {
+  
+  int i, Npix, ZeroImaginary;
+  float *t1, *t2, *temp;
+  Vector *Ire, *Iim, *Ore, *Oim;
+
+  if ((argc != 6) || (strcmp (argv[3], "to"))) {
+    gprint (GP_ERR, "USAGE: fft1d (real) (imag) to (real) (imag)\n");
+    return (FALSE);
+  }
+
+  Iim = NULL;
+  ZeroImaginary = TRUE;
+  if (strcmp (argv[2], "0")) {
+    ZeroImaginary = FALSE;
+    if ((Iim = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  }    
+  if ((Ire = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Ore = SelectVector (argv[4], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Oim = SelectVector (argv[5], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  Npix = Ire[0].Nelements;
+  if (!ZeroImaginary && (Npix != Iim[0].Nelements)) {
+    gprint (GP_ERR, "vector mismatch in size\n");
+    return (FALSE);
+  }
+
+  if (!IsBinaryOld (Npix)) {
+    gprint (GP_ERR, "Npix is not a binary number!\n");
+    return (FALSE);
+  }
+  
+  ALLOCATE (temp, float, 2*Npix);
+  if (ZeroImaginary) {
+    t1 = Ire[0].elements;
+    for (i = 0; i < Npix; i++, t1++) {
+      temp[2*i  ] = *t1;
+      temp[2*i+1] = 0;
+    }
+  } else {
+    t1 = Ire[0].elements;
+    t2 = Iim[0].elements;
+    for (i = 0; i < Npix; i++, t1++, t2++) {
+      temp[2*i  ] = *t1;
+      temp[2*i+1] = *t2;
+    }
+  }    
+    
+  fftold (temp, Npix, 1); 
+
+  Ore[0].Nelements = Npix;
+  Oim[0].Nelements = Npix;
+  REALLOCATE (Ore[0].elements, float, Npix);
+  REALLOCATE (Oim[0].elements, float, Npix);
+ 
+  t1 = Ore[0].elements;
+  t2 = Oim[0].elements;
+  for (i = 0; i < Npix; i++, t1++, t2++) {
+    *t1 = temp[2*i  ] / Npix;
+    *t2 = temp[2*i+1] / Npix;
+  }    
+  
+  free (temp);
+  
+  return (TRUE);
+}
+
+  
Index: /trunk/Ohana/src/opihi/cmd.data/fft2d.c
===================================================================
--- /trunk/Ohana/src/opihi/cmd.data/fft2d.c	(revision 16098)
+++ /trunk/Ohana/src/opihi/cmd.data/fft2d.c	(revision 16099)
@@ -42,5 +42,5 @@
   Ny = Ire[0].header.Naxis[1];
   Naxis[0] = Ny; Naxis[1] = Nx;
-  if (!IsBinary (Npix)) {
+  if (!IsBinaryOld (Npix)) {
     gprint (GP_ERR, "dimensions are not binary!\n");
     return (FALSE);
Index: /trunk/Ohana/src/opihi/cmd.data/init.c
===================================================================
--- /trunk/Ohana/src/opihi/cmd.data/init.c	(revision 16098)
+++ /trunk/Ohana/src/opihi/cmd.data/init.c	(revision 16099)
@@ -26,4 +26,5 @@
 int erase            PROTO((int, char **));
 int extract          PROTO((int, char **));
+int fft1dold         PROTO((int, char **));
 int fft1d            PROTO((int, char **));
 int fft2d            PROTO((int, char **));
@@ -145,4 +146,5 @@
   {"erase",        erase,	     "erase objects on an overlay"},
   {"extract",      extract,	     "extract a portion of a buffer into another buffer"},
+  {"fft1dold",     fft1dold,	     "fft on the pixel-stream in an image"},
   {"fft1d",        fft1d,	     "fft on the pixel-stream in an image"},
   {"fft2d",        fft2d,	     "fft on an image"},
Index: /trunk/Ohana/src/opihi/cmd.data/test/fft1d.sh
===================================================================
--- /trunk/Ohana/src/opihi/cmd.data/test/fft1d.sh	(revision 16099)
+++ /trunk/Ohana/src/opihi/cmd.data/test/fft1d.sh	(revision 16099)
@@ -0,0 +1,25 @@
+
+list tests
+ test1
+ test2
+ test3
+end
+
+macro test1
+ $PASS = 1
+ break -auto off
+
+ create t 0 4096 1.0
+ set f = dsin(20*t)
+
+ fft1d f 0 to Frn Fin
+
+ fft1dold f 0 to Fro Fio
+
+ clear
+ section a 0.0 0.0 1.0 0.5
+ lim t Fro; box; plot -pt 7 -c blue t Fro; plot -pt 2 -c red t Frn
+
+ section b 0.0 0.5 1.0 0.5
+ lim t Fio; box; plot -pt 7 -c blue t Fio; plot -pt 2 -c red t Fin
+end
Index: /trunk/Ohana/src/opihi/include/data.h
===================================================================
--- /trunk/Ohana/src/opihi/include/data.h	(revision 16098)
+++ /trunk/Ohana/src/opihi/include/data.h	(revision 16099)
@@ -84,8 +84,12 @@
 
 /* in fft.c */
-void fft (float *dataRe, float *dataIm, int N, int Nbit);
+void fft (float *dataRe, float *dataIm, int N, int Nbit, int forward);
 void fftN (float *data, int *nn, int ndim, int isign);
-int IsBinary (int N);
-void fourn (float *data, int *nn, int ndim, int isign);
+int IsBinary (int N, int *Nbit);
+
+// fftold.c
+void fftold (float *Data, int N, int isign);
+void fftNold (float *data, int *nn, int ndim, int isign);
+int IsBinaryOld (int N);
 
 /* in gaussj.c */
Index: /trunk/Ohana/src/opihi/lib.data/Makefile
===================================================================
--- /trunk/Ohana/src/opihi/lib.data/Makefile	(revision 16098)
+++ /trunk/Ohana/src/opihi/lib.data/Makefile	(revision 16099)
@@ -19,4 +19,5 @@
 $(SDIR)/page.$(ARCH).o                  \
 $(SDIR)/fft.$(ARCH).o			\
+$(SDIR)/fftold.$(ARCH).o			\
 $(SDIR)/svdcmp.$(ARCH).o		\
 $(SDIR)/convert.$(ARCH).o		\
Index: /trunk/Ohana/src/opihi/lib.data/fft.c
===================================================================
--- /trunk/Ohana/src/opihi/lib.data/fft.c	(revision 16098)
+++ /trunk/Ohana/src/opihi/lib.data/fft.c	(revision 16099)
@@ -8,8 +8,9 @@
 }
 
-void fft (double *x, double *y, int n, int Nbit) {
+void fft (float *x, float *y, int n, int Nbit, int forward) {
 
   int i,j,k,n1,n2;
-  double c,s,e,a,t1,t2;        
+  float c,s,e,a,t1,t2;        
+  double factor;
          
   // bit-reverse
@@ -19,8 +20,8 @@
     n1 = n2;
     while ( j >= n1 ) {
-      j = j - n1;
-      n1 = n1/2;
+      j -= n1;
+      n1 /= 2;
     }
-    j = j + n1;
+    j += n1;
                
     if (i < j) {
@@ -33,13 +34,18 @@
     }
   }
-  
                                           
   n1 = 0; /* FFT */
   n2 = 1;
                                              
+  if (forward) {
+    factor = +2.0*M_PI;
+  } else {
+    factor = -2.0*M_PI;
+  }
+
   for (i=0; i < Nbit; i++) {
     n1 = n2;
     n2 = n2 + n2;
-    e = -6.283185307179586/n2;
+    e = factor/n2;
     a = 0.0;
                                              
@@ -60,21 +66,27 @@
   }
                                       
+  // re-normalize
+  for (i = 0; i < n; i++) {
+    x[i] /= n;
+    y[i] /= n;
+  }
+
   return;
 }                          
 
-int IsBinary (int N, int *nbin) {
+// check that a number is binary (2^Nbit).  returns int(log_2(N)) in Nbit
+int IsBinary (int N, int *Nbit) {
 
-  // check that a number is binary (2^nbit)
+  int i, Nset;
 
-  Nbit = 0;
+  *Nbit = 0;
   Nset = 0;
   for (i = 0; i < 8*sizeof(int); i++) {
     if (N & 0x01) {
       Nset ++;
-      Nbin = i;
+      *Nbit = i;
     }
     N >>= 1;
   }
-  *nbin = Nbit;
   if (Nset > 1) return (FALSE);
 }
@@ -92,6 +104,6 @@
 /* m: n = 2**m                                            */
 /*   input/output                                         */
-/* x: double array of length n with real part of data     */
-/* y: double array of length n with imag part of data     */
+/* x: float array of length n with real part of data     */
+/* y: float array of length n with imag part of data     */
 /*                                                        */
 /*   Permission to copy and use this program is granted   */
Index: /trunk/Ohana/src/opihi/lib.data/fftold.c
===================================================================
--- /trunk/Ohana/src/opihi/lib.data/fftold.c	(revision 16098)
+++ /trunk/Ohana/src/opihi/lib.data/fftold.c	(revision 16099)
@@ -3,5 +3,5 @@
 #define FSWAP(a,b) tempr=(a);(a)=(b);(b)=tempr
 
-void fft (float *Data, int N, int isign) {
+void fftold (float *Data, int N, int isign) {
 
   int n,mmax,m,j,istep,i;
@@ -51,30 +51,33 @@
 }
 
-void fftold (float *Data, int N, int isign) {
-
-  int n,mmax,m,j,istep,i;
+void fftold_Nindices (float *Data, int N, int isign) {
+
+  int n,ifp1,m,j,ipf2,i;
   double wtemp,wr,wpr,wpi,wi,theta;
   float tempr,tempi, *data;
 
-  data = Data - 1;
+  data = Data;
   n = N << 1;
-  j = 1;
-
-  for (i = 1; i < n; i+=2) {
-    if (j > i) {
-      FSWAP (data[j], data[i]);
-      FSWAP (data[j+1], data[i+1]);
-    }
-    m = n >> 1;
-    while (m >= 2 && j > m) {
-      j -= m;
-      m >>= 1;
-    }
-    j += m;
-  }
-  mmax = 2;
-  while (n > mmax) {
-    istep = 2*mmax;
-    theta = 6.28318530717959 / (isign*mmax);
+
+  ip1 = 2;
+  ip2 = n;
+  i2rev = 0;
+  for (i2 = 0; i2 < ip2; i2 += ip1) {
+    if (i2 < i2ref) {
+      FSWAP (data[i2rev], data[i2]);
+      FSWAP (data[i2rev+1], data[i2+1]);
+    }
+    ibit = ip2 >> 1;
+    while (ibit >= ip1 && i2rev >= ibit) {
+      i2rev -= ibit;
+      ibit >>= 1;
+    }
+    i2rev += ibit;
+  }
+
+  ifp1 = 2;
+  while (ifp1 < ip2) {
+    ifp2 = ifp1 << 1;
+    theta = isign*6.28318530717959 / (2*ifp1/ip1);
     wtemp = sin(0.5*theta);
     wpr = -2.0*wtemp*wtemp;
@@ -82,23 +85,26 @@
     wr = 1.0;
     wi = 0.0;
-    for (m = 1; m < mmax; m+=2) {
-      for (i = m; i <= n; i+=istep) {
-	j = i + mmax;
-	tempr = wr*data[j] - wi*data[j+1];
-	tempi = wr*data[j+1] + wi*data[j];
-	data[j] = data[i] - tempr;
-	data[j+1] = data[i+1] - tempi;
-	data[i] += tempr;
-	data[i+1] += tempi;
+    for (i3 = 0; i3 < ifp1; i3+=ip1) {
+      // in the N-d version, this loop is broken into two loops
+      for (i1 = i3; i1 < n; i1+=ipf2) {
+	k1 = i1;
+	k2 = k1 + ifp1;
+	tempr = wr*data[k2] - wi*data[k2+1];
+	tempi = wr*data[k2+1] + wi*data[k2];
+	data[k2] = data[k1] - tempr;
+	data[k2+1] = data[k1+1] - tempi;
+	data[k1] += tempr;
+	data[k1+1] += tempi;
       }
       wr = (wtemp = wr)*wpr - wi*wpi+wr;
       wi = wi*wpr + wtemp*wpi + wi;
     }
-    mmax = istep;
-  }
+    ifp1 = ifp2;
+  }
+
 }
 
 /* convert indices to zero reference */
-void fftN (float *data, int *nn, int ndim, int isign) {
+void fftNold (float *data, int *nn, int ndim, int isign) {
 
   int i1,i2,i3,i2rev,i3rev,ip1,ip2,ip3,ifp1,ifp2;
@@ -138,5 +144,5 @@
     while (ifp1 < ip2) {
       ifp2 = ifp1 << 1;
-      theta = isign*6.28318530717959/(ifp2/ip1);
+      theta = isign*6.28318530717959/(2*ifp1/ip1);
       wtemp = sin(0.5*theta);
       wpr = -2.0*wtemp*wtemp;
@@ -182,5 +188,5 @@
 */ 
 
-int IsBinary (int N) {
+int IsBinaryOld (int N) {
 
   int i, nbit;
@@ -199,72 +205,2 @@
 
 }  
-
-#define FSWAP(a,b) tempr=(a);(a)=(b);(b)=tempr
-
-void fourn (float *data, int *nn, int ndim, int isign) {
-
-  int i1,i2,i3,i2rev,i3rev,ip1,ip2,ip3,ifp1,ifp2;
-  int ibit,idim,k1,k2,n,nprev,nrem,ntot;
-  float tempi,tempr;
-  double theta,wi,wpi,wpr,wr,wtemp;
-
-  ntot = 1;
-  for (idim = 1; idim <= ndim; idim++)
-    ntot *= nn[idim];
-  nprev = 1;
-  for (idim = ndim; idim >= 1; idim--) {
-    n = nn[idim];
-    nrem = ntot/(n*nprev);
-    ip1 = nprev << 1;
-    ip2 = ip1*n;
-    ip3 = ip2*nrem;
-    i2rev = 1;
-    for (i2 = 1; i2 <= ip2; i2+=ip1) {
-      if (i2 < i2rev) {
-	for (i1 = i2; i1 <= i2+ip1-2; i1+=2) {
-	  for (i3 = i1; i3 <= ip3;i3+=ip2) {
-	    i3rev = i2rev+i3-i2;
-	    FSWAP(data[i3],data[i3rev]);
-	    FSWAP(data[i3+1],data[i3rev+1]);
-	  }
-	}
-      }
-      ibit = ip2 >> 1;
-      while (ibit >= ip1 && i2rev > ibit) {
-	i2rev -=  ibit;
-	ibit >>= 1;
-      }
-      i2rev += ibit;
-    }
-    ifp1 = ip1;
-    while (ifp1 < ip2) {
-      ifp2 = ifp1 << 1;
-      theta = isign*6.28318530717959/(ifp2/ip1);
-      wtemp = sin(0.5*theta);
-      wpr  =  -2.0*wtemp*wtemp;
-      wpi = sin(theta);
-      wr = 1.0;
-      wi = 0.0;
-      for (i3 = 1;i3<=ifp1;i3+=ip1) {
-	for (i1 = i3;i1<=i3+ip1-2;i1+=2) {
-	  for (i2 = i1;i2<=ip3;i2+=ifp2) {
-	    k1 = i2;
-	    k2 = k1+ifp1;
-	    tempr = wr*data[k2]-wi*data[k2+1];
-	    tempi = wr*data[k2+1]+wi*data[k2];
-	    data[k2] = data[k1]-tempr;
-	    data[k2+1] = data[k1+1]-tempi;
-	    data[k1] += tempr;
-	    data[k1+1] += tempi;
-	  }
-	}
-	wr = (wtemp = wr)*wpr-wi*wpi+wr;
-	wi = wi*wpr+wtemp*wpi+wi;
-      }
-      ifp1 = ifp2;
-    }
-    nprev *= n;
-  }
-}
-
-#undef FSWAP