18 #ifndef itkFFTWCommon_h
19 #define itkFFTWCommon_h
21 #if defined(ITK_USE_FFTWF) || defined(ITK_USE_FFTWD)
22 # if defined(ITK_USE_CUFFTW)
28 # if !defined(FFTW_WISDOM_ONLY)
31 # error "FFTW 3.3.2 or later is required so that FFTW_WISDOM_ONLY is defined."
54 template <
typename TPixel>
64 #if defined(ITK_USE_FFTWF)
76 # ifdef ITK_USE_CUFFTW
88 bool canDestroyInput =
false)
90 return Plan_dft_c2r(1, &n, in, out, flags, threads, canDestroyInput);
100 bool canDestroyInput =
false)
102 auto * sizes =
new int[2];
105 PlanType plan = Plan_dft_c2r(2, sizes, in, out, flags, threads, canDestroyInput);
118 bool canDestroyInput =
false)
120 auto * sizes =
new int[3];
124 PlanType plan = Plan_dft_c2r(3, sizes, in, out, flags, threads, canDestroyInput);
136 bool canDestroyInput =
false)
138 # ifndef ITK_USE_CUFFTW
140 fftwf_plan_with_nthreads(threads);
146 unsigned roflags = flags;
147 if (!(flags & FFTW_ESTIMATE))
149 roflags = flags | FFTW_WISDOM_ONLY;
151 PlanType plan = fftwf_plan_dft_c2r(rank, n, in, out, roflags);
158 plan = fftwf_plan_dft_c2r(rank, n, in, out, flags);
164 for (
int i = 0; i < rank; i++)
169 fftwf_plan_dft_c2r(rank, n, din, out, flags);
172 plan = fftwf_plan_dft_c2r(rank, n, in, out, roflags);
174 # ifndef ITK_USE_CUFFTW
178 itkAssertOrThrowMacro(plan !=
nullptr,
"PLAN_CREATION_FAILED ");
189 bool canDestroyInput =
false)
191 return Plan_dft_r2c(1, &n, in, out, flags, threads, canDestroyInput);
201 bool canDestroyInput =
false)
203 auto * sizes =
new int[2];
206 PlanType plan = Plan_dft_r2c(2, sizes, in, out, flags, threads, canDestroyInput);
219 bool canDestroyInput =
false)
221 auto * sizes =
new int[3];
225 PlanType plan = Plan_dft_r2c(3, sizes, in, out, flags, threads, canDestroyInput);
237 bool canDestroyInput =
false)
240 # ifndef ITK_USE_CUFFTW
242 fftwf_plan_with_nthreads(threads);
248 unsigned roflags = flags;
249 if (!(flags & FFTW_ESTIMATE))
251 roflags = flags | FFTW_WISDOM_ONLY;
253 PlanType plan = fftwf_plan_dft_r2c(rank, n, in, out, roflags);
260 plan = fftwf_plan_dft_r2c(rank, n, in, out, flags);
266 for (
int i = 0; i < rank; i++)
271 fftwf_plan_dft_r2c(rank, n, din, out, flags);
274 plan = fftwf_plan_dft_r2c(rank, n, in, out, roflags);
276 # ifndef ITK_USE_CUFFTW
280 itkAssertOrThrowMacro(plan !=
nullptr,
"PLAN_CREATION_FAILED ");
291 bool canDestroyInput =
false)
293 return Plan_dft(1, &n, in, out, sign, flags, threads, canDestroyInput);
304 bool canDestroyInput =
false)
306 auto * sizes =
new int[2];
309 PlanType plan = Plan_dft(2, sizes, in, out, sign, flags, threads, canDestroyInput);
323 bool canDestroyInput =
false)
325 auto * sizes =
new int[3];
329 PlanType plan = Plan_dft(3, sizes, in, out, sign, flags, threads, canDestroyInput);
342 bool canDestroyInput =
false)
344 # ifndef ITK_USE_CUFFTW
346 fftwf_plan_with_nthreads(threads);
352 unsigned roflags = flags;
353 if (!(flags & FFTW_ESTIMATE))
355 roflags = flags | FFTW_WISDOM_ONLY;
357 PlanType plan = fftwf_plan_dft(rank, n, in, out, sign, roflags);
364 plan = fftwf_plan_dft(rank, n, in, out, sign, flags);
370 for (
int i = 0; i < rank; i++)
375 fftwf_plan_dft(rank, n, din, out, sign, flags);
378 plan = fftwf_plan_dft(rank, n, in, out, sign, roflags);
380 # ifndef ITK_USE_CUFFTW
384 itkAssertOrThrowMacro(plan !=
nullptr,
"PLAN_CREATION_FAILED ");
397 # ifndef ITK_USE_CUFFTW
400 fftwf_destroy_plan(p);
404 #endif // ITK_USE_FFTWF
407 #if defined(ITK_USE_FFTWD)
418 # ifdef ITK_USE_CUFFTW
430 bool canDestroyInput =
false)
432 return Plan_dft_c2r(1, &n, in, out, flags, threads, canDestroyInput);
442 bool canDestroyInput =
false)
444 auto * sizes =
new int[2];
447 PlanType plan = Plan_dft_c2r(2, sizes, in, out, flags, threads, canDestroyInput);
460 bool canDestroyInput =
false)
462 auto * sizes =
new int[3];
466 PlanType plan = Plan_dft_c2r(3, sizes, in, out, flags, threads, canDestroyInput);
478 bool canDestroyInput =
false)
480 # ifndef ITK_USE_CUFFTW
482 fftw_plan_with_nthreads(threads);
488 unsigned roflags = flags;
489 if (!(flags & FFTW_ESTIMATE))
491 roflags = flags | FFTW_WISDOM_ONLY;
493 PlanType plan = fftw_plan_dft_c2r(rank, n, in, out, roflags);
500 plan = fftw_plan_dft_c2r(rank, n, in, out, flags);
506 for (
int i = 0; i < rank; i++)
511 fftw_plan_dft_c2r(rank, n, din, out, flags);
514 plan = fftw_plan_dft_c2r(rank, n, in, out, roflags);
516 # ifndef ITK_USE_CUFFTW
520 itkAssertOrThrowMacro(plan !=
nullptr,
"PLAN_CREATION_FAILED ");
531 bool canDestroyInput =
false)
533 return Plan_dft_r2c(1, &n, in, out, flags, threads, canDestroyInput);
543 bool canDestroyInput =
false)
545 auto * sizes =
new int[2];
548 PlanType plan = Plan_dft_r2c(2, sizes, in, out, flags, threads, canDestroyInput);
561 bool canDestroyInput =
false)
563 auto * sizes =
new int[3];
567 PlanType plan = Plan_dft_r2c(3, sizes, in, out, flags, threads, canDestroyInput);
579 bool canDestroyInput =
false)
581 # ifndef ITK_USE_CUFFTW
583 fftw_plan_with_nthreads(threads);
589 unsigned roflags = flags;
590 if (!(flags & FFTW_ESTIMATE))
592 roflags = flags | FFTW_WISDOM_ONLY;
594 PlanType plan = fftw_plan_dft_r2c(rank, n, in, out, roflags);
601 plan = fftw_plan_dft_r2c(rank, n, in, out, flags);
607 for (
int i = 0; i < rank; i++)
612 fftw_plan_dft_r2c(rank, n, din, out, flags);
615 plan = fftw_plan_dft_r2c(rank, n, in, out, roflags);
617 # ifndef ITK_USE_CUFFTW
621 itkAssertOrThrowMacro(plan !=
nullptr,
"PLAN_CREATION_FAILED ");
632 bool canDestroyInput =
false)
634 return Plan_dft(1, &n, in, out, sign, flags, threads, canDestroyInput);
645 bool canDestroyInput =
false)
647 auto * sizes =
new int[2];
650 PlanType plan = Plan_dft(2, sizes, in, out, sign, flags, threads, canDestroyInput);
664 bool canDestroyInput =
false)
666 auto * sizes =
new int[3];
670 PlanType plan = Plan_dft(3, sizes, in, out, sign, flags, threads, canDestroyInput);
683 bool canDestroyInput =
false)
685 # ifndef ITK_USE_CUFFTW
687 fftw_plan_with_nthreads(threads);
693 unsigned roflags = flags;
694 if (!(flags & FFTW_ESTIMATE))
696 roflags = flags | FFTW_WISDOM_ONLY;
698 PlanType plan = fftw_plan_dft(rank, n, in, out, sign, roflags);
705 plan = fftw_plan_dft(rank, n, in, out, sign, flags);
711 for (
int i = 0; i < rank; i++)
716 fftw_plan_dft(rank, n, din, out, sign, flags);
719 plan = fftw_plan_dft(rank, n, in, out, sign, roflags);
721 # ifndef ITK_USE_CUFFTW
725 itkAssertOrThrowMacro(plan !=
nullptr,
"PLAN_CREATION_FAILED ");
738 # ifndef ITK_USE_CUFFTW
741 fftw_destroy_plan(p);