36 #ifndef VIGRA_MULTI_FFT_HXX
37 #define VIGRA_MULTI_FFT_HXX
40 #include "multi_array.hxx"
41 #include "navigator.hxx"
42 #include "copyimage.hxx"
58 template <
unsigned int N,
class T,
class C>
59 void moveDCToCenterImpl(MultiArrayView<N, T, C> a,
unsigned int startDimension)
61 typedef typename MultiArrayView<N, T, C>::traverser Traverser;
62 typedef MultiArrayNavigator<Traverser, N> Navigator;
63 typedef typename Navigator::iterator Iterator;
65 for(
unsigned int d = startDimension; d < N; ++d)
67 Navigator nav(a.traverser_begin(), a.shape(), d);
69 for( ; nav.hasMore(); nav++ )
71 Iterator i = nav.begin();
72 int s = nav.end() - i;
77 for(
int k=0; k<s2; ++k)
79 std::swap(i[k], i[k+s2]);
85 for(
int k=0; k<s2; ++k)
96 template <
unsigned int N,
class T,
class C>
97 void moveDCToUpperLeftImpl(MultiArrayView<N, T, C> a,
unsigned int startDimension)
99 typedef typename MultiArrayView<N, T, C>::traverser Traverser;
100 typedef MultiArrayNavigator<Traverser, N> Navigator;
101 typedef typename Navigator::iterator Iterator;
103 for(
unsigned int d = startDimension; d < N; ++d)
105 Navigator nav(a.traverser_begin(), a.shape(), d);
107 for( ; nav.hasMore(); nav++ )
109 Iterator i = nav.begin();
110 int s = nav.end() - i;
115 for(
int k=0; k<s2; ++k)
117 std::swap(i[k], i[k+s2]);
123 for(
int k=s2; k>0; --k)
142 template <
unsigned int N,
class T,
class C>
145 detail::moveDCToCenterImpl(a, 0);
148 template <
unsigned int N,
class T,
class C>
151 detail::moveDCToUpperLeftImpl(a, 0);
154 template <
unsigned int N,
class T,
class C>
155 inline void moveDCToHalfspaceCenter(MultiArrayView<N, T, C> a)
157 detail::moveDCToCenterImpl(a, 1);
160 template <
unsigned int N,
class T,
class C>
161 inline void moveDCToHalfspaceUpperLeft(MultiArrayView<N, T, C> a)
163 detail::moveDCToUpperLeftImpl(a, 1);
170 fftwPlanCreate(
unsigned int N,
int* shape,
171 FFTWComplex<double> * in,
int* instrides,
int instep,
172 FFTWComplex<double> * out,
int* outstrides,
int outstep,
173 int sign,
unsigned int planner_flags)
175 return fftw_plan_many_dft(N, shape, 1,
176 (fftw_complex *)in, instrides, instep, 0,
177 (fftw_complex *)out, outstrides, outstep, 0,
178 sign, planner_flags);
182 fftwPlanCreate(
unsigned int N,
int* shape,
183 double * in,
int* instrides,
int instep,
184 FFTWComplex<double> * out,
int* outstrides,
int outstep,
185 int ,
unsigned int planner_flags)
187 return fftw_plan_many_dft_r2c(N, shape, 1,
188 in, instrides, instep, 0,
189 (fftw_complex *)out, outstrides, outstep, 0,
194 fftwPlanCreate(
unsigned int N,
int* shape,
195 FFTWComplex<double> * in,
int* instrides,
int instep,
196 double * out,
int* outstrides,
int outstep,
197 int ,
unsigned int planner_flags)
199 return fftw_plan_many_dft_c2r(N, shape, 1,
200 (fftw_complex *)in, instrides, instep, 0,
201 out, outstrides, outstep, 0,
206 fftwPlanCreate(
unsigned int N,
int* shape,
207 FFTWComplex<float> * in,
int* instrides,
int instep,
208 FFTWComplex<float> * out,
int* outstrides,
int outstep,
209 int sign,
unsigned int planner_flags)
211 return fftwf_plan_many_dft(N, shape, 1,
212 (fftwf_complex *)in, instrides, instep, 0,
213 (fftwf_complex *)out, outstrides, outstep, 0,
214 sign, planner_flags);
218 fftwPlanCreate(
unsigned int N,
int* shape,
219 float * in,
int* instrides,
int instep,
220 FFTWComplex<float> * out,
int* outstrides,
int outstep,
221 int ,
unsigned int planner_flags)
223 return fftwf_plan_many_dft_r2c(N, shape, 1,
224 in, instrides, instep, 0,
225 (fftwf_complex *)out, outstrides, outstep, 0,
230 fftwPlanCreate(
unsigned int N,
int* shape,
231 FFTWComplex<float> * in,
int* instrides,
int instep,
232 float * out,
int* outstrides,
int outstep,
233 int ,
unsigned int planner_flags)
235 return fftwf_plan_many_dft_c2r(N, shape, 1,
236 (fftwf_complex *)in, instrides, instep, 0,
237 out, outstrides, outstep, 0,
242 fftwPlanCreate(
unsigned int N,
int* shape,
243 FFTWComplex<long double> * in,
int* instrides,
int instep,
244 FFTWComplex<long double> * out,
int* outstrides,
int outstep,
245 int sign,
unsigned int planner_flags)
247 return fftwl_plan_many_dft(N, shape, 1,
248 (fftwl_complex *)in, instrides, instep, 0,
249 (fftwl_complex *)out, outstrides, outstep, 0,
250 sign, planner_flags);
254 fftwPlanCreate(
unsigned int N,
int* shape,
255 long double * in,
int* instrides,
int instep,
256 FFTWComplex<long double> * out,
int* outstrides,
int outstep,
257 int ,
unsigned int planner_flags)
259 return fftwl_plan_many_dft_r2c(N, shape, 1,
260 in, instrides, instep, 0,
261 (fftwl_complex *)out, outstrides, outstep, 0,
266 fftwPlanCreate(
unsigned int N,
int* shape,
267 FFTWComplex<long double> * in,
int* instrides,
int instep,
268 long double * out,
int* outstrides,
int outstep,
269 int ,
unsigned int planner_flags)
271 return fftwl_plan_many_dft_c2r(N, shape, 1,
272 (fftwl_complex *)in, instrides, instep, 0,
273 out, outstrides, outstep, 0,
277 inline void fftwPlanDestroy(fftw_plan plan)
280 fftw_destroy_plan(plan);
283 inline void fftwPlanDestroy(fftwf_plan plan)
286 fftwf_destroy_plan(plan);
289 inline void fftwPlanDestroy(fftwl_plan plan)
292 fftwl_destroy_plan(plan);
296 fftwPlanExecute(fftw_plan plan)
302 fftwPlanExecute(fftwf_plan plan)
308 fftwPlanExecute(fftwl_plan plan)
314 fftwPlanExecute(fftw_plan plan, FFTWComplex<double> * in, FFTWComplex<double> * out)
316 fftw_execute_dft(plan, (fftw_complex *)in, (fftw_complex *)out);
320 fftwPlanExecute(fftw_plan plan,
double * in, FFTWComplex<double> * out)
322 fftw_execute_dft_r2c(plan, in, (fftw_complex *)out);
326 fftwPlanExecute(fftw_plan plan, FFTWComplex<double> * in,
double * out)
328 fftw_execute_dft_c2r(plan, (fftw_complex *)in, out);
332 fftwPlanExecute(fftwf_plan plan, FFTWComplex<float> * in, FFTWComplex<float> * out)
334 fftwf_execute_dft(plan, (fftwf_complex *)in, (fftwf_complex *)out);
338 fftwPlanExecute(fftwf_plan plan,
float * in, FFTWComplex<float> * out)
340 fftwf_execute_dft_r2c(plan, in, (fftwf_complex *)out);
344 fftwPlanExecute(fftwf_plan plan, FFTWComplex<float> * in,
float * out)
346 fftwf_execute_dft_c2r(plan, (fftwf_complex *)in, out);
350 fftwPlanExecute(fftwl_plan plan, FFTWComplex<long double> * in, FFTWComplex<long double> * out)
352 fftwl_execute_dft(plan, (fftwl_complex *)in, (fftwl_complex *)out);
356 fftwPlanExecute(fftwl_plan plan,
long double * in, FFTWComplex<long double> * out)
358 fftwl_execute_dft_r2c(plan, in, (fftwl_complex *)out);
362 fftwPlanExecute(fftwl_plan plan, FFTWComplex<long double> * in,
long double * out)
364 fftwl_execute_dft_c2r(plan, (fftwl_complex *)in, out);
368 int fftwPaddingSize(
int s)
374 static const int size = 1330;
375 static int goodSizes[size] = { 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16,
376 18, 20, 21, 22, 24, 25, 26, 27, 28, 30, 32, 33, 35, 36, 39, 40, 42, 44, 45, 48,
377 49, 50, 52, 54, 55, 56, 60, 63, 64, 65, 66, 70, 72, 75, 77, 78, 80, 81,
378 84, 88, 90, 91, 96, 98, 99, 100, 104, 105, 108, 110, 112, 117, 120, 125,
379 126, 128, 130, 132, 135, 140, 144, 147, 150, 154, 156, 160, 162, 165,
380 168, 175, 176, 180, 182, 189, 192, 195, 196, 198, 200, 208, 210, 216,
381 220, 224, 225, 231, 234, 240, 243, 245, 250, 252, 256, 260, 264, 270,
382 273, 275, 280, 288, 294, 297, 300, 308, 312, 315, 320, 324, 325, 330,
383 336, 343, 350, 351, 352, 360, 364, 375, 378, 384, 385, 390, 392, 396,
384 400, 405, 416, 420, 432, 440, 441, 448, 450, 455, 462, 468, 480, 486,
385 490, 495, 500, 504, 512, 520, 525, 528, 539, 540, 546, 550, 560, 567,
386 576, 585, 588, 594, 600, 616, 624, 625, 630, 637, 640, 648, 650, 660,
387 672, 675, 686, 693, 700, 702, 704, 720, 728, 729, 735, 750, 756, 768,
388 770, 780, 784, 792, 800, 810, 819, 825, 832, 840, 864, 875, 880, 882,
389 891, 896, 900, 910, 924, 936, 945, 960, 972, 975, 980, 990, 1000, 1008,
390 1024, 1029, 1040, 1050, 1053, 1056, 1078, 1080, 1092, 1100, 1120, 1125,
391 1134, 1152, 1155, 1170, 1176, 1188, 1200, 1215, 1225, 1232, 1248, 1250,
392 1260, 1274, 1280, 1296, 1300, 1320, 1323, 1344, 1350, 1365, 1372, 1375,
393 1386, 1400, 1404, 1408, 1440, 1456, 1458, 1470, 1485, 1500, 1512, 1536,
394 1540, 1560, 1568, 1575, 1584, 1600, 1617, 1620, 1625, 1638, 1650, 1664,
395 1680, 1701, 1715, 1728, 1750, 1755, 1760, 1764, 1782, 1792, 1800, 1820,
396 1848, 1872, 1875, 1890, 1911, 1920, 1925, 1944, 1950, 1960, 1980, 2000,
397 2016, 2025, 2048, 2058, 2079, 2080, 2100, 2106, 2112, 2156, 2160, 2184,
398 2187, 2200, 2205, 2240, 2250, 2268, 2275, 2304, 2310, 2340, 2352, 2376,
399 2400, 2401, 2430, 2450, 2457, 2464, 2475, 2496, 2500, 2520, 2548, 2560,
400 2592, 2600, 2625, 2640, 2646, 2673, 2688, 2695, 2700, 2730, 2744, 2750,
401 2772, 2800, 2808, 2816, 2835, 2880, 2912, 2916, 2925, 2940, 2970, 3000,
402 3024, 3072, 3080, 3087, 3120, 3125, 3136, 3150, 3159, 3168, 3185, 3200,
403 3234, 3240, 3250, 3276, 3300, 3328, 3360, 3375, 3402, 3430, 3456, 3465,
404 3500, 3510, 3520, 3528, 3564, 3584, 3600, 3640, 3645, 3675, 3696, 3744,
405 3750, 3773, 3780, 3822, 3840, 3850, 3888, 3900, 3920, 3960, 3969, 4000,
406 4032, 4050, 4095, 4096, 4116, 4125, 4158, 4160, 4200, 4212, 4224, 4312,
407 4320, 4368, 4374, 4375, 4400, 4410, 4455, 4459, 4480, 4500, 4536, 4550,
408 4608, 4620, 4680, 4704, 4725, 4752, 4800, 4802, 4851, 4860, 4875, 4900,
409 4914, 4928, 4950, 4992, 5000, 5040, 5096, 5103, 5120, 5145, 5184, 5200,
410 5250, 5265, 5280, 5292, 5346, 5376, 5390, 5400, 5460, 5488, 5500, 5544,
411 5600, 5616, 5625, 5632, 5670, 5733, 5760, 5775, 5824, 5832, 5850, 5880,
412 5940, 6000, 6048, 6075, 6125, 6144, 6160, 6174, 6237, 6240, 6250, 6272,
413 6300, 6318, 6336, 6370, 6400, 6468, 6480, 6500, 6552, 6561, 6600, 6615,
414 6656, 6720, 6750, 6804, 6825, 6860, 6875, 6912, 6930, 7000, 7020, 7040,
415 7056, 7128, 7168, 7200, 7203, 7280, 7290, 7350, 7371, 7392, 7425, 7488,
416 7500, 7546, 7560, 7644, 7680, 7700, 7776, 7800, 7840, 7875, 7920, 7938,
417 8000, 8019, 8064, 8085, 8100, 8125, 8190, 8192, 8232, 8250, 8316, 8320,
418 8400, 8424, 8448, 8505, 8575, 8624, 8640, 8736, 8748, 8750, 8775, 8800,
419 8820, 8910, 8918, 8960, 9000, 9072, 9100, 9216, 9240, 9261, 9360, 9375,
420 9408, 9450, 9477, 9504, 9555, 9600, 9604, 9625, 9702, 9720, 9750, 9800,
421 9828, 9856, 9900, 9984, 10000, 10080, 10125, 10192, 10206, 10240, 10290,
422 10368, 10395, 10400, 10500, 10530, 10560, 10584, 10692, 10752, 10780,
423 10800, 10920, 10935, 10976, 11000, 11025, 11088, 11200, 11232, 11250,
424 11264, 11319, 11340, 11375, 11466, 11520, 11550, 11648, 11664, 11700,
425 11760, 11880, 11907, 12000, 12005, 12096, 12150, 12250, 12285, 12288,
426 12320, 12348, 12375, 12474, 12480, 12500, 12544, 12600, 12636, 12672,
427 12740, 12800, 12936, 12960, 13000, 13104, 13122, 13125, 13200, 13230,
428 13312, 13365, 13377, 13440, 13475, 13500, 13608, 13650, 13720, 13750,
429 13824, 13860, 14000, 14040, 14080, 14112, 14175, 14256, 14336, 14400,
430 14406, 14553, 14560, 14580, 14625, 14700, 14742, 14784, 14850, 14976,
431 15000, 15092, 15120, 15288, 15309, 15360, 15400, 15435, 15552, 15600,
432 15625, 15680, 15750, 15795, 15840, 15876, 15925, 16000, 16038, 16128,
433 16170, 16200, 16250, 16380, 16384, 16464, 16500, 16632, 16640, 16800,
434 16807, 16848, 16875, 16896, 17010, 17150, 17199, 17248, 17280, 17325,
435 17472, 17496, 17500, 17550, 17600, 17640, 17820, 17836, 17920, 18000,
436 18144, 18200, 18225, 18375, 18432, 18480, 18522, 18711, 18720, 18750,
437 18816, 18865, 18900, 18954, 19008, 19110, 19200, 19208, 19250, 19404,
438 19440, 19500, 19600, 19656, 19683, 19712, 19800, 19845, 19968, 20000,
439 20160, 20250, 20384, 20412, 20475, 20480, 20580, 20625, 20736, 20790,
440 20800, 21000, 21060, 21120, 21168, 21384, 21504, 21560, 21600, 21609,
441 21840, 21870, 21875, 21952, 22000, 22050, 22113, 22176, 22275, 22295,
442 22400, 22464, 22500, 22528, 22638, 22680, 22750, 22932, 23040, 23100,
443 23296, 23328, 23400, 23520, 23625, 23760, 23814, 24000, 24010, 24057,
444 24192, 24255, 24300, 24375, 24500, 24570, 24576, 24640, 24696, 24750,
445 24948, 24960, 25000, 25088, 25200, 25272, 25344, 25480, 25515, 25600,
446 25725, 25872, 25920, 26000, 26208, 26244, 26250, 26325, 26400, 26411,
447 26460, 26624, 26730, 26754, 26880, 26950, 27000, 27216, 27300, 27440,
448 27500, 27648, 27720, 27783, 28000, 28080, 28125, 28160, 28224, 28350,
449 28431, 28512, 28665, 28672, 28800, 28812, 28875, 29106, 29120, 29160,
450 29250, 29400, 29484, 29568, 29700, 29952, 30000, 30184, 30240, 30375,
451 30576, 30618, 30625, 30720, 30800, 30870, 31104, 31185, 31200, 31213,
452 31250, 31360, 31500, 31590, 31680, 31752, 31850, 32000, 32076, 32256,
453 32340, 32400, 32500, 32760, 32768, 32805, 32928, 33000, 33075, 33264,
454 33280, 33600, 33614, 33696, 33750, 33792, 33957, 34020, 34125, 34300,
455 34375, 34398, 34496, 34560, 34650, 34944, 34992, 35000, 35100, 35200,
456 35280, 35640, 35672, 35721, 35840, 36000, 36015, 36288, 36400, 36450,
457 36750, 36855, 36864, 36960, 37044, 37125, 37422, 37440, 37500, 37632,
458 37730, 37800, 37908, 38016, 38220, 38400, 38416, 38500, 38808, 38880,
459 39000, 39200, 39312, 39366, 39375, 39424, 39600, 39690, 39936, 40000,
460 40095, 40131, 40320, 40425, 40500, 40625, 40768, 40824, 40950, 40960,
461 41160, 41250, 41472, 41580, 41600, 42000, 42120, 42240, 42336, 42525,
462 42768, 42875, 43008, 43120, 43200, 43218, 43659, 43680, 43740, 43750,
463 43875, 43904, 44000, 44100, 44226, 44352, 44550, 44590, 44800, 44928,
464 45000, 45056, 45276, 45360, 45500, 45864, 45927, 46080, 46200, 46305,
465 46592, 46656, 46800, 46875, 47040, 47250, 47385, 47520, 47628, 47775,
466 48000, 48020, 48114, 48125, 48384, 48510, 48600, 48750, 49000, 49140,
467 49152, 49280, 49392, 49500, 49896, 49920, 50000, 50176, 50400, 50421,
468 50544, 50625, 50688, 50960, 51030, 51200, 51450, 51597, 51744, 51840,
469 51975, 52000, 52416, 52488, 52500, 52650, 52800, 52822, 52920, 53248,
470 53460, 53508, 53760, 53900, 54000, 54432, 54600, 54675, 54880, 55000,
471 55125, 55296, 55440, 55566, 56000, 56133, 56160, 56250, 56320, 56448,
472 56595, 56700, 56862, 56875, 57024, 57330, 57344, 57600, 57624, 57750,
473 58212, 58240, 58320, 58500, 58800, 58968, 59049, 59136, 59400, 59535,
474 59904, 60000, 60025, 60368, 60480, 60750, 61152, 61236, 61250, 61425,
475 61440, 61600, 61740, 61875, 62208, 62370, 62400, 62426, 62500, 62720,
476 63000, 63180, 63360, 63504, 63700, 64000, 64152, 64512, 64680, 64800,
477 64827, 65000, 65520, 65536, 65610, 65625, 65856, 66000, 66150, 66339,
478 66528, 66560, 66825, 66885, 67200, 67228, 67375, 67392, 67500, 67584,
479 67914, 68040, 68250, 68600, 68750, 68796, 68992, 69120, 69300, 69888,
480 69984, 70000, 70200, 70400, 70560, 70875, 71280, 71344, 71442, 71680,
481 72000, 72030, 72171, 72576, 72765, 72800, 72900, 73125, 73500, 73710,
482 73728, 73920, 74088, 74250, 74844, 74880, 75000, 75264, 75460, 75600,
483 75816, 76032, 76440, 76545, 76800, 76832, 77000, 77175, 77616, 77760,
484 78000, 78125, 78400, 78624, 78732, 78750, 78848, 78975, 79200, 79233,
485 79380, 79625, 79872, 80000, 80190, 80262, 80640, 80850, 81000, 81250,
486 81536, 81648, 81900, 81920, 82320, 82500, 82944, 83160, 83200, 83349,
487 84000, 84035, 84240, 84375, 84480, 84672, 85050, 85293, 85536, 85750,
488 85995, 86016, 86240, 86400, 86436, 86625, 87318, 87360, 87480, 87500,
489 87750, 87808, 88000, 88200, 88452, 88704, 89100, 89180, 89600, 89856,
490 90000, 90112, 90552, 90720, 91000, 91125, 91728, 91854, 91875, 92160,
491 92400, 92610, 93184, 93312, 93555, 93600, 93639, 93750, 94080, 94325,
492 94500, 94770, 95040, 95256, 95550, 96000, 96040, 96228, 96250, 96768,
493 97020, 97200, 97500, 98000, 98280, 98304, 98415, 98560, 98784, 99000,
494 99225, 99792, 99840 };
496 if(s <= 0 || s >= goodSizes[size-1])
499 int * upperBound = std::upper_bound(goodSizes, goodSizes+size, s);
500 if(upperBound > goodSizes && upperBound[-1] == s)
507 int fftwEvenPaddingSize(
int s)
513 static const int size = 1063;
514 static int goodSizes[size] = { 2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22,
515 24, 26, 28, 30, 32, 36, 40, 42, 44, 48, 50, 52, 54, 56, 60, 64, 66, 70,
516 72, 78, 80, 84, 88, 90, 96, 98, 100, 104, 108, 110, 112, 120, 126, 128,
517 130, 132, 140, 144, 150, 154, 156, 160, 162, 168, 176, 180, 182, 192,
518 196, 198, 200, 208, 210, 216, 220, 224, 234, 240, 250, 252, 256, 260,
519 264, 270, 280, 288, 294, 300, 308, 312, 320, 324, 330, 336, 350, 352,
520 360, 364, 378, 384, 390, 392, 396, 400, 416, 420, 432, 440, 448, 450,
521 462, 468, 480, 486, 490, 500, 504, 512, 520, 528, 540, 546, 550, 560,
522 576, 588, 594, 600, 616, 624, 630, 640, 648, 650, 660, 672, 686, 700,
523 702, 704, 720, 728, 750, 756, 768, 770, 780, 784, 792, 800, 810, 832,
524 840, 864, 880, 882, 896, 900, 910, 924, 936, 960, 972, 980, 990, 1000,
525 1008, 1024, 1040, 1050, 1056, 1078, 1080, 1092, 1100, 1120, 1134, 1152,
526 1170, 1176, 1188, 1200, 1232, 1248, 1250, 1260, 1274, 1280, 1296, 1300,
527 1320, 1344, 1350, 1372, 1386, 1400, 1404, 1408, 1440, 1456, 1458, 1470,
528 1500, 1512, 1536, 1540, 1560, 1568, 1584, 1600, 1620, 1638, 1650, 1664,
529 1680, 1728, 1750, 1760, 1764, 1782, 1792, 1800, 1820, 1848, 1872, 1890,
530 1920, 1944, 1950, 1960, 1980, 2000, 2016, 2048, 2058, 2080, 2100, 2106,
531 2112, 2156, 2160, 2184, 2200, 2240, 2250, 2268, 2304, 2310, 2340, 2352,
532 2376, 2400, 2430, 2450, 2464, 2496, 2500, 2520, 2548, 2560, 2592, 2600,
533 2640, 2646, 2688, 2700, 2730, 2744, 2750, 2772, 2800, 2808, 2816, 2880,
534 2912, 2916, 2940, 2970, 3000, 3024, 3072, 3080, 3120, 3136, 3150, 3168,
535 3200, 3234, 3240, 3250, 3276, 3300, 3328, 3360, 3402, 3430, 3456, 3500,
536 3510, 3520, 3528, 3564, 3584, 3600, 3640, 3696, 3744, 3750, 3780, 3822,
537 3840, 3850, 3888, 3900, 3920, 3960, 4000, 4032, 4050, 4096, 4116, 4158,
538 4160, 4200, 4212, 4224, 4312, 4320, 4368, 4374, 4400, 4410, 4480, 4500,
539 4536, 4550, 4608, 4620, 4680, 4704, 4752, 4800, 4802, 4860, 4900, 4914,
540 4928, 4950, 4992, 5000, 5040, 5096, 5120, 5184, 5200, 5250, 5280, 5292,
541 5346, 5376, 5390, 5400, 5460, 5488, 5500, 5544, 5600, 5616, 5632, 5670,
542 5760, 5824, 5832, 5850, 5880, 5940, 6000, 6048, 6144, 6160, 6174, 6240,
543 6250, 6272, 6300, 6318, 6336, 6370, 6400, 6468, 6480, 6500, 6552, 6600,
544 6656, 6720, 6750, 6804, 6860, 6912, 6930, 7000, 7020, 7040, 7056, 7128,
545 7168, 7200, 7280, 7290, 7350, 7392, 7488, 7500, 7546, 7560, 7644, 7680,
546 7700, 7776, 7800, 7840, 7920, 7938, 8000, 8064, 8100, 8190, 8192, 8232,
547 8250, 8316, 8320, 8400, 8424, 8448, 8624, 8640, 8736, 8748, 8750, 8800,
548 8820, 8910, 8918, 8960, 9000, 9072, 9100, 9216, 9240, 9360, 9408, 9450,
549 9504, 9600, 9604, 9702, 9720, 9750, 9800, 9828, 9856, 9900, 9984, 10000,
550 10080, 10192, 10206, 10240, 10290, 10368, 10400, 10500, 10530, 10560,
551 10584, 10692, 10752, 10780, 10800, 10920, 10976, 11000, 11088, 11200,
552 11232, 11250, 11264, 11340, 11466, 11520, 11550, 11648, 11664, 11700,
553 11760, 11880, 12000, 12096, 12150, 12250, 12288, 12320, 12348, 12474,
554 12480, 12500, 12544, 12600, 12636, 12672, 12740, 12800, 12936, 12960,
555 13000, 13104, 13122, 13200, 13230, 13312, 13440, 13500, 13608, 13650,
556 13720, 13750, 13824, 13860, 14000, 14040, 14080, 14112, 14256, 14336,
557 14400, 14406, 14560, 14580, 14700, 14742, 14784, 14850, 14976, 15000,
558 15092, 15120, 15288, 15360, 15400, 15552, 15600, 15680, 15750, 15840,
559 15876, 16000, 16038, 16128, 16170, 16200, 16250, 16380, 16384, 16464,
560 16500, 16632, 16640, 16800, 16848, 16896, 17010, 17150, 17248, 17280,
561 17472, 17496, 17500, 17550, 17600, 17640, 17820, 17836, 17920, 18000,
562 18144, 18200, 18432, 18480, 18522, 18720, 18750, 18816, 18900, 18954,
563 19008, 19110, 19200, 19208, 19250, 19404, 19440, 19500, 19600, 19656,
564 19712, 19800, 19968, 20000, 20160, 20250, 20384, 20412, 20480, 20580,
565 20736, 20790, 20800, 21000, 21060, 21120, 21168, 21384, 21504, 21560,
566 21600, 21840, 21870, 21952, 22000, 22050, 22176, 22400, 22464, 22500,
567 22528, 22638, 22680, 22750, 22932, 23040, 23100, 23296, 23328, 23400,
568 23520, 23760, 23814, 24000, 24010, 24192, 24300, 24500, 24570, 24576,
569 24640, 24696, 24750, 24948, 24960, 25000, 25088, 25200, 25272, 25344,
570 25480, 25600, 25872, 25920, 26000, 26208, 26244, 26250, 26400, 26460,
571 26624, 26730, 26754, 26880, 26950, 27000, 27216, 27300, 27440, 27500,
572 27648, 27720, 28000, 28080, 28160, 28224, 28350, 28512, 28672, 28800,
573 28812, 29106, 29120, 29160, 29250, 29400, 29484, 29568, 29700, 29952,
574 30000, 30184, 30240, 30576, 30618, 30720, 30800, 30870, 31104, 31200,
575 31250, 31360, 31500, 31590, 31680, 31752, 31850, 32000, 32076, 32256,
576 32340, 32400, 32500, 32760, 32768, 32928, 33000, 33264, 33280, 33600,
577 33614, 33696, 33750, 33792, 34020, 34300, 34398, 34496, 34560, 34650,
578 34944, 34992, 35000, 35100, 35200, 35280, 35640, 35672, 35840, 36000,
579 36288, 36400, 36450, 36750, 36864, 36960, 37044, 37422, 37440, 37500,
580 37632, 37730, 37800, 37908, 38016, 38220, 38400, 38416, 38500, 38808,
581 38880, 39000, 39200, 39312, 39366, 39424, 39600, 39690, 39936, 40000,
582 40320, 40500, 40768, 40824, 40950, 40960, 41160, 41250, 41472, 41580,
583 41600, 42000, 42120, 42240, 42336, 42768, 43008, 43120, 43200, 43218,
584 43680, 43740, 43750, 43904, 44000, 44100, 44226, 44352, 44550, 44590,
585 44800, 44928, 45000, 45056, 45276, 45360, 45500, 45864, 46080, 46200,
586 46592, 46656, 46800, 47040, 47250, 47520, 47628, 48000, 48020, 48114,
587 48384, 48510, 48600, 48750, 49000, 49140, 49152, 49280, 49392, 49500,
588 49896, 49920, 50000, 50176, 50400, 50544, 50688, 50960, 51030, 51200,
589 51450, 51744, 51840, 52000, 52416, 52488, 52500, 52650, 52800, 52822,
590 52920, 53248, 53460, 53508, 53760, 53900, 54000, 54432, 54600, 54880,
591 55000, 55296, 55440, 55566, 56000, 56160, 56250, 56320, 56448, 56700,
592 56862, 57024, 57330, 57344, 57600, 57624, 57750, 58212, 58240, 58320,
593 58500, 58800, 58968, 59136, 59400, 59904, 60000, 60368, 60480, 60750,
594 61152, 61236, 61250, 61440, 61600, 61740, 62208, 62370, 62400, 62426,
595 62500, 62720, 63000, 63180, 63360, 63504, 63700, 64000, 64152, 64512,
596 64680, 64800, 65000, 65520, 65536, 65610, 65856, 66000, 66150, 66528,
597 66560, 67200, 67228, 67392, 67500, 67584, 67914, 68040, 68250, 68600,
598 68750, 68796, 68992, 69120, 69300, 69888, 69984, 70000, 70200, 70400,
599 70560, 71280, 71344, 71442, 71680, 72000, 72030, 72576, 72800, 72900,
600 73500, 73710, 73728, 73920, 74088, 74250, 74844, 74880, 75000, 75264,
601 75460, 75600, 75816, 76032, 76440, 76800, 76832, 77000, 77616, 77760,
602 78000, 78400, 78624, 78732, 78750, 78848, 79200, 79380, 79872, 80000,
603 80190, 80262, 80640, 80850, 81000, 81250, 81536, 81648, 81900, 81920,
604 82320, 82500, 82944, 83160, 83200, 84000, 84240, 84480, 84672, 85050,
605 85536, 85750, 86016, 86240, 86400, 86436, 87318, 87360, 87480, 87500,
606 87750, 87808, 88000, 88200, 88452, 88704, 89100, 89180, 89600, 89856,
607 90000, 90112, 90552, 90720, 91000, 91728, 91854, 92160, 92400, 92610,
608 93184, 93312, 93600, 93750, 94080, 94500, 94770, 95040, 95256, 95550,
609 96000, 96040, 96228, 96250, 96768, 97020, 97200, 97500, 98000, 98280,
610 98304, 98560, 98784, 99000, 99792, 99840 };
612 if(s <= 0 || s >= goodSizes[size-1])
615 int * upperBound = std::upper_bound(goodSizes, goodSizes+size, s);
616 if(upperBound > goodSizes && upperBound[-1] == s)
623 struct FFTEmbedKernel
625 template <
unsigned int N,
class Real,
class C,
class Shape>
627 exec(MultiArrayView<N, Real, C> & out, Shape
const & kernelShape,
628 Shape & srcPoint, Shape & destPoint,
bool copyIt)
630 for(srcPoint[M]=0; srcPoint[M]<kernelShape[M]; ++srcPoint[M])
632 if(srcPoint[M] < (kernelShape[M] + 1) / 2)
634 destPoint[M] = srcPoint[M];
638 destPoint[M] = srcPoint[M] + out.shape(M) - kernelShape[M];
641 FFTEmbedKernel<M-1>::exec(out, kernelShape, srcPoint, destPoint, copyIt);
647 struct FFTEmbedKernel<0>
649 template <
unsigned int N,
class Real,
class C,
class Shape>
651 exec(MultiArrayView<N, Real, C> & out, Shape
const & kernelShape,
652 Shape & srcPoint, Shape & destPoint,
bool copyIt)
654 for(srcPoint[0]=0; srcPoint[0]<kernelShape[0]; ++srcPoint[0])
656 if(srcPoint[0] < (kernelShape[0] + 1) / 2)
658 destPoint[0] = srcPoint[0];
662 destPoint[0] = srcPoint[0] + out.shape(0) - kernelShape[0];
667 out[destPoint] = out[srcPoint];
674 template <
unsigned int N,
class Real,
class C1,
class C2>
676 fftEmbedKernel(MultiArrayView<N, Real, C1> kernel,
677 MultiArrayView<N, Real, C2> out,
680 typedef typename MultiArrayShape<N>::type Shape;
682 MultiArrayView<N, Real, C2> kout = out.subarray(Shape(), kernel.shape());
690 Shape srcPoint, destPoint;
691 FFTEmbedKernel<(int)N-1>::exec(out, kernel.shape(), srcPoint, destPoint,
false);
694 template <
unsigned int N,
class Real,
class C1,
class C2>
696 fftEmbedArray(MultiArrayView<N, Real, C1> in,
697 MultiArrayView<N, Real, C2> out)
699 typedef typename MultiArrayShape<N>::type Shape;
701 Shape diff = out.shape() - in.shape(),
703 rightDiff = diff - leftDiff,
704 right = in.shape() + leftDiff;
706 out.subarray(leftDiff, right) = in;
708 typedef typename MultiArrayView<N, Real, C2>::traverser Traverser;
709 typedef MultiArrayNavigator<Traverser, N> Navigator;
710 typedef typename Navigator::iterator Iterator;
712 for(
unsigned int d = 0; d < N; ++d)
714 Navigator nav(out.traverser_begin(), out.shape(), d);
716 for( ; nav.hasMore(); nav++ )
718 Iterator i = nav.begin();
719 for(
int k=1; k<=leftDiff[d]; ++k)
720 i[leftDiff[d] - k] = i[leftDiff[d] + k];
721 for(
int k=0; k<rightDiff[d]; ++k)
722 i[right[d] + k] = i[right[d] - k - 2];
729 template <
class T,
int N>
731 fftwBestPaddedShape(TinyVector<T, N> shape)
733 for(
unsigned int k=0; k<N; ++k)
734 shape[k] = detail::fftwPaddingSize(shape[k]);
738 template <
class T,
int N>
740 fftwBestPaddedShapeR2C(TinyVector<T, N> shape)
742 shape[0] = detail::fftwEvenPaddingSize(shape[0]);
743 for(
unsigned int k=1; k<N; ++k)
744 shape[k] = detail::fftwPaddingSize(shape[k]);
758 template <
class T,
int N>
762 shape[0] = shape[0] / 2 + 1;
766 template <
class T,
int N>
768 fftwCorrespondingShapeC2R(TinyVector<T, N> shape,
bool oddDimension0 =
false)
770 shape[0] = oddDimension0
771 ? (shape[0] - 1) * 2 + 1
772 : (shape[0] - 1) * 2;
809 template <
unsigned int N,
class Real =
double>
813 typedef typename FFTWReal2Complex<Real>::plan_type PlanType;
817 Shape shape, instrides, outstrides;
837 template <
class C1,
class C2>
840 int SIGN,
unsigned int planner_flags = FFTW_ESTIMATE)
843 init(in, out, SIGN, planner_flags);
857 template <
class C1,
class C2>
860 unsigned int planner_flags = FFTW_ESTIMATE)
863 init(in, out, planner_flags);
877 template <
class C1,
class C2>
880 unsigned int planner_flags = FFTW_ESTIMATE)
883 init(in, out, planner_flags);
894 instrides.swap(o.instrides);
895 outstrides.swap(o.outstrides);
908 instrides.swap(o.instrides);
909 outstrides.swap(o.outstrides);
920 detail::fftwPlanDestroy(plan);
927 template <
class C1,
class C2>
930 int SIGN,
unsigned int planner_flags = FFTW_ESTIMATE)
932 vigra_precondition(in.strideOrdering() == out.strideOrdering(),
933 "FFTWPlan.init(): input and output must have the same stride ordering.");
935 initImpl(in.permuteStridesDescending(), out.permuteStridesDescending(),
936 SIGN, planner_flags);
943 template <
class C1,
class C2>
946 unsigned int planner_flags = FFTW_ESTIMATE)
949 "FFTWPlan.init(): input and output must have the same stride ordering.");
952 FFTW_FORWARD, planner_flags);
959 template <
class C1,
class C2>
962 unsigned int planner_flags = FFTW_ESTIMATE)
965 "FFTWPlan.init(): input and output must have the same stride ordering.");
968 FFTW_BACKWARD, planner_flags);
978 template <
class C1,
class C2>
982 executeImpl(in.permuteStridesDescending(), out.permuteStridesDescending());
992 template <
class C1,
class C2>
1006 template <
class C1,
class C2>
1015 template <
class MI,
class MO>
1016 void initImpl(MI ins, MO outs,
int SIGN,
unsigned int planner_flags);
1018 template <
class MI,
class MO>
1019 void executeImpl(MI ins, MO outs)
const;
1024 vigra_precondition(in.shape() == out.shape(),
1025 "FFTWPlan.init(): input and output must have the same shape.");
1028 void checkShapes(MultiArrayView<N, Real, StridedArrayTag> ins,
1029 MultiArrayView<N, FFTWComplex<Real>, StridedArrayTag> outs)
const
1031 for(
int k=0; k<(int)N-1; ++k)
1032 vigra_precondition(ins.shape(k) == outs.shape(k),
1033 "FFTWPlan.init(): input and output must have matching shapes.");
1034 vigra_precondition(ins.shape(N-1) / 2 + 1 == outs.shape(N-1),
1035 "FFTWPlan.init(): input and output must have matching shapes.");
1038 void checkShapes(MultiArrayView<N, FFTWComplex<Real>, StridedArrayTag> ins,
1039 MultiArrayView<N, Real, StridedArrayTag> outs)
const
1041 for(
int k=0; k<(int)N-1; ++k)
1042 vigra_precondition(ins.shape(k) == outs.shape(k),
1043 "FFTWPlan.init(): input and output must have matching shapes.");
1044 vigra_precondition(outs.shape(N-1) / 2 + 1 == ins.shape(N-1),
1045 "FFTWPlan.init(): input and output must have matching shapes.");
1049 template <
unsigned int N,
class Real>
1050 template <
class MI,
class MO>
1052 FFTWPlan<N, Real>::initImpl(MI ins, MO outs,
int SIGN,
unsigned int planner_flags)
1054 checkShapes(ins, outs);
1056 typename MultiArrayShape<N>::type logicalShape(SIGN == FFTW_FORWARD
1060 Shape newShape(logicalShape.begin(), logicalShape.end()),
1061 newIStrides(ins.stride().begin(), ins.stride().end()),
1062 newOStrides(outs.stride().begin(), outs.stride().end()),
1063 itotal(ins.shape().begin(), ins.shape().end()),
1064 ototal(outs.shape().begin(), outs.shape().end());
1066 for(
unsigned int j=1; j<N; ++j)
1068 itotal[j] = ins.stride(j-1) / ins.stride(j);
1069 ototal[j] = outs.stride(j-1) / outs.stride(j);
1072 PlanType newPlan = detail::fftwPlanCreate(N, newShape.begin(),
1073 ins.data(), itotal.begin(), ins.stride(N-1),
1074 outs.data(), ototal.begin(), outs.stride(N-1),
1075 SIGN, planner_flags);
1076 detail::fftwPlanDestroy(plan);
1078 shape.swap(newShape);
1079 instrides.swap(newIStrides);
1080 outstrides.swap(newOStrides);
1084 template <
unsigned int N,
class Real>
1085 template <
class MI,
class MO>
1086 void FFTWPlan<N, Real>::executeImpl(MI ins, MO outs)
const
1088 vigra_precondition(plan != 0,
"FFTWPlan::execute(): plan is NULL.");
1090 typename MultiArrayShape<N>::type lshape(sign == FFTW_FORWARD
1094 vigra_precondition((lshape == TinyVectorView<int, N>(shape.data())),
1095 "FFTWPlan::execute(): shape mismatch between plan and data.");
1096 vigra_precondition((ins.stride() == TinyVectorView<int, N>(instrides.data())),
1097 "FFTWPlan::execute(): strides mismatch between plan and input data.");
1098 vigra_precondition((outs.stride() == TinyVectorView<int, N>(outstrides.data())),
1099 "FFTWPlan::execute(): strides mismatch between plan and output data.");
1101 detail::fftwPlanExecute(plan, ins.data(), outs.data());
1103 typedef typename MO::value_type V;
1104 if(sign == FFTW_BACKWARD)
1105 outs *= V(1.0) / Real(outs.size());
1147 template <
unsigned int N,
class Real>
1155 RArray realArray, realKernel;
1156 CArray fourierArray, fourierKernel;
1157 bool useFourierKernel;
1168 : useFourierKernel(false)
1182 template <
class C1,
class C2,
class C3>
1186 unsigned int planner_flags = FFTW_ESTIMATE)
1187 : useFourierKernel(false)
1189 init(in, kernel, out, planner_flags);
1203 template <
class C1,
class C2,
class C3>
1207 unsigned int planner_flags = FFTW_ESTIMATE)
1208 : useFourierKernel(true)
1210 init(in, kernel, out, planner_flags);
1225 template <
class C1,
class C2,
class C3>
1229 bool fourierDomainKernel,
1230 unsigned int planner_flags = FFTW_ESTIMATE)
1232 init(in, kernel, out, fourierDomainKernel, planner_flags);
1248 template <
class C1,
class C2,
class C3>
1250 bool useFourierKernel =
false,
1251 unsigned int planner_flags = FFTW_ESTIMATE)
1253 if(useFourierKernel)
1254 init(inOut, kernel, planner_flags);
1256 initFourierKernel(inOut, kernel, planner_flags);
1263 template <
class C1,
class C2,
class C3>
1267 unsigned int planner_flags = FFTW_ESTIMATE)
1269 vigra_precondition(in.
shape() == out.
shape(),
1270 "FFTWConvolvePlan::init(): input and output must have the same shape.");
1278 template <
class C1,
class C2,
class C3>
1282 unsigned int planner_flags = FFTW_ESTIMATE)
1284 vigra_precondition(in.
shape() == out.
shape(),
1285 "FFTWConvolvePlan::init(): input and output must have the same shape.");
1286 initFourierKernel(in.
shape(), kernel.shape(), planner_flags);
1293 template <
class C1,
class C2,
class C3>
1297 bool fourierDomainKernel,
1298 unsigned int planner_flags = FFTW_ESTIMATE)
1300 vigra_precondition(in.shape() == out.shape(),
1301 "FFTWConvolvePlan::init(): input and output must have the same shape.");
1302 useFourierKernel = fourierDomainKernel;
1303 initComplex(in.shape(), kernel.shape(), planner_flags);
1312 template <
class C1,
class KernelIterator,
class OutIterator>
1314 KernelIterator kernels, KernelIterator kernelsEnd,
1315 OutIterator outs,
unsigned int planner_flags = FFTW_ESTIMATE)
1317 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1318 typedef typename KernelArray::value_type KernelValue;
1319 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1320 typedef typename OutArray::value_type OutValue;
1322 bool realKernel = IsSameType<KernelValue, Real>::value;
1323 bool fourierKernel = IsSameType<KernelValue, Complex>::value;
1325 vigra_precondition(realKernel || fourierKernel,
1326 "FFTWConvolvePlan::initMany(): kernels have unsuitable value_type.");
1327 vigra_precondition((IsSameType<OutValue, Real>::value),
1328 "FFTWConvolvePlan::initMany(): outputs have unsuitable value_type.");
1337 initFourierKernelMany(in.
shape(),
1338 checkShapesFourier(in.
shape(), kernels, kernelsEnd, outs),
1349 template <
class C1,
class KernelIterator,
class OutIterator>
1351 KernelIterator kernels, KernelIterator kernelsEnd,
1353 bool fourierDomainKernels,
1354 unsigned int planner_flags = FFTW_ESTIMATE)
1356 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1357 typedef typename KernelArray::value_type KernelValue;
1358 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1359 typedef typename OutArray::value_type OutValue;
1361 vigra_precondition((IsSameType<KernelValue, Complex>::value),
1362 "FFTWConvolvePlan::initMany(): kernels have unsuitable value_type.");
1363 vigra_precondition((IsSameType<OutValue, Complex>::value),
1364 "FFTWConvolvePlan::initMany(): outputs have unsuitable value_type.");
1366 useFourierKernel = fourierDomainKernels;
1368 Shape paddedShape = checkShapesComplex(in.shape(), kernels, kernelsEnd, outs);
1370 CArray newFourierArray(paddedShape), newFourierKernel(paddedShape);
1372 FFTWPlan<N, Real> fplan(newFourierArray, newFourierArray, FFTW_FORWARD, planner_flags);
1373 FFTWPlan<N, Real> bplan(newFourierArray, newFourierArray, FFTW_BACKWARD, planner_flags);
1375 forward_plan = fplan;
1376 backward_plan = bplan;
1377 fourierArray.
swap(newFourierArray);
1378 fourierKernel.
swap(newFourierKernel);
1381 void init(Shape inOut, Shape kernel,
1382 unsigned int planner_flags = FFTW_ESTIMATE);
1384 void initFourierKernel(Shape inOut, Shape kernel,
1385 unsigned int planner_flags = FFTW_ESTIMATE);
1387 void initComplex(Shape inOut, Shape kernel,
1388 unsigned int planner_flags = FFTW_ESTIMATE);
1390 void initMany(Shape inOut, Shape maxKernel,
1391 unsigned int planner_flags = FFTW_ESTIMATE)
1393 init(inOut, maxKernel, planner_flags);
1396 void initFourierKernelMany(Shape inOut, Shape kernels,
1397 unsigned int planner_flags = FFTW_ESTIMATE)
1399 initFourierKernel(inOut, kernels, planner_flags);
1409 template <
class C1,
class C2,
class C3>
1410 void execute(MultiArrayView<N, Real, C1> in,
1411 MultiArrayView<N, Real, C2> kernel,
1412 MultiArrayView<N, Real, C3> out);
1421 template <
class C1,
class C2,
class C3>
1422 void execute(MultiArrayView<N, Real, C1> in,
1423 MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
1424 MultiArrayView<N, Real, C3> out);
1433 template <
class C1,
class C2,
class C3>
1434 void execute(MultiArrayView<N, FFTWComplex<Real>, C1> in,
1435 MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
1436 MultiArrayView<N, FFTWComplex<Real>, C3> out);
1446 template <
class C1,
class KernelIterator,
class OutIterator>
1447 void executeMany(MultiArrayView<N, FFTWComplex<Real>, C1> in,
1448 KernelIterator kernels, KernelIterator kernelsEnd,
1458 template <
class C1,
class KernelIterator,
class OutIterator>
1460 KernelIterator kernels, KernelIterator kernelsEnd,
1463 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1464 typedef typename KernelArray::value_type KernelValue;
1465 typedef typename IsSameType<KernelValue, Complex>::type UseFourierKernel;
1466 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1467 typedef typename OutArray::value_type OutValue;
1469 bool realKernel = IsSameType<KernelValue, Real>::value;
1470 bool fourierKernel = IsSameType<KernelValue, Complex>::value;
1472 vigra_precondition(realKernel || fourierKernel,
1473 "FFTWConvolvePlan::executeMany(): kernels have unsuitable value_type.");
1474 vigra_precondition((IsSameType<OutValue, Real>::value),
1475 "FFTWConvolvePlan::executeMany(): outputs have unsuitable value_type.");
1477 executeManyImpl(in, kernels, kernelsEnd, outs, UseFourierKernel());
1482 template <
class KernelIterator,
class OutIterator>
1483 Shape checkShapes(Shape in,
1484 KernelIterator kernels, KernelIterator kernelsEnd,
1487 template <
class KernelIterator,
class OutIterator>
1488 Shape checkShapesFourier(Shape in,
1489 KernelIterator kernels, KernelIterator kernelsEnd,
1492 template <
class KernelIterator,
class OutIterator>
1493 Shape checkShapesComplex(Shape in,
1494 KernelIterator kernels, KernelIterator kernelsEnd,
1497 template <
class C1,
class KernelIterator,
class OutIterator>
1500 KernelIterator kernels, KernelIterator kernelsEnd,
1501 OutIterator outs, VigraFalseType );
1503 template <
class C1,
class KernelIterator,
class OutIterator>
1506 KernelIterator kernels, KernelIterator kernelsEnd,
1507 OutIterator outs, VigraTrueType );
1511 template <
unsigned int N,
class Real>
1514 unsigned int planner_flags)
1516 Shape paddedShape = fftwBestPaddedShapeR2C(in + kernel - Shape(1)),
1519 CArray newFourierArray(complexShape), newFourierKernel(complexShape);
1521 Shape realStrides = 2*newFourierArray.stride();
1523 RArray newRealArray(paddedShape, realStrides, (Real*)newFourierArray.data());
1524 RArray newRealKernel(paddedShape, realStrides, (Real*)newFourierKernel.data());
1526 FFTWPlan<N, Real> fplan(newRealArray, newFourierArray, planner_flags);
1527 FFTWPlan<N, Real> bplan(newFourierArray, newRealArray, planner_flags);
1529 forward_plan = fplan;
1530 backward_plan = bplan;
1531 realArray = newRealArray;
1532 realKernel = newRealKernel;
1533 fourierArray.swap(newFourierArray);
1534 fourierKernel.swap(newFourierKernel);
1535 useFourierKernel =
false;
1538 template <
unsigned int N,
class Real>
1540 FFTWConvolvePlan<N, Real>::initFourierKernel(Shape in, Shape kernel,
1541 unsigned int planner_flags)
1543 Shape complexShape = kernel,
1544 paddedShape = fftwCorrespondingShapeC2R(complexShape);
1546 for(
unsigned int k=0; k<N; ++k)
1547 vigra_precondition(in[k] <= paddedShape[k],
1548 "FFTWConvolvePlan::init(): kernel too small for given input.");
1550 CArray newFourierArray(complexShape), newFourierKernel(complexShape);
1552 Shape realStrides = 2*newFourierArray.stride();
1554 RArray newRealArray(paddedShape, realStrides, (Real*)newFourierArray.data());
1555 RArray newRealKernel(paddedShape, realStrides, (Real*)newFourierKernel.data());
1557 FFTWPlan<N, Real> fplan(newRealArray, newFourierArray, planner_flags);
1558 FFTWPlan<N, Real> bplan(newFourierArray, newRealArray, planner_flags);
1560 forward_plan = fplan;
1561 backward_plan = bplan;
1562 realArray = newRealArray;
1563 realKernel = newRealKernel;
1564 fourierArray.swap(newFourierArray);
1565 fourierKernel.swap(newFourierKernel);
1566 useFourierKernel =
true;
1569 template <
unsigned int N,
class Real>
1571 FFTWConvolvePlan<N, Real>::initComplex(Shape in, Shape kernel,
1572 unsigned int planner_flags)
1576 if(useFourierKernel)
1578 for(
unsigned int k=0; k<N; ++k)
1579 vigra_precondition(in[k] <= kernel[k],
1580 "FFTWConvolvePlan::init(): kernel too small for given input.");
1582 paddedShape = kernel;
1586 paddedShape = fftwBestPaddedShape(in + kernel - Shape(1));
1589 CArray newFourierArray(paddedShape), newFourierKernel(paddedShape);
1591 FFTWPlan<N, Real> fplan(newFourierArray, newFourierArray, FFTW_FORWARD, planner_flags);
1592 FFTWPlan<N, Real> bplan(newFourierArray, newFourierArray, FFTW_BACKWARD, planner_flags);
1594 forward_plan = fplan;
1595 backward_plan = bplan;
1596 fourierArray.swap(newFourierArray);
1597 fourierKernel.swap(newFourierKernel);
1600 #ifndef DOXYGEN // doxygen documents these functions as free functions
1602 template <
unsigned int N,
class Real>
1603 template <
class C1,
class C2,
class C3>
1606 MultiArrayView<N, Real, C2> kernel,
1607 MultiArrayView<N, Real, C3> out)
1609 vigra_precondition(!useFourierKernel,
1610 "FFTWConvolvePlan::execute(): plan was generated for Fourier kernel, got spatial kernel.");
1612 vigra_precondition(in.shape() == out.shape(),
1613 "FFTWConvolvePlan::execute(): input and output must have the same shape.");
1615 Shape paddedShape = fftwBestPaddedShapeR2C(in.shape() + kernel.shape() - Shape(1)),
1616 diff = paddedShape - in.shape(),
1618 right = in.shape() + left;
1620 vigra_precondition(paddedShape == realArray.shape(),
1621 "FFTWConvolvePlan::execute(): shape mismatch between input and plan.");
1623 detail::fftEmbedArray(in, realArray);
1624 forward_plan.execute(realArray, fourierArray);
1626 detail::fftEmbedKernel(kernel, realKernel);
1627 forward_plan.execute(realKernel, fourierKernel);
1629 fourierArray *= fourierKernel;
1631 backward_plan.execute(fourierArray, realArray);
1633 out = realArray.subarray(left, right);
1636 template <
unsigned int N,
class Real>
1637 template <
class C1,
class C2,
class C3>
1640 MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
1641 MultiArrayView<N, Real, C3> out)
1643 vigra_precondition(useFourierKernel,
1644 "FFTWConvolvePlan::execute(): plan was generated for spatial kernel, got Fourier kernel.");
1646 vigra_precondition(in.shape() == out.shape(),
1647 "FFTWConvolvePlan::execute(): input and output must have the same shape.");
1649 vigra_precondition(kernel.shape() == fourierArray.shape(),
1650 "FFTWConvolvePlan::execute(): shape mismatch between kernel and plan.");
1652 Shape paddedShape = fftwCorrespondingShapeC2R(kernel.shape(),
odd(in.shape(0))),
1653 diff = paddedShape - in.shape(),
1655 right = in.shape() + left;
1657 vigra_precondition(paddedShape == realArray.shape(),
1658 "FFTWConvolvePlan::execute(): shape mismatch between input and plan.");
1660 detail::fftEmbedArray(in, realArray);
1661 forward_plan.execute(realArray, fourierArray);
1663 fourierKernel = kernel;
1664 moveDCToHalfspaceUpperLeft(fourierKernel);
1666 fourierArray *= fourierKernel;
1668 backward_plan.execute(fourierArray, realArray);
1670 out = realArray.subarray(left, right);
1673 template <
unsigned int N,
class Real>
1674 template <
class C1,
class C2,
class C3>
1677 MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
1678 MultiArrayView<N, FFTWComplex<Real>, C3> out)
1680 vigra_precondition(in.shape() == out.shape(),
1681 "FFTWConvolvePlan::execute(): input and output must have the same shape.");
1683 Shape paddedShape = fourierArray.shape(),
1684 diff = paddedShape - in.shape(),
1686 right = in.shape() + left;
1688 if(useFourierKernel)
1690 vigra_precondition(kernel.shape() == fourierArray.shape(),
1691 "FFTWConvolvePlan::execute(): shape mismatch between kernel and plan.");
1693 fourierKernel = kernel;
1698 detail::fftEmbedKernel(kernel, fourierKernel);
1699 forward_plan.execute(fourierKernel, fourierKernel);
1702 detail::fftEmbedArray(in, fourierArray);
1703 forward_plan.execute(fourierArray, fourierArray);
1705 fourierArray *= fourierKernel;
1707 backward_plan.execute(fourierArray, fourierArray);
1709 out = fourierArray.subarray(left, right);
1712 template <
unsigned int N,
class Real>
1713 template <
class C1,
class KernelIterator,
class OutIterator>
1715 FFTWConvolvePlan<N, Real>::executeManyImpl(MultiArrayView<N, Real, C1> in,
1716 KernelIterator kernels, KernelIterator kernelsEnd,
1717 OutIterator outs, VigraFalseType )
1719 vigra_precondition(!useFourierKernel,
1720 "FFTWConvolvePlan::execute(): plan was generated for Fourier kernel, got spatial kernel.");
1722 Shape kernelMax = checkShapes(in.shape(), kernels, kernelsEnd, outs),
1723 paddedShape = fftwBestPaddedShapeR2C(in.shape() + kernelMax - Shape(1)),
1724 diff = paddedShape - in.shape(),
1726 right = in.shape() + left;
1728 vigra_precondition(paddedShape == realArray.shape(),
1729 "FFTWConvolvePlan::executeMany(): shape mismatch between input and plan.");
1731 detail::fftEmbedArray(in, realArray);
1732 forward_plan.execute(realArray, fourierArray);
1734 for(; kernels != kernelsEnd; ++kernels, ++outs)
1736 detail::fftEmbedKernel(*kernels, realKernel);
1737 forward_plan.execute(realKernel, fourierKernel);
1739 fourierKernel *= fourierArray;
1741 backward_plan.execute(fourierKernel, realKernel);
1743 *outs = realKernel.subarray(left, right);
1747 template <
unsigned int N,
class Real>
1748 template <
class C1,
class KernelIterator,
class OutIterator>
1750 FFTWConvolvePlan<N, Real>::executeManyImpl(MultiArrayView<N, Real, C1> in,
1751 KernelIterator kernels, KernelIterator kernelsEnd,
1752 OutIterator outs, VigraTrueType )
1754 vigra_precondition(useFourierKernel,
1755 "FFTWConvolvePlan::execute(): plan was generated for spatial kernel, got Fourier kernel.");
1757 Shape complexShape = checkShapesFourier(in.shape(), kernels, kernelsEnd, outs),
1758 paddedShape = fftwCorrespondingShapeC2R(complexShape,
odd(in.shape(0))),
1759 diff = paddedShape - in.shape(),
1761 right = in.shape() + left;
1763 vigra_precondition(complexShape == fourierArray.shape(),
1764 "FFTWConvolvePlan::executeFourierKernelMany(): shape mismatch between kernels and plan.");
1766 vigra_precondition(paddedShape == realArray.shape(),
1767 "FFTWConvolvePlan::executeFourierKernelMany(): shape mismatch between input and plan.");
1769 detail::fftEmbedArray(in, realArray);
1770 forward_plan.execute(realArray, fourierArray);
1772 for(; kernels != kernelsEnd; ++kernels, ++outs)
1774 fourierKernel = *kernels;
1775 moveDCToHalfspaceUpperLeft(fourierKernel);
1776 fourierKernel *= fourierArray;
1778 backward_plan.execute(fourierKernel, realKernel);
1780 *outs = realKernel.subarray(left, right);
1784 template <
unsigned int N,
class Real>
1785 template <
class C1,
class KernelIterator,
class OutIterator>
1788 KernelIterator kernels, KernelIterator kernelsEnd,
1791 typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
1792 typedef typename KernelArray::value_type KernelValue;
1793 typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
1794 typedef typename OutArray::value_type OutValue;
1796 vigra_precondition((IsSameType<KernelValue, Complex>::value),
1797 "FFTWConvolvePlan::executeMany(): kernels have unsuitable value_type.");
1798 vigra_precondition((IsSameType<OutValue, Complex>::value),
1799 "FFTWConvolvePlan::executeMany(): outputs have unsuitable value_type.");
1801 Shape paddedShape = checkShapesComplex(in.shape(), kernels, kernelsEnd, outs),
1802 diff = paddedShape - in.shape(),
1804 right = in.shape() + left;
1806 detail::fftEmbedArray(in, fourierArray);
1807 forward_plan.execute(fourierArray, fourierArray);
1809 for(; kernels != kernelsEnd; ++kernels, ++outs)
1811 if(useFourierKernel)
1813 fourierKernel = *kernels;
1818 detail::fftEmbedKernel(*kernels, fourierKernel);
1819 forward_plan.execute(fourierKernel, fourierKernel);
1822 fourierKernel *= fourierArray;
1824 backward_plan.execute(fourierKernel, fourierKernel);
1826 *outs = fourierKernel.subarray(left, right);
1832 template <
unsigned int N,
class Real>
1833 template <
class KernelIterator,
class OutIterator>
1834 typename FFTWConvolvePlan<N, Real>::Shape
1835 FFTWConvolvePlan<N, Real>::checkShapes(Shape in,
1836 KernelIterator kernels, KernelIterator kernelsEnd,
1839 vigra_precondition(kernels != kernelsEnd,
1840 "FFTWConvolvePlan::checkShapes(): empty kernel sequence.");
1843 for(; kernels != kernelsEnd; ++kernels, ++outs)
1845 vigra_precondition(in == outs->shape(),
1846 "FFTWConvolvePlan::checkShapes(): shape mismatch between input and (one) output.");
1847 kernelMax = max(kernelMax, kernels->shape());
1849 vigra_precondition(
prod(kernelMax) > 0,
1850 "FFTWConvolvePlan::checkShapes(): all kernels have size 0.");
1854 template <
unsigned int N,
class Real>
1855 template <
class KernelIterator,
class OutIterator>
1856 typename FFTWConvolvePlan<N, Real>::Shape
1857 FFTWConvolvePlan<N, Real>::checkShapesFourier(Shape in,
1858 KernelIterator kernels, KernelIterator kernelsEnd,
1861 vigra_precondition(kernels != kernelsEnd,
1862 "FFTWConvolvePlan::checkShapesFourier(): empty kernel sequence.");
1864 Shape complexShape = kernels->shape(),
1865 paddedShape = fftwCorrespondingShapeC2R(complexShape);
1867 for(
unsigned int k=0; k<N; ++k)
1868 vigra_precondition(in[k] <= paddedShape[k],
1869 "FFTWConvolvePlan::checkShapesFourier(): kernels too small for given input.");
1871 for(; kernels != kernelsEnd; ++kernels, ++outs)
1873 vigra_precondition(in == outs->shape(),
1874 "FFTWConvolvePlan::checkShapesFourier(): shape mismatch between input and (one) output.");
1875 vigra_precondition(complexShape == kernels->shape(),
1876 "FFTWConvolvePlan::checkShapesFourier(): all kernels must have the same size.");
1878 return complexShape;
1881 template <
unsigned int N,
class Real>
1882 template <
class KernelIterator,
class OutIterator>
1883 typename FFTWConvolvePlan<N, Real>::Shape
1884 FFTWConvolvePlan<N, Real>::checkShapesComplex(Shape in,
1885 KernelIterator kernels, KernelIterator kernelsEnd,
1888 vigra_precondition(kernels != kernelsEnd,
1889 "FFTWConvolvePlan::checkShapesComplex(): empty kernel sequence.");
1891 Shape kernelShape = kernels->shape();
1892 for(; kernels != kernelsEnd; ++kernels, ++outs)
1894 vigra_precondition(in == outs->shape(),
1895 "FFTWConvolvePlan::checkShapesComplex(): shape mismatch between input and (one) output.");
1896 if(useFourierKernel)
1898 vigra_precondition(kernelShape == kernels->shape(),
1899 "FFTWConvolvePlan::checkShapesComplex(): Fourier domain kernels must have identical size.");
1903 kernelShape = max(kernelShape, kernels->shape());
1906 vigra_precondition(
prod(kernelShape) > 0,
1907 "FFTWConvolvePlan::checkShapesComplex(): all kernels have size 0.");
1909 if(useFourierKernel)
1911 for(
unsigned int k=0; k<N; ++k)
1912 vigra_precondition(in[k] <= kernelShape[k],
1913 "FFTWConvolvePlan::checkShapesComplex(): kernels too small for given input.");
1918 return fftwBestPaddedShape(in + kernelShape - Shape(1));
1929 template <
unsigned int N,
class Real,
class C1,
class C2>
1932 MultiArrayView<N, FFTWComplex<Real>, C2> out)
1934 FFTWPlan<N, Real>(in, out, FFTW_FORWARD).execute(in, out);
1937 template <
unsigned int N,
class Real,
class C1,
class C2>
1940 MultiArrayView<N, FFTWComplex<Real>, C2> out)
1942 FFTWPlan<N, Real>(in, out, FFTW_BACKWARD).execute(in, out);
1945 template <
unsigned int N,
class Real,
class C1,
class C2>
1948 MultiArrayView<N, FFTWComplex<Real>, C2> out)
1950 if(in.shape() == out.shape())
1954 FFTWPlan<N, Real>(out, out, FFTW_FORWARD).execute(out, out);
1958 FFTWPlan<N, Real>(in, out).execute(in, out);
1961 vigra_precondition(
false,
1962 "fourierTransform(): shape mismatch between input and output.");
1965 template <
unsigned int N,
class Real,
class C1,
class C2>
1968 MultiArrayView<N, Real, C2> out)
1971 "fourierTransformInverse(): shape mismatch between input and output.");
1972 FFTWPlan<N, Real>(in, out).execute(in, out);
2164 doxygen_overloaded_function(template <...>
void convolveFFT)
2166 template <
unsigned int N,
class Real,
class C1,
class C2,
class C3>
2169 MultiArrayView<N, Real, C2> kernel,
2170 MultiArrayView<N, Real, C3> out)
2172 FFTWConvolvePlan<N, Real>(in, kernel, out).execute(in, kernel, out);
2175 template <
unsigned int N,
class Real,
class C1,
class C2,
class C3>
2178 MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
2179 MultiArrayView<N, Real, C3> out)
2181 FFTWConvolvePlan<N, Real>(in, kernel, out).execute(in, kernel, out);
2190 template <
unsigned int N,
class Real,
class C1,
class C2,
class C3>
2193 MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
2194 MultiArrayView<N, FFTWComplex<Real>, C3> out,
2195 bool fourierDomainKernel)
2197 FFTWConvolvePlan<N, Real>(in, kernel, out, fourierDomainKernel).execute(in, kernel, out);
2206 template <
unsigned int N,
class Real,
class C1,
2207 class KernelIterator,
class OutIterator>
2210 KernelIterator kernels, KernelIterator kernelsEnd,
2213 FFTWConvolvePlan<N, Real> plan;
2214 plan.initMany(in, kernels, kernelsEnd, outs);
2215 plan.executeMany(in, kernels, kernelsEnd, outs);
2224 template <
unsigned int N,
class Real,
class C1,
2225 class KernelIterator,
class OutIterator>
2228 KernelIterator kernels, KernelIterator kernelsEnd,
2230 bool fourierDomainKernel)
2232 FFTWConvolvePlan<N, Real> plan;
2233 plan.initMany(in, kernels, kernelsEnd, outs, fourierDomainKernel);
2234 plan.executeMany(in, kernels, kernelsEnd, outs);
2241 #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:979
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.
void convolveFFT(...)
Convolve an array with a kernel by means of the Fourier transform.
~FFTWPlan()
Destructor.
Definition: multi_fft.hxx:918
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:960
const difference_type & shape() const
Definition: multi_array.hxx:1602
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:1279
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:1226
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:1264
Main MultiArray class containing the memory management.
Definition: multi_array.hxx:595
void execute(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out) const
Execute a real-to-complex transform.
Definition: multi_fft.hxx:993
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:1313
std::ptrdiff_t MultiArrayIndex
Definition: multi_iterator.hxx:348
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:928
void convolveFFTComplexMany(...)
Convolve a complex-valued array with a sequence of kernels by means of the Fourier transform...
FFTWComplex< R >::NormType norm(const FFTWComplex< R > &a)
norm (= magnitude)
Definition: fftw3.hxx:1037
NumericTraits< V >::Promote prod(TinyVectorBase< V, SIZE, D1, D2 > const &l)
product of the vector's elements
Definition: tinyvector.hxx:1707
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:838
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:1294
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:858
Definition: multi_fft.hxx:810
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:1350
MultiArrayView< N, T, StridedArrayTag > permuteStridesDescending() const
Definition: multi_array.hxx:2090
Definition: metaprogramming.hxx:102
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:1459
Class for fixed size vectors.This class contains an array of size SIZE of the specified VALUETYPE...
Definition: accessor.hxx:939
FixedPoint16< IntBits3, OverflowHandling > & div(FixedPoint16< IntBits1, OverflowHandling > l, FixedPoint16< IntBits2, OverflowHandling > r, FixedPoint16< IntBits3, OverflowHandling > &result)
division with enforced result type.
Definition: fixedpoint.hxx:1602
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:1204
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:878
Definition: multi_fft.hxx:1148
FFTWPlan(FFTWPlan const &other)
Copy constructor.
Definition: multi_fft.hxx:888
Base class for, and view to, vigra::MultiArray.
Definition: multi_array.hxx:593
void convolveFFTComplex(...)
Convolve a complex-valued array by means of the Fourier transform.
FFTWPlan()
Create an empty plan.
Definition: multi_fft.hxx:825
FFTWPlan & operator=(FFTWPlan const &other)
Copy assigment.
Definition: multi_fft.hxx:901
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:1249
T sign(T t)
Definition: mathutil.hxx:551
void execute(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, Real, C2 > out) const
Execute a complex-to-real transform.
Definition: multi_fft.hxx:1007
void swap(MultiArray &other)
Definition: multi_array.hxx:2924
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:944
FFTWConvolvePlan()
Create an empty plan.
Definition: multi_fft.hxx:1167
Wrapper class for the FFTW complex types 'fftw_complex'.
Definition: fftw3.hxx:131
difference_type strideOrdering() const
Definition: multi_array.hxx:1571
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:1183