31# include <itpp/config.h>
33# include <itpp/config_msvc.h>
37#if defined(HAVE_FFT_MKL)
45# include <mkl_service.h>
46# undef DftiCreateDescriptor
49#elif defined(HAVE_FFT_ACML)
56#elif defined(HAVE_FFTW3)
67enum MultithreadingTag {SingleThreaded = 1, OmpThreaded};
73static const MultithreadingTag ThreadingTag = OmpThreaded;
76static const int contexts_per_transform_type = 10;
83 Mutex& operator=(
const Mutex&);
85 Mutex() {omp_init_lock(&_lck);}
86 ~Mutex() {omp_destroy_lock(&_lck);}
88 void lock() {omp_set_lock(&_lck);}
90 bool try_lock() {
return (omp_test_lock(&_lck)) != 0;}
92 void unlock() {omp_unset_lock(&_lck);}
98static const MultithreadingTag ThreadingTag = SingleThreaded;
101static const int contexts_per_transform_type = 1;
108 Mutex& operator=(
const Mutex&);
113 bool try_lock() {
return true;}
126 Lock& operator=(
const Lock&);
128 Lock(Mutex& m): _m(m) {_m.lock();}
129 ~Lock() {_m.unlock();}
138struct FFTCplx_Traits {
140 typedef cvec OutType;
143struct IFFTCplx_Traits {
145 typedef cvec OutType;
148struct FFTReal_Traits {
150 typedef cvec OutType;
153struct IFFTReal_Traits {
169template<
typename TransformTraits>
class Transform;
171template<MultithreadingTag>
inline void init_fft_library();
173#if defined(HAVE_FFT_MKL)
184template<>
inline void init_fft_library<SingleThreaded>() {}
185template<>
inline void init_fft_library<OmpThreaded>()
189 mkl::mkl_domain_set_num_threads(1, MKL_FFT);
196inline void release_descriptor(mkl::DFTI_DESCRIPTOR* h)
199 MKL_LONG status = mkl::DftiFreeDescriptor(&h);
201 it_info(mkl::DftiErrorMessage(status));
202 it_error(
"MKL library release_descriptor() failed on DftiFreeDescriptor.");
207template<>
class Transform<FFTCplx_Traits>
209 mkl::DFTI_DESCRIPTOR* _h;
210 int _transform_length;
212 Transform(): _h(NULL), _transform_length(0) {}
214 void compute_transform(
const cvec &in, cvec &out) {
215 out.set_size(in.size(),
false);
216 if(_transform_length != in.size()) {
217 release_descriptor(_h);
218 _transform_length = in.size();
220 MKL_LONG status = mkl::DftiCreateDescriptor(&_h, mkl::DFTI_DOUBLE, mkl::DFTI_COMPLEX, 1, _transform_length);
222 it_info(mkl::DftiErrorMessage(status));
223 it_error(
"MKL library compute_transform() failed on DftiCreateDescriptor.");
226 mkl::DftiSetValue(_h, mkl::DFTI_PLACEMENT, mkl::DFTI_NOT_INPLACE);
228 status = mkl::DftiCommitDescriptor(_h);
230 it_info(mkl::DftiErrorMessage(status));
231 it_error(
"MKL library compute_transform() failed on DftiCommitDescriptor.");
235 mkl::DftiComputeForward(_h, (
void *)in._data(), out._data());
238 void reset() {release_descriptor(_h); *
this = Transform();}
241template<>
class Transform<IFFTCplx_Traits>
243 mkl::DFTI_DESCRIPTOR* _h;
244 int _transform_length;
246 Transform(): _h(NULL), _transform_length(0) {}
248 void compute_transform(
const cvec &in, cvec &out) {
249 out.set_size(in.size(),
false);
250 if(_transform_length != in.size()) {
251 release_descriptor(_h);
252 _transform_length = in.size();
253 MKL_LONG status = mkl::DftiCreateDescriptor(&_h, mkl::DFTI_DOUBLE, mkl::DFTI_COMPLEX, 1, _transform_length);
255 it_info(mkl::DftiErrorMessage(status));
256 it_error(
"MKL library compute_transform() failed on DftiCreateDescriptor.");
259 mkl::DftiSetValue(_h, mkl::DFTI_PLACEMENT, mkl::DFTI_NOT_INPLACE);
260 mkl::DftiSetValue(_h, mkl::DFTI_BACKWARD_SCALE, 1.0 / _transform_length);
262 status = mkl::DftiCommitDescriptor(_h);
264 it_info(mkl::DftiErrorMessage(status));
265 it_error(
"MKL library compute_transform() failed on DftiCommitDescriptor.");
269 mkl::DftiComputeBackward(_h, (
void *)in._data(), out._data());
272 void reset() {release_descriptor(_h); *
this = Transform();}
275template<>
class Transform<FFTReal_Traits>
277 mkl::DFTI_DESCRIPTOR* _h;
278 int _transform_length;
280 Transform(): _h(NULL), _transform_length(0) {}
282 void compute_transform(
const vec &in, cvec &out) {
283 out.set_size(in.size(),
false);
284 if(_transform_length != in.size()) {
285 release_descriptor(_h);
286 _transform_length = in.size();
288 MKL_LONG status = mkl::DftiCreateDescriptor(&_h, mkl::DFTI_DOUBLE, mkl::DFTI_REAL, 1, _transform_length);
290 it_info(mkl::DftiErrorMessage(status));
291 it_error(
"MKL library compute_transform() failed on DftiCreateDescriptor.");
294 mkl::DftiSetValue(_h, mkl::DFTI_PLACEMENT, mkl::DFTI_NOT_INPLACE);
296 status = mkl::DftiCommitDescriptor(_h);
298 it_info(mkl::DftiErrorMessage(status));
299 it_error(
"MKL library compute_transform() failed on DftiCommitDescriptor.");
303 mkl::DftiComputeForward(_h, (
void *)in._data(), out._data());
307 int istart =
ceil_i(in.size() / 2.0);
308 int idelta = in.size() - istart;
309 out.set_subvector(istart,
reverse(
conj(out(1, idelta))));
312 void reset() {release_descriptor(_h); *
this = Transform();}
315template<>
class Transform<IFFTReal_Traits>
317 mkl::DFTI_DESCRIPTOR* _h;
318 int _transform_length;
320 Transform(): _h(NULL), _transform_length(0) {}
322 void compute_transform(
const cvec &in, vec &out) {
323 out.set_size(in.size(),
false);
324 if(_transform_length != in.size()) {
325 release_descriptor(_h);
326 _transform_length = in.size();
328 MKL_LONG status = mkl::DftiCreateDescriptor(&_h, mkl::DFTI_DOUBLE, mkl::DFTI_REAL, 1, _transform_length);
330 it_info(mkl::DftiErrorMessage(status));
331 it_error(
"MKL library compute_transform() failed on DftiCreateDescriptor.");
334 mkl::DftiSetValue(_h, mkl::DFTI_PLACEMENT, mkl::DFTI_NOT_INPLACE);
335 mkl::DftiSetValue(_h, mkl::DFTI_BACKWARD_SCALE, 1.0 / _transform_length);
337 status = mkl::DftiCommitDescriptor(_h);
339 it_info(mkl::DftiErrorMessage(status));
340 it_error(
"MKL library compute_transform() failed on DftiCommitDescriptor.");
344 mkl::DftiComputeBackward(_h, (
void *)in._data(), out._data());
347 void reset() {release_descriptor(_h); *
this = Transform();}
353#if defined(HAVE_FFT_ACML)
363template<>
inline void init_fft_library<SingleThreaded>() {}
364template<>
inline void init_fft_library<OmpThreaded>() {}
369template<>
class Transform<FFTCplx_Traits>
372 int _transform_length;
374 Transform(): _transform_length(0) {}
376 void compute_transform(
const cvec &in, cvec &out) {
378 out.set_size(in.size(),
false);
379 if(_transform_length != in.size()) {
380 _transform_length = in.size();
382 int min_required_size = 5 * _transform_length + 100;
383 if(_scratchpad.size() < min_required_size) _scratchpad.set_size(min_required_size);
385 acml::zfft1dx(0, 1.0,
false, _transform_length, (acml::doublecomplex *)in._data(), 1,
386 (acml::doublecomplex *)out._data(), 1,
387 (acml::doublecomplex *)_scratchpad._data(), &info);
389 acml::zfft1dx(-1, 1.0,
false, _transform_length, (acml::doublecomplex *)in._data(), 1,
390 (acml::doublecomplex *)out._data(), 1,
391 (acml::doublecomplex *)_scratchpad._data(), &info);
394 void reset() {*
this = Transform();}
397template<>
class Transform<IFFTCplx_Traits>
400 int _transform_length;
402 Transform(): _transform_length(0) {}
404 void compute_transform(
const cvec &in, cvec &out) {
406 out.set_size(in.size(),
false);
407 if(_transform_length != in.size()) {
408 _transform_length = in.size();
410 int min_required_size = 5 * _transform_length + 100;
411 if(_scratchpad.size() < min_required_size) _scratchpad.set_size(min_required_size);
413 acml::zfft1dx(0, 1.0 / _transform_length,
false, _transform_length, (acml::doublecomplex *)in._data(), 1,
414 (acml::doublecomplex *)out._data(), 1,
415 (acml::doublecomplex *)_scratchpad._data(), &info);
417 acml::zfft1dx(1, 1.0 / _transform_length,
false, _transform_length, (acml::doublecomplex *)in._data(), 1,
418 (acml::doublecomplex *)out._data(), 1,
419 (acml::doublecomplex *)_scratchpad._data(), &info);
422 void reset() {*
this = Transform();}
425template<>
class Transform<FFTReal_Traits>
428 int _transform_length;
430 Transform(): _transform_length(0) {}
432 void compute_transform(
const vec &in, cvec &out) {
436 if(_transform_length != in.size()) {
437 _transform_length = in.size();
439 int min_required_size = 5 * _transform_length + 100;
440 if(_scratchpad.size() < min_required_size) _scratchpad.set_size(min_required_size);
442 acml::dzfft(0, _transform_length, out_re._data(), _scratchpad._data(), &info);
444 acml::dzfft(1, _transform_length, out_re._data(), _scratchpad._data(), &info);
447 double factor = std::sqrt(
static_cast<double>(_transform_length));
451 vec out_im(_transform_length);
453 if(!(_transform_length % 2)) out_im(_transform_length / 2) = 0.0;
454 out_im.set_subvector(1,
reverse(out_re(_transform_length / 2 + 1, _transform_length - 1)));
455 out_im.set_subvector(_transform_length / 2 + 1, -out_re(_transform_length / 2 + 1, _transform_length - 1));
456 out_re.set_subvector(_transform_length / 2 + 1,
reverse(out_re(1, (_transform_length - 1) / 2)));
458 out.set_size(_transform_length,
false);
462 void reset() {*
this = Transform();}
465template<>
class Transform<IFFTReal_Traits>
468 int _transform_length;
470 Transform(): _transform_length(0) {}
472 void compute_transform(
const cvec &in, vec &out) {
474 out.set_size(in.size());
475 out.set_subvector(0,
real(in(0, in.size() / 2)));
476 out.set_subvector(in.size() / 2 + 1, -
imag(in(in.size() / 2 + 1, in.size() - 1)));
479 if(_transform_length != in.size()) {
480 _transform_length = in.size();
482 int min_required_size = 5 * _transform_length + 100;
483 if(_scratchpad.size() < min_required_size) _scratchpad.set_size(min_required_size);
485 acml::zdfft(0, _transform_length, out._data(), _scratchpad._data(), &info);
487 acml::zdfft(1, _transform_length, out._data(), _scratchpad._data(), &info);
488 out.set_subvector(1,
reverse(out(1, _transform_length - 1)));
491 double factor = 1.0 / std::sqrt(
static_cast<double>(_transform_length));
495 void reset() {*
this = Transform();}
503#if defined(HAVE_FFTW3)
511template<>
inline void init_fft_library<SingleThreaded>() {}
512template<>
inline void init_fft_library<OmpThreaded>() {}
514Mutex& get_library_lock()
516 static Mutex FFTW3LibraryLock;
517 return FFTW3LibraryLock;
522inline void destroy_plan(fftw_plan p)
524 if(p != NULL) fftw_destroy_plan(p);
527template<>
class Transform<FFTCplx_Traits>
530 int _transform_length;
532 Transform(): _p(NULL), _transform_length(0) {}
534 void compute_transform(
const cvec &in, cvec &out) {
535 out.set_size(in.size(),
false);
536 if(_transform_length != in.size()) {
537 Lock l(get_library_lock());
539 _transform_length = in.size();
542 _p = fftw_plan_dft_1d(_transform_length, (fftw_complex *)in._data(),
543 (fftw_complex *)out._data(),
544 FFTW_FORWARD, FFTW_ESTIMATE);
547 fftw_execute_dft(_p, (fftw_complex *)in._data(),
548 (fftw_complex *)out._data());
551 void reset() {destroy_plan(_p); *
this = Transform();}
554template<>
class Transform<IFFTCplx_Traits>
557 int _transform_length;
559 Transform(): _p(NULL), _transform_length(0) {}
561 void compute_transform(
const cvec &in, cvec &out) {
562 out.set_size(in.size(),
false);
563 if(_transform_length != in.size()) {
564 Lock l(get_library_lock());
566 _transform_length = in.size();
570 _p = fftw_plan_dft_1d(_transform_length, (fftw_complex *)in._data(),
571 (fftw_complex *)out._data(),
572 FFTW_BACKWARD, FFTW_ESTIMATE);
575 fftw_execute_dft(_p, (fftw_complex *)in._data(),
576 (fftw_complex *)out._data());
578 double inv_N = 1.0 / _transform_length;
583 void reset() {destroy_plan(_p); *
this = Transform();}
586template<>
class Transform<FFTReal_Traits>
589 int _transform_length;
591 Transform(): _p(NULL), _transform_length(0) {}
593 void compute_transform(
const vec &in, cvec &out) {
594 out.set_size(in.size(),
false);
595 if(_transform_length != in.size()) {
596 Lock l(get_library_lock());
598 _transform_length = in.size();
602 _p = fftw_plan_dft_r2c_1d(_transform_length, (
double *)in._data(),
603 (fftw_complex *)out._data(),
607 fftw_execute_dft_r2c(_p, (
double *)in._data(),
608 (fftw_complex *)out._data());
612 int offset =
ceil_i(_transform_length / 2.0);
613 int n_elem = _transform_length - offset;
614 for(
int i = 0; i < n_elem; ++i) {
615 out(offset + i) = std::conj(out(n_elem - i));
619 void reset() {destroy_plan(_p); *
this = Transform();}
622template<>
class Transform<IFFTReal_Traits>
625 int _transform_length;
627 Transform(): _p(NULL), _transform_length(0) {}
629 void compute_transform(
const cvec &in, vec &out) {
630 out.set_size(in.size(),
false);
631 if(_transform_length != in.size()) {
632 Lock l(get_library_lock());
634 _transform_length = in.size();
638 _p = fftw_plan_dft_c2r_1d(_transform_length, (fftw_complex *)in._data(),
639 (
double *)out._data(),
640 FFTW_ESTIMATE | FFTW_PRESERVE_INPUT);
643 fftw_execute_dft_c2r(_p, (fftw_complex *)in._data(),
644 (
double *)out._data());
646 double inv_N = 1.0 / _transform_length;
651 void reset() {destroy_plan(_p); *
this = Transform();}
657template<>
class Transform<DCT_Traits>
660 int _transform_length;
662 Transform(): _p(NULL), _transform_length(0) {}
664 void compute_transform(
const vec &in, vec &out) {
665 out.set_size(in.size(),
false);
666 if(_transform_length != in.size()) {
667 Lock l(get_library_lock());
669 _transform_length = in.size();
673 _p = fftw_plan_r2r_1d(_transform_length, (
double *)in._data(),
674 (
double *)out._data(),
675 FFTW_REDFT10, FFTW_ESTIMATE);
678 fftw_execute_r2r(_p, (
double *)in._data(), (
double *)out._data());
681 out /= std::sqrt(2.0 * _transform_length);
682 out(0) /= std::sqrt(2.0);
685 void reset() {destroy_plan(_p); *
this = Transform();}
688template<>
class Transform<IDCT_Traits>
691 int _transform_length;
693 Transform(): _p(NULL), _transform_length(0) {}
695 void compute_transform(
const vec &in, vec &out) {
699 out(0) *= std::sqrt(2.0);
700 out /= std::sqrt(2.0 * in.size());
702 if(_transform_length != in.size()) {
703 Lock l(get_library_lock());
705 _transform_length = in.size();
709 _p = fftw_plan_r2r_1d(_transform_length, (
double *)in._data(),
710 (
double *)out._data(),
711 FFTW_REDFT01, FFTW_ESTIMATE);
714 fftw_execute_r2r(_p, (
double *)out._data(), (
double *)out._data());
717 void reset() {destroy_plan(_p); *
this = Transform();}
722#if defined(HAVE_FFT_MKL) || defined(HAVE_FFT_ACML)
729template<>
class Transform<DCT_Traits>
731 Transform<FFTReal_Traits> _tr;
735 void compute_transform(
const vec &in, vec &out) {
743 for(
int i = 0; i < N; i++) {
744 c(i) *= std::complex<double>(std::cos(
pi * i / N / 2), std::sin(-
pi * i / N / 2))
745 / std::sqrt(2.0 * N);
748 out(0) /= std::sqrt(2.0);
752 void reset() {_tr.reset();}
756template<>
class Transform<IDCT_Traits>
758 Transform<IFFTReal_Traits> _tr;
762 void compute_transform(
const vec &in, vec &out) {
768 c.set_size(2 * N,
true);
769 c(0) *= std::sqrt(2.0);
770 for(
int i = 0; i < N; i++) {
771 c(i) *= std::complex<double>(std::cos(
pi * i / N / 2), std::sin(
pi * i / N / 2))
772 * std::sqrt(2.0 * N);
774 for(
int i = N - 1; i >= 1; i--) {
775 c(c.size() - i) = c(i) * std::complex<double>(std::cos(
pi * i / N),
776 std::sin(-
pi * i / N));
778 _tr.compute_transform(c, out);
779 out.set_size(N,
true);
783 void reset() {_tr.reset();}
790template<
typename TransformTraits>
class Locked_Transform :
private Transform<TransformTraits>
792 typedef Transform<TransformTraits> Base;
795 Locked_Transform() {}
797 void release_context() {Lock l(_m); Base::reset();}
798 void run_transform(
const typename TransformTraits::InType& in,
typename TransformTraits::OutType& out) {Lock l(_m); Base::compute_transform(in, out);}
817static bool is_library_initialized =
false;
819template<
typename TransformTraits>
class Transform_Provider
821 typedef Locked_Transform<TransformTraits> Transform;
822 Transform _transforms[contexts_per_transform_type];
825 Transform_Provider(): _id(0) {
826 if(!is_library_initialized) {
828 init_fft_library<ThreadingTag>();
829 is_library_initialized =
true;
832 int get_context_id() {
835 if(ret == contexts_per_transform_type)
841 void run_transform(
int id,
const typename TransformTraits::InType& in,
typename TransformTraits::OutType& out) {
842 _transforms[
id - 1].run_transform(in, out);
846 ~Transform_Provider() {
847 for(
int i = 0; i < contexts_per_transform_type; ++i)
848 _transforms[i].release_context();
853template<
typename TransformTraits> Transform_Provider<TransformTraits>& get_transform_provider()
855 static Transform_Provider<TransformTraits> p;
859void fft(
const cvec &in, cvec &out)
861 static int context_id = 0;
862 #pragma omp threadprivate(context_id)
864 if(context_id == 0) {
869 context_id = get_transform_provider<FFTCplx_Traits>().get_context_id();
872 it_assert(in.size() > 0,
"fft(): zero-sized input detected");
874 get_transform_provider<FFTCplx_Traits>().run_transform(context_id, in, out);
877void ifft(
const cvec &in, cvec &out)
879 static int context_id = 0;
880 #pragma omp threadprivate(context_id)
882 if(context_id == 0) {
887 context_id = get_transform_provider<IFFTCplx_Traits>().get_context_id();
890 it_assert(in.size() > 0,
"ifft(): zero-sized input detected");
892 get_transform_provider<IFFTCplx_Traits>().run_transform(context_id, in, out);
895void fft_real(
const vec &in, cvec &out)
897 static int context_id = 0;
898 #pragma omp threadprivate(context_id)
900 if(context_id == 0) {
905 context_id = get_transform_provider<FFTReal_Traits>().get_context_id();
908 it_assert(in.size() > 0,
"fft_real(): zero-sized input detected");
910 get_transform_provider<FFTReal_Traits>().run_transform(context_id, in, out);
915 static int context_id = 0;
916 #pragma omp threadprivate(context_id)
918 if(context_id == 0) {
923 context_id = get_transform_provider<IFFTReal_Traits>().get_context_id();
926 it_assert(in.size() > 0,
"ifft_real(): zero-sized input detected");
928 get_transform_provider<IFFTReal_Traits>().run_transform(context_id, in, out);
931void dct(
const vec &in, vec &out)
933 static int context_id = 0;
934 #pragma omp threadprivate(context_id)
936 if(context_id == 0) {
941 context_id = get_transform_provider<DCT_Traits>().get_context_id();
944 it_assert(in.size() > 0,
"dct(): zero-sized input detected");
946 get_transform_provider<DCT_Traits>().run_transform(context_id, in, out);
949void idct(
const vec &in, vec &out)
951 static int context_id = 0;
952 #pragma omp threadprivate(context_id)
954 if(context_id == 0) {
959 context_id = get_transform_provider<IDCT_Traits>().get_context_id();
962 it_assert(in.size() > 0,
"dct(): zero-sized input detected");
964 get_transform_provider<IDCT_Traits>().run_transform(context_id, in, out);
971void fft(
const cvec &in, cvec &out)
973 it_error(
"FFT library is needed to use fft() function");
976void ifft(
const cvec &in, cvec &out)
978 it_error(
"FFT library is needed to use ifft() function");
981void fft_real(
const vec &in, cvec &out)
983 it_error(
"FFT library is needed to use fft_real() function");
988 it_error(
"FFT library is needed to use ifft_real() function");
991void dct(
const vec &in, vec &out)
993 it_error(
"FFT library is needed to use dct() function");
996void idct(
const vec &in, vec &out)
998 it_error(
"FFT library is needed to use idct() function");
1007cvec
fft(
const cvec &in)
1014cvec
fft(
const cvec &in,
const int N)
1023cvec
ifft(
const cvec &in)
1030cvec
ifft(
const cvec &in,
const int N)
1046cvec
fft_real(
const vec& in,
const int N)
1062vec
ifft_real(
const cvec &in,
const int N)
1071vec
dct(
const vec &in)
1078vec
dct(
const vec &in,
const int N)
1088vec
idct(
const vec &in)
1095vec
idct(
const vec &in,
const int N)
1107template ITPP_EXPORT vec
dht(
const vec &v);
1108template ITPP_EXPORT cvec
dht(
const cvec &v);
1110template ITPP_EXPORT
void dht(
const vec &vin, vec &vout);
1111template ITPP_EXPORT
void dht(
const cvec &vin, cvec &vout);
1113template ITPP_EXPORT
void self_dht(vec &v);
1114template ITPP_EXPORT
void self_dht(cvec &v);
1116template ITPP_EXPORT vec
dwht(
const vec &v);
1117template ITPP_EXPORT cvec
dwht(
const cvec &v);
1119template ITPP_EXPORT
void dwht(
const vec &vin, vec &vout);
1120template ITPP_EXPORT
void dwht(
const cvec &vin, cvec &vout);
1122template ITPP_EXPORT
void self_dwht(vec &v);
1123template ITPP_EXPORT
void self_dwht(cvec &v);
1125template ITPP_EXPORT mat
dht2(
const mat &m);
1126template ITPP_EXPORT cmat
dht2(
const cmat &m);
1128template ITPP_EXPORT mat
dwht2(
const mat &m);
1129template ITPP_EXPORT cmat
dwht2(
const cmat &m);
void set_size(int n, bool copy=false)
Resizing an Array<T>.
friend const Array< T > concat(const Array< T > &a1, const T &e)
Append element e to the end of the Array a.
ITPP_EXPORT void idct(const vec &in, vec &out)
Inverse Discrete Cosine Transform (IDCT)
ITPP_EXPORT void dct(const vec &in, vec &out)
Discrete Cosine Transform (DCT)
ITPP_EXPORT bool have_cosine_transforms()
Run-time test if library is built with cosine transforms enabled.
#define it_error(s)
Abort unconditionally.
#define it_info(s)
Print information message.
#define it_assert(t, s)
Abort if t is not true.
ITPP_EXPORT void ifft(const cvec &in, cvec &out)
Inverse Fast Fourier Transform.
ITPP_EXPORT void fft_real(const vec &in, cvec &out)
Real Fast Fourier Transform.
ITPP_EXPORT void ifft_real(const cvec &in, vec &out)
Inverse Real Fast Fourier Transform.
ITPP_EXPORT bool have_fourier_transforms()
Run-time test if library is built with Fast Fourier Transforms enabled.
ITPP_EXPORT void fft(const cvec &in, cvec &out)
Fast Fourier Transform.
Mat< T > dwht2(const Mat< T > &m)
Fast 2D Walsh Hadamard Transform.
Vec< T > dwht(const Vec< T > &v)
Fast Walsh Hadamard Transform.
Vec< T > dht(const Vec< T > &v)
Fast Hadamard Transform.
void self_dht(Vec< T > &v)
Fast Hadamard Transform - memory efficient. Stores the result in v.
Mat< T > dht2(const Mat< T > &m)
Fast 2D Hadamard Transform.
void self_dwht(Vec< T > &v)
Fast Walsh Hadamard Transform - memory efficient (result in v)
vec imag(const cvec &data)
Imaginary part of complex values.
vec real(const cvec &data)
Real part of complex values.
cvec conj(const cvec &x)
Conjugate of complex value.
Vec< T > reverse(const Vec< T > &in)
Reverse the input vector.
cvec to_cvec(const Vec< T > &v)
Converts a Vec<T> to cvec.
const double pi
Constant Pi.
int ceil_i(double x)
The nearest larger integer.