46 #include <blitz/array.h> 49 #include <blitz/tinyvec-et.h> 56 template <
typename DataT>
74 template<
typename DataT>
137 void saveWisdom()
const;
157 template<
int Dim,
typename IndexT>
158 blitz::TinyVector<IndexT,Dim>
159 getPaddedShape(
const blitz::TinyVector<IndexT,Dim>& minimumExtents)
const;
192 template<
int Dim,
typename IndexT>
193 void pad(
const blitz::Array<DataT,Dim>& in,
194 blitz::Array<DataT,Dim>& out,
195 blitz::TinyVector<IndexT,Dim>& lb,
196 blitz::TinyVector<IndexT,Dim>& ub,
197 const blitz::TinyVector<IndexT,Dim>& minimumExtents = IndexT(0),
199 const DataT paddingValue = static_cast<DataT>(0))
const;
233 template<
int Dim,
typename IndexT>
234 void pad(
const blitz::Array<std::complex<DataT>,Dim>& in,
235 blitz::Array<std::complex<DataT>,Dim>& out,
236 blitz::TinyVector<IndexT,Dim>& lb,
237 blitz::TinyVector<IndexT,Dim>& ub,
238 const blitz::TinyVector<IndexT,Dim>& minimumExtents = IndexT(0),
240 const std::complex<DataT> paddingValue =
241 static_cast< std::complex<DataT>
>(0))
const;
259 template<
int Dim,
typename IndexT>
260 void unpad(
const blitz::Array<DataT,Dim>& in,
261 blitz::Array<DataT,Dim>& out,
262 const blitz::TinyVector<IndexT,Dim>& lb,
263 const blitz::TinyVector<IndexT,Dim>& ub)
const;
281 template<
int Dim,
typename IndexT>
282 void unpad(
const blitz::Array<std::complex<DataT>,Dim>& in,
283 blitz::Array<std::complex<DataT>,Dim>& out,
284 const blitz::TinyVector<IndexT,Dim>& lb,
285 const blitz::TinyVector<IndexT,Dim>& ub)
const;
319 void plan_forward(blitz::Array<DataT,Dim>& in,
320 blitz::Array<std::complex<DataT>,Dim>& out,
322 const unsigned int plan_flags = FFTW_MEASURE)
const;
356 void plan_backward(blitz::Array<std::complex<DataT>,Dim>& in,
357 blitz::Array<DataT,Dim>& out,
359 const unsigned int plan_flags = FFTW_MEASURE)
const;
393 void plan_forward(blitz::Array<std::complex<DataT>,Dim>& in,
394 blitz::Array<std::complex<DataT>,Dim>& out,
396 const unsigned int plan_flags = FFTW_MEASURE)
const;
404 get_plan_forward(blitz::Array<DataT,Dim>& in,
405 blitz::Array<std::complex<DataT>,Dim>& out,
407 const unsigned int plan_flags = FFTW_MEASURE)
const;
411 get_plan_forward(blitz::Array<std::complex<DataT>,Dim>& in,
412 blitz::Array<std::complex<DataT>,Dim>& out,
414 const unsigned int plan_flags = FFTW_MEASURE)
const;
448 void plan_backward(blitz::Array<std::complex<DataT>,Dim>& in,
449 blitz::Array<std::complex<DataT>,Dim>& out,
451 const unsigned int plan_flags = FFTW_MEASURE)
const;
460 get_plan_backward(blitz::Array<std::complex<DataT>,Dim>& in,
461 blitz::Array<DataT,Dim>& out,
463 const unsigned int plan_flags = FFTW_MEASURE)
const;
467 get_plan_backward(blitz::Array<std::complex<DataT>,Dim>& in,
468 blitz::Array<std::complex<DataT>,Dim>& out,
470 const unsigned int plan_flags = FFTW_MEASURE)
const;
488 void forward(
const blitz::Array<DataT,Dim>& in,
489 blitz::Array<std::complex<DataT>,Dim>& out)
const;
507 void forward(blitz::Array<std::complex<DataT>,Dim>& in,
508 blitz::Array<std::complex<DataT>,Dim>& out)
const;
527 void backward(blitz::Array<std::complex<DataT>,Dim>& in,
528 blitz::Array<DataT,Dim>& out,
548 void backward(blitz::Array<std::complex<DataT>,Dim>& in,
549 blitz::Array<std::complex<DataT>,Dim>& out,
559 void exec_guru_plan(blitz::Array<std::complex<DataT>,Dim>& in,
560 blitz::Array<std::complex<DataT>,Dim>& out,
561 blitz_fftw_plan plan,
566 void exec_guru_plan_r2c(blitz::Array<DataT,Dim>& in,
567 blitz::Array<std::complex<DataT>,Dim>& out,
568 blitz_fftw_plan plan,
572 void exec_guru_plan_c2r(blitz::Array<std::complex<DataT>,Dim>& in,
573 blitz::Array<DataT,Dim>& out,
574 blitz_fftw_plan plan,
589 void unShuffle(
const blitz::Array<DataT,Dim>& in,
590 blitz::Array<DataT,Dim>& out)
const;
602 static size_t nextBestFFTSize(
const size_t size);
603 static size_t prevBestFFTSize(
const size_t size);
604 static const std::set<size_t>& nextBestFFTSizes();
607 void translate(
const blitz::Array<DataT,Dim>& data,
608 blitz::Array<DataT,Dim>& dataTrans,
613 void translate(
const blitz::Array<DataT,Dim>& data,
614 blitz::Array<DataT,Dim>& dataTrans,
615 const blitz::TinyVector<BlitzIndexT,Dim>& t)
const;
623 class BlitzFFTWDestructor
626 ~BlitzFFTWDestructor()
628 if (BlitzFFTW::p_instance != 0) {
629 delete BlitzFFTW::p_instance;
630 BlitzFFTW::p_instance = 0;
636 friend class BlitzFFTWDestructor;
638 static const size_t maxPrepareFFTSize = 65535;
639 static void prepareFFTSizes();
640 static std::set<size_t> _bestFFTSizes;
648 typedef DataT blitz_fftw_complex[2];
650 void (*blitz_fftw_execute)(
const blitz_fftw_plan p);
651 blitz_fftw_plan (*blitz_fftw_plan_dft)(
int rank,
const int *n, blitz_fftw_complex *in, blitz_fftw_complex *out,
int sign,
unsigned flags);
652 void (*blitz_fftw_execute_dft)(
const blitz_fftw_plan p, blitz_fftw_complex *in, blitz_fftw_complex *out);
653 blitz_fftw_plan (*blitz_fftw_plan_dft_r2c)(
int rank,
const int *n, DataT *in, blitz_fftw_complex *out,
unsigned flags);
654 blitz_fftw_plan (*blitz_fftw_plan_dft_c2r)(
int rank,
const int *n, blitz_fftw_complex *in, DataT *out,
unsigned flags);
655 void (*blitz_fftw_execute_dft_r2c)(
const blitz_fftw_plan p, DataT *in, blitz_fftw_complex *out);
656 void (*blitz_fftw_execute_dft_c2r)(
const blitz_fftw_plan p, blitz_fftw_complex *in, DataT *out);
657 void (*blitz_fftw_destroy_plan)(blitz_fftw_plan p);
658 void (*blitz_fftw_cleanup)(void);
659 void *(*blitz_fftw_malloc)(
size_t n);
673 #ifndef BLITZFFTW_HH_NO_INCLUDE_ICC 674 #include "BlitzFFTW.icc" static BlitzFFTW * instance()
Get a handle to the Singleton BlitzFFTW Object.
fftwf_plan real_fftw_plan
BlitzFFTWPlan< DataT >::real_fftw_plan blitz_fftw_plan