36 #ifndef VIGRA_MULTI_FFT_HXX 37 #define VIGRA_MULTI_FFT_HXX 40 #include "multi_array.hxx" 41 #include "multi_math.hxx" 42 #include "navigator.hxx" 43 #include "copyimage.hxx" 44 #include "threading.hxx" 60 template <
unsigned int N,
class T,
class C>
61 void moveDCToCenterImpl(MultiArrayView<N, T, C> a,
unsigned int startDimension)
64 typedef MultiArrayNavigator<Traverser, N> Navigator;
65 typedef typename Navigator::iterator Iterator;
67 for(
unsigned int d = startDimension; d < N; ++d)
69 Navigator nav(a.traverser_begin(), a.shape(), d);
71 for( ; nav.hasMore(); nav++ )
73 Iterator i = nav.begin();
74 int s = nav.end() - i;
79 for(
int k=0; k<s2; ++k)
81 std::swap(i[k], i[k+s2]);
87 for(
int k=0; k<s2; ++k)
98 template <
unsigned int N,
class T,
class C>
99 void moveDCToUpperLeftImpl(MultiArrayView<N, T, C> a,
unsigned int startDimension)
102 typedef MultiArrayNavigator<Traverser, N> Navigator;
103 typedef typename Navigator::iterator Iterator;
105 for(
unsigned int d = startDimension; d < N; ++d)
107 Navigator nav(a.traverser_begin(), a.shape(), d);
109 for( ; nav.hasMore(); nav++ )
111 Iterator i = nav.begin();
112 int s = nav.end() - i;
117 for(
int k=0; k<s2; ++k)
119 std::swap(i[k], i[k+s2]);
125 for(
int k=s2; k>0; --k)
144 template <
unsigned int N,
class T,
class C>
147 detail::moveDCToCenterImpl(a, 0);
150 template <
unsigned int N,
class T,
class C>
153 detail::moveDCToUpperLeftImpl(a, 0);
156 template <
unsigned int N,
class T,
class C>
157 inline void moveDCToHalfspaceCenter(MultiArrayView<N, T, C> a)
159 detail::moveDCToCenterImpl(a, 1);
162 template <
unsigned int N,
class T,
class C>
163 inline void moveDCToHalfspaceUpperLeft(MultiArrayView<N, T, C> a)
165 detail::moveDCToUpperLeftImpl(a, 1);
171 #ifndef VIGRA_SINGLE_THREADED 173 template <
int DUMMY=0>
177 threading::lock_guard<threading::mutex> guard_;
180 : guard_(plan_mutex_)
183 static threading::mutex plan_mutex_;
187 threading::mutex FFTWLock<DUMMY>::plan_mutex_;
189 #else // VIGRA_SINGLE_THREADED 191 template <
int DUMMY=0>
200 #endif // not VIGRA_SINGLE_THREADED 203 fftwPlanCreate(
unsigned int N,
int* shape,
204 FFTWComplex<double> * in,
int* instrides,
int instep,
205 FFTWComplex<double> * out,
int* outstrides,
int outstep,
206 int sign,
unsigned int planner_flags)
208 return fftw_plan_many_dft(N, shape, 1,
209 (fftw_complex *)in, instrides, instep, 0,
210 (fftw_complex *)out, outstrides, outstep, 0,
211 sign, planner_flags);
215 fftwPlanCreate(
unsigned int N,
int* shape,
216 double * in,
int* instrides,
int instep,
217 FFTWComplex<double> * out,
int* outstrides,
int outstep,
218 int ,
unsigned int planner_flags)
220 return fftw_plan_many_dft_r2c(N, shape, 1,
221 in, instrides, instep, 0,
222 (fftw_complex *)out, outstrides, outstep, 0,
227 fftwPlanCreate(
unsigned int N,
int* shape,
228 FFTWComplex<double> * in,
int* instrides,
int instep,
229 double * out,
int* outstrides,
int outstep,
230 int ,
unsigned int planner_flags)
232 return fftw_plan_many_dft_c2r(N, shape, 1,
233 (fftw_complex *)in, instrides, instep, 0,
234 out, outstrides, outstep, 0,
239 fftwPlanCreate(
unsigned int N,
int* shape,
240 FFTWComplex<float> * in,
int* instrides,
int instep,
241 FFTWComplex<float> * out,
int* outstrides,
int outstep,
242 int sign,
unsigned int planner_flags)
244 return fftwf_plan_many_dft(N, shape, 1,
245 (fftwf_complex *)in, instrides, instep, 0,
246 (fftwf_complex *)out, outstrides, outstep, 0,
247 sign, planner_flags);
251 fftwPlanCreate(
unsigned int N,
int* shape,
252 float * in,
int* instrides,
int instep,
253 FFTWComplex<float> * out,
int* outstrides,
int outstep,
254 int ,
unsigned int planner_flags)
256 return fftwf_plan_many_dft_r2c(N, shape, 1,
257 in, instrides, instep, 0,
258 (fftwf_complex *)out, outstrides, outstep, 0,
263 fftwPlanCreate(
unsigned int N,
int* shape,
264 FFTWComplex<float> * in,
int* instrides,
int instep,
265 float * out,
int* outstrides,
int outstep,
266 int ,
unsigned int planner_flags)
268 return fftwf_plan_many_dft_c2r(N, shape, 1,
269 (fftwf_complex *)in, instrides, instep, 0,
270 out, outstrides, outstep, 0,
275 fftwPlanCreate(
unsigned int N,
int* shape,
276 FFTWComplex<long double> * in,
int* instrides,
int instep,
277 FFTWComplex<long double> * out,
int* outstrides,
int outstep,
278 int sign,
unsigned int planner_flags)
280 return fftwl_plan_many_dft(N, shape, 1,
281 (fftwl_complex *)in, instrides, instep, 0,
282 (fftwl_complex *)out, outstrides, outstep, 0,
283 sign, planner_flags);
287 fftwPlanCreate(
unsigned int N,
int* shape,
288 long double * in,
int* instrides,
int instep,
289 FFTWComplex<long double> * out,
int* outstrides,
int outstep,
290 int ,
unsigned int planner_flags)
292 return fftwl_plan_many_dft_r2c(N, shape, 1,
293 in, instrides, instep, 0,
294 (fftwl_complex *)out, outstrides, outstep, 0,
299 fftwPlanCreate(
unsigned int N,
int* shape,
300 FFTWComplex<long double> * in,
int* instrides,
int instep,
301 long double * out,
int* outstrides,
int outstep,
302 int ,
unsigned int planner_flags)
304 return fftwl_plan_many_dft_c2r(N, shape, 1,
305 (fftwl_complex *)in, instrides, instep, 0,
306 out, outstrides, outstep, 0,
310 inline void fftwPlanDestroy(fftw_plan plan)
313 fftw_destroy_plan(plan);
316 inline void fftwPlanDestroy(fftwf_plan plan)
319 fftwf_destroy_plan(plan);
322 inline void fftwPlanDestroy(fftwl_plan plan)
325 fftwl_destroy_plan(plan);
329 fftwPlanExecute(fftw_plan plan)
335 fftwPlanExecute(fftwf_plan plan)
341 fftwPlanExecute(fftwl_plan plan)
347 fftwPlanExecute(fftw_plan plan, FFTWComplex<double> * in, FFTWComplex<double> * out)
349 fftw_execute_dft(plan, (fftw_complex *)in, (fftw_complex *)out);
353 fftwPlanExecute(fftw_plan plan,
double * in, FFTWComplex<double> * out)
355 fftw_execute_dft_r2c(plan, in, (fftw_complex *)out);
359 fftwPlanExecute(fftw_plan plan, FFTWComplex<double> * in,
double * out)
361 fftw_execute_dft_c2r(plan, (fftw_complex *)in, out);
365 fftwPlanExecute(fftwf_plan plan, FFTWComplex<float> * in, FFTWComplex<float> * out)
367 fftwf_execute_dft(plan, (fftwf_complex *)in, (fftwf_complex *)out);
371 fftwPlanExecute(fftwf_plan plan,
float * in, FFTWComplex<float> * out)
373 fftwf_execute_dft_r2c(plan, in, (fftwf_complex *)out);
377 fftwPlanExecute(fftwf_plan plan, FFTWComplex<float> * in,
float * out)
379 fftwf_execute_dft_c2r(plan, (fftwf_complex *)in, out);
383 fftwPlanExecute(fftwl_plan plan, FFTWComplex<long double> * in, FFTWComplex<long double> * out)
385 fftwl_execute_dft(plan, (fftwl_complex *)in, (fftwl_complex *)out);
389 fftwPlanExecute(fftwl_plan plan,
long double * in, FFTWComplex<long double> * out)
391 fftwl_execute_dft_r2c(plan, in, (fftwl_complex *)out);
395 fftwPlanExecute(fftwl_plan plan, FFTWComplex<long double> * in,
long double * out)
397 fftwl_execute_dft_c2r(plan, (fftwl_complex *)in, out);
401 struct FFTWPaddingSize
403 static const int size = 1330;
404 static const int evenSize = 1063;
405 static int goodSizes[size];
406 static int goodEvenSizes[evenSize];
408 static inline int find(
int s)
410 if(s <= 0 || s >= goodSizes[size-1])
413 int * upperBound = std::upper_bound(goodSizes, goodSizes+size, s);
414 if(upperBound > goodSizes && upperBound[-1] == s)
420 static inline int findEven(
int s)
422 if(s <= 0 || s >= goodEvenSizes[evenSize-1])
425 int * upperBound = std::upper_bound(goodEvenSizes, goodEvenSizes+evenSize, s);
426 if(upperBound > goodEvenSizes && upperBound[-1] == s)
438 int FFTWPaddingSize<DUMMY>::goodSizes[size] = {
439 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16,
440 18, 20, 21, 22, 24, 25, 26, 27, 28, 30, 32, 33, 35, 36, 39, 40, 42, 44, 45, 48,
441 49, 50, 52, 54, 55, 56, 60, 63, 64, 65, 66, 70, 72, 75, 77, 78, 80, 81,
442 84, 88, 90, 91, 96, 98, 99, 100, 104, 105, 108, 110, 112, 117, 120, 125,
443 126, 128, 130, 132, 135, 140, 144, 147, 150, 154, 156, 160, 162, 165,
444 168, 175, 176, 180, 182, 189, 192, 195, 196, 198, 200, 208, 210, 216,
445 220, 224, 225, 231, 234, 240, 243, 245, 250, 252, 256, 260, 264, 270,
446 273, 275, 280, 288, 294, 297, 300, 308, 312, 315, 320, 324, 325, 330,
447 336, 343, 350, 351, 352, 360, 364, 375, 378, 384, 385, 390, 392, 396,
448 400, 405, 416, 420, 432, 440, 441, 448, 450, 455, 462, 468, 480, 486,
449 490, 495, 500, 504, 512, 520, 525, 528, 539, 540, 546, 550, 560, 567,
450 576, 585, 588, 594, 600, 616, 624, 625, 630, 637, 640, 648, 650, 660,
451 672, 675, 686, 693, 700, 702, 704, 720, 728, 729, 735, 750, 756, 768,
452 770, 780, 784, 792, 800, 810, 819, 825, 832, 840, 864, 875, 880, 882,
453 891, 896, 900, 910, 924, 936, 945, 960, 972, 975, 980, 990, 1000, 1008,
454 1024, 1029, 1040, 1050, 1053, 1056, 1078, 1080, 1092, 1100, 1120, 1125,
455 1134, 1152, 1155, 1170, 1176, 1188, 1200, 1215, 1225, 1232, 1248, 1250,
456 1260, 1274, 1280, 1296, 1300, 1320, 1323, 1344, 1350, 1365, 1372, 1375,
457 1386, 1400, 1404, 1408, 1440, 1456, 1458, 1470, 1485, 1500, 1512, 1536,
458 1540, 1560, 1568, 1575, 1584, 1600, 1617, 1620, 1625, 1638, 1650, 1664,
459 1680, 1701, 1715, 1728, 1750, 1755, 1760, 1764, 1782, 1792, 1800, 1820,
460 1848, 1872, 1875, 1890, 1911, 1920, 1925, 1944, 1950, 1960, 1980, 2000,
461 2016, 2025, 2048, 2058, 2079, 2080, 2100, 2106, 2112, 2156, 2160, 2184,
462 2187, 2200, 2205, 2240, 2250, 2268, 2275, 2304, 2310, 2340, 2352, 2376,
463 2400, 2401, 2430, 2450, 2457, 2464, 2475, 2496, 2500, 2520, 2548, 2560,
464 2592, 2600, 2625, 2640, 2646, 2673, 2688, 2695, 2700, 2730, 2744, 2750,
465 2772, 2800, 2808, 2816, 2835, 2880, 2912, 2916, 2925, 2940, 2970, 3000,
466 3024, 3072, 3080, 3087, 3120, 3125, 3136, 3150, 3159, 3168, 3185, 3200,
467 3234, 3240, 3250, 3276, 3300, 3328, 3360, 3375, 3402, 3430, 3456, 3465,
468 3500, 3510, 3520, 3528, 3564, 3584, 3600, 3640, 3645, 3675, 3696, 3744,
469 3750, 3773, 3780, 3822, 3840, 3850, 3888, 3900, 3920, 3960, 3969, 4000,
470 4032, 4050, 4095, 4096, 4116, 4125, 4158, 4160, 4200, 4212, 4224, 4312,
471 4320, 4368, 4374, 4375, 4400, 4410, 4455, 4459, 4480, 4500, 4536, 4550,
472 4608, 4620, 4680, 4704, 4725, 4752, 4800, 4802, 4851, 4860, 4875, 4900,
473 4914, 4928, 4950, 4992, 5000, 5040, 5096, 5103, 5120, 5145, 5184, 5200,
474 5250, 5265, 5280, 5292, 5346, 5376, 5390, 5400, 5460, 5488, 5500, 5544,
475 5600, 5616, 5625, 5632, 5670, 5733, 5760, 5775, 5824, 5832, 5850, 5880,
476 5940, 6000, 6048, 6075, 6125, 6144, 6160, 6174, 6237, 6240, 6250, 6272,
477 6300, 6318, 6336, 6370, 6400, 6468, 6480, 6500, 6552, 6561, 6600, 6615,
478 6656, 6720, 6750, 6804, 6825, 6860, 6875, 6912, 6930, 7000, 7020, 7040,
479 7056, 7128, 7168, 7200, 7203, 7280, 7290, 7350, 7371, 7392, 7425, 7488,
480 7500, 7546, 7560, 7644, 7680, 7700, 7776, 7800, 7840, 7875, 7920, 7938,
481 8000, 8019, 8064, 8085, 8100, 8125, 8190, 8192, 8232, 8250, 8316, 8320,
482 8400, 8424, 8448, 8505, 8575, 8624, 8640, 8736, 8748, 8750, 8775, 8800,
483 8820, 8910, 8918, 8960, 9000, 9072, 9100, 9216, 9240, 9261, 9360, 9375,
484 9408, 9450, 9477, 9504, 9555, 9600, 9604, 9625, 9702, 9720, 9750, 9800,
485 9828, 9856, 9900, 9984, 10000, 10080, 10125, 10192, 10206, 10240, 10290,
486 10368, 10395, 10400, 10500, 10530, 10560, 10584, 10692, 10752, 10780,
487 10800, 10920, 10935, 10976, 11000, 11025, 11088, 11200, 11232, 11250,
488 11264, 11319, 11340, 11375, 11466, 11520, 11550, 11648, 11664, 11700,
489 11760, 11880, 11907, 12000, 12005, 12096, 12150, 12250, 12285, 12288,
490 12320, 12348, 12375, 12474, 12480, 12500, 12544, 12600, 12636, 12672,
491 12740, 12800, 12936, 12960, 13000, 13104, 13122, 13125, 13200, 13230,
492 13312, 13365, 13377, 13440, 13475, 13500, 13608, 13650, 13720, 13750,
493 13824, 13860, 14000, 14040, 14080, 14112, 14175, 14256, 14336, 14400,
494 14406, 14553, 14560, 14580, 14625, 14700, 14742, 14784, 14850, 14976,
495 15000, 15092, 15120, 15288, 15309, 15360, 15400, 15435, 15552, 15600,
496 15625, 15680, 15750, 15795, 15840, 15876, 15925, 16000, 16038, 16128,
497 16170, 16200, 16250, 16380, 16384, 16464, 16500, 16632, 16640, 16800,
498 16807, 16848, 16875, 16896, 17010, 17150, 17199, 17248, 17280, 17325,
499 17472, 17496, 17500, 17550, 17600, 17640, 17820, 17836, 17920, 18000,
500 18144, 18200, 18225, 18375, 18432, 18480, 18522, 18711, 18720, 18750,
501 18816, 18865, 18900, 18954, 19008, 19110, 19200, 19208, 19250, 19404,
502 19440, 19500, 19600, 19656, 19683, 19712, 19800, 19845, 19968, 20000,
503 20160, 20250, 20384, 20412, 20475, 20480, 20580, 20625, 20736, 20790,
504 20800, 21000, 21060, 21120, 21168, 21384, 21504, 21560, 21600, 21609,
505 21840, 21870, 21875, 21952, 22000, 22050, 22113, 22176, 22275, 22295,
506 22400, 22464, 22500, 22528, 22638, 22680, 22750, 22932, 23040, 23100,
507 23296, 23328, 23400, 23520, 23625, 23760, 23814, 24000, 24010, 24057,
508 24192, 24255, 24300, 24375, 24500, 24570, 24576, 24640, 24696, 24750,
509 24948, 24960, 25000, 25088, 25200, 25272, 25344, 25480, 25515, 25600,
510 25725, 25872, 25920, 26000, 26208, 26244, 26250, 26325, 26400, 26411,
511 26460, 26624, 26730, 26754, 26880, 26950, 27000, 27216, 27300, 27440,
512 27500, 27648, 27720, 27783, 28000, 28080, 28125, 28160, 28224, 28350,
513 28431, 28512, 28665, 28672, 28800, 28812, 28875, 29106, 29120, 29160,
514 29250, 29400, 29484, 29568, 29700, 29952, 30000, 30184, 30240, 30375,
515 30576, 30618, 30625, 30720, 30800, 30870, 31104, 31185, 31200, 31213,
516 31250, 31360, 31500, 31590, 31680, 31752, 31850, 32000, 32076, 32256,
517 32340, 32400, 32500, 32760, 32768, 32805, 32928, 33000, 33075, 33264,
518 33280, 33600, 33614, 33696, 33750, 33792, 33957, 34020, 34125, 34300,
519 34375, 34398, 34496, 34560, 34650, 34944, 34992, 35000, 35100, 35200,
520 35280, 35640, 35672, 35721, 35840, 36000, 36015, 36288, 36400, 36450,
521 36750, 36855, 36864, 36960, 37044, 37125, 37422, 37440, 37500, 37632,
522 37730, 37800, 37908, 38016, 38220, 38400, 38416, 38500, 38808, 38880,
523 39000, 39200, 39312, 39366, 39375, 39424, 39600, 39690, 39936, 40000,
524 40095, 40131, 40320, 40425, 40500, 40625, 40768, 40824, 40950, 40960,
525 41160, 41250, 41472, 41580, 41600, 42000, 42120, 42240, 42336, 42525,
526 42768, 42875, 43008, 43120, 43200, 43218, 43659, 43680, 43740, 43750,
527 43875, 43904, 44000, 44100, 44226, 44352, 44550, 44590, 44800, 44928,
528 45000, 45056, 45276, 45360, 45500, 45864, 45927, 46080, 46200, 46305,
529 46592, 46656, 46800, 46875, 47040, 47250, 47385, 47520, 47628, 47775,
530 48000, 48020, 48114, 48125, 48384, 48510, 48600, 48750, 49000, 49140,
531 49152, 49280, 49392, 49500, 49896, 49920, 50000, 50176, 50400, 50421,
532 50544, 50625, 50688, 50960, 51030, 51200, 51450, 51597, 51744, 51840,
533 51975, 52000, 52416, 52488, 52500, 52650, 52800, 52822, 52920, 53248,
534 53460, 53508, 53760, 53900, 54000, 54432, 54600, 54675, 54880, 55000,
535 55125, 55296, 55440, 55566, 56000, 56133, 56160, 56250, 56320, 56448,
536 56595, 56700, 56862, 56875, 57024, 57330, 57344, 57600, 57624, 57750,
537 58212, 58240, 58320, 58500, 58800, 58968, 59049, 59136, 59400, 59535,
538 59904, 60000, 60025, 60368, 60480, 60750, 61152, 61236, 61250, 61425,
539 61440, 61600, 61740, 61875, 62208, 62370, 62400, 62426, 62500, 62720,
540 63000, 63180, 63360, 63504, 63700, 64000, 64152, 64512, 64680, 64800,
541 64827, 65000, 65520, 65536, 65610, 65625, 65856, 66000, 66150, 66339,
542 66528, 66560, 66825, 66885, 67200, 67228, 67375, 67392, 67500, 67584,
543 67914, 68040, 68250, 68600, 68750, 68796, 68992, 69120, 69300, 69888,
544 69984, 70000, 70200, 70400, 70560, 70875, 71280, 71344, 71442, 71680,
545 72000, 72030, 72171, 72576, 72765, 72800, 72900, 73125, 73500, 73710,
546 73728, 73920, 74088, 74250, 74844, 74880, 75000, 75264, 75460, 75600,
547 75816, 76032, 76440, 76545, 76800, 76832, 77000, 77175, 77616, 77760,
548 78000, 78125, 78400, 78624, 78732, 78750, 78848, 78975, 79200, 79233,
549 79380, 79625, 79872, 80000, 80190, 80262, 80640, 80850, 81000, 81250,
550 81536, 81648, 81900, 81920, 82320, 82500, 82944, 83160, 83200, 83349,
551 84000, 84035, 84240, 84375, 84480, 84672, 85050, 85293, 85536, 85750,
552 85995, 86016, 86240, 86400, 86436, 86625, 87318, 87360, 87480, 87500,
553 87750, 87808, 88000, 88200, 88452, 88704, 89100, 89180, 89600, 89856,
554 90000, 90112, 90552, 90720, 91000, 91125, 91728, 91854, 91875, 92160,
555 92400, 92610, 93184, 93312, 93555, 93600, 93639, 93750, 94080, 94325,
556 94500, 94770, 95040, 95256, 95550, 96000, 96040, 96228, 96250, 96768,
557 97020, 97200, 97500, 98000, 98280, 98304, 98415, 98560, 98784, 99000,
562 int FFTWPaddingSize<DUMMY>::goodEvenSizes[evenSize] = {
563 2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22,
564 24, 26, 28, 30, 32, 36, 40, 42, 44, 48, 50, 52, 54, 56, 60, 64, 66, 70,
565 72, 78, 80, 84, 88, 90, 96, 98, 100, 104, 108, 110, 112, 120, 126, 128,
566 130, 132, 140, 144, 150, 154, 156, 160, 162, 168, 176, 180, 182, 192,
567 196, 198, 200, 208, 210, 216, 220, 224, 234, 240, 250, 252, 256, 260,
568 264, 270, 280, 288, 294, 300, 308, 312, 320, 324, 330, 336, 350, 352,
569 360, 364, 378, 384, 390, 392, 396, 400, 416, 420, 432, 440, 448, 450,
570 462, 468, 480, 486, 490, 500, 504, 512, 520, 528, 540, 546, 550, 560,
571 576, 588, 594, 600, 616, 624, 630, 640, 648, 650, 660, 672, 686, 700,
572 702, 704, 720, 728, 750, 756, 768, 770, 780, 784, 792, 800, 810, 832,
573 840, 864, 880, 882, 896, 900, 910, 924, 936, 960, 972, 980, 990, 1000,
574 1008, 1024, 1040, 1050, 1056, 1078, 1080, 1092, 1100, 1120, 1134, 1152,
575 1170, 1176, 1188, 1200, 1232, 1248, 1250, 1260, 1274, 1280, 1296, 1300,
576 1320, 1344, 1350, 1372, 1386, 1400, 1404, 1408, 1440, 1456, 1458, 1470,
577 1500, 1512, 1536, 1540, 1560, 1568, 1584, 1600, 1620, 1638, 1650, 1664,
578 1680, 1728, 1750, 1760, 1764, 1782, 1792, 1800, 1820, 1848, 1872, 1890,
579 1920, 1944, 1950, 1960, 1980, 2000, 2016, 2048, 2058, 2080, 2100, 2106,
580 2112, 2156, 2160, 2184, 2200, 2240, 2250, 2268, 2304, 2310, 2340, 2352,
581 2376, 2400, 2430, 2450, 2464, 2496, 2500, 2520, 2548, 2560, 2592, 2600,
582 2640, 2646, 2688, 2700, 2730, 2744, 2750, 2772, 2800, 2808, 2816, 2880,
583 2912, 2916, 2940, 2970, 3000, 3024, 3072, 3080, 3120, 3136, 3150, 3168,
584 3200, 3234, 3240, 3250, 3276, 3300, 3328, 3360, 3402, 3430, 3456, 3500,
585 3510, 3520, 3528, 3564, 3584, 3600, 3640, 3696, 3744, 3750, 3780, 3822,
586 3840, 3850, 3888, 3900, 3920, 3960, 4000, 4032, 4050, 4096, 4116, 4158,
587 4160, 4200, 4212, 4224, 4312, 4320, 4368, 4374, 4400, 4410, 4480, 4500,
588 4536, 4550, 4608, 4620, 4680, 4704, 4752, 4800, 4802, 4860, 4900, 4914,
589 4928, 4950, 4992, 5000, 5040, 5096, 5120, 5184, 5200, 5250, 5280, 5292,
590 5346, 5376, 5390, 5400, 5460, 5488, 5500, 5544, 5600, 5616, 5632, 5670,
591 5760, 5824, 5832, 5850, 5880, 5940, 6000, 6048, 6144, 6160, 6174, 6240,
592 6250, 6272, 6300, 6318, 6336, 6370, 6400, 6468, 6480, 6500, 6552, 6600,
593 6656, 6720, 6750, 6804, 6860, 6912, 6930, 7000, 7020, 7040, 7056, 7128,
594 7168, 7200, 7280, 7290, 7350, 7392, 7488, 7500, 7546, 7560, 7644, 7680,
595 7700, 7776, 7800, 7840, 7920, 7938, 8000, 8064, 8100, 8190, 8192, 8232,
596 8250, 8316, 8320, 8400, 8424, 8448, 8624, 8640, 8736, 8748, 8750, 8800,
597 8820, 8910, 8918, 8960, 9000, 9072, 9100, 9216, 9240, 9360, 9408, 9450,
598 9504, 9600, 9604, 9702, 9720, 9750, 9800, 9828, 9856, 9900, 9984, 10000,
599 10080, 10192, 10206, 10240, 10290, 10368, 10400, 10500, 10530, 10560,
600 10584, 10692, 10752, 10780, 10800, 10920, 10976, 11000, 11088, 11200,
601 11232, 11250, 11264, 11340, 11466, 11520, 11550, 11648, 11664, 11700,
602 11760, 11880, 12000, 12096, 12150, 12250, 12288, 12320, 12348, 12474,
603 12480, 12500, 12544, 12600, 12636, 12672, 12740, 12800, 12936, 12960,
604 13000, 13104, 13122, 13200, 13230, 13312, 13440, 13500, 13608, 13650,
605 13720, 13750, 13824, 13860, 14000, 14040, 14080, 14112, 14256, 14336,
606 14400, 14406, 14560, 14580, 14700, 14742, 14784, 14850, 14976, 15000,
607 15092, 15120, 15288, 15360, 15400, 15552, 15600, 15680, 15750, 15840,
608 15876, 16000, 16038, 16128, 16170, 16200, 16250, 16380, 16384, 16464,
609 16500, 16632, 16640, 16800, 16848, 16896, 17010, 17150, 17248, 17280,
610 17472, 17496, 17500, 17550, 17600, 17640, 17820, 17836, 17920, 18000,
611 18144, 18200, 18432, 18480, 18522, 18720, 18750, 18816, 18900, 18954,
612 19008, 19110, 19200, 19208, 19250, 19404, 19440, 19500, 19600, 19656,
613 19712, 19800, 19968, 20000, 20160, 20250, 20384, 20412, 20480, 20580,
614 20736, 20790, 20800, 21000, 21060, 21120, 21168, 21384, 21504, 21560,
615 21600, 21840, 21870, 21952, 22000, 22050, 22176, 22400, 22464, 22500,
616 22528, 22638, 22680, 22750, 22932, 23040, 23100, 23296, 23328, 23400,
617 23520, 23760, 23814, 24000, 24010, 24192, 24300, 24500, 24570, 24576,
618 24640, 24696, 24750, 24948, 24960, 25000, 25088, 25200, 25272, 25344,
619 25480, 25600, 25872, 25920, 26000, 26208, 26244, 26250, 26400, 26460,
620 26624, 26730, 26754, 26880, 26950, 27000, 27216, 27300, 27440, 27500,
621 27648, 27720, 28000, 28080, 28160, 28224, 28350, 28512, 28672, 28800,
622 28812, 29106, 29120, 29160, 29250, 29400, 29484, 29568, 29700, 29952,
623 30000, 30184, 30240, 30576, 30618, 30720, 30800, 30870, 31104, 31200,
624 31250, 31360, 31500, 31590, 31680, 31752, 31850, 32000, 32076, 32256,
625 32340, 32400, 32500, 32760, 32768, 32928, 33000, 33264, 33280, 33600,
626 33614, 33696, 33750, 33792, 34020, 34300, 34398, 34496, 34560, 34650,
627 34944, 34992, 35000, 35100, 35200, 35280, 35640, 35672, 35840, 36000,
628 36288, 36400, 36450, 36750, 36864, 36960, 37044, 37422, 37440, 37500,
629 37632, 37730, 37800, 37908, 38016, 38220, 38400, 38416, 38500, 38808,
630 38880, 39000, 39200, 39312, 39366, 39424, 39600, 39690, 39936, 40000,
631 40320, 40500, 40768, 40824, 40950, 40960, 41160, 41250, 41472, 41580,
632 41600, 42000, 42120, 42240, 42336, 42768, 43008, 43120, 43200, 43218,
633 43680, 43740, 43750, 43904, 44000, 44100, 44226, 44352, 44550, 44590,
634 44800, 44928, 45000, 45056, 45276, 45360, 45500, 45864, 46080, 46200,
635 46592, 46656, 46800, 47040, 47250, 47520, 47628, 48000, 48020, 48114,
636 48384, 48510, 48600, 48750, 49000, 49140, 49152, 49280, 49392, 49500,
637 49896, 49920, 50000, 50176, 50400, 50544, 50688, 50960, 51030, 51200,
638 51450, 51744, 51840, 52000, 52416, 52488, 52500, 52650, 52800, 52822,
639 52920, 53248, 53460, 53508, 53760, 53900, 54000, 54432, 54600, 54880,
640 55000, 55296, 55440, 55566, 56000, 56160, 56250, 56320, 56448, 56700,
641 56862, 57024, 57330, 57344, 57600, 57624, 57750, 58212, 58240, 58320,
642 58500, 58800, 58968, 59136, 59400, 59904, 60000, 60368, 60480, 60750,
643 61152, 61236, 61250, 61440, 61600, 61740, 62208, 62370, 62400, 62426,
644 62500, 62720, 63000, 63180, 63360, 63504, 63700, 64000, 64152, 64512,
645 64680, 64800, 65000, 65520, 65536, 65610, 65856, 66000, 66150, 66528,
646 66560, 67200, 67228, 67392, 67500, 67584, 67914, 68040, 68250, 68600,
647 68750, 68796, 68992, 69120, 69300, 69888, 69984, 70000, 70200, 70400,
648 70560, 71280, 71344, 71442, 71680, 72000, 72030, 72576, 72800, 72900,
649 73500, 73710, 73728, 73920, 74088, 74250, 74844, 74880, 75000, 75264,
650 75460, 75600, 75816, 76032, 76440, 76800, 76832, 77000, 77616, 77760,
651 78000, 78400, 78624, 78732, 78750, 78848, 79200, 79380, 79872, 80000,
652 80190, 80262, 80640, 80850, 81000, 81250, 81536, 81648, 81900, 81920,
653 82320, 82500, 82944, 83160, 83200, 84000, 84240, 84480, 84672, 85050,
654 85536, 85750, 86016, 86240, 86400, 86436, 87318, 87360, 87480, 87500,
655 87750, 87808, 88000, 88200, 88452, 88704, 89100, 89180, 89600, 89856,
656 90000, 90112, 90552, 90720, 91000, 91728, 91854, 92160, 92400, 92610,
657 93184, 93312, 93600, 93750, 94080, 94500, 94770, 95040, 95256, 95550,
658 96000, 96040, 96228, 96250, 96768, 97020, 97200, 97500, 98000, 98280,
659 98304, 98560, 98784, 99000, 99792, 99840
663 struct FFTEmbedKernel
665 template <
unsigned int N,
class Real,
class C,
class Shape>
667 exec(MultiArrayView<N, Real, C> & out, Shape
const & kernelShape,
668 Shape & srcPoint, Shape & destPoint,
bool copyIt)
670 for(srcPoint[M]=0; srcPoint[M]<kernelShape[M]; ++srcPoint[M])
672 if(srcPoint[M] < (kernelShape[M] + 1) / 2)
674 destPoint[M] = srcPoint[M];
678 destPoint[M] = srcPoint[M] + out.shape(M) - kernelShape[M];
681 FFTEmbedKernel<M-1>::exec(out, kernelShape, srcPoint, destPoint, copyIt);
687 struct FFTEmbedKernel<0>
689 template <
unsigned int N,
class Real,
class C,
class Shape>
691 exec(MultiArrayView<N, Real, C> & out, Shape
const & kernelShape,
692 Shape & srcPoint, Shape & destPoint,
bool copyIt)
694 for(srcPoint[0]=0; srcPoint[0]<kernelShape[0]; ++srcPoint[0])
696 if(srcPoint[0] < (kernelShape[0] + 1) / 2)
698 destPoint[0] = srcPoint[0];
702 destPoint[0] = srcPoint[0] + out.shape(0) - kernelShape[0];
707 out[destPoint] = out[srcPoint];
714 template <
unsigned int N,
class Real,
class C1,
class C2>
716 fftEmbedKernel(MultiArrayView<N, Real, C1> kernel,
717 MultiArrayView<N, Real, C2> out,
722 MultiArrayView<N, Real, C2> kout = out.subarray(Shape(), kernel.shape());
730 Shape srcPoint, destPoint;
731 FFTEmbedKernel<(int)N-1>::exec(out, kernel.shape(), srcPoint, destPoint,
false);
734 template <
unsigned int N,
class Real,
class C1,
class C2>
736 fftEmbedArray(MultiArrayView<N, Real, C1> in,
737 MultiArrayView<N, Real, C2> out)
741 Shape diff = out.shape() - in.shape(),
743 rightDiff = diff - leftDiff,
744 right = in.shape() + leftDiff;
746 out.subarray(leftDiff, right) = in;
749 typedef MultiArrayNavigator<Traverser, N> Navigator;
750 typedef typename Navigator::iterator Iterator;
752 for(
unsigned int d = 0; d < N; ++d)
754 Navigator nav(out.traverser_begin(), out.shape(), d);
756 for( ; nav.hasMore(); nav++ )
758 Iterator i = nav.begin();
759 for(
int k=1; k<=leftDiff[d]; ++k)
760 i[leftDiff[d] - k] = i[leftDiff[d] + k];
761 for(
int k=0; k<rightDiff[d]; ++k)
762 i[right[d] + k] = i[right[d] - k - 2];
769 template <
class T,
int N>
771 fftwBestPaddedShape(TinyVector<T, N> shape)
773 for(
unsigned int k=0; k<N; ++k)
774 shape[k] = detail::FFTWPaddingSize<0>::find(shape[k]);
778 template <
class T,
int N>
780 fftwBestPaddedShapeR2C(TinyVector<T, N> shape)
782 shape[0] = detail::FFTWPaddingSize<0>::findEven(shape[0]);
783 for(
unsigned int k=1; k<N; ++k)
784 shape[k] = detail::FFTWPaddingSize<0>::find(shape[k]);
799 template <
class T,
int N>
803 shape[0] = shape[0] / 2 + 1;
807 template <
class T,
int N>
809 fftwCorrespondingShapeC2R(
TinyVector<T, N> shape,
bool oddDimension0 =
false)
811 shape[0] = oddDimension0
812 ? (shape[0] - 1) * 2 + 1
813 : (shape[0] - 1) * 2;
850 template <
unsigned int N,
class Real =
double>
854 typedef typename FFTWReal2Complex<Real>::plan_type PlanType;
858 Shape shape, instrides, outstrides;
878 template <
class C1,
class C2>
881 int SIGN,
unsigned int planner_flags = FFTW_ESTIMATE)
884 init(in, out, SIGN, planner_flags);
898 template <
class C1,
class C2>
901 unsigned int planner_flags = FFTW_ESTIMATE)
904 init(in, out, planner_flags);
918 template <
class C1,
class C2>
921 unsigned int planner_flags = FFTW_ESTIMATE)
924 init(in, out, planner_flags);
935 instrides.swap(o.instrides);
936 outstrides.swap(o.outstrides);
949 instrides.swap(o.instrides);
950 outstrides.swap(o.outstrides);
961 detail::FFTWLock<> lock;
962 detail::fftwPlanDestroy(plan);
969 template <
class C1,
class C2>
972 int SIGN,
unsigned int planner_flags = FFTW_ESTIMATE)
974 vigra_precondition(in.strideOrdering() == out.strideOrdering(),
975 "FFTWPlan.init(): input and output must have the same stride ordering.");
977 initImpl(in.permuteStridesDescending(), out.permuteStridesDescending(),
978 SIGN, planner_flags);
985 template <
class C1,
class C2>
988 unsigned int planner_flags = FFTW_ESTIMATE)
991 "FFTWPlan.init(): input and output must have the same stride ordering.");
994 FFTW_FORWARD, planner_flags);
1001 template <
class C1,
class C2>
1004 unsigned int planner_flags = FFTW_ESTIMATE)
1007 "FFTWPlan.init(): input and output must have the same stride ordering.");
1010 FFTW_BACKWARD, planner_flags);
1020 template <
class C1,
class C2>
1024 executeImpl(in.permuteStridesDescending(), out.permuteStridesDescending());
1034 template <
class C1,
class C2>
1048 template <
class C1,
class C2>
1057 template <
class MI,
class MO>
1058 void initImpl(MI ins, MO outs,
int SIGN,
unsigned int planner_flags);
1060 template <
class MI,
class MO>
1061 void executeImpl(MI ins, MO outs)
const;
1066 vigra_precondition(in.shape() == out.shape(),
1067 "FFTWPlan.init(): input and output must have the same shape.");
1073 for(
int k=0; k<(int)N-1; ++k)
1074 vigra_precondition(ins.
shape(k) == outs.shape(k),
1075 "FFTWPlan.init(): input and output must have matching shapes.");
1076 vigra_precondition(ins.
shape(N-1) / 2 + 1 == outs.shape(N-1),
1077 "FFTWPlan.init(): input and output must have matching shapes.");
1083 for(
int k=0; k<(int)N-1; ++k)
1084 vigra_precondition(ins.
shape(k) == outs.
shape(k),
1085 "FFTWPlan.init(): input and output must have matching shapes.");
1086 vigra_precondition(outs.
shape(N-1) / 2 + 1 == ins.
shape(N-1),
1087 "FFTWPlan.init(): input and output must have matching shapes.");
1091 template <
unsigned int N,
class Real>
1092 template <
class MI,
class MO>
1096 checkShapes(ins, outs);
1102 Shape newShape(logicalShape.begin(), logicalShape.end()),
1103 newIStrides(ins.stride().begin(), ins.stride().end()),
1104 newOStrides(outs.stride().begin(), outs.stride().end()),
1105 itotal(ins.shape().begin(), ins.shape().end()),
1106 ototal(outs.shape().begin(), outs.shape().end());
1108 for(
unsigned int j=1; j<N; ++j)
1110 itotal[j] = ins.stride(j-1) / ins.stride(j);
1111 ototal[j] = outs.stride(j-1) / outs.stride(j);
1115 detail::FFTWLock<> lock;
1116 PlanType newPlan = detail::fftwPlanCreate(N, newShape.begin(),
1117 ins.data(), itotal.begin(), ins.stride(N-1),
1118 outs.data(), ototal.begin(), outs.stride(N-1),
1119 SIGN, planner_flags);
1120 detail::fftwPlanDestroy(plan);
1124 shape.swap(newShape);
1125 instrides.swap(newIStrides);
1126 outstrides.swap(newOStrides);
1130 template <
unsigned int N,
class Real>
1131 template <
class MI,
class MO>
1134 vigra_precondition(plan != 0,
"FFTWPlan::execute(): plan is NULL.");
1141 "FFTWPlan::execute(): shape mismatch between plan and data.");
1143 "FFTWPlan::execute(): strides mismatch between plan and input data.");
1145 "FFTWPlan::execute(): strides mismatch between plan and output data.");
1147 detail::fftwPlanExecute(plan, ins.data(), outs.data());
1149 typedef typename MO::value_type V;
1150 if(sign == FFTW_BACKWARD)
1151 outs *= V(1.0) / Real(outs.size());
1193 template <
unsigned int N,
class Real>
1201 RArray realArray, realKernel;
1202 CArray fourierArray, fourierKernel;
1203 bool useFourierKernel;
1214 : useFourierKernel(false)
1228 template <
class C1,
class C2,
class C3>
1232 unsigned int planner_flags = FFTW_ESTIMATE)
1233 : useFourierKernel(false)
1235 init(in, kernel, out, planner_flags);
1249 template <
class C1,
class C2,
class C3>
1253 unsigned int planner_flags = FFTW_ESTIMATE)
1254 : useFourierKernel(true)
1256 init(in, kernel, out, planner_flags);
1271 template <
class C1,
class C2,
class C3>
1275 bool fourierDomainKernel,
1276 unsigned int planner_flags = FFTW_ESTIMATE)
1278 init(in, kernel, out, fourierDomainKernel, planner_flags);
1294 template <
class C1,
class C2,
class C3>
1296 bool useFourierKernel =
false,
1297 unsigned int planner_flags = FFTW_ESTIMATE)
1299 if(useFourierKernel)
1300 init(inOut, kernel, planner_flags);
1302 initFourierKernel(inOut, kernel, planner_flags);
1309 template <
class C1,
class C2,
class C3>
1313 unsigned int planner_flags = FFTW_ESTIMATE)
1315 vigra_precondition(in.
shape() == out.
shape(),
1316 "FFTWConvolvePlan::init(): input and output must have the same shape.");
1317 init(in.
shape(), kernel.
shape(), planner_flags);
1324 template <
class C1,
class C2,
class C3>
1328 unsigned int planner_flags = FFTW_ESTIMATE)
1330 vigra_precondition(in.
shape() == out.
shape(),
1331 "FFTWConvolvePlan::init(): input and output must have the same shape.");
1332 initFourierKernel(in.
shape(), kernel.shape(), planner_flags);
1339 template <
class C1,
class C2,
class C3>
1343 bool fourierDomainKernel,
1344 unsigned int planner_flags = FFTW_ESTIMATE)
1346 vigra_precondition(in.shape() == out.shape(),
1347 "FFTWConvolvePlan::init(): input and output must have the same shape.");
1348 useFourierKernel = fourierDomainKernel;
1349 initComplex(in.shape(), kernel.shape(), planner_flags);
1358 template <
class C1,
class KernelIterator,
class OutIterator>
1360 KernelIterator kernels, KernelIterator kernelsEnd,
1361 OutIterator outs,
unsigned int planner_flags = FFTW_ESTIMATE)
1363 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1364 typedef typename KernelArray::value_type KernelValue;
1365 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1366 typedef typename OutArray::value_type OutValue;
1368 bool realKernel = IsSameType<KernelValue, Real>::value;
1369 bool fourierKernel = IsSameType<KernelValue, Complex>::value;
1371 vigra_precondition(realKernel || fourierKernel,
1372 "FFTWConvolvePlan::initMany(): kernels have unsuitable value_type.");
1373 vigra_precondition((IsSameType<OutValue, Real>::value),
1374 "FFTWConvolvePlan::initMany(): outputs have unsuitable value_type.");
1378 initMany(in.
shape(), checkShapes(in.
shape(), kernels, kernelsEnd, outs),
1383 initFourierKernelMany(in.
shape(),
1384 checkShapesFourier(in.
shape(), kernels, kernelsEnd, outs),
1395 template <
class C1,
class KernelIterator,
class OutIterator>
1397 KernelIterator kernels, KernelIterator kernelsEnd,
1399 bool fourierDomainKernels,
1400 unsigned int planner_flags = FFTW_ESTIMATE)
1402 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1403 typedef typename KernelArray::value_type KernelValue;
1404 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1405 typedef typename OutArray::value_type OutValue;
1407 vigra_precondition((IsSameType<KernelValue, Complex>::value),
1408 "FFTWConvolvePlan::initMany(): kernels have unsuitable value_type.");
1409 vigra_precondition((IsSameType<OutValue, Complex>::value),
1410 "FFTWConvolvePlan::initMany(): outputs have unsuitable value_type.");
1412 useFourierKernel = fourierDomainKernels;
1414 Shape paddedShape = checkShapesComplex(in.shape(), kernels, kernelsEnd, outs);
1416 CArray newFourierArray(paddedShape), newFourierKernel(paddedShape);
1418 FFTWPlan<N, Real> fplan(newFourierArray, newFourierArray, FFTW_FORWARD, planner_flags);
1419 FFTWPlan<N, Real> bplan(newFourierArray, newFourierArray, FFTW_BACKWARD, planner_flags);
1421 forward_plan = fplan;
1422 backward_plan = bplan;
1423 fourierArray.
swap(newFourierArray);
1424 fourierKernel.
swap(newFourierKernel);
1427 void init(Shape inOut, Shape kernel,
1428 unsigned int planner_flags = FFTW_ESTIMATE);
1430 void initFourierKernel(Shape inOut, Shape kernel,
1431 unsigned int planner_flags = FFTW_ESTIMATE);
1433 void initComplex(Shape inOut, Shape kernel,
1434 unsigned int planner_flags = FFTW_ESTIMATE);
1436 void initMany(Shape inOut, Shape maxKernel,
1437 unsigned int planner_flags = FFTW_ESTIMATE)
1439 init(inOut, maxKernel, planner_flags);
1442 void initFourierKernelMany(Shape inOut, Shape kernels,
1443 unsigned int planner_flags = FFTW_ESTIMATE)
1445 initFourierKernel(inOut, kernels, planner_flags);
1455 template <
class C1,
class C2,
class C3>
1460 executeImpl(in, kernel, out);
1470 template <
class C1,
class C2,
class C3>
1482 template <
class C1,
class C2,
class C3>
1495 template <
class C1,
class KernelIterator,
class OutIterator>
1497 KernelIterator kernels, KernelIterator kernelsEnd,
1507 template <
class C1,
class KernelIterator,
class OutIterator>
1509 KernelIterator kernels, KernelIterator kernelsEnd,
1512 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1513 typedef typename KernelArray::value_type KernelValue;
1514 typedef typename IsSameType<KernelValue, Complex>::type UseFourierKernel;
1515 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1516 typedef typename OutArray::value_type OutValue;
1518 bool realKernel = IsSameType<KernelValue, Real>::value;
1519 bool fourierKernel = IsSameType<KernelValue, Complex>::value;
1521 vigra_precondition(realKernel || fourierKernel,
1522 "FFTWConvolvePlan::executeMany(): kernels have unsuitable value_type.");
1523 vigra_precondition((IsSameType<OutValue, Real>::value),
1524 "FFTWConvolvePlan::executeMany(): outputs have unsuitable value_type.");
1526 executeManyImpl(in, kernels, kernelsEnd, outs, UseFourierKernel());
1531 template <
class KernelIterator,
class OutIterator>
1532 Shape checkShapes(Shape in,
1533 KernelIterator kernels, KernelIterator kernelsEnd,
1536 template <
class KernelIterator,
class OutIterator>
1537 Shape checkShapesFourier(Shape in,
1538 KernelIterator kernels, KernelIterator kernelsEnd,
1541 template <
class KernelIterator,
class OutIterator>
1542 Shape checkShapesComplex(Shape in,
1543 KernelIterator kernels, KernelIterator kernelsEnd,
1546 template <
class C1,
class C2,
class C3>
1550 bool do_correlation=
false);
1552 template <
class C1,
class KernelIterator,
class OutIterator>
1555 KernelIterator kernels, KernelIterator kernelsEnd,
1556 OutIterator outs, VigraFalseType );
1558 template <
class C1,
class KernelIterator,
class OutIterator>
1561 KernelIterator kernels, KernelIterator kernelsEnd,
1562 OutIterator outs, VigraTrueType );
1566 template <
unsigned int N,
class Real>
1569 unsigned int planner_flags)
1571 Shape paddedShape = fftwBestPaddedShapeR2C(in + kernel - Shape(1)),
1574 CArray newFourierArray(complexShape), newFourierKernel(complexShape);
1576 Shape realStrides = 2*newFourierArray.stride();
1578 RArray newRealArray(paddedShape, realStrides, (Real*)newFourierArray.data());
1579 RArray newRealKernel(paddedShape, realStrides, (Real*)newFourierKernel.
data());
1584 forward_plan = fplan;
1585 backward_plan = bplan;
1586 realArray = newRealArray;
1587 realKernel = newRealKernel;
1588 fourierArray.swap(newFourierArray);
1589 fourierKernel.swap(newFourierKernel);
1590 useFourierKernel =
false;
1593 template <
unsigned int N,
class Real>
1596 unsigned int planner_flags)
1598 Shape complexShape = kernel,
1599 paddedShape = fftwCorrespondingShapeC2R(complexShape);
1601 for(
unsigned int k=0; k<N; ++k)
1602 vigra_precondition(in[k] <= paddedShape[k],
1603 "FFTWConvolvePlan::init(): kernel too small for given input.");
1605 CArray newFourierArray(complexShape), newFourierKernel(complexShape);
1607 Shape realStrides = 2*newFourierArray.stride();
1609 RArray newRealArray(paddedShape, realStrides, (Real*)newFourierArray.data());
1610 RArray newRealKernel(paddedShape, realStrides, (Real*)newFourierKernel.
data());
1615 forward_plan = fplan;
1616 backward_plan = bplan;
1617 realArray = newRealArray;
1618 realKernel = newRealKernel;
1619 fourierArray.swap(newFourierArray);
1620 fourierKernel.swap(newFourierKernel);
1621 useFourierKernel =
true;
1624 template <
unsigned int N,
class Real>
1627 unsigned int planner_flags)
1631 if(useFourierKernel)
1633 for(
unsigned int k=0; k<N; ++k)
1634 vigra_precondition(in[k] <= kernel[k],
1635 "FFTWConvolvePlan::init(): kernel too small for given input.");
1637 paddedShape = kernel;
1641 paddedShape = fftwBestPaddedShape(in + kernel - Shape(1));
1644 CArray newFourierArray(paddedShape), newFourierKernel(paddedShape);
1646 FFTWPlan<N, Real> fplan(newFourierArray, newFourierArray, FFTW_FORWARD, planner_flags);
1647 FFTWPlan<N, Real> bplan(newFourierArray, newFourierArray, FFTW_BACKWARD, planner_flags);
1649 forward_plan = fplan;
1650 backward_plan = bplan;
1651 fourierArray.swap(newFourierArray);
1652 fourierKernel.swap(newFourierKernel);
1655 #ifndef DOXYGEN // doxygen documents these functions as free functions 1657 template <
unsigned int N,
class Real>
1658 template <
class C1,
class C2,
class C3>
1663 bool do_correlation)
1665 vigra_precondition(!useFourierKernel,
1666 "FFTWConvolvePlan::execute(): plan was generated for Fourier kernel, got spatial kernel.");
1668 vigra_precondition(in.
shape() == out.
shape(),
1669 "FFTWConvolvePlan::execute(): input and output must have the same shape.");
1671 Shape paddedShape = fftwBestPaddedShapeR2C(in.
shape() + kernel.
shape() - Shape(1)),
1672 diff = paddedShape - in.
shape(),
1674 right = in.
shape() + left;
1676 vigra_precondition(paddedShape == realArray.shape(),
1677 "FFTWConvolvePlan::execute(): shape mismatch between input and plan.");
1679 detail::fftEmbedArray(in, realArray);
1680 forward_plan.execute(realArray, fourierArray);
1682 detail::fftEmbedKernel(kernel, realKernel);
1683 forward_plan.execute(realKernel, fourierKernel);
1687 using namespace multi_math;
1688 fourierArray *=
conj(fourierKernel);
1692 fourierArray *= fourierKernel;
1695 backward_plan.execute(fourierArray, realArray);
1697 out = realArray.
subarray(left, right);
1700 template <
unsigned int N,
class Real>
1701 template <
class C1,
class C2,
class C3>
1707 vigra_precondition(useFourierKernel,
1708 "FFTWConvolvePlan::execute(): plan was generated for spatial kernel, got Fourier kernel.");
1710 vigra_precondition(in.
shape() == out.
shape(),
1711 "FFTWConvolvePlan::execute(): input and output must have the same shape.");
1713 vigra_precondition(kernel.
shape() == fourierArray.shape(),
1714 "FFTWConvolvePlan::execute(): shape mismatch between kernel and plan.");
1716 Shape paddedShape = fftwCorrespondingShapeC2R(kernel.
shape(),
odd(in.
shape(0))),
1717 diff = paddedShape - in.
shape(),
1719 right = in.
shape() + left;
1721 vigra_precondition(paddedShape == realArray.shape(),
1722 "FFTWConvolvePlan::execute(): shape mismatch between input and plan.");
1724 detail::fftEmbedArray(in, realArray);
1725 forward_plan.execute(realArray, fourierArray);
1727 fourierKernel = kernel;
1728 moveDCToHalfspaceUpperLeft(fourierKernel);
1730 fourierArray *= fourierKernel;
1732 backward_plan.execute(fourierArray, realArray);
1734 out = realArray.
subarray(left, right);
1737 template <
unsigned int N,
class Real>
1738 template <
class C1,
class C2,
class C3>
1744 vigra_precondition(in.
shape() == out.
shape(),
1745 "FFTWConvolvePlan::execute(): input and output must have the same shape.");
1747 Shape paddedShape = fourierArray.shape(),
1748 diff = paddedShape - in.
shape(),
1750 right = in.
shape() + left;
1752 if(useFourierKernel)
1754 vigra_precondition(kernel.
shape() == fourierArray.shape(),
1755 "FFTWConvolvePlan::execute(): shape mismatch between kernel and plan.");
1757 fourierKernel = kernel;
1762 detail::fftEmbedKernel(kernel, fourierKernel);
1763 forward_plan.execute(fourierKernel, fourierKernel);
1766 detail::fftEmbedArray(in, fourierArray);
1767 forward_plan.execute(fourierArray, fourierArray);
1769 fourierArray *= fourierKernel;
1771 backward_plan.execute(fourierArray, fourierArray);
1773 out = fourierArray.
subarray(left, right);
1776 template <
unsigned int N,
class Real>
1777 template <
class C1,
class KernelIterator,
class OutIterator>
1780 KernelIterator kernels, KernelIterator kernelsEnd,
1781 OutIterator outs, VigraFalseType )
1783 vigra_precondition(!useFourierKernel,
1784 "FFTWConvolvePlan::execute(): plan was generated for Fourier kernel, got spatial kernel.");
1786 Shape kernelMax = checkShapes(in.
shape(), kernels, kernelsEnd, outs),
1787 paddedShape = fftwBestPaddedShapeR2C(in.
shape() + kernelMax - Shape(1)),
1788 diff = paddedShape - in.
shape(),
1790 right = in.
shape() + left;
1792 vigra_precondition(paddedShape == realArray.shape(),
1793 "FFTWConvolvePlan::executeMany(): shape mismatch between input and plan.");
1795 detail::fftEmbedArray(in, realArray);
1796 forward_plan.execute(realArray, fourierArray);
1798 for(; kernels != kernelsEnd; ++kernels, ++outs)
1800 detail::fftEmbedKernel(*kernels, realKernel);
1801 forward_plan.execute(realKernel, fourierKernel);
1803 fourierKernel *= fourierArray;
1805 backward_plan.execute(fourierKernel, realKernel);
1807 *outs = realKernel.subarray(left, right);
1811 template <
unsigned int N,
class Real>
1812 template <
class C1,
class KernelIterator,
class OutIterator>
1815 KernelIterator kernels, KernelIterator kernelsEnd,
1816 OutIterator outs, VigraTrueType )
1818 vigra_precondition(useFourierKernel,
1819 "FFTWConvolvePlan::execute(): plan was generated for spatial kernel, got Fourier kernel.");
1821 Shape complexShape = checkShapesFourier(in.
shape(), kernels, kernelsEnd, outs),
1822 paddedShape = fftwCorrespondingShapeC2R(complexShape,
odd(in.
shape(0))),
1823 diff = paddedShape - in.
shape(),
1825 right = in.
shape() + left;
1827 vigra_precondition(complexShape == fourierArray.shape(),
1828 "FFTWConvolvePlan::executeFourierKernelMany(): shape mismatch between kernels and plan.");
1830 vigra_precondition(paddedShape == realArray.shape(),
1831 "FFTWConvolvePlan::executeFourierKernelMany(): shape mismatch between input and plan.");
1833 detail::fftEmbedArray(in, realArray);
1834 forward_plan.execute(realArray, fourierArray);
1836 for(; kernels != kernelsEnd; ++kernels, ++outs)
1838 fourierKernel = *kernels;
1839 moveDCToHalfspaceUpperLeft(fourierKernel);
1840 fourierKernel *= fourierArray;
1842 backward_plan.execute(fourierKernel, realKernel);
1844 *outs = realKernel.subarray(left, right);
1848 template <
unsigned int N,
class Real>
1849 template <
class C1,
class KernelIterator,
class OutIterator>
1852 KernelIterator kernels, KernelIterator kernelsEnd,
1855 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1856 typedef typename KernelArray::value_type KernelValue;
1857 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1858 typedef typename OutArray::value_type OutValue;
1860 vigra_precondition((IsSameType<KernelValue, Complex>::value),
1861 "FFTWConvolvePlan::executeMany(): kernels have unsuitable value_type.");
1862 vigra_precondition((IsSameType<OutValue, Complex>::value),
1863 "FFTWConvolvePlan::executeMany(): outputs have unsuitable value_type.");
1865 Shape paddedShape = checkShapesComplex(in.
shape(), kernels, kernelsEnd, outs),
1866 diff = paddedShape - in.
shape(),
1868 right = in.
shape() + left;
1870 detail::fftEmbedArray(in, fourierArray);
1871 forward_plan.execute(fourierArray, fourierArray);
1873 for(; kernels != kernelsEnd; ++kernels, ++outs)
1875 if(useFourierKernel)
1877 fourierKernel = *kernels;
1882 detail::fftEmbedKernel(*kernels, fourierKernel);
1883 forward_plan.execute(fourierKernel, fourierKernel);
1886 fourierKernel *= fourierArray;
1888 backward_plan.execute(fourierKernel, fourierKernel);
1890 *outs = fourierKernel.subarray(left, right);
1896 template <
unsigned int N,
class Real>
1897 template <
class KernelIterator,
class OutIterator>
1900 KernelIterator kernels, KernelIterator kernelsEnd,
1903 vigra_precondition(kernels != kernelsEnd,
1904 "FFTWConvolvePlan::checkShapes(): empty kernel sequence.");
1907 for(; kernels != kernelsEnd; ++kernels, ++outs)
1909 vigra_precondition(in == outs->shape(),
1910 "FFTWConvolvePlan::checkShapes(): shape mismatch between input and (one) output.");
1911 kernelMax = max(kernelMax, kernels->shape());
1913 vigra_precondition(
prod(kernelMax) > 0,
1914 "FFTWConvolvePlan::checkShapes(): all kernels have size 0.");
1918 template <
unsigned int N,
class Real>
1919 template <
class KernelIterator,
class OutIterator>
1922 KernelIterator kernels, KernelIterator kernelsEnd,
1925 vigra_precondition(kernels != kernelsEnd,
1926 "FFTWConvolvePlan::checkShapesFourier(): empty kernel sequence.");
1928 Shape complexShape = kernels->shape(),
1929 paddedShape = fftwCorrespondingShapeC2R(complexShape);
1931 for(
unsigned int k=0; k<N; ++k)
1932 vigra_precondition(in[k] <= paddedShape[k],
1933 "FFTWConvolvePlan::checkShapesFourier(): kernels too small for given input.");
1935 for(; kernels != kernelsEnd; ++kernels, ++outs)
1937 vigra_precondition(in == outs->shape(),
1938 "FFTWConvolvePlan::checkShapesFourier(): shape mismatch between input and (one) output.");
1939 vigra_precondition(complexShape == kernels->shape(),
1940 "FFTWConvolvePlan::checkShapesFourier(): all kernels must have the same size.");
1942 return complexShape;
1945 template <
unsigned int N,
class Real>
1946 template <
class KernelIterator,
class OutIterator>
1949 KernelIterator kernels, KernelIterator kernelsEnd,
1952 vigra_precondition(kernels != kernelsEnd,
1953 "FFTWConvolvePlan::checkShapesComplex(): empty kernel sequence.");
1955 Shape kernelShape = kernels->shape();
1956 for(; kernels != kernelsEnd; ++kernels, ++outs)
1958 vigra_precondition(in == outs->shape(),
1959 "FFTWConvolvePlan::checkShapesComplex(): shape mismatch between input and (one) output.");
1960 if(useFourierKernel)
1962 vigra_precondition(kernelShape == kernels->shape(),
1963 "FFTWConvolvePlan::checkShapesComplex(): Fourier domain kernels must have identical size.");
1967 kernelShape = max(kernelShape, kernels->shape());
1970 vigra_precondition(
prod(kernelShape) > 0,
1971 "FFTWConvolvePlan::checkShapesComplex(): all kernels have size 0.");
1973 if(useFourierKernel)
1975 for(
unsigned int k=0; k<N; ++k)
1976 vigra_precondition(in[k] <= kernelShape[k],
1977 "FFTWConvolvePlan::checkShapesComplex(): kernels too small for given input.");
1982 return fftwBestPaddedShape(in + kernelShape - Shape(1));
2018 template <
unsigned int N,
class Real>
2046 template <
class C1,
class C2,
class C3>
2050 unsigned int planner_flags = FFTW_ESTIMATE)
2051 : BaseType(in, kernel, out, planner_flags)
2066 template <
class C1,
class C2,
class C3>
2068 bool useFourierKernel =
false,
2069 unsigned int planner_flags = FFTW_ESTIMATE)
2070 : BaseType(inOut, kernel, false, planner_flags)
2077 template <
class C1,
class C2,
class C3>
2081 unsigned int planner_flags = FFTW_ESTIMATE)
2083 vigra_precondition(in.
shape() == out.
shape(),
2084 "FFTWCorrelatePlan::init(): input and output must have the same shape.");
2085 BaseType::init(in.
shape(), kernel.
shape(), planner_flags);
2095 template <
class C1,
class C2,
class C3>
2100 BaseType::executeImpl(in, kernel, out,
true);
2110 template <
unsigned int N,
class Real,
class C1,
class C2>
2118 template <
unsigned int N,
class Real,
class C1,
class C2>
2126 template <
unsigned int N,
class Real,
class C1,
class C2>
2142 vigra_precondition(
false,
2143 "fourierTransform(): shape mismatch between input and output.");
2146 template <
unsigned int N,
class Real,
class C1,
class C2>
2152 "fourierTransformInverse(): shape mismatch between input and output.");
2347 template <
unsigned int N,
class Real,
class C1,
class C2,
class C3>
2356 template <
unsigned int N,
class Real,
class C1,
class C2,
class C3>
2371 template <
unsigned int N,
class Real,
class C1,
class C2,
class C3>
2376 bool fourierDomainKernel)
2387 template <
unsigned int N,
class Real,
class C1,
2388 class KernelIterator,
class OutIterator>
2391 KernelIterator kernels, KernelIterator kernelsEnd,
2395 plan.
initMany(in, kernels, kernelsEnd, outs);
2405 template <
unsigned int N,
class Real,
class C1,
2406 class KernelIterator,
class OutIterator>
2409 KernelIterator kernels, KernelIterator kernelsEnd,
2411 bool fourierDomainKernel)
2414 plan.
initMany(in, kernels, kernelsEnd, outs, fourierDomainKernel);
2468 template <
unsigned int N,
class Real,
class C1,
class C2,
class C3>
2481 #endif // VIGRA_MULTI_FFT_HXX void execute(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out) const
Execute a complex-to-complex transform.
Definition: multi_fft.hxx:1021
FFTWCorrelatePlan()
Create an empty plan.
Definition: multi_fft.hxx:2031
void convolveFFTMany(...)
Convolve a real-valued array with a sequence of kernels by means of the Fourier transform.
void execute(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, Real, C2 > kernel, MultiArrayView< N, Real, C3 > out)
Execute a plan to convolve a real array with a real kernel.
Definition: multi_fft.hxx:1456
void convolveFFT(...)
Convolve an array with a kernel by means of the Fourier transform.
~FFTWPlan()
Destructor.
Definition: multi_fft.hxx:959
FFTWComplex< R > conj(const FFTWComplex< R > &a)
complex conjugate
Definition: fftw3.hxx:1030
FFTWReal2Complex< Real >::type complex_type
Definition: fftw3.hxx:136
void init(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, Real, C2 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Init a complex-to-real transform.
Definition: multi_fft.hxx:1002
const difference_type & shape() const
Definition: multi_array.hxx:1596
void init(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > kernel, MultiArrayView< N, Real, C3 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Init a plan to convolve a real array with a complex kernel.
Definition: multi_fft.hxx:1325
FFTWConvolvePlan(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > kernel, MultiArrayView< N, FFTWComplex< Real >, C3 > out, bool fourierDomainKernel, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan to convolve a complex array with a complex kernel.
Definition: multi_fft.hxx:1272
pointer data() const
Definition: multi_array.hxx:1846
void init(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, Real, C2 > kernel, MultiArrayView< N, Real, C3 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Init a plan to convolve a real array with a real kernel.
Definition: multi_fft.hxx:1310
Main MultiArray class containing the memory management.
Definition: multi_array.hxx:2422
void execute(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out) const
Execute a real-to-complex transform.
Definition: multi_fft.hxx:1035
void initMany(MultiArrayView< N, Real, C1 > in, KernelIterator kernels, KernelIterator kernelsEnd, OutIterator outs, unsigned int planner_flags=FFTW_ESTIMATE)
Init a plan to convolve a real array with a sequence of kernels.
Definition: multi_fft.hxx:1359
std::ptrdiff_t MultiArrayIndex
Definition: multi_fwd.hxx:60
void init(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out, int SIGN, unsigned int planner_flags=FFTW_ESTIMATE)
Init a complex-to-complex transform.
Definition: multi_fft.hxx:970
void convolveFFTComplexMany(...)
Convolve a complex-valued array with a sequence of kernels by means of the Fourier transform...
Definition: accessor.hxx:43
bool odd(int t)
Check if an integer is odd.
doxygen_overloaded_function(template<... > void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays.
FFTWComplex< R >::NormType norm(const FFTWComplex< R > &a)
norm (= magnitude)
Definition: fftw3.hxx:1037
vigra::detail::MultiIteratorChooser< StrideTag >::template Traverser< actual_dimension, T, T &, T * >::type traverser
Definition: multi_array.hxx:712
NumericTraits< V >::Promote prod(TinyVectorBase< V, SIZE, D1, D2 > const &l)
product of the vector's elements
Definition: tinyvector.hxx:2097
void executeMany(MultiArrayView< N, FFTWComplex< Real >, C1 > in, KernelIterator kernels, KernelIterator kernelsEnd, OutIterator outs)
Execute a plan to convolve a complex array with a sequence of kernels.
FFTWPlan(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out, int SIGN, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan for a complex-to-complex transform.
Definition: multi_fft.hxx:879
void init(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > kernel, MultiArrayView< N, FFTWComplex< Real >, C3 > out, bool fourierDomainKernel, unsigned int planner_flags=FFTW_ESTIMATE)
Init a plan to convolve a complex array with a complex kernel.
Definition: multi_fft.hxx:1340
FFTWPlan(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan for a real-to-complex transform.
Definition: multi_fft.hxx:899
Definition: multi_fft.hxx:851
void initMany(MultiArrayView< N, FFTWComplex< Real >, C1 > in, KernelIterator kernels, KernelIterator kernelsEnd, OutIterator outs, bool fourierDomainKernels, unsigned int planner_flags=FFTW_ESTIMATE)
Init a plan to convolve a complex array with a sequence of kernels.
Definition: multi_fft.hxx:1396
void init(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, Real, C2 > kernel, MultiArrayView< N, Real, C3 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Init a plan to convolve a real array with a real kernel.
Definition: multi_fft.hxx:2078
FFTWCorrelatePlan(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, Real, C2 > kernel, MultiArrayView< N, Real, C3 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan to correlate a real array with a real kernel.
Definition: multi_fft.hxx:2047
MultiArrayView< N, T, StridedArrayTag > permuteStridesDescending() const
Definition: multi_array.hxx:2121
Definition: metaprogramming.hxx:105
TinyVector< MultiArrayIndex, N > type
Definition: multi_shape.hxx:250
bool even(int t)
Check if an integer is even.
Wrapper for fixed size vectors.
Definition: tinyvector.hxx:621
void executeMany(MultiArrayView< N, Real, C1 > in, KernelIterator kernels, KernelIterator kernelsEnd, OutIterator outs)
Execute a plan to convolve a real array with a sequence of kernels.
Definition: multi_fft.hxx:1508
void execute(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, Real, C2 > kernel, MultiArrayView< N, Real, C3 > out)
Execute a plan to correlate a real array with a real kernel.
Definition: multi_fft.hxx:2096
Definition: multi_fft.hxx:2019
FixedPoint16< IntBits3, OverflowHandling > & div(FixedPoint16< IntBits1, OverflowHandling > l, FixedPoint16< IntBits2, OverflowHandling > r, FixedPoint16< IntBits3, OverflowHandling > &result)
division with enforced result type.
Definition: fixedpoint.hxx:1616
FFTWConvolvePlan(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > kernel, MultiArrayView< N, Real, C3 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan to convolve a real array with a complex kernel.
Definition: multi_fft.hxx:1250
void correlateFFT(...)
Correlate an array with a kernel by means of the Fourier transform.
FFTWPlan(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, Real, C2 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan for a complex-to-real transform.
Definition: multi_fft.hxx:919
Definition: multi_fft.hxx:1194
FFTWPlan(FFTWPlan const &other)
Copy constructor.
Definition: multi_fft.hxx:929
Base class for, and view to, vigra::MultiArray.
Definition: multi_array.hxx:652
void convolveFFTComplex(...)
Convolve a complex-valued array by means of the Fourier transform.
FFTWPlan()
Create an empty plan.
Definition: multi_fft.hxx:866
FFTWPlan & operator=(FFTWPlan const &other)
Copy assigment.
Definition: multi_fft.hxx:942
MultiArrayView subarray(difference_type p, difference_type q) const
Definition: multi_array.hxx:1476
FFTWConvolvePlan(Shape inOut, Shape kernel, bool useFourierKernel=false, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan from just the shape information.
Definition: multi_fft.hxx:1295
T sign(T t)
The sign function.
Definition: mathutil.hxx:574
void execute(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, Real, C2 > out) const
Execute a complex-to-real transform.
Definition: multi_fft.hxx:1049
void swap(MultiArray &other)
Definition: multi_array.hxx:3018
void init(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Init a real-to-complex transform.
Definition: multi_fft.hxx:986
FFTWConvolvePlan()
Create an empty plan.
Definition: multi_fft.hxx:1213
Wrapper class for the FFTW complex types 'fftw_complex'.
Definition: fftw3.hxx:131
FFTWCorrelatePlan(Shape inOut, Shape kernel, bool useFourierKernel=false, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan from just the shape information.
Definition: multi_fft.hxx:2067
difference_type strideOrdering() const
Definition: multi_array.hxx:1565
FFTWConvolvePlan(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, Real, C2 > kernel, MultiArrayView< N, Real, C3 > out, unsigned int planner_flags=FFTW_ESTIMATE)
Create a plan to convolve a real array with a real kernel.
Definition: multi_fft.hxx:1229