[Scipy-svn] r3213 - in trunk/Lib/fftpack: . src

scipy-svn@scip... scipy-svn@scip...
Wed Aug 1 08:10:15 CDT 2007


Author: cdavid
Date: 2007-08-01 08:09:45 -0500 (Wed, 01 Aug 2007)
New Revision: 3213

Added:
   trunk/Lib/fftpack/src/drfft_djbfft.c
   trunk/Lib/fftpack/src/drfft_fftpack.c
   trunk/Lib/fftpack/src/drfft_fftw.c
   trunk/Lib/fftpack/src/drfft_fftw3.c
Modified:
   trunk/Lib/fftpack/setup.py
   trunk/Lib/fftpack/src/drfft.c
Log:
Clean fft code for real input

Modified: trunk/Lib/fftpack/setup.py
===================================================================
--- trunk/Lib/fftpack/setup.py	2007-08-01 11:39:07 UTC (rev 3212)
+++ trunk/Lib/fftpack/setup.py	2007-08-01 13:09:45 UTC (rev 3213)
@@ -31,7 +31,9 @@
         libraries=['dfftpack'],
         extra_info=[fft_opt_info, djbfft_info],
         depends=['src/zfft_djbfft.c', 'src/zfft_fftpack.c', 'src/zfft_fftw.c',
-            'src/zfft_fftw3.c', 'src/zfft_mkl.c'],
+            'src/zfft_fftw3.c', 'src/zfft_mkl.c', 
+            'src/drfft_djbfft.c', 'src/drfft_fftpack.c',
+            'src/drfft_fftw3.c', 'src/drfft_fftw.c'],
     )
 
     config.add_extension('convolve',

Modified: trunk/Lib/fftpack/src/drfft.c
===================================================================
--- trunk/Lib/fftpack/src/drfft.c	2007-08-01 11:39:07 UTC (rev 3212)
+++ trunk/Lib/fftpack/src/drfft.c	2007-08-01 13:09:45 UTC (rev 3213)
@@ -6,211 +6,66 @@
 
 #include "fftpack.h"
 
-/**************** DJBFFT *****************************/
-#ifdef WITH_DJBFFT
-GEN_CACHE(ddjbfft,(int n)
-	  ,unsigned int* f;
-	   double* ptr;
-	  ,caches_ddjbfft[i].n==n
-	  ,caches_ddjbfft[id].f = (unsigned int*)malloc(sizeof(unsigned int)*(n));
-	   caches_ddjbfft[id].ptr = (double*)malloc(sizeof(double)*n);
-	   fftfreq_rtable(caches_ddjbfft[id].f,n);
-	  ,free(caches_ddjbfft[id].f);
-	   free(caches_ddjbfft[id].ptr);
-	  ,10)
-#endif
-
-#if defined WITH_FFTW3
-/**************** FFTW3 *****************************/
-GEN_CACHE(drfftw,(int n,int d,int flags)
-          ,int direction;
-          int flags;
-          fftw_plan plan;
-          double *ptr;
-          ,((caches_drfftw[i].n==n) && 
-            (caches_drfftw[i].direction==d) &&
-            (caches_drfftw[i].flags==flags))
-          ,caches_drfftw[id].direction = d;
-          caches_drfftw[id].flags = flags;
-          caches_drfftw[id].ptr = (double*)fftw_malloc(sizeof(double)*(n));
-          caches_drfftw[id].plan = fftw_plan_r2r_1d(n,caches_drfftw[id].ptr,
-                       caches_drfftw[id].ptr,(d>0?FFTW_R2HC:FFTW_HC2R),flags);
-          ,fftw_destroy_plan(caches_drfftw[id].plan);
-          fftw_free(caches_drfftw[id].ptr);
-          ,10)
-#elif defined WITH_FFTW
-/**************** FFTW2 *****************************/
-GEN_CACHE(drfftw,(int n,int d,int flags)
-	  ,int direction;
-	   int flags;
-	   rfftw_plan plan;
-	   double *ptr;
-	  ,((caches_drfftw[i].n==n) && 
-	   (caches_drfftw[i].direction==d) &&
-	   (caches_drfftw[i].flags==flags))
-	  ,caches_drfftw[id].direction = d;
-	   caches_drfftw[id].flags = flags;
-	   caches_drfftw[id].plan = rfftw_create_plan(n,
-		(d>0?FFTW_REAL_TO_COMPLEX:FFTW_COMPLEX_TO_REAL),flags);
-	   caches_drfftw[id].ptr = (double*)malloc(sizeof(double)*(n));
-	  ,rfftw_destroy_plan(caches_drfftw[id].plan);
-	   free(caches_drfftw[id].ptr);
-	  ,10)
-#else
-/**************** FFTPACK ZFFT **********************/
-extern void F_FUNC(dfftf,DFFTF)(int*,double*,double*);
-extern void F_FUNC(dfftb,DFFTB)(int*,double*,double*);
-extern void F_FUNC(dffti,DFFTI)(int*,double*);
-GEN_CACHE(dfftpack,(int n)
-	  ,double* wsave;
-	  ,(caches_dfftpack[i].n==n)
-	  ,caches_dfftpack[id].wsave = (double*)malloc(sizeof(double)*(2*n+15));
-	   F_FUNC(dffti,DFFTI)(&n,caches_dfftpack[id].wsave);
-	  ,free(caches_dfftpack[id].wsave);
-	  ,10)
-#endif
-
-extern void destroy_drfft_cache(void) {
-#ifdef WITH_DJBFFT
-  destroy_ddjbfft_caches();
-#endif
-#if defined(WITH_FFTW3) || defined(WITH_FFTW)
-  destroy_drfftw_caches();
-#else
-  destroy_dfftpack_caches();
-#endif
+/* The following macro convert private backend specific function to the public
+ * functions exported by the module  */
+#define GEN_PUBLIC_API(name) \
+void destroy_drfft_cache(void)\
+{\
+        destroy_dr##name##_caches();\
+}\
+\
+void drfft(double *inout, int n, \
+        int direction, int howmany, int normalize)\
+{\
+        drfft_##name(inout, n, direction, howmany, normalize);\
 }
 
-/**************** DRFFT function **********************/
+/* ************** Definition of backend specific functions ********* */
 
+/*
+ * To add a backend :
+ *  - create a file drfft_name.c, where you define a function drfft_name where
+ *  name is the name of your backend. If you do not use the GEN_CACHE macro,
+ *  you will need to define a function void destroy_drname_caches(void), 
+ *  which can do nothing
+ *  - in drfft.c, include the drfft_name.c file, and add the 3 following lines
+ *  just after it:
+ *  #ifndef WITH_DJBFFT
+ *      GEN_PUBLIC_API(name)
+ *  #endif
+ */
 
-extern void drfft(double *inout,
-		  int n,int direction,int howmany,int normalize) {
-  int i;
-  double *ptr = inout;
-#if defined(WITH_FFTW3) || defined(WITH_FFTW) || defined(WITH_DJBFFT)
-  double *ptrc = NULL;
-#endif
-#if defined WITH_FFTW3
-  fftw_plan plan = NULL;
+#ifdef WITH_FFTW3
+    #include "drfft_fftw3.c"
+    #ifndef WITH_DJBFFT
+        GEN_PUBLIC_API(fftw3)
+    #endif
 #elif defined WITH_FFTW
-  rfftw_plan plan = NULL;
-#else
-  double* wsave = NULL;
+    #include "drfft_fftw.c"
+    #ifndef WITH_DJBFFT
+        GEN_PUBLIC_API(fftw)
+    #endif
+#else /* Use fftpack by default */
+    #include "drfft_fftpack.c"
+    #ifndef WITH_DJBFFT
+        GEN_PUBLIC_API(fftpack)
+    #endif 
 #endif
-#ifdef WITH_DJBFFT
-  unsigned int *f = NULL;
-#endif
 
+/* 
+ * djbfft must be used at the end, because it needs another backend (defined
+ * above) for non 2^n * size 
+ */
 #ifdef WITH_DJBFFT
-  switch (n) {
-  case 2:;case 4:;case 8:;case 16:;case 32:;case 64:;case 128:;case 256:;
-  case 512:;case 1024:;case 2048:;case 4096:;case 8192:
-    i = get_cache_id_ddjbfft(n);
-    f = caches_ddjbfft[i].f;
-    ptrc = caches_ddjbfft[i].ptr;
-  }
-  if (f==NULL)
-#endif
-#ifdef WITH_FFTW3
+    #include "drfft_djbfft.c"
+    void destroy_drfft_cache(void)
     {
-      i = get_cache_id_drfftw(n,direction,FFTW_ESTIMATE);
-      plan = caches_drfftw[i].plan;
-      ptrc = caches_drfftw[i].ptr;
+        destroy_drdjbfft_caches();
+        drfft_def_destroy_cache();
     }
-#elif defined WITH_FFTW
+    void drfft(double *inout, int n, 
+            int direction, int howmany, int normalize)
     {
-      i = get_cache_id_drfftw(n,direction,FFTW_IN_PLACE|FFTW_ESTIMATE);
-      plan = caches_drfftw[i].plan;
-      ptrc = caches_drfftw[i].ptr;
+        drfft_djbfft(inout, n, direction, howmany, normalize);
     }
-#else
-    wsave = caches_dfftpack[get_cache_id_dfftpack(n)].wsave;
 #endif
-
-  switch (direction) {
-
-  case 1:
-    for (i=0;i<howmany;++i,ptr+=n) {
-#ifdef WITH_DJBFFT
-      if (f!=NULL) {
-	COPYSTD2DJB(ptr,ptrc,n);
-	switch (n) {
-#define TMPCASE(N) case N: fftr8_##N(ptrc); break
-	  TMPCASE(2);TMPCASE(4);TMPCASE(8);TMPCASE(16);TMPCASE(32);
-	  TMPCASE(64);TMPCASE(128);TMPCASE(256);TMPCASE(512);
-	  TMPCASE(1024);TMPCASE(2048);TMPCASE(4096);TMPCASE(8192);
-#undef TMPCASE
-	}
-	COPYDJB2STD(ptrc,ptr,f,n);
-      } else
-#endif
-#if defined WITH_FFTW3
-	{
-	  memcpy(ptrc,ptr,sizeof(double)*n);
-	  fftw_execute(plan);
-	  COPYRFFTW2STD(ptrc,ptr,n);
-	}
-#elif defined WITH_FFTW
-	{
-	  memcpy(ptrc,ptr,sizeof(double)*n);
-	  rfftw(plan,1,(fftw_real*)ptrc,1,1,NULL,1,1);
-	  COPYRFFTW2STD(ptrc,ptr,n);
-	}
-#else
-	F_FUNC(dfftf,DFFTF)(&n,ptr,wsave);
-#endif
-    }
-    break;
-
-  case -1:
-    for (i=0;i<howmany;++i,ptr+=n) {
-#ifdef WITH_DJBFFT
-      if (f!=NULL) {
-	COPYINVSTD2DJB(ptr,ptrc,normalize,f,n);
-	switch (n) {
-#define TMPCASE(N)case N:if(normalize)fftr8_scale##N(ptrc);fftr8_un##N(ptrc);break
-	  TMPCASE(2);TMPCASE(4);TMPCASE(8);TMPCASE(16);TMPCASE(32);
-	  TMPCASE(64);TMPCASE(128);TMPCASE(256);TMPCASE(512);
-	  TMPCASE(1024);TMPCASE(2048);TMPCASE(4096);TMPCASE(8192);
-#undef TMPCASE
-	}
-	COPYINVDJB2STD(ptrc,ptr,n);
-      } else
-#endif
-#if defined WITH_FFTW3
-	{
-	  COPYINVRFFTW2STD(ptr,ptrc,n);
-	  fftw_execute(plan);
-	  memcpy(ptr,ptrc,sizeof(double)*n);
-	}
-#elif defined WITH_FFTW
-	{
-	  COPYINVRFFTW2STD(ptr,ptrc,n);
-	  rfftw(plan,1,(fftw_real*)ptrc,1,1,NULL,1,1);
-	  memcpy(ptr,ptrc,sizeof(double)*n);
-	}
-#else
-	F_FUNC(dfftb,DFFTB)(&n,ptr,wsave);
-#endif
-    }
-    break;
-
-  default:
-    fprintf(stderr,"drfft: invalid direction=%d\n",direction);
-  }
-
-  if (normalize 
-#ifdef WITH_DJBFFT
-      && (f==NULL||direction==1)
-#endif
-      ) {
-    double d = 1.0/n;
-    ptr = inout;
-    for (i=n*howmany-1;i>=0;--i)
-      (*(ptr++)) *= d;
-  }
-}
-
-
-

Added: trunk/Lib/fftpack/src/drfft_djbfft.c
===================================================================
--- trunk/Lib/fftpack/src/drfft_djbfft.c	2007-08-01 11:39:07 UTC (rev 3212)
+++ trunk/Lib/fftpack/src/drfft_djbfft.c	2007-08-01 13:09:45 UTC (rev 3213)
@@ -0,0 +1,131 @@
+/*
+ * Last Change: Wed Aug 01 08:00 PM 2007 J
+ *
+ * Original code by Pearu Peterson.
+ */
+
+/*
+ * DJBFFT only implements size 2^N !
+ *
+ * drfft_def and drfft_def_destroy_cache are the functions used for size different
+ * than 2^N
+ */
+#ifdef WITH_FFTW3
+#define drfft_def drfft_fftw3
+#define drfft_def_destroy_cache destroy_drfftw3_caches
+#elif defined WITH_FFTW
+#define drfft_def drfft_fftw
+#define drfft_def_destroy_cache destroy_drfftw_caches
+#else
+#define drfft_def drfft_fftpack
+#define drfft_def_destroy_cache destroy_drfftpack_caches
+#endif
+
+GEN_CACHE(drdjbfft, (int n)
+	    , unsigned int *f;
+	    double *ptr;, 
+        caches_drdjbfft[i].n == n, 
+        caches_drdjbfft[id].f = (unsigned int *) malloc(sizeof(unsigned int) * (n));
+	    caches_drdjbfft[id].ptr = (double *) malloc(sizeof(double) * n);
+	    fftfreq_rtable(caches_drdjbfft[id].f, n);,
+	    free(caches_drdjbfft[id].f); 
+        free(caches_drdjbfft[id].ptr);, 
+        10)
+
+/**************** ZFFT function **********************/
+static void drfft_djbfft(double * inout,
+			 int n, int direction, int howmany, int normalize)
+{
+    int i;
+    double *ptr = inout;
+    double *ptrc = NULL;
+    unsigned int *f = NULL;
+
+    switch (n) {
+        case 2:;
+        case 4:;
+        case 8:;
+        case 16:;
+        case 32:;
+        case 64:;
+        case 128:;
+        case 256:;
+        case 512:;
+        case 1024:;
+        case 2048:;
+        case 4096:;
+        case 8192:
+        i = get_cache_id_drdjbfft(n);
+        f = caches_drdjbfft[i].f;
+        ptrc = caches_drdjbfft[i].ptr;
+    }
+    if (f == NULL) {
+        drfft_def(inout, n, direction, howmany, normalize);
+    }
+
+    switch (direction) {
+    case 1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            if (f != NULL) {
+                COPYSTD2DJB(ptr, ptrc, n);
+                switch (n) {
+#define TMPCASE(N) case N: fftr8_##N(ptrc); break
+                    TMPCASE(2);
+                    TMPCASE(4);
+                    TMPCASE(8);
+                    TMPCASE(16);
+                    TMPCASE(32);
+                    TMPCASE(64);
+                    TMPCASE(128);
+                    TMPCASE(256);
+                    TMPCASE(512);
+                    TMPCASE(1024);
+                    TMPCASE(2048);
+                    TMPCASE(4096);
+                    TMPCASE(8192);
+#undef TMPCASE
+                }
+                COPYDJB2STD(ptrc, ptr, f, n);
+            } 
+        }
+        break;
+
+    case -1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            if (f != NULL) {
+                COPYINVSTD2DJB(ptr, ptrc, normalize, f, n);
+                switch (n) {
+
+#define TMPCASE(N)case N:if(normalize)fftr8_scale##N(ptrc);fftr8_un##N(ptrc);break
+                    TMPCASE(2);
+                    TMPCASE(4);
+                    TMPCASE(8);
+                    TMPCASE(16);
+                    TMPCASE(32);
+                    TMPCASE(64);
+                    TMPCASE(128);
+                    TMPCASE(256);
+                    TMPCASE(512);
+                    TMPCASE(1024);
+                    TMPCASE(2048);
+                    TMPCASE(4096);
+                    TMPCASE(8192);
+#undef TMPCASE
+                }
+                COPYINVDJB2STD(ptrc, ptr, n);
+            } 
+        }
+        break;
+
+    default:
+        fprintf(stderr, "drfft: invalid direction=%d\n", direction);
+    }
+
+    if (normalize && f != NULL && direction == 1) {
+        double d = 1.0 / n;
+        ptr = inout;
+        for (i = n * howmany - 1; i >= 0; --i) {
+            (*(ptr++)) *= d;
+        }
+    }
+}

Added: trunk/Lib/fftpack/src/drfft_fftpack.c
===================================================================
--- trunk/Lib/fftpack/src/drfft_fftpack.c	2007-08-01 11:39:07 UTC (rev 3212)
+++ trunk/Lib/fftpack/src/drfft_fftpack.c	2007-08-01 13:09:45 UTC (rev 3213)
@@ -0,0 +1,54 @@
+/*
+ * Last Change: Wed Aug 01 07:00 PM 2007 J
+ *
+ * FFTPACK implementation
+ *
+ * Original code by Pearu Peterson.
+ */
+
+extern void F_FUNC(dfftf, DFFTF) (int *, double *, double *);
+extern void F_FUNC(dfftb, DFFTB) (int *, double *, double *);
+extern void F_FUNC(dffti, DFFTI) (int *, double *);
+GEN_CACHE(drfftpack, (int n)
+	  , double *wsave;
+	  , (caches_drfftpack[i].n == n)
+	  , caches_drfftpack[id].wsave =
+	  (double *) malloc(sizeof(double) * (2 * n + 15));
+	  F_FUNC(dffti, DFFTI) (&n, caches_drfftpack[id].wsave);
+	  , free(caches_drfftpack[id].wsave);
+	  , 10)
+
+static void drfft_fftpack(double *inout, int n, int direction, int howmany,
+			  int normalize)
+{
+    int i;
+    double *ptr = inout;
+    double *wsave = NULL;
+    wsave = caches_drfftpack[get_cache_id_drfftpack(n)].wsave;
+
+
+    switch (direction) {
+        case 1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            dfftf_(&n, ptr, wsave);
+        }
+        break;
+
+    case -1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            dfftb_(&n, ptr, wsave);
+        }
+        break;
+
+    default:
+        fprintf(stderr, "drfft: invalid direction=%d\n", direction);
+    }
+
+    if (normalize) {
+        double d = 1.0 / n;
+        ptr = inout;
+        for (i = n * howmany - 1; i >= 0; --i) {
+            (*(ptr++)) *= d;
+        }
+    }
+}

Added: trunk/Lib/fftpack/src/drfft_fftw.c
===================================================================
--- trunk/Lib/fftpack/src/drfft_fftw.c	2007-08-01 11:39:07 UTC (rev 3212)
+++ trunk/Lib/fftpack/src/drfft_fftw.c	2007-08-01 13:09:45 UTC (rev 3213)
@@ -0,0 +1,70 @@
+/*
+ * Last Change: Wed Aug 01 07:00 PM 2007 J
+ *
+ * FFTW2 implementation
+ *
+ * Original code by Pearu Peterson.
+ */
+
+GEN_CACHE(drfftw, (int n, int d, int flags)
+	  , int direction;
+	  int flags;
+	  rfftw_plan plan;
+	  double *ptr;, ((caches_drfftw[i].n == n) &&
+			 (caches_drfftw[i].direction == d) &&
+			 (caches_drfftw[i].flags == flags))
+	  , caches_drfftw[id].direction = d;
+	  caches_drfftw[id].flags = flags;
+	  caches_drfftw[id].plan = rfftw_create_plan(n,
+						     (d >
+						      0 ?
+						      FFTW_REAL_TO_COMPLEX
+						      :
+						      FFTW_COMPLEX_TO_REAL),
+						     flags);
+	  caches_drfftw[id].ptr =
+	  (double *) malloc(sizeof(double) * (n));,
+	  rfftw_destroy_plan(caches_drfftw[id].plan);
+	  free(caches_drfftw[id].ptr);, 10)
+
+static void drfft_fftw(double *inout, int n, int dir, int
+			howmany, int normalize)
+{
+    int i;
+    double *ptr = inout;
+    double *ptrc = NULL;
+    rfftw_plan plan = NULL;
+
+    i = get_cache_id_drfftw(n, dir, FFTW_IN_PLACE | FFTW_ESTIMATE);
+    plan = caches_drfftw[i].plan;
+    ptrc = caches_drfftw[i].ptr;
+
+    switch (dir) {
+    case 1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            memcpy(ptrc, ptr, sizeof(double) * n);
+            rfftw(plan, 1, (fftw_real *) ptrc, 1, 1, NULL, 1, 1);
+            COPYRFFTW2STD(ptrc, ptr, n);
+        }
+        break;
+
+    case -1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            COPYINVRFFTW2STD(ptr, ptrc, n);
+            rfftw(plan, 1, (fftw_real *) ptrc, 1, 1, NULL, 1, 1);
+            memcpy(ptr, ptrc, sizeof(double) * n);
+        }
+        break;
+
+    default:
+        fprintf(stderr, "drfft: invalid direction=%d\n", dir);
+    }
+
+    if (normalize) {
+        double d = 1.0 / n;
+        ptr = inout;
+        for (i = n * howmany - 1; i >= 0; --i) {
+            (*(ptr++)) *= d;
+        }
+    }
+}

Added: trunk/Lib/fftpack/src/drfft_fftw3.c
===================================================================
--- trunk/Lib/fftpack/src/drfft_fftw3.c	2007-08-01 11:39:07 UTC (rev 3212)
+++ trunk/Lib/fftpack/src/drfft_fftw3.c	2007-08-01 13:09:45 UTC (rev 3213)
@@ -0,0 +1,65 @@
+/*
+ * Last Change: Wed Aug 01 07:00 PM 2007 J
+ *
+ * FFTW3 implementation
+ *
+ * Original code by Pearu Peterson.
+ */
+
+GEN_CACHE(drfftw3, (int n, int d, int flags)
+	  , int direction;
+	  int flags;
+	  fftw_plan plan;
+	  double *ptr;, ((caches_drfftw3[i].n == n) &&
+			 (caches_drfftw3[i].direction == d) &&
+			 (caches_drfftw3[i].flags == flags))
+	  , caches_drfftw3[id].direction = d;
+	  caches_drfftw3[id].flags = flags;
+	  caches_drfftw3[id].ptr =
+	  (double *) fftw_malloc(sizeof(double) * (n));
+	  caches_drfftw3[id].plan =
+	  fftw_plan_r2r_1d(n, caches_drfftw3[id].ptr, caches_drfftw3[id].ptr,
+			   (d > 0 ? FFTW_R2HC : FFTW_HC2R), flags);,
+	  fftw_destroy_plan(caches_drfftw3[id].plan);
+	  fftw_free(caches_drfftw3[id].ptr);, 10)
+
+static void drfft_fftw3(double *inout, int n, int direction, int
+			howmany, int normalize)
+{
+    int i;
+    double *ptr = inout;
+
+    double *ptrc = NULL;
+    fftw_plan plan = NULL;
+
+    i = get_cache_id_drfftw3(n, direction, (1U << 6));
+    plan = caches_drfftw3[i].plan;
+    ptrc = caches_drfftw3[i].ptr;
+    switch (direction) {
+    case 1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            memcpy(ptrc, ptr, sizeof(double) * n);
+            fftw_execute(plan);
+            COPYRFFTW2STD(ptrc, ptr, n);
+        }
+        break;
+
+    case -1:
+        for (i = 0; i < howmany; ++i, ptr += n) {
+            COPYINVRFFTW2STD(ptr, ptrc, n);
+            fftw_execute(plan);
+            memcpy(ptr, ptrc, sizeof(double) * n);
+        }
+        break;
+    default:
+        fprintf(stderr, "drfft: invalid direction=%d\n", direction);
+    }
+
+    if (normalize) {
+        double d = 1.0 / n;
+        ptr = inout;
+        for (i = n * howmany - 1; i >= 0; --i) {
+            (*(ptr++)) *= d;
+        }
+    }
+}



More information about the Scipy-svn mailing list