[Scipy-svn] r4286 - branches/refactor_fft/scipy/fftpack/src/djbfft
scipy-svn at scipy.org
scipy-svn at scipy.org
Mon May 12 09:09:23 EDT 2008
Author: cdavid
Date: 2008-05-12 08:09:20 -0500 (Mon, 12 May 2008)
New Revision: 4286
Modified:
branches/refactor_fft/scipy/fftpack/src/djbfft/drfft.cxx
Log:
Fix drfft implementation by djbfft backend.
Modified: branches/refactor_fft/scipy/fftpack/src/djbfft/drfft.cxx
===================================================================
--- branches/refactor_fft/scipy/fftpack/src/djbfft/drfft.cxx 2008-05-12 13:08:51 UTC (rev 4285)
+++ branches/refactor_fft/scipy/fftpack/src/djbfft/drfft.cxx 2008-05-12 13:09:20 UTC (rev 4286)
@@ -1,13 +1,14 @@
/*
- * Last Change: Sun May 11 09:00 PM 2008 J
+ * Last Change: Wed Aug 01 08:00 PM 2007 J
*
* Original code by Pearu Peterson.
*/
/*
- * RDJBFFT only implements size 2^N !
+ * DJBFFT only implements size 2^N !
*
- * drfft_def is the functions used for size different * than 2^N
+ * drfft_def and drfft_def_destroy_cache are the functions used for size different
+ * than 2^N
*/
#include <new>
#include <cassert>
@@ -67,117 +68,130 @@
free(m_f);
}
-int RDJBFFTCache::compute_forward(double *inout) const
+int RDJBFFTCache::compute_forward(double *inout) const
{
- const int n = m_id.m_n;
+ double *ptr = inout;
+ double *ptrc = m_ptr;
- COPYSTD2DJB(inout, m_ptr, n);
- switch (n) {
-#define TMPCASE(N) case N: fftr8_##N(m_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);
+ int n = m_id.m_n;
+ 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(m_ptr, inout, m_f, n);
-
- return 0;
+ }
+ COPYDJB2STD(ptrc, ptr, m_f, n);
+ return 0;
}
-int RDJBFFTCache::compute_backward(double *inout, int normalize) const
+int RDJBFFTCache::compute_backward(double *inout, int normalize) const
{
- const int n = m_id.m_n;
+ double *ptr = inout;
+ double *ptrc = m_ptr;
- COPYINVSTD2DJB(inout, m_ptr, normalize, m_f, n);
- switch (n) {
-#define TMPCASE(N)case N:if(normalize)fftr8_scale##N(m_ptr);fftr8_un##N(m_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);
+ int n = m_id.m_n;
+
+ COPYINVSTD2DJB(ptr, ptrc, normalize, m_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(m_ptr, inout, n);
-
- return 0;
+ }
+ COPYINVDJB2STD(ptrc, ptr, n);
+ return 0;
}
static CacheManager<DJBFFTCacheId, RDJBFFTCache> rdjbfft_cmgr(10);
+#if 0
+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)
+#endif
+
/**************** ZFFT function **********************/
static void drfft_djbfft(double * inout,
int n, int direction, int howmany, int normalize)
{
- int i;
- double *ptr = inout;
- RDJBFFTCache *cache;
- unsigned int *f = NULL;
+ int i;
+ double *ptr = inout;
+ RDJBFFTCache *cache;
- 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:
- cache = rdjbfft_cmgr.get_cache(DJBFFTCacheId(n));
- break;
- default:
- /* For sizes not handled by djbfft, use default
- * implementation and returns */
- drfft_def(inout, n, direction, howmany, normalize);
- return;
- }
+ 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:
+ cache = rdjbfft_cmgr.get_cache(DJBFFTCacheId(n));
+ break;
+ default:
+ drfft_def(inout, n, direction, howmany, normalize);
+ return;
+ }
- switch (direction) {
- case 1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- cache->compute_forward(ptr);
- }
- break;
+ switch (direction) {
+ case 1:
+ for (i = 0; i < howmany; ++i, ptr += n) {
+ cache->compute_forward(ptr);
+ }
+ break;
- case -1:
- for (i = 0; i < howmany; ++i, ptr += n) {
- cache->compute_backward(ptr, normalize);
- }
- break;
+ case -1:
+ for (i = 0; i < howmany; ++i, ptr += n) {
+ cache->compute_backward(ptr, normalize);
+ }
+ break;
- default:
- fprintf(stderr, "drfft: invalid direction=%d\n",
- direction);
- }
+ 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;
- }
- }
+ if (normalize && direction == 1) {
+ 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