[Scipy-svn] r4214 - in branches/refactor_fft/scipy/fftpack: . src
scipy-svn at scipy.org
scipy-svn at scipy.org
Mon May 5 03:24:59 EDT 2008
Author: cdavid
Date: 2008-05-05 02:24:48 -0500 (Mon, 05 May 2008)
New Revision: 4214
Added:
branches/refactor_fft/scipy/fftpack/src/zfft.cxx
branches/refactor_fft/scipy/fftpack/src/zfft_djbfft.cxx
branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.cxx
Removed:
branches/refactor_fft/scipy/fftpack/src/zfft.c
branches/refactor_fft/scipy/fftpack/src/zfft_djbfft.c
branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.c
Modified:
branches/refactor_fft/scipy/fftpack/SConstruct
branches/refactor_fft/scipy/fftpack/setup.py
branches/refactor_fft/scipy/fftpack/src/fftpack.h
Log:
Put zfft API + djbfft/fftw3 implementation into c++.
Modified: branches/refactor_fft/scipy/fftpack/SConstruct
===================================================================
--- branches/refactor_fft/scipy/fftpack/SConstruct 2008-05-05 06:52:28 UTC (rev 4213)
+++ branches/refactor_fft/scipy/fftpack/SConstruct 2008-05-05 07:24:48 UTC (rev 4214)
@@ -40,7 +40,7 @@
env.PrependUnique(LIBPATH = env['build_dir'])
# Build _fftpack
-src = ['src/zfft.c','src/drfft.c','src/zrfft.c', 'src/zfftnd.c', 'fftpack.pyf']
+src = ['src/zfft.cxx','src/drfft.c','src/zrfft.c', 'src/zfftnd.c', 'fftpack.pyf']
env.NumpyPythonExtension('_fftpack', src)
# Build convolve
Modified: branches/refactor_fft/scipy/fftpack/setup.py
===================================================================
--- branches/refactor_fft/scipy/fftpack/setup.py 2008-05-05 06:52:28 UTC (rev 4213)
+++ branches/refactor_fft/scipy/fftpack/setup.py 2008-05-05 07:24:48 UTC (rev 4214)
@@ -23,7 +23,7 @@
config.add_library('dfftpack',
sources=[join('dfftpack','*.f')])
- sources = ['fftpack.pyf','src/zfft.c','src/drfft.c','src/zrfft.c',
+ sources = ['fftpack.pyf','src/zfft.cxx','src/drfft.c','src/zrfft.c',
'src/zfftnd.c']
config.add_extension('_fftpack',
@@ -31,8 +31,8 @@
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/drfft_djbfft.c', 'src/drfft_fftpack.c',
+ 'src/zfft_fftw3.cxx', 'src/zfft_mkl.c',
+ 'src/drfft_djbfft.cxx', 'src/drfft_fftpack.c',
'src/drfft_fftw3.c', 'src/drfft_fftw.c',
'src/zfftnd_fftpack.c', 'src/zfftnd_fftw.c',
'src/zfftnd_fftw3.c', 'src/zfftnd_mkl.c',
Modified: branches/refactor_fft/scipy/fftpack/src/fftpack.h
===================================================================
--- branches/refactor_fft/scipy/fftpack/src/fftpack.h 2008-05-05 06:52:28 UTC (rev 4213)
+++ branches/refactor_fft/scipy/fftpack/src/fftpack.h 2008-05-05 07:24:48 UTC (rev 4214)
@@ -31,13 +31,19 @@
#endif
#ifdef SCIPY_DJBFFT_H
+#ifdef __cplusplus
+extern "C" {
+#endif
#define WITH_DJBFFT
#define complex8 complex_double
#define COMPLEX8_H
#include <fftfreq.h>
#include <fftc8.h>
#include <fftr8.h>
+#ifdef __cplusplus
+}
#endif
+#endif
#ifdef SCIPY_MKL_H
#define WITH_MKL
Deleted: branches/refactor_fft/scipy/fftpack/src/zfft.c
===================================================================
--- branches/refactor_fft/scipy/fftpack/src/zfft.c 2008-05-05 06:52:28 UTC (rev 4213)
+++ branches/refactor_fft/scipy/fftpack/src/zfft.c 2008-05-05 07:24:48 UTC (rev 4214)
@@ -1,76 +0,0 @@
-/*
- Interface to various FFT libraries.
- Double complex FFT and IFFT.
- Author: Pearu Peterson, August 2002
- */
-
-#include "fftpack.h"
-
-/* The following macro convert private backend specific function to the public
- * functions exported by the module */
-#define GEN_PUBLIC_API(name) \
-void destroy_zfft_cache(void)\
-{\
- destroy_z##name##_caches();\
-}\
-\
-void zfft(complex_double *inout, int n, \
- int direction, int howmany, int normalize)\
-{\
- zfft_##name(inout, n, direction, howmany, normalize);\
-}
-
-/* ************** Definition of backend specific functions ********* */
-
-/*
- * To add a backend :
- * - create a file zfft_name.c, where you define a function zfft_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_zname_caches(void),
- * which can do nothing
- * - in zfft.c, include the zfft_name.c file, and add the 3 following lines
- * just after it:
- * #ifndef WITH_DJBFFT
- * GEN_PUBLIC_API(name)
- * #endif
- */
-
-#ifdef WITH_FFTW3
- #include "zfft_fftw3.c"
- #ifndef WITH_DJBFFT
- GEN_PUBLIC_API(fftw3)
- #endif
-#elif defined WITH_FFTW
- #include "zfft_fftw.c"
- #ifndef WITH_DJBFFT
- GEN_PUBLIC_API(fftw)
- #endif
-#elif defined WITH_MKL
- #include "zfft_mkl.c"
- #ifndef WITH_DJBFFT
- GEN_PUBLIC_API(mkl)
- #endif
-#else /* Use fftpack by default */
- #include "zfft_fftpack.c"
- #ifndef WITH_DJBFFT
- GEN_PUBLIC_API(fftpack)
- #endif
-#endif
-
-/*
- * djbfft must be used at the end, because it needs another backend (defined
- * above) for non 2^n * size
- */
-#ifdef WITH_DJBFFT
- #include "zfft_djbfft.c"
- void destroy_zfft_cache(void)
- {
- destroy_zdjbfft_caches();
- zfft_def_destroy_cache();
- }
- void zfft(complex_double *inout, int n,
- int direction, int howmany, int normalize)
- {
- zfft_djbfft(inout, n, direction, howmany, normalize);
- }
-#endif
Copied: branches/refactor_fft/scipy/fftpack/src/zfft.cxx (from rev 4212, branches/refactor_fft/scipy/fftpack/src/zfft.c)
===================================================================
--- branches/refactor_fft/scipy/fftpack/src/zfft.c 2008-05-04 13:25:43 UTC (rev 4212)
+++ branches/refactor_fft/scipy/fftpack/src/zfft.cxx 2008-05-05 07:24:48 UTC (rev 4214)
@@ -0,0 +1,76 @@
+/*
+ Interface to various FFT libraries.
+ Double complex FFT and IFFT.
+ Author: Pearu Peterson, August 2002
+ */
+
+#include "fftpack.h"
+
+/* The following macro convert private backend specific function to the public
+ * functions exported by the module */
+#define GEN_PUBLIC_API(name) \
+extern "C" void destroy_zfft_cache(void)\
+{\
+ destroy_z##name##_caches();\
+}\
+\
+extern "C" void zfft(complex_double *inout, int n, \
+ int direction, int howmany, int normalize)\
+{\
+ zfft_##name(inout, n, direction, howmany, normalize);\
+}
+
+/* ************** Definition of backend specific functions ********* */
+
+/*
+ * To add a backend :
+ * - create a file zfft_name.c, where you define a function zfft_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_zname_caches(void),
+ * which can do nothing
+ * - in zfft.c, include the zfft_name.c file, and add the 3 following lines
+ * just after it:
+ * #ifndef WITH_DJBFFT
+ * GEN_PUBLIC_API(name)
+ * #endif
+ */
+
+#ifdef WITH_FFTW3
+ #include "zfft_fftw3.cxx"
+ #ifndef WITH_DJBFFT
+ GEN_PUBLIC_API(fftw3)
+ #endif
+#elif defined WITH_FFTW
+ #include "zfft_fftw.c"
+ #ifndef WITH_DJBFFT
+ GEN_PUBLIC_API(fftw)
+ #endif
+#elif defined WITH_MKL
+ #include "zfft_mkl.c"
+ #ifndef WITH_DJBFFT
+ GEN_PUBLIC_API(mkl)
+ #endif
+#else /* Use fftpack by default */
+ #include "zfft_fftpack.c"
+ #ifndef WITH_DJBFFT
+ GEN_PUBLIC_API(fftpack)
+ #endif
+#endif
+
+/*
+ * djbfft must be used at the end, because it needs another backend (defined
+ * above) for non 2^n * size
+ */
+#ifdef WITH_DJBFFT
+ #include "zfft_djbfft.cxx"
+ extern "C" void destroy_zfft_cache(void)
+ {
+ destroy_zdjbfft_caches();
+ zfft_def_destroy_cache();
+ }
+ extern "C" void zfft(complex_double *inout, int n,
+ int direction, int howmany, int normalize)
+ {
+ zfft_djbfft(inout, n, direction, howmany, normalize);
+ }
+#endif
Deleted: branches/refactor_fft/scipy/fftpack/src/zfft_djbfft.c
===================================================================
--- branches/refactor_fft/scipy/fftpack/src/zfft_djbfft.c 2008-05-05 06:52:28 UTC (rev 4213)
+++ branches/refactor_fft/scipy/fftpack/src/zfft_djbfft.c 2008-05-05 07:24:48 UTC (rev 4214)
@@ -1,151 +0,0 @@
-/*
-* DJBFFT only implements size 2^N !
-*
-* zfft_def and zfft_def_destroy_cache are the functions
-* used for size different than 2^N
-*/
-#ifdef WITH_FFTWORK
-#define zfft_def zfft_fftwork
-#define zfft_def_destroy_cache destroy_zfftwork_cache
-#elif defined WITH_FFTW3
-#define zfft_def zfft_fftw3
-#define zfft_def_destroy_cache destroy_zfftw3_caches
-#elif defined WITH_FFTW
-#define zfft_def zfft_fftw
-#define zfft_def_destroy_cache destroy_zfftw_caches
-#else
-#define zfft_def zfft_fftpack
-#define zfft_def_destroy_cache destroy_zfftpack_caches
-#endif
-
-GEN_CACHE(zdjbfft,(int n)
- ,unsigned int* f;
- double* ptr;
- ,caches_zdjbfft[i].n==n
- ,caches_zdjbfft[id].f = (unsigned int*)malloc(sizeof(unsigned int)*(n));
- caches_zdjbfft[id].ptr = (double*)malloc(sizeof(double)*(2*n));
- fftfreq_ctable(caches_zdjbfft[id].f,n);
- for(i=0;i<n;++i)
- caches_zdjbfft[id].f[i] = (n-caches_zdjbfft[id].f[i])%n;
- ,free(caches_zdjbfft[id].f);
- free(caches_zdjbfft[id].ptr);
- ,10)
-
-/**************** ZFFT function **********************/
-static void zfft_djbfft(complex_double * inout,
- int n, int direction, int howmany, int normalize)
-{
- int i;
- complex_double *ptr = inout;
- int j;
- complex_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_zdjbfft(n);
- f = caches_zdjbfft[i].f;
- ptrc = (complex_double *) caches_zdjbfft[i].ptr;
- }
- if (f == 0) {
- zfft_def(inout, n, direction, howmany, normalize);
- }
-
- switch (direction) {
- case 1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- if (f != NULL) {
- memcpy(ptrc, ptr, 2 * n * sizeof(double));
- switch (n) {
-#define TMPCASE(N) case N: fftc8_##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
- }
- for (j = 0; j < n; ++j) {
- *(ptr + f[j]) = *(ptrc + j);
- }
- }
-
- }
- break;
-
- case -1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- if (f != NULL) {
- for (j = 0; j < n; ++j) {
- *(ptrc + j) = *(ptr + f[j]);
- }
- switch (n) {
-#define TMPCASE(N) case N: fftc8_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
- }
- memcpy(ptr, ptrc, 2 * n * sizeof(double));
- }
- }
- break;
- default:
- fprintf(stderr, "zfft: invalid direction=%d\n", direction);
- }
-
- if (normalize) {
- ptr = inout;
- if (f != NULL) {
- for (i = 0; i < howmany; ++i, ptr += n) {
- switch (n) {
-#define TMPCASE(N) case N: fftc8_scale##N(ptr); 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
- }
- }
- }
- }
-}
Copied: branches/refactor_fft/scipy/fftpack/src/zfft_djbfft.cxx (from rev 4212, branches/refactor_fft/scipy/fftpack/src/zfft_djbfft.c)
Deleted: branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.c
===================================================================
--- branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.c 2008-05-05 06:52:28 UTC (rev 4213)
+++ branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.c 2008-05-05 07:24:48 UTC (rev 4214)
@@ -1,121 +0,0 @@
-
-/* This cache uses FFTW_MEASURE for the plans, and do not copy the data. */
-GEN_CACHE(zfftw3,(int n,int d)
- ,int direction;
- fftw_plan plan;
- fftw_complex *wrk;
- ,((caches_zfftw3[i].n==n) &&
- (caches_zfftw3[i].direction==d))
- ,caches_zfftw3[id].direction = d;
- /* This working buffer is only used to compute the plan: we need it
- since FFTW_MEASURE destroys its input when computing a plan */
- caches_zfftw3[id].wrk = fftw_malloc(n * sizeof(double) * 2);
- caches_zfftw3[id].plan = fftw_plan_dft_1d(n,
- caches_zfftw3[id].wrk,
- caches_zfftw3[id].wrk,
- (d>0?FFTW_FORWARD:FFTW_BACKWARD),
- FFTW_ESTIMATE | FFTW_UNALIGNED);
- ,
- fftw_destroy_plan(caches_zfftw3[id].plan);
- fftw_free(caches_zfftw3[id].wrk);
- ,10)
-
-static void zfft_fftw3(complex_double * inout, int n, int dir, int
-howmany, int normalize)
-{
- fftw_complex *ptr = (fftw_complex*)inout;
- fftw_complex *ptrm;
- fftw_plan plan = NULL;
- double factor = 1./n;
-
- int i;
-
- plan = caches_zfftw3[get_cache_id_zfftw3(n, dir)].plan;
-
- switch (dir) {
- case 1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- fftw_execute_dft(plan, ptr, ptr);
- }
- break;
-
- case -1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- fftw_execute_dft(plan, ptr, ptr);
- }
- break;
-
- default:
- fprintf(stderr, "zfft: invalid dir=%d\n", dir);
- }
-
- if (normalize) {
- ptr =(fftw_complex*)inout;
- for (i = n * howmany - 1; i >= 0; --i) {
- *((double *) (ptr)) *= factor;
- *((double *) (ptr++) + 1) *= factor;
- }
- }
-}
-#if 0
-GEN_CACHE(zfftw3,(int n,int d)
- ,int direction;
- fftw_plan plan;
- fftw_complex* ptr;
- ,((caches_zfftw3[i].n==n) &&
- (caches_zfftw3[i].direction==d))
- ,caches_zfftw3[id].direction = d;
- caches_zfftw3[id].ptr = fftw_malloc(sizeof(fftw_complex)*(n));
- caches_zfftw3[id].plan = fftw_plan_dft_1d(n, caches_zfftw3[id].ptr,
- caches_zfftw3[id].ptr,
- (d>0?FFTW_FORWARD:FFTW_BACKWARD),
- FFTW_ESTIMATE);
- ,fftw_destroy_plan(caches_zfftw3[id].plan);
- fftw_free(caches_zfftw3[id].ptr);
- ,10)
-
-static void zfft_fftw3(complex_double * inout, int n, int dir, int
-howmany, int normalize)
-{
- complex_double *ptr = inout;
- fftw_complex *ptrm = NULL;
- fftw_plan plan = NULL;
-
- int i;
-
- plan = caches_zfftw3[get_cache_id_zfftw3(n, dir)].plan;
-
- switch (dir) {
- case 1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- ptrm =
- caches_zfftw3[get_cache_id_zfftw3(n, dir)].ptr;
- memcpy(ptrm, ptr, sizeof(double) * 2 * n);
- fftw_execute(plan);
- memcpy(ptr, ptrm, sizeof(double) * 2 * n);
- }
- break;
-
- case -1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- ptrm =
- caches_zfftw3[get_cache_id_zfftw3(n, dir)].ptr;
- memcpy(ptrm, ptr, sizeof(double) * 2 * n);
- fftw_execute(plan);
- memcpy(ptr, ptrm, sizeof(double) * 2 * n);
- }
- break;
-
- default:
- fprintf(stderr, "zfft: invalid dir=%d\n", dir);
- }
-
- if (normalize) {
- ptr = inout;
- for (i = n * howmany - 1; i >= 0; --i) {
- *((double *) (ptr)) /= n;
- *((double *) (ptr++) + 1) /= n;
- }
- }
-}
-#endif
Copied: branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.cxx (from rev 4212, branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.c)
===================================================================
--- branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.c 2008-05-04 13:25:43 UTC (rev 4212)
+++ branches/refactor_fft/scipy/fftpack/src/zfft_fftw3.cxx 2008-05-05 07:24:48 UTC (rev 4214)
@@ -0,0 +1,120 @@
+/* This cache uses FFTW_MEASURE for the plans, and do not copy the data. */
+GEN_CACHE(zfftw3,(int n,int d)
+ ,int direction;
+ fftw_plan plan;
+ fftw_complex *wrk;
+ ,((caches_zfftw3[i].n==n) &&
+ (caches_zfftw3[i].direction==d))
+ ,caches_zfftw3[id].direction = d;
+ /* This working buffer is only used to compute the plan: we need it
+ since FFTW_MEASURE destroys its input when computing a plan */
+ caches_zfftw3[id].wrk = (fftw_complex*)fftw_malloc(n * sizeof(double) * 2);
+ caches_zfftw3[id].plan = fftw_plan_dft_1d(n,
+ caches_zfftw3[id].wrk,
+ caches_zfftw3[id].wrk,
+ (d>0?FFTW_FORWARD:FFTW_BACKWARD),
+ FFTW_ESTIMATE | FFTW_UNALIGNED);
+ ,
+ fftw_destroy_plan(caches_zfftw3[id].plan);
+ fftw_free(caches_zfftw3[id].wrk);
+ ,10)
+
+static void zfft_fftw3(complex_double * inout, int n, int dir, int
+howmany, int normalize)
+{
+ fftw_complex *ptr = (fftw_complex*)inout;
+ fftw_complex *ptrm;
+ fftw_plan plan = NULL;
+ double factor = 1./n;
+
+ int i;
+
+ plan = caches_zfftw3[get_cache_id_zfftw3(n, dir)].plan;
+
+ switch (dir) {
+ case 1:
+ for (i = 0; i < howmany; ++i, ptr += n) {
+ fftw_execute_dft(plan, ptr, ptr);
+ }
+ break;
+
+ case -1:
+ for (i = 0; i < howmany; ++i, ptr += n) {
+ fftw_execute_dft(plan, ptr, ptr);
+ }
+ break;
+
+ default:
+ fprintf(stderr, "zfft: invalid dir=%d\n", dir);
+ }
+
+ if (normalize) {
+ ptr =(fftw_complex*)inout;
+ for (i = n * howmany - 1; i >= 0; --i) {
+ *((double *) (ptr)) *= factor;
+ *((double *) (ptr++) + 1) *= factor;
+ }
+ }
+}
+#if 0
+GEN_CACHE(zfftw3,(int n,int d)
+ ,int direction;
+ fftw_plan plan;
+ fftw_complex* ptr;
+ ,((caches_zfftw3[i].n==n) &&
+ (caches_zfftw3[i].direction==d))
+ ,caches_zfftw3[id].direction = d;
+ caches_zfftw3[id].ptr = fftw_malloc(sizeof(fftw_complex)*(n));
+ caches_zfftw3[id].plan = fftw_plan_dft_1d(n, caches_zfftw3[id].ptr,
+ caches_zfftw3[id].ptr,
+ (d>0?FFTW_FORWARD:FFTW_BACKWARD),
+ FFTW_ESTIMATE);
+ ,fftw_destroy_plan(caches_zfftw3[id].plan);
+ fftw_free(caches_zfftw3[id].ptr);
+ ,10)
+
+static void zfft_fftw3(complex_double * inout, int n, int dir, int
+howmany, int normalize)
+{
+ complex_double *ptr = inout;
+ fftw_complex *ptrm = NULL;
+ fftw_plan plan = NULL;
+
+ int i;
+
+ plan = caches_zfftw3[get_cache_id_zfftw3(n, dir)].plan;
+
+ switch (dir) {
+ case 1:
+ for (i = 0; i < howmany; ++i, ptr += n) {
+ ptrm =
+ caches_zfftw3[get_cache_id_zfftw3(n, dir)].ptr;
+ memcpy(ptrm, ptr, sizeof(double) * 2 * n);
+ fftw_execute(plan);
+ memcpy(ptr, ptrm, sizeof(double) * 2 * n);
+ }
+ break;
+
+ case -1:
+ for (i = 0; i < howmany; ++i, ptr += n) {
+ ptrm =
+ caches_zfftw3[get_cache_id_zfftw3(n, dir)].ptr;
+ memcpy(ptrm, ptr, sizeof(double) * 2 * n);
+ fftw_execute(plan);
+ memcpy(ptr, ptrm, sizeof(double) * 2 * n);
+ }
+ break;
+
+ default:
+ fprintf(stderr, "zfft: invalid dir=%d\n", dir);
+ }
+
+ if (normalize) {
+ ptr = inout;
+ for (i = n * howmany - 1; i >= 0; --i) {
+ *((double *) (ptr)) /= n;
+ *((double *) (ptr++) + 1) /= n;
+ }
+ }
+}
+#endif
More information about the Scipy-svn
mailing list