36 #ifndef VIGRA_MULTI_FFT_HXX 
   37 #define VIGRA_MULTI_FFT_HXX 
   40 #include "metaprogramming.hxx" 
   41 #include "multi_array.hxx" 
   42 #include "multi_math.hxx" 
   43 #include "navigator.hxx" 
   44 #include "copyimage.hxx" 
   45 #include "threading.hxx" 
   61 template <
unsigned int N, 
class T, 
class C>
 
   62 void moveDCToCenterImpl(MultiArrayView<N, T, C> a, 
unsigned int startDimension)
 
   64     typedef typename MultiArrayView<N, T, C>::traverser Traverser;
 
   65     typedef MultiArrayNavigator<Traverser, N> Navigator;
 
   66     typedef typename Navigator::iterator Iterator;
 
   68     for(
unsigned int d = startDimension; d < N; ++d)
 
   70         Navigator nav(a.traverser_begin(), a.shape(), d);
 
   72         for( ; nav.hasMore(); nav++ )
 
   74             Iterator i = nav.begin();
 
   75             int s  = nav.end() - i;
 
   80                 for(
int k=0; k<s2; ++k)
 
   82                     std::swap(i[k], i[k+s2]);
 
   88                 for(
int k=0; k<s2; ++k)
 
   99 template <
unsigned int N, 
class T, 
class C>
 
  100 void moveDCToUpperLeftImpl(MultiArrayView<N, T, C> a, 
unsigned int startDimension)
 
  102     typedef typename MultiArrayView<N, T, C>::traverser Traverser;
 
  103     typedef MultiArrayNavigator<Traverser, N> Navigator;
 
  104     typedef typename Navigator::iterator Iterator;
 
  106     for(
unsigned int d = startDimension; d < N; ++d)
 
  108         Navigator nav(a.traverser_begin(), a.shape(), d);
 
  110         for( ; nav.hasMore(); nav++ )
 
  112             Iterator i = nav.begin();
 
  113             int s  = nav.end() - i;
 
  118                 for(
int k=0; k<s2; ++k)
 
  120                     std::swap(i[k], i[k+s2]);
 
  126                 for(
int k=s2; k>0; --k)
 
  145 template <
unsigned int N, 
class T, 
class C>
 
  148     detail::moveDCToCenterImpl(a, 0);
 
  151 template <
unsigned int N, 
class T, 
class C>
 
  154     detail::moveDCToUpperLeftImpl(a, 0);
 
  157 template <
unsigned int N, 
class T, 
class C>
 
  158 inline void moveDCToHalfspaceCenter(MultiArrayView<N, T, C> a)
 
  160     detail::moveDCToCenterImpl(a, 1);
 
  163 template <
unsigned int N, 
class T, 
class C>
 
  164 inline void moveDCToHalfspaceUpperLeft(MultiArrayView<N, T, C> a)
 
  166     detail::moveDCToUpperLeftImpl(a, 1);
 
  172 #ifndef VIGRA_SINGLE_THREADED 
  174 template <
int DUMMY=0>
 
  178     threading::lock_guard<threading::mutex> guard_;
 
  181     : guard_(plan_mutex_)
 
  184     static threading::mutex plan_mutex_;
 
  188 threading::mutex FFTWLock<DUMMY>::plan_mutex_;
 
  190 #else // VIGRA_SINGLE_THREADED 
  192 template <
int DUMMY=0>
 
  201 #endif // not VIGRA_SINGLE_THREADED 
  204 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  205                FFTWComplex<double> * in,  
int* instrides,  
int instep,
 
  206                FFTWComplex<double> * out, 
int* outstrides, 
int outstep,
 
  207                int sign, 
unsigned int planner_flags)
 
  209     return fftw_plan_many_dft(N, shape, 1,
 
  210                               (fftw_complex *)in, instrides, instep, 0,
 
  211                               (fftw_complex *)out, outstrides, outstep, 0,
 
  212                               sign, planner_flags);
 
  216 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  217                double * in,  
int* instrides,  
int instep,
 
  218                FFTWComplex<double> * out, 
int* outstrides, 
int outstep,
 
  219                int , 
unsigned int planner_flags)
 
  221     return fftw_plan_many_dft_r2c(N, shape, 1,
 
  222                                    in, instrides, instep, 0,
 
  223                                    (fftw_complex *)out, outstrides, outstep, 0,
 
  228 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  229                FFTWComplex<double> * in,  
int* instrides,  
int instep,
 
  230                double * out, 
int* outstrides, 
int outstep,
 
  231                int , 
unsigned int planner_flags)
 
  233     return fftw_plan_many_dft_c2r(N, shape, 1,
 
  234                                   (fftw_complex *)in, instrides, instep, 0,
 
  235                                   out, outstrides, outstep, 0,
 
  240 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  241                FFTWComplex<float> * in,  
int* instrides,  
int instep,
 
  242                FFTWComplex<float> * out, 
int* outstrides, 
int outstep,
 
  243                int sign, 
unsigned int planner_flags)
 
  245     return fftwf_plan_many_dft(N, shape, 1,
 
  246                                (fftwf_complex *)in, instrides, instep, 0,
 
  247                                (fftwf_complex *)out, outstrides, outstep, 0,
 
  248                                sign, planner_flags);
 
  252 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  253                float * in,  
int* instrides,  
int instep,
 
  254                FFTWComplex<float> * out, 
int* outstrides, 
int outstep,
 
  255                int , 
unsigned int planner_flags)
 
  257     return fftwf_plan_many_dft_r2c(N, shape, 1,
 
  258                                     in, instrides, instep, 0,
 
  259                                     (fftwf_complex *)out, outstrides, outstep, 0,
 
  264 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  265                FFTWComplex<float> * in,  
int* instrides,  
int instep,
 
  266                float * out, 
int* outstrides, 
int outstep,
 
  267                int , 
unsigned int planner_flags)
 
  269     return fftwf_plan_many_dft_c2r(N, shape, 1,
 
  270                                    (fftwf_complex *)in, instrides, instep, 0,
 
  271                                    out, outstrides, outstep, 0,
 
  276 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  277                FFTWComplex<long double> * in,  
int* instrides,  
int instep,
 
  278                FFTWComplex<long double> * out, 
int* outstrides, 
int outstep,
 
  279                int sign, 
unsigned int planner_flags)
 
  281     return fftwl_plan_many_dft(N, shape, 1,
 
  282                                (fftwl_complex *)in, instrides, instep, 0,
 
  283                                (fftwl_complex *)out, outstrides, outstep, 0,
 
  284                                sign, planner_flags);
 
  288 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  289                long double * in,  
int* instrides,  
int instep,
 
  290                FFTWComplex<long double> * out, 
int* outstrides, 
int outstep,
 
  291                int , 
unsigned int planner_flags)
 
  293     return fftwl_plan_many_dft_r2c(N, shape, 1,
 
  294                                     in, instrides, instep, 0,
 
  295                                     (fftwl_complex *)out, outstrides, outstep, 0,
 
  300 fftwPlanCreate(
unsigned int N, 
int* shape,
 
  301                FFTWComplex<long double> * in,  
int* instrides,  
int instep,
 
  302                long double * out, 
int* outstrides, 
int outstep,
 
  303                int , 
unsigned int planner_flags)
 
  305     return fftwl_plan_many_dft_c2r(N, shape, 1,
 
  306                                    (fftwl_complex *)in, instrides, instep, 0,
 
  307                                    out, outstrides, outstep, 0,
 
  311 inline void fftwPlanDestroy(fftw_plan plan)
 
  314         fftw_destroy_plan(plan);
 
  317 inline void fftwPlanDestroy(fftwf_plan plan)
 
  320         fftwf_destroy_plan(plan);
 
  323 inline void fftwPlanDestroy(fftwl_plan plan)
 
  326         fftwl_destroy_plan(plan);
 
  330 fftwPlanExecute(fftw_plan plan)
 
  336 fftwPlanExecute(fftwf_plan plan)
 
  342 fftwPlanExecute(fftwl_plan plan)
 
  348 fftwPlanExecute(fftw_plan plan, FFTWComplex<double> * in,  FFTWComplex<double> * out)
 
  350     fftw_execute_dft(plan, (fftw_complex *)in, (fftw_complex *)out);
 
  354 fftwPlanExecute(fftw_plan plan, 
double * in,  FFTWComplex<double> * out)
 
  356     fftw_execute_dft_r2c(plan, in, (fftw_complex *)out);
 
  360 fftwPlanExecute(fftw_plan plan, FFTWComplex<double> * in,  
double * out)
 
  362     fftw_execute_dft_c2r(plan, (fftw_complex *)in, out);
 
  366 fftwPlanExecute(fftwf_plan plan, FFTWComplex<float> * in,  FFTWComplex<float> * out)
 
  368     fftwf_execute_dft(plan, (fftwf_complex *)in, (fftwf_complex *)out);
 
  372 fftwPlanExecute(fftwf_plan plan, 
float * in,  FFTWComplex<float> * out)
 
  374     fftwf_execute_dft_r2c(plan, in, (fftwf_complex *)out);
 
  378 fftwPlanExecute(fftwf_plan plan, FFTWComplex<float> * in,  
float * out)
 
  380     fftwf_execute_dft_c2r(plan, (fftwf_complex *)in, out);
 
  384 fftwPlanExecute(fftwl_plan plan, FFTWComplex<long double> * in,  FFTWComplex<long double> * out)
 
  386     fftwl_execute_dft(plan, (fftwl_complex *)in, (fftwl_complex *)out);
 
  390 fftwPlanExecute(fftwl_plan plan, 
long double * in,  FFTWComplex<long double> * out)
 
  392     fftwl_execute_dft_r2c(plan, in, (fftwl_complex *)out);
 
  396 fftwPlanExecute(fftwl_plan plan, FFTWComplex<long double> * in,  
long double * out)
 
  398     fftwl_execute_dft_c2r(plan, (fftwl_complex *)in, out);
 
  402 struct FFTWPaddingSize
 
  404     static const int size = 1330;
 
  405     static const int evenSize = 1063;
 
  406     static int goodSizes[size];
 
  407     static int goodEvenSizes[evenSize];
 
  409     static inline int find(
int s)
 
  411         if(s <= 0 || s >= goodSizes[size-1])
 
  414         int * upperBound = std::upper_bound(goodSizes, goodSizes+size, s);
 
  415         if(upperBound > goodSizes && upperBound[-1] == s)
 
  421     static inline int findEven(
int s)
 
  423         if(s <= 0 || s >= goodEvenSizes[evenSize-1])
 
  426         int * upperBound = std::upper_bound(goodEvenSizes, goodEvenSizes+evenSize, s);
 
  427         if(upperBound > goodEvenSizes && upperBound[-1] == s)
 
  439 int FFTWPaddingSize<DUMMY>::goodSizes[size] = {
 
  440         1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16,
 
  441         18, 20, 21, 22, 24, 25, 26, 27, 28, 30, 32, 33, 35, 36, 39, 40, 42, 44, 45, 48,
 
  442         49, 50, 52, 54, 55, 56, 60, 63, 64, 65, 66, 70, 72, 75, 77, 78, 80, 81,
 
  443         84, 88, 90, 91, 96, 98, 99, 100, 104, 105, 108, 110, 112, 117, 120, 125,
 
  444         126, 128, 130, 132, 135, 140, 144, 147, 150, 154, 156, 160, 162, 165,
 
  445         168, 175, 176, 180, 182, 189, 192, 195, 196, 198, 200, 208, 210, 216,
 
  446         220, 224, 225, 231, 234, 240, 243, 245, 250, 252, 256, 260, 264, 270,
 
  447         273, 275, 280, 288, 294, 297, 300, 308, 312, 315, 320, 324, 325, 330,
 
  448         336, 343, 350, 351, 352, 360, 364, 375, 378, 384, 385, 390, 392, 396,
 
  449         400, 405, 416, 420, 432, 440, 441, 448, 450, 455, 462, 468, 480, 486,
 
  450         490, 495, 500, 504, 512, 520, 525, 528, 539, 540, 546, 550, 560, 567,
 
  451         576, 585, 588, 594, 600, 616, 624, 625, 630, 637, 640, 648, 650, 660,
 
  452         672, 675, 686, 693, 700, 702, 704, 720, 728, 729, 735, 750, 756, 768,
 
  453         770, 780, 784, 792, 800, 810, 819, 825, 832, 840, 864, 875, 880, 882,
 
  454         891, 896, 900, 910, 924, 936, 945, 960, 972, 975, 980, 990, 1000, 1008,
 
  455         1024, 1029, 1040, 1050, 1053, 1056, 1078, 1080, 1092, 1100, 1120, 1125,
 
  456         1134, 1152, 1155, 1170, 1176, 1188, 1200, 1215, 1225, 1232, 1248, 1250,
 
  457         1260, 1274, 1280, 1296, 1300, 1320, 1323, 1344, 1350, 1365, 1372, 1375,
 
  458         1386, 1400, 1404, 1408, 1440, 1456, 1458, 1470, 1485, 1500, 1512, 1536,
 
  459         1540, 1560, 1568, 1575, 1584, 1600, 1617, 1620, 1625, 1638, 1650, 1664,
 
  460         1680, 1701, 1715, 1728, 1750, 1755, 1760, 1764, 1782, 1792, 1800, 1820,
 
  461         1848, 1872, 1875, 1890, 1911, 1920, 1925, 1944, 1950, 1960, 1980, 2000,
 
  462         2016, 2025, 2048, 2058, 2079, 2080, 2100, 2106, 2112, 2156, 2160, 2184,
 
  463         2187, 2200, 2205, 2240, 2250, 2268, 2275, 2304, 2310, 2340, 2352, 2376,
 
  464         2400, 2401, 2430, 2450, 2457, 2464, 2475, 2496, 2500, 2520, 2548, 2560,
 
  465         2592, 2600, 2625, 2640, 2646, 2673, 2688, 2695, 2700, 2730, 2744, 2750,
 
  466         2772, 2800, 2808, 2816, 2835, 2880, 2912, 2916, 2925, 2940, 2970, 3000,
 
  467         3024, 3072, 3080, 3087, 3120, 3125, 3136, 3150, 3159, 3168, 3185, 3200,
 
  468         3234, 3240, 3250, 3276, 3300, 3328, 3360, 3375, 3402, 3430, 3456, 3465,
 
  469         3500, 3510, 3520, 3528, 3564, 3584, 3600, 3640, 3645, 3675, 3696, 3744,
 
  470         3750, 3773, 3780, 3822, 3840, 3850, 3888, 3900, 3920, 3960, 3969, 4000,
 
  471         4032, 4050, 4095, 4096, 4116, 4125, 4158, 4160, 4200, 4212, 4224, 4312,
 
  472         4320, 4368, 4374, 4375, 4400, 4410, 4455, 4459, 4480, 4500, 4536, 4550,
 
  473         4608, 4620, 4680, 4704, 4725, 4752, 4800, 4802, 4851, 4860, 4875, 4900,
 
  474         4914, 4928, 4950, 4992, 5000, 5040, 5096, 5103, 5120, 5145, 5184, 5200,
 
  475         5250, 5265, 5280, 5292, 5346, 5376, 5390, 5400, 5460, 5488, 5500, 5544,
 
  476         5600, 5616, 5625, 5632, 5670, 5733, 5760, 5775, 5824, 5832, 5850, 5880,
 
  477         5940, 6000, 6048, 6075, 6125, 6144, 6160, 6174, 6237, 6240, 6250, 6272,
 
  478         6300, 6318, 6336, 6370, 6400, 6468, 6480, 6500, 6552, 6561, 6600, 6615,
 
  479         6656, 6720, 6750, 6804, 6825, 6860, 6875, 6912, 6930, 7000, 7020, 7040,
 
  480         7056, 7128, 7168, 7200, 7203, 7280, 7290, 7350, 7371, 7392, 7425, 7488,
 
  481         7500, 7546, 7560, 7644, 7680, 7700, 7776, 7800, 7840, 7875, 7920, 7938,
 
  482         8000, 8019, 8064, 8085, 8100, 8125, 8190, 8192, 8232, 8250, 8316, 8320,
 
  483         8400, 8424, 8448, 8505, 8575, 8624, 8640, 8736, 8748, 8750, 8775, 8800,
 
  484         8820, 8910, 8918, 8960, 9000, 9072, 9100, 9216, 9240, 9261, 9360, 9375,
 
  485         9408, 9450, 9477, 9504, 9555, 9600, 9604, 9625, 9702, 9720, 9750, 9800,
 
  486         9828, 9856, 9900, 9984, 10000, 10080, 10125, 10192, 10206, 10240, 10290,
 
  487         10368, 10395, 10400, 10500, 10530, 10560, 10584, 10692, 10752, 10780,
 
  488         10800, 10920, 10935, 10976, 11000, 11025, 11088, 11200, 11232, 11250,
 
  489         11264, 11319, 11340, 11375, 11466, 11520, 11550, 11648, 11664, 11700,
 
  490         11760, 11880, 11907, 12000, 12005, 12096, 12150, 12250, 12285, 12288,
 
  491         12320, 12348, 12375, 12474, 12480, 12500, 12544, 12600, 12636, 12672,
 
  492         12740, 12800, 12936, 12960, 13000, 13104, 13122, 13125, 13200, 13230,
 
  493         13312, 13365, 13377, 13440, 13475, 13500, 13608, 13650, 13720, 13750,
 
  494         13824, 13860, 14000, 14040, 14080, 14112, 14175, 14256, 14336, 14400,
 
  495         14406, 14553, 14560, 14580, 14625, 14700, 14742, 14784, 14850, 14976,
 
  496         15000, 15092, 15120, 15288, 15309, 15360, 15400, 15435, 15552, 15600,
 
  497         15625, 15680, 15750, 15795, 15840, 15876, 15925, 16000, 16038, 16128,
 
  498         16170, 16200, 16250, 16380, 16384, 16464, 16500, 16632, 16640, 16800,
 
  499         16807, 16848, 16875, 16896, 17010, 17150, 17199, 17248, 17280, 17325,
 
  500         17472, 17496, 17500, 17550, 17600, 17640, 17820, 17836, 17920, 18000,
 
  501         18144, 18200, 18225, 18375, 18432, 18480, 18522, 18711, 18720, 18750,
 
  502         18816, 18865, 18900, 18954, 19008, 19110, 19200, 19208, 19250, 19404,
 
  503         19440, 19500, 19600, 19656, 19683, 19712, 19800, 19845, 19968, 20000,
 
  504         20160, 20250, 20384, 20412, 20475, 20480, 20580, 20625, 20736, 20790,
 
  505         20800, 21000, 21060, 21120, 21168, 21384, 21504, 21560, 21600, 21609,
 
  506         21840, 21870, 21875, 21952, 22000, 22050, 22113, 22176, 22275, 22295,
 
  507         22400, 22464, 22500, 22528, 22638, 22680, 22750, 22932, 23040, 23100,
 
  508         23296, 23328, 23400, 23520, 23625, 23760, 23814, 24000, 24010, 24057,
 
  509         24192, 24255, 24300, 24375, 24500, 24570, 24576, 24640, 24696, 24750,
 
  510         24948, 24960, 25000, 25088, 25200, 25272, 25344, 25480, 25515, 25600,
 
  511         25725, 25872, 25920, 26000, 26208, 26244, 26250, 26325, 26400, 26411,
 
  512         26460, 26624, 26730, 26754, 26880, 26950, 27000, 27216, 27300, 27440,
 
  513         27500, 27648, 27720, 27783, 28000, 28080, 28125, 28160, 28224, 28350,
 
  514         28431, 28512, 28665, 28672, 28800, 28812, 28875, 29106, 29120, 29160,
 
  515         29250, 29400, 29484, 29568, 29700, 29952, 30000, 30184, 30240, 30375,
 
  516         30576, 30618, 30625, 30720, 30800, 30870, 31104, 31185, 31200, 31213,
 
  517         31250, 31360, 31500, 31590, 31680, 31752, 31850, 32000, 32076, 32256,
 
  518         32340, 32400, 32500, 32760, 32768, 32805, 32928, 33000, 33075, 33264,
 
  519         33280, 33600, 33614, 33696, 33750, 33792, 33957, 34020, 34125, 34300,
 
  520         34375, 34398, 34496, 34560, 34650, 34944, 34992, 35000, 35100, 35200,
 
  521         35280, 35640, 35672, 35721, 35840, 36000, 36015, 36288, 36400, 36450,
 
  522         36750, 36855, 36864, 36960, 37044, 37125, 37422, 37440, 37500, 37632,
 
  523         37730, 37800, 37908, 38016, 38220, 38400, 38416, 38500, 38808, 38880,
 
  524         39000, 39200, 39312, 39366, 39375, 39424, 39600, 39690, 39936, 40000,
 
  525         40095, 40131, 40320, 40425, 40500, 40625, 40768, 40824, 40950, 40960,
 
  526         41160, 41250, 41472, 41580, 41600, 42000, 42120, 42240, 42336, 42525,
 
  527         42768, 42875, 43008, 43120, 43200, 43218, 43659, 43680, 43740, 43750,
 
  528         43875, 43904, 44000, 44100, 44226, 44352, 44550, 44590, 44800, 44928,
 
  529         45000, 45056, 45276, 45360, 45500, 45864, 45927, 46080, 46200, 46305,
 
  530         46592, 46656, 46800, 46875, 47040, 47250, 47385, 47520, 47628, 47775,
 
  531         48000, 48020, 48114, 48125, 48384, 48510, 48600, 48750, 49000, 49140,
 
  532         49152, 49280, 49392, 49500, 49896, 49920, 50000, 50176, 50400, 50421,
 
  533         50544, 50625, 50688, 50960, 51030, 51200, 51450, 51597, 51744, 51840,
 
  534         51975, 52000, 52416, 52488, 52500, 52650, 52800, 52822, 52920, 53248,
 
  535         53460, 53508, 53760, 53900, 54000, 54432, 54600, 54675, 54880, 55000,
 
  536         55125, 55296, 55440, 55566, 56000, 56133, 56160, 56250, 56320, 56448,
 
  537         56595, 56700, 56862, 56875, 57024, 57330, 57344, 57600, 57624, 57750,
 
  538         58212, 58240, 58320, 58500, 58800, 58968, 59049, 59136, 59400, 59535,
 
  539         59904, 60000, 60025, 60368, 60480, 60750, 61152, 61236, 61250, 61425,
 
  540         61440, 61600, 61740, 61875, 62208, 62370, 62400, 62426, 62500, 62720,
 
  541         63000, 63180, 63360, 63504, 63700, 64000, 64152, 64512, 64680, 64800,
 
  542         64827, 65000, 65520, 65536, 65610, 65625, 65856, 66000, 66150, 66339,
 
  543         66528, 66560, 66825, 66885, 67200, 67228, 67375, 67392, 67500, 67584,
 
  544         67914, 68040, 68250, 68600, 68750, 68796, 68992, 69120, 69300, 69888,
 
  545         69984, 70000, 70200, 70400, 70560, 70875, 71280, 71344, 71442, 71680,
 
  546         72000, 72030, 72171, 72576, 72765, 72800, 72900, 73125, 73500, 73710,
 
  547         73728, 73920, 74088, 74250, 74844, 74880, 75000, 75264, 75460, 75600,
 
  548         75816, 76032, 76440, 76545, 76800, 76832, 77000, 77175, 77616, 77760,
 
  549         78000, 78125, 78400, 78624, 78732, 78750, 78848, 78975, 79200, 79233,
 
  550         79380, 79625, 79872, 80000, 80190, 80262, 80640, 80850, 81000, 81250,
 
  551         81536, 81648, 81900, 81920, 82320, 82500, 82944, 83160, 83200, 83349,
 
  552         84000, 84035, 84240, 84375, 84480, 84672, 85050, 85293, 85536, 85750,
 
  553         85995, 86016, 86240, 86400, 86436, 86625, 87318, 87360, 87480, 87500,
 
  554         87750, 87808, 88000, 88200, 88452, 88704, 89100, 89180, 89600, 89856,
 
  555         90000, 90112, 90552, 90720, 91000, 91125, 91728, 91854, 91875, 92160,
 
  556         92400, 92610, 93184, 93312, 93555, 93600, 93639, 93750, 94080, 94325,
 
  557         94500, 94770, 95040, 95256, 95550, 96000, 96040, 96228, 96250, 96768,
 
  558         97020, 97200, 97500, 98000, 98280, 98304, 98415, 98560, 98784, 99000,
 
  563 int FFTWPaddingSize<DUMMY>::goodEvenSizes[evenSize] = {
 
  564         2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22,
 
  565         24, 26, 28, 30, 32, 36, 40, 42, 44, 48, 50, 52, 54, 56, 60, 64, 66, 70,
 
  566         72, 78, 80, 84, 88, 90, 96, 98, 100, 104, 108, 110, 112, 120, 126, 128,
 
  567         130, 132, 140, 144, 150, 154, 156, 160, 162, 168, 176, 180, 182, 192,
 
  568         196, 198, 200, 208, 210, 216, 220, 224, 234, 240, 250, 252, 256, 260,
 
  569         264, 270, 280, 288, 294, 300, 308, 312, 320, 324, 330, 336, 350, 352,
 
  570         360, 364, 378, 384, 390, 392, 396, 400, 416, 420, 432, 440, 448, 450,
 
  571         462, 468, 480, 486, 490, 500, 504, 512, 520, 528, 540, 546, 550, 560,
 
  572         576, 588, 594, 600, 616, 624, 630, 640, 648, 650, 660, 672, 686, 700,
 
  573         702, 704, 720, 728, 750, 756, 768, 770, 780, 784, 792, 800, 810, 832,
 
  574         840, 864, 880, 882, 896, 900, 910, 924, 936, 960, 972, 980, 990, 1000,
 
  575         1008, 1024, 1040, 1050, 1056, 1078, 1080, 1092, 1100, 1120, 1134, 1152,
 
  576         1170, 1176, 1188, 1200, 1232, 1248, 1250, 1260, 1274, 1280, 1296, 1300,
 
  577         1320, 1344, 1350, 1372, 1386, 1400, 1404, 1408, 1440, 1456, 1458, 1470,
 
  578         1500, 1512, 1536, 1540, 1560, 1568, 1584, 1600, 1620, 1638, 1650, 1664,
 
  579         1680, 1728, 1750, 1760, 1764, 1782, 1792, 1800, 1820, 1848, 1872, 1890,
 
  580         1920, 1944, 1950, 1960, 1980, 2000, 2016, 2048, 2058, 2080, 2100, 2106,
 
  581         2112, 2156, 2160, 2184, 2200, 2240, 2250, 2268, 2304, 2310, 2340, 2352,
 
  582         2376, 2400, 2430, 2450, 2464, 2496, 2500, 2520, 2548, 2560, 2592, 2600,
 
  583         2640, 2646, 2688, 2700, 2730, 2744, 2750, 2772, 2800, 2808, 2816, 2880,
 
  584         2912, 2916, 2940, 2970, 3000, 3024, 3072, 3080, 3120, 3136, 3150, 3168,
 
  585         3200, 3234, 3240, 3250, 3276, 3300, 3328, 3360, 3402, 3430, 3456, 3500,
 
  586         3510, 3520, 3528, 3564, 3584, 3600, 3640, 3696, 3744, 3750, 3780, 3822,
 
  587         3840, 3850, 3888, 3900, 3920, 3960, 4000, 4032, 4050, 4096, 4116, 4158,
 
  588         4160, 4200, 4212, 4224, 4312, 4320, 4368, 4374, 4400, 4410, 4480, 4500,
 
  589         4536, 4550, 4608, 4620, 4680, 4704, 4752, 4800, 4802, 4860, 4900, 4914,
 
  590         4928, 4950, 4992, 5000, 5040, 5096, 5120, 5184, 5200, 5250, 5280, 5292,
 
  591         5346, 5376, 5390, 5400, 5460, 5488, 5500, 5544, 5600, 5616, 5632, 5670,
 
  592         5760, 5824, 5832, 5850, 5880, 5940, 6000, 6048, 6144, 6160, 6174, 6240,
 
  593         6250, 6272, 6300, 6318, 6336, 6370, 6400, 6468, 6480, 6500, 6552, 6600,
 
  594         6656, 6720, 6750, 6804, 6860, 6912, 6930, 7000, 7020, 7040, 7056, 7128,
 
  595         7168, 7200, 7280, 7290, 7350, 7392, 7488, 7500, 7546, 7560, 7644, 7680,
 
  596         7700, 7776, 7800, 7840, 7920, 7938, 8000, 8064, 8100, 8190, 8192, 8232,
 
  597         8250, 8316, 8320, 8400, 8424, 8448, 8624, 8640, 8736, 8748, 8750, 8800,
 
  598         8820, 8910, 8918, 8960, 9000, 9072, 9100, 9216, 9240, 9360, 9408, 9450,
 
  599         9504, 9600, 9604, 9702, 9720, 9750, 9800, 9828, 9856, 9900, 9984, 10000,
 
  600         10080, 10192, 10206, 10240, 10290, 10368, 10400, 10500, 10530, 10560,
 
  601         10584, 10692, 10752, 10780, 10800, 10920, 10976, 11000, 11088, 11200,
 
  602         11232, 11250, 11264, 11340, 11466, 11520, 11550, 11648, 11664, 11700,
 
  603         11760, 11880, 12000, 12096, 12150, 12250, 12288, 12320, 12348, 12474,
 
  604         12480, 12500, 12544, 12600, 12636, 12672, 12740, 12800, 12936, 12960,
 
  605         13000, 13104, 13122, 13200, 13230, 13312, 13440, 13500, 13608, 13650,
 
  606         13720, 13750, 13824, 13860, 14000, 14040, 14080, 14112, 14256, 14336,
 
  607         14400, 14406, 14560, 14580, 14700, 14742, 14784, 14850, 14976, 15000,
 
  608         15092, 15120, 15288, 15360, 15400, 15552, 15600, 15680, 15750, 15840,
 
  609         15876, 16000, 16038, 16128, 16170, 16200, 16250, 16380, 16384, 16464,
 
  610         16500, 16632, 16640, 16800, 16848, 16896, 17010, 17150, 17248, 17280,
 
  611         17472, 17496, 17500, 17550, 17600, 17640, 17820, 17836, 17920, 18000,
 
  612         18144, 18200, 18432, 18480, 18522, 18720, 18750, 18816, 18900, 18954,
 
  613         19008, 19110, 19200, 19208, 19250, 19404, 19440, 19500, 19600, 19656,
 
  614         19712, 19800, 19968, 20000, 20160, 20250, 20384, 20412, 20480, 20580,
 
  615         20736, 20790, 20800, 21000, 21060, 21120, 21168, 21384, 21504, 21560,
 
  616         21600, 21840, 21870, 21952, 22000, 22050, 22176, 22400, 22464, 22500,
 
  617         22528, 22638, 22680, 22750, 22932, 23040, 23100, 23296, 23328, 23400,
 
  618         23520, 23760, 23814, 24000, 24010, 24192, 24300, 24500, 24570, 24576,
 
  619         24640, 24696, 24750, 24948, 24960, 25000, 25088, 25200, 25272, 25344,
 
  620         25480, 25600, 25872, 25920, 26000, 26208, 26244, 26250, 26400, 26460,
 
  621         26624, 26730, 26754, 26880, 26950, 27000, 27216, 27300, 27440, 27500,
 
  622         27648, 27720, 28000, 28080, 28160, 28224, 28350, 28512, 28672, 28800,
 
  623         28812, 29106, 29120, 29160, 29250, 29400, 29484, 29568, 29700, 29952,
 
  624         30000, 30184, 30240, 30576, 30618, 30720, 30800, 30870, 31104, 31200,
 
  625         31250, 31360, 31500, 31590, 31680, 31752, 31850, 32000, 32076, 32256,
 
  626         32340, 32400, 32500, 32760, 32768, 32928, 33000, 33264, 33280, 33600,
 
  627         33614, 33696, 33750, 33792, 34020, 34300, 34398, 34496, 34560, 34650,
 
  628         34944, 34992, 35000, 35100, 35200, 35280, 35640, 35672, 35840, 36000,
 
  629         36288, 36400, 36450, 36750, 36864, 36960, 37044, 37422, 37440, 37500,
 
  630         37632, 37730, 37800, 37908, 38016, 38220, 38400, 38416, 38500, 38808,
 
  631         38880, 39000, 39200, 39312, 39366, 39424, 39600, 39690, 39936, 40000,
 
  632         40320, 40500, 40768, 40824, 40950, 40960, 41160, 41250, 41472, 41580,
 
  633         41600, 42000, 42120, 42240, 42336, 42768, 43008, 43120, 43200, 43218,
 
  634         43680, 43740, 43750, 43904, 44000, 44100, 44226, 44352, 44550, 44590,
 
  635         44800, 44928, 45000, 45056, 45276, 45360, 45500, 45864, 46080, 46200,
 
  636         46592, 46656, 46800, 47040, 47250, 47520, 47628, 48000, 48020, 48114,
 
  637         48384, 48510, 48600, 48750, 49000, 49140, 49152, 49280, 49392, 49500,
 
  638         49896, 49920, 50000, 50176, 50400, 50544, 50688, 50960, 51030, 51200,
 
  639         51450, 51744, 51840, 52000, 52416, 52488, 52500, 52650, 52800, 52822,
 
  640         52920, 53248, 53460, 53508, 53760, 53900, 54000, 54432, 54600, 54880,
 
  641         55000, 55296, 55440, 55566, 56000, 56160, 56250, 56320, 56448, 56700,
 
  642         56862, 57024, 57330, 57344, 57600, 57624, 57750, 58212, 58240, 58320,
 
  643         58500, 58800, 58968, 59136, 59400, 59904, 60000, 60368, 60480, 60750,
 
  644         61152, 61236, 61250, 61440, 61600, 61740, 62208, 62370, 62400, 62426,
 
  645         62500, 62720, 63000, 63180, 63360, 63504, 63700, 64000, 64152, 64512,
 
  646         64680, 64800, 65000, 65520, 65536, 65610, 65856, 66000, 66150, 66528,
 
  647         66560, 67200, 67228, 67392, 67500, 67584, 67914, 68040, 68250, 68600,
 
  648         68750, 68796, 68992, 69120, 69300, 69888, 69984, 70000, 70200, 70400,
 
  649         70560, 71280, 71344, 71442, 71680, 72000, 72030, 72576, 72800, 72900,
 
  650         73500, 73710, 73728, 73920, 74088, 74250, 74844, 74880, 75000, 75264,
 
  651         75460, 75600, 75816, 76032, 76440, 76800, 76832, 77000, 77616, 77760,
 
  652         78000, 78400, 78624, 78732, 78750, 78848, 79200, 79380, 79872, 80000,
 
  653         80190, 80262, 80640, 80850, 81000, 81250, 81536, 81648, 81900, 81920,
 
  654         82320, 82500, 82944, 83160, 83200, 84000, 84240, 84480, 84672, 85050,
 
  655         85536, 85750, 86016, 86240, 86400, 86436, 87318, 87360, 87480, 87500,
 
  656         87750, 87808, 88000, 88200, 88452, 88704, 89100, 89180, 89600, 89856,
 
  657         90000, 90112, 90552, 90720, 91000, 91728, 91854, 92160, 92400, 92610,
 
  658         93184, 93312, 93600, 93750, 94080, 94500, 94770, 95040, 95256, 95550,
 
  659         96000, 96040, 96228, 96250, 96768, 97020, 97200, 97500, 98000, 98280,
 
  660         98304, 98560, 98784, 99000, 99792, 99840
 
  664 struct FFTEmbedKernel
 
  666     template <
unsigned int N, 
class Real, 
class C, 
class Shape>
 
  668     exec(MultiArrayView<N, Real, C> & out, Shape 
const & kernelShape,
 
  669          Shape & srcPoint, Shape & destPoint, 
bool copyIt)
 
  671         for(srcPoint[M]=0; srcPoint[M]<kernelShape[M]; ++srcPoint[M])
 
  673             if(srcPoint[M] < (kernelShape[M] + 1) / 2)
 
  675                 destPoint[M] = srcPoint[M];
 
  679                 destPoint[M] = srcPoint[M] + out.shape(M) - kernelShape[M];
 
  682             FFTEmbedKernel<M-1>::exec(out, kernelShape, srcPoint, destPoint, copyIt);
 
  688 struct FFTEmbedKernel<0>
 
  690     template <
unsigned int N, 
class Real, 
class C, 
class Shape>
 
  692     exec(MultiArrayView<N, Real, C> & out, Shape 
const & kernelShape,
 
  693          Shape & srcPoint, Shape & destPoint, 
bool copyIt)
 
  695         for(srcPoint[0]=0; srcPoint[0]<kernelShape[0]; ++srcPoint[0])
 
  697             if(srcPoint[0] < (kernelShape[0] + 1) / 2)
 
  699                 destPoint[0] = srcPoint[0];
 
  703                 destPoint[0] = srcPoint[0] + out.shape(0) - kernelShape[0];
 
  708                 out[destPoint] = out[srcPoint];
 
  715 template <
unsigned int N, 
class Real, 
class C1, 
class C2>
 
  717 fftEmbedKernel(MultiArrayView<N, Real, C1> kernel,
 
  718                MultiArrayView<N, Real, C2> out,
 
  721     typedef typename MultiArrayShape<N>::type Shape;
 
  723     MultiArrayView<N, Real, C2> kout = out.subarray(Shape(), kernel.shape());
 
  731     Shape srcPoint, destPoint;
 
  732     FFTEmbedKernel<(int)N-1>::exec(out, kernel.shape(), srcPoint, destPoint, 
false);
 
  735 template <
unsigned int N, 
class Real, 
class C1, 
class C2>
 
  737 fftEmbedArray(MultiArrayView<N, Real, C1> in,
 
  738               MultiArrayView<N, Real, C2> out)
 
  740     typedef typename MultiArrayShape<N>::type Shape;
 
  742     Shape diff = out.shape() - in.shape(),
 
  744           rightDiff = diff - leftDiff,
 
  745           right = in.shape() + leftDiff;
 
  747     out.subarray(leftDiff, right) = in;
 
  749     typedef typename MultiArrayView<N, Real, C2>::traverser Traverser;
 
  750     typedef MultiArrayNavigator<Traverser, N> Navigator;
 
  751     typedef typename Navigator::iterator Iterator;
 
  753     for(
unsigned int d = 0; d < N; ++d)
 
  755         Navigator nav(out.traverser_begin(), out.shape(), d);
 
  757         for( ; nav.hasMore(); nav++ )
 
  759             Iterator i = nav.begin();
 
  760             for(
int k=1; k<=leftDiff[d]; ++k)
 
  761                 i[leftDiff[d] - k] = i[leftDiff[d] + k];
 
  762             for(
int k=0; k<rightDiff[d]; ++k)
 
  763                 i[right[d] + k] = i[right[d] - k - 2];
 
  770 template <
class T, 
int N>
 
  772 fftwBestPaddedShape(TinyVector<T, N> shape)
 
  774     for(
unsigned int k=0; k<N; ++k)
 
  775         shape[k] = detail::FFTWPaddingSize<0>::find(shape[k]);
 
  779 template <
class T, 
int N>
 
  781 fftwBestPaddedShapeR2C(TinyVector<T, N> shape)
 
  783     shape[0] = detail::FFTWPaddingSize<0>::findEven(shape[0]);
 
  784     for(
unsigned int k=1; k<N; ++k)
 
  785         shape[k] = detail::FFTWPaddingSize<0>::find(shape[k]);
 
  800 template <
class T, 
int N>
 
  804     shape[0] = shape[0] / 2 + 1;
 
  808 template <
class T, 
int N>
 
  810 fftwCorrespondingShapeC2R(TinyVector<T, N> shape, 
bool oddDimension0 = 
false)
 
  812     shape[0] = oddDimension0
 
  813                   ? (shape[0] - 1) * 2 + 1
 
  814                   : (shape[0] - 1) * 2;
 
  851 template <
unsigned int N, 
class Real = 
double>
 
  855     typedef typename FFTWReal2Complex<Real>::plan_type PlanType;
 
  859     Shape shape, instrides, outstrides;
 
  879     template <
class C1, 
class C2>
 
  882              int SIGN, 
unsigned int planner_flags = FFTW_ESTIMATE)
 
  885         init(in, out, SIGN, planner_flags);
 
  899     template <
class C1, 
class C2>
 
  902              unsigned int planner_flags = FFTW_ESTIMATE)
 
  905         init(in, out, planner_flags);
 
  919     template <
class C1, 
class C2>
 
  922              unsigned int planner_flags = FFTW_ESTIMATE)
 
  925         init(in, out, planner_flags);
 
  936         instrides.swap(o.instrides);
 
  937         outstrides.swap(o.outstrides);
 
  950             instrides.swap(o.instrides);
 
  951             outstrides.swap(o.outstrides);
 
  962         detail::FFTWLock<> lock;
 
  963         detail::fftwPlanDestroy(plan);
 
  970     template <
class C1, 
class C2>
 
  973               int SIGN, 
unsigned int planner_flags = FFTW_ESTIMATE)
 
  975         vigra_precondition(in.strideOrdering() == out.strideOrdering(),
 
  976             "FFTWPlan.init(): input and output must have the same stride ordering.");
 
  978         initImpl(in.permuteStridesDescending(), out.permuteStridesDescending(),
 
  979                  SIGN, planner_flags);
 
  986     template <
class C1, 
class C2>
 
  989               unsigned int planner_flags = FFTW_ESTIMATE)
 
  992             "FFTWPlan.init(): input and output must have the same stride ordering.");
 
  995                  FFTW_FORWARD, planner_flags);
 
 1002     template <
class C1, 
class C2>
 
 1005               unsigned int planner_flags = FFTW_ESTIMATE)
 
 1008             "FFTWPlan.init(): input and output must have the same stride ordering.");
 
 1011                  FFTW_BACKWARD, planner_flags);
 
 1021     template <
class C1, 
class C2>
 
 1025         executeImpl(in.permuteStridesDescending(), out.permuteStridesDescending());
 
 1035     template <
class C1, 
class C2>
 
 1049     template <
class C1, 
class C2>
 
 1058     template <
class MI, 
class MO>
 
 1059     void initImpl(MI ins, MO outs, 
int SIGN, 
unsigned int planner_flags);
 
 1061     template <
class MI, 
class MO>
 
 1062     void executeImpl(MI ins, MO outs) 
const;
 
 1067         vigra_precondition(in.shape() == out.shape(),
 
 1068             "FFTWPlan.init(): input and output must have the same shape.");
 
 1071     void checkShapes(MultiArrayView<N, Real, StridedArrayTag> ins,
 
 1072                      MultiArrayView<N, FFTWComplex<Real>, StridedArrayTag> outs)
 const 
 1074         for(
int k=0; k<(int)N-1; ++k)
 
 1075             vigra_precondition(ins.shape(k) == outs.shape(k),
 
 1076                 "FFTWPlan.init(): input and output must have matching shapes.");
 
 1077         vigra_precondition(ins.shape(N-1) / 2 + 1 == outs.shape(N-1),
 
 1078             "FFTWPlan.init(): input and output must have matching shapes.");
 
 1081     void checkShapes(MultiArrayView<N, FFTWComplex<Real>, StridedArrayTag> ins,
 
 1082                      MultiArrayView<N, Real, StridedArrayTag> outs)
 const 
 1084         for(
int k=0; k<(int)N-1; ++k)
 
 1085             vigra_precondition(ins.shape(k) == outs.shape(k),
 
 1086                 "FFTWPlan.init(): input and output must have matching shapes.");
 
 1087         vigra_precondition(outs.shape(N-1) / 2 + 1 == ins.shape(N-1),
 
 1088             "FFTWPlan.init(): input and output must have matching shapes.");
 
 1092 template <
unsigned int N, 
class Real>
 
 1093 template <
class MI, 
class MO>
 
 1095 FFTWPlan<N, Real>::initImpl(MI ins, MO outs, 
int SIGN, 
unsigned int planner_flags)
 
 1097     checkShapes(ins, outs);
 
 1099     typename MultiArrayShape<N>::type logicalShape(SIGN == FFTW_FORWARD
 
 1103     Shape newShape(logicalShape.begin(), logicalShape.end()),
 
 1104           newIStrides(ins.stride().begin(), ins.stride().end()),
 
 1105           newOStrides(outs.stride().begin(), outs.stride().end()),
 
 1106           itotal(ins.shape().begin(), ins.shape().end()),
 
 1107           ototal(outs.shape().begin(), outs.shape().end());
 
 1109     for(
unsigned int j=1; j<N; ++j)
 
 1111         itotal[j] = ins.stride(j-1) / ins.stride(j);
 
 1112         ototal[j] = outs.stride(j-1) / outs.stride(j);
 
 1116         detail::FFTWLock<> lock;
 
 1117         PlanType newPlan = detail::fftwPlanCreate(N, newShape.begin(),
 
 1118                                       ins.data(), itotal.begin(), ins.stride(N-1),
 
 1119                                       outs.data(), ototal.begin(), outs.stride(N-1),
 
 1120                                       SIGN, planner_flags);
 
 1121         detail::fftwPlanDestroy(plan);
 
 1125     shape.swap(newShape);
 
 1126     instrides.swap(newIStrides);
 
 1127     outstrides.swap(newOStrides);
 
 1131 template <
unsigned int N, 
class Real>
 
 1132 template <
class MI, 
class MO>
 
 1133 void FFTWPlan<N, Real>::executeImpl(MI ins, MO outs)
 const 
 1135     vigra_precondition(plan != 0, 
"FFTWPlan::execute(): plan is NULL.");
 
 1137     typename MultiArrayShape<N>::type lshape(sign == FFTW_FORWARD
 
 1141     vigra_precondition((lshape == TinyVectorView<int, N>(shape.data())),
 
 1142         "FFTWPlan::execute(): shape mismatch between plan and data.");
 
 1143     vigra_precondition((ins.stride() == TinyVectorView<int, N>(instrides.data())),
 
 1144         "FFTWPlan::execute(): strides mismatch between plan and input data.");
 
 1145     vigra_precondition((outs.stride() == TinyVectorView<int, N>(outstrides.data())),
 
 1146         "FFTWPlan::execute(): strides mismatch between plan and output data.");
 
 1148     detail::fftwPlanExecute(plan, ins.data(), outs.data());
 
 1150     typedef typename MO::value_type V;
 
 1151     if(sign == FFTW_BACKWARD)
 
 1152         outs *= V(1.0) / Real(outs.size());
 
 1194 template <
unsigned int N, 
class Real>
 
 1202     RArray realArray, realKernel;
 
 1203     CArray fourierArray, fourierKernel;
 
 1204     bool useFourierKernel;
 
 1215     : useFourierKernel(false)
 
 1229     template <
class C1, 
class C2, 
class C3>
 
 1233                      unsigned int planner_flags = FFTW_ESTIMATE)
 
 1234     : useFourierKernel(false)
 
 1236         init(in, kernel, out, planner_flags);
 
 1250     template <
class C1, 
class C2, 
class C3>
 
 1254                      unsigned int planner_flags = FFTW_ESTIMATE)
 
 1255     : useFourierKernel(true)
 
 1257         init(in, kernel, out, planner_flags);
 
 1272     template <
class C1, 
class C2, 
class C3>
 
 1276                      bool fourierDomainKernel,
 
 1277                      unsigned int planner_flags = FFTW_ESTIMATE)
 
 1279         init(in, kernel, out, fourierDomainKernel, planner_flags);
 
 1295     template <
class C1, 
class C2, 
class C3>
 
 1297                      bool useFourierKernel = 
false,
 
 1298                      unsigned int planner_flags = FFTW_ESTIMATE)
 
 1300         if(useFourierKernel)
 
 1301             init(inOut, kernel, planner_flags);
 
 1303             initFourierKernel(inOut, kernel, planner_flags);
 
 1310     template <
class C1, 
class C2, 
class C3>
 
 1314               unsigned int planner_flags = FFTW_ESTIMATE)
 
 1316         vigra_precondition(in.
shape() == out.
shape(),
 
 1317             "FFTWConvolvePlan::init(): input and output must have the same shape.");
 
 1325     template <
class C1, 
class C2, 
class C3>
 
 1329               unsigned int planner_flags = FFTW_ESTIMATE)
 
 1331         vigra_precondition(in.
shape() == out.
shape(),
 
 1332             "FFTWConvolvePlan::init(): input and output must have the same shape.");
 
 1333         initFourierKernel(in.
shape(), kernel.shape(), planner_flags);
 
 1340     template <
class C1, 
class C2, 
class C3>
 
 1344               bool fourierDomainKernel,
 
 1345               unsigned int planner_flags = FFTW_ESTIMATE)
 
 1347         vigra_precondition(in.shape() == out.shape(),
 
 1348             "FFTWConvolvePlan::init(): input and output must have the same shape.");
 
 1349         useFourierKernel = fourierDomainKernel;
 
 1350         initComplex(in.shape(), kernel.shape(), planner_flags);
 
 1359     template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1361                   KernelIterator kernels, KernelIterator kernelsEnd,
 
 1362                   OutIterator outs, 
unsigned int planner_flags = FFTW_ESTIMATE)
 
 1364         typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
 
 1365         typedef typename KernelArray::value_type KernelValue;
 
 1366         typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
 
 1367         typedef typename OutArray::value_type OutValue;
 
 1369         bool realKernel = IsSameType<KernelValue, Real>::value;
 
 1370         bool fourierKernel = IsSameType<KernelValue, Complex>::value;
 
 1372         vigra_precondition(realKernel || fourierKernel,
 
 1373              "FFTWConvolvePlan::initMany(): kernels have unsuitable value_type.");
 
 1374         vigra_precondition((IsSameType<OutValue, Real>::value),
 
 1375              "FFTWConvolvePlan::initMany(): outputs have unsuitable value_type.");
 
 1384             initFourierKernelMany(in.
shape(),
 
 1385                                   checkShapesFourier(in.
shape(), kernels, kernelsEnd, outs),
 
 1396     template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1398                   KernelIterator kernels, KernelIterator kernelsEnd,
 
 1400                   bool fourierDomainKernels,
 
 1401                   unsigned int planner_flags = FFTW_ESTIMATE)
 
 1403         typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
 
 1404         typedef typename KernelArray::value_type KernelValue;
 
 1405         typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
 
 1406         typedef typename OutArray::value_type OutValue;
 
 1408         vigra_precondition((IsSameType<KernelValue, Complex>::value),
 
 1409              "FFTWConvolvePlan::initMany(): kernels have unsuitable value_type.");
 
 1410         vigra_precondition((IsSameType<OutValue, Complex>::value),
 
 1411              "FFTWConvolvePlan::initMany(): outputs have unsuitable value_type.");
 
 1413         useFourierKernel = fourierDomainKernels;
 
 1415         Shape paddedShape = checkShapesComplex(in.shape(), kernels, kernelsEnd, outs);
 
 1417         CArray newFourierArray(paddedShape), newFourierKernel(paddedShape);
 
 1419         FFTWPlan<N, Real> fplan(newFourierArray, newFourierArray, FFTW_FORWARD, planner_flags);
 
 1420         FFTWPlan<N, Real> bplan(newFourierArray, newFourierArray, FFTW_BACKWARD, planner_flags);
 
 1422         forward_plan = fplan;
 
 1423         backward_plan = bplan;
 
 1424         fourierArray.
swap(newFourierArray);
 
 1425         fourierKernel.
swap(newFourierKernel);
 
 1428     void init(Shape inOut, Shape kernel,
 
 1429               unsigned int planner_flags = FFTW_ESTIMATE);
 
 1431     void initFourierKernel(Shape inOut, Shape kernel,
 
 1432                            unsigned int planner_flags = FFTW_ESTIMATE);
 
 1434     void initComplex(Shape inOut, Shape kernel,
 
 1435                      unsigned int planner_flags = FFTW_ESTIMATE);
 
 1437     void initMany(Shape inOut, Shape maxKernel,
 
 1438                   unsigned int planner_flags = FFTW_ESTIMATE)
 
 1440         init(inOut, maxKernel, planner_flags);
 
 1443     void initFourierKernelMany(Shape inOut, Shape kernels,
 
 1444                                unsigned int planner_flags = FFTW_ESTIMATE)
 
 1446         initFourierKernel(inOut, kernels, planner_flags);
 
 1456     template <
class C1, 
class C2, 
class C3>
 
 1461         executeImpl(in, kernel, out);
 
 1471     template <
class C1, 
class C2, 
class C3>
 
 1483     template <
class C1, 
class C2, 
class C3>
 
 1496     template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1498                      KernelIterator kernels, KernelIterator kernelsEnd,
 
 1508     template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1510                      KernelIterator kernels, KernelIterator kernelsEnd,
 
 1513         typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
 
 1514         typedef typename KernelArray::value_type KernelValue;
 
 1515         typedef typename IsSameType<KernelValue, Complex>::type UseFourierKernel;
 
 1516         typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
 
 1517         typedef typename OutArray::value_type OutValue;
 
 1519         bool realKernel = IsSameType<KernelValue, Real>::value;
 
 1520         bool fourierKernel = IsSameType<KernelValue, Complex>::value;
 
 1522         vigra_precondition(realKernel || fourierKernel,
 
 1523              "FFTWConvolvePlan::executeMany(): kernels have unsuitable value_type.");
 
 1524         vigra_precondition((IsSameType<OutValue, Real>::value),
 
 1525              "FFTWConvolvePlan::executeMany(): outputs have unsuitable value_type.");
 
 1527         executeManyImpl(in, kernels, kernelsEnd, outs, UseFourierKernel());
 
 1532     template <
class KernelIterator, 
class OutIterator>
 
 1533     Shape checkShapes(Shape in,
 
 1534                       KernelIterator kernels, KernelIterator kernelsEnd,
 
 1537     template <
class KernelIterator, 
class OutIterator>
 
 1538     Shape checkShapesFourier(Shape in,
 
 1539                              KernelIterator kernels, KernelIterator kernelsEnd,
 
 1542     template <
class KernelIterator, 
class OutIterator>
 
 1543     Shape checkShapesComplex(Shape in,
 
 1544                              KernelIterator kernels, KernelIterator kernelsEnd,
 
 1547     template <
class C1, 
class C2, 
class C3>
 
 1551                      bool do_correlation=
false);
 
 1553     template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1556                     KernelIterator kernels, KernelIterator kernelsEnd,
 
 1557                     OutIterator outs, VigraFalseType );
 
 1559     template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1562                     KernelIterator kernels, KernelIterator kernelsEnd,
 
 1563                     OutIterator outs, VigraTrueType );
 
 1567 template <
unsigned int N, 
class Real>
 
 1570                                 unsigned int planner_flags)
 
 1572     Shape paddedShape = fftwBestPaddedShapeR2C(in + kernel - Shape(1)),
 
 1575     CArray newFourierArray(complexShape), newFourierKernel(complexShape);
 
 1577     Shape realStrides = 2*newFourierArray.stride();
 
 1579     RArray newRealArray(paddedShape, realStrides, (Real*)newFourierArray.data());
 
 1580     RArray newRealKernel(paddedShape, realStrides, (Real*)newFourierKernel.data());
 
 1582     FFTWPlan<N, Real> fplan(newRealArray, newFourierArray, planner_flags);
 
 1583     FFTWPlan<N, Real> bplan(newFourierArray, newRealArray, planner_flags);
 
 1585     forward_plan = fplan;
 
 1586     backward_plan = bplan;
 
 1587     realArray = newRealArray;
 
 1588     realKernel = newRealKernel;
 
 1589     fourierArray.swap(newFourierArray);
 
 1590     fourierKernel.swap(newFourierKernel);
 
 1591     useFourierKernel = 
false;
 
 1594 template <
unsigned int N, 
class Real>
 
 1596 FFTWConvolvePlan<N, Real>::initFourierKernel(Shape in, Shape kernel,
 
 1597                                              unsigned int planner_flags)
 
 1599     Shape complexShape = kernel,
 
 1600           paddedShape  = fftwCorrespondingShapeC2R(complexShape);
 
 1602     for(
unsigned int k=0; k<N; ++k)
 
 1603         vigra_precondition(in[k] <= paddedShape[k],
 
 1604              "FFTWConvolvePlan::init(): kernel too small for given input.");
 
 1606     CArray newFourierArray(complexShape), newFourierKernel(complexShape);
 
 1608     Shape realStrides = 2*newFourierArray.stride();
 
 1610     RArray newRealArray(paddedShape, realStrides, (Real*)newFourierArray.data());
 
 1611     RArray newRealKernel(paddedShape, realStrides, (Real*)newFourierKernel.data());
 
 1613     FFTWPlan<N, Real> fplan(newRealArray, newFourierArray, planner_flags);
 
 1614     FFTWPlan<N, Real> bplan(newFourierArray, newRealArray, planner_flags);
 
 1616     forward_plan = fplan;
 
 1617     backward_plan = bplan;
 
 1618     realArray = newRealArray;
 
 1619     realKernel = newRealKernel;
 
 1620     fourierArray.swap(newFourierArray);
 
 1621     fourierKernel.swap(newFourierKernel);
 
 1622     useFourierKernel = 
true;
 
 1625 template <
unsigned int N, 
class Real>
 
 1627 FFTWConvolvePlan<N, Real>::initComplex(Shape in, Shape kernel,
 
 1628                                         unsigned int planner_flags)
 
 1632     if(useFourierKernel)
 
 1634         for(
unsigned int k=0; k<N; ++k)
 
 1635             vigra_precondition(in[k] <= kernel[k],
 
 1636                  "FFTWConvolvePlan::init(): kernel too small for given input.");
 
 1638         paddedShape = kernel;
 
 1642         paddedShape  = fftwBestPaddedShape(in + kernel - Shape(1));
 
 1645     CArray newFourierArray(paddedShape), newFourierKernel(paddedShape);
 
 1647     FFTWPlan<N, Real> fplan(newFourierArray, newFourierArray, FFTW_FORWARD, planner_flags);
 
 1648     FFTWPlan<N, Real> bplan(newFourierArray, newFourierArray, FFTW_BACKWARD, planner_flags);
 
 1650     forward_plan = fplan;
 
 1651     backward_plan = bplan;
 
 1652     fourierArray.swap(newFourierArray);
 
 1653     fourierKernel.swap(newFourierKernel);
 
 1656 #ifndef DOXYGEN // doxygen documents these functions as free functions 
 1658 template <
unsigned int N, 
class Real>
 
 1659 template <
class C1, 
class C2, 
class C3>
 
 1661 FFTWConvolvePlan<N, Real>::executeImpl(MultiArrayView<N, Real, C1> in,
 
 1662                                        MultiArrayView<N, Real, C2> kernel,
 
 1663                                        MultiArrayView<N, Real, C3> out,
 
 1664                                        bool do_correlation)
 
 1666     vigra_precondition(!useFourierKernel,
 
 1667        "FFTWConvolvePlan::execute(): plan was generated for Fourier kernel, got spatial kernel.");
 
 1669     vigra_precondition(in.shape() == out.shape(),
 
 1670         "FFTWConvolvePlan::execute(): input and output must have the same shape.");
 
 1672     Shape paddedShape = fftwBestPaddedShapeR2C(in.shape() + kernel.shape() - Shape(1)),
 
 1673           diff = paddedShape - in.shape(),
 
 1675           right = in.shape() + left;
 
 1677     vigra_precondition(paddedShape == realArray.shape(),
 
 1678        "FFTWConvolvePlan::execute(): shape mismatch between input and plan.");
 
 1680     detail::fftEmbedArray(in, realArray);
 
 1681     forward_plan.execute(realArray, fourierArray);
 
 1683     detail::fftEmbedKernel(kernel, realKernel);
 
 1684     forward_plan.execute(realKernel, fourierKernel);
 
 1688         using namespace multi_math;
 
 1689         fourierArray *= 
conj(fourierKernel);
 
 1693         fourierArray *= fourierKernel;
 
 1696     backward_plan.execute(fourierArray, realArray);
 
 1698     out = realArray.subarray(left, right);
 
 1701 template <
unsigned int N, 
class Real>
 
 1702 template <
class C1, 
class C2, 
class C3>
 
 1705                                    MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
 
 1706                                    MultiArrayView<N, Real, C3> out)
 
 1708     vigra_precondition(useFourierKernel,
 
 1709        "FFTWConvolvePlan::execute(): plan was generated for spatial kernel, got Fourier kernel.");
 
 1711     vigra_precondition(in.shape() == out.shape(),
 
 1712         "FFTWConvolvePlan::execute(): input and output must have the same shape.");
 
 1714     vigra_precondition(kernel.shape() == fourierArray.shape(),
 
 1715        "FFTWConvolvePlan::execute(): shape mismatch between kernel and plan.");
 
 1717     Shape paddedShape = fftwCorrespondingShapeC2R(kernel.shape(), 
odd(in.shape(0))),
 
 1718           diff = paddedShape - in.shape(),
 
 1720           right = in.shape() + left;
 
 1722     vigra_precondition(paddedShape == realArray.shape(),
 
 1723        "FFTWConvolvePlan::execute(): shape mismatch between input and plan.");
 
 1725     detail::fftEmbedArray(in, realArray);
 
 1726     forward_plan.execute(realArray, fourierArray);
 
 1728     fourierKernel = kernel;
 
 1729     moveDCToHalfspaceUpperLeft(fourierKernel);
 
 1731     fourierArray *= fourierKernel;
 
 1733     backward_plan.execute(fourierArray, realArray);
 
 1735     out = realArray.subarray(left, right);
 
 1738 template <
unsigned int N, 
class Real>
 
 1739 template <
class C1, 
class C2, 
class C3>
 
 1742                                     MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
 
 1743                                     MultiArrayView<N, FFTWComplex<Real>, C3> out)
 
 1745     vigra_precondition(in.shape() == out.shape(),
 
 1746         "FFTWConvolvePlan::execute(): input and output must have the same shape.");
 
 1748     Shape paddedShape = fourierArray.shape(),
 
 1749           diff = paddedShape - in.shape(),
 
 1751           right = in.shape() + left;
 
 1753     if(useFourierKernel)
 
 1755         vigra_precondition(kernel.shape() == fourierArray.shape(),
 
 1756            "FFTWConvolvePlan::execute(): shape mismatch between kernel and plan.");
 
 1758         fourierKernel = kernel;
 
 1763         detail::fftEmbedKernel(kernel, fourierKernel);
 
 1764         forward_plan.execute(fourierKernel, fourierKernel);
 
 1767     detail::fftEmbedArray(in, fourierArray);
 
 1768     forward_plan.execute(fourierArray, fourierArray);
 
 1770     fourierArray *= fourierKernel;
 
 1772     backward_plan.execute(fourierArray, fourierArray);
 
 1774     out = fourierArray.subarray(left, right);
 
 1777 template <
unsigned int N, 
class Real>
 
 1778 template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1780 FFTWConvolvePlan<N, Real>::executeManyImpl(MultiArrayView<N, Real, C1> in,
 
 1781                                            KernelIterator kernels, KernelIterator kernelsEnd,
 
 1782                                            OutIterator outs, VigraFalseType )
 
 1784     vigra_precondition(!useFourierKernel,
 
 1785        "FFTWConvolvePlan::execute(): plan was generated for Fourier kernel, got spatial kernel.");
 
 1787     Shape kernelMax = checkShapes(in.shape(), kernels, kernelsEnd, outs),
 
 1788           paddedShape = fftwBestPaddedShapeR2C(in.shape() + kernelMax - Shape(1)),
 
 1789           diff = paddedShape - in.shape(),
 
 1791           right = in.shape() + left;
 
 1793     vigra_precondition(paddedShape == realArray.shape(),
 
 1794        "FFTWConvolvePlan::executeMany(): shape mismatch between input and plan.");
 
 1796     detail::fftEmbedArray(in, realArray);
 
 1797     forward_plan.execute(realArray, fourierArray);
 
 1799     for(; kernels != kernelsEnd; ++kernels, ++outs)
 
 1801         detail::fftEmbedKernel(*kernels, realKernel);
 
 1802         forward_plan.execute(realKernel, fourierKernel);
 
 1804         fourierKernel *= fourierArray;
 
 1806         backward_plan.execute(fourierKernel, realKernel);
 
 1808         *outs = realKernel.subarray(left, right);
 
 1812 template <
unsigned int N, 
class Real>
 
 1813 template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1815 FFTWConvolvePlan<N, Real>::executeManyImpl(MultiArrayView<N, Real, C1> in,
 
 1816                                            KernelIterator kernels, KernelIterator kernelsEnd,
 
 1817                                            OutIterator outs, VigraTrueType )
 
 1819     vigra_precondition(useFourierKernel,
 
 1820        "FFTWConvolvePlan::execute(): plan was generated for spatial kernel, got Fourier kernel.");
 
 1822     Shape complexShape = checkShapesFourier(in.shape(), kernels, kernelsEnd, outs),
 
 1823           paddedShape = fftwCorrespondingShapeC2R(complexShape, 
odd(in.shape(0))),
 
 1824           diff = paddedShape - in.shape(),
 
 1826           right = in.shape() + left;
 
 1828     vigra_precondition(complexShape == fourierArray.shape(),
 
 1829        "FFTWConvolvePlan::executeFourierKernelMany(): shape mismatch between kernels and plan.");
 
 1831     vigra_precondition(paddedShape == realArray.shape(),
 
 1832        "FFTWConvolvePlan::executeFourierKernelMany(): shape mismatch between input and plan.");
 
 1834     detail::fftEmbedArray(in, realArray);
 
 1835     forward_plan.execute(realArray, fourierArray);
 
 1837     for(; kernels != kernelsEnd; ++kernels, ++outs)
 
 1839         fourierKernel = *kernels;
 
 1840         moveDCToHalfspaceUpperLeft(fourierKernel);
 
 1841         fourierKernel *= fourierArray;
 
 1843         backward_plan.execute(fourierKernel, realKernel);
 
 1845         *outs = realKernel.subarray(left, right);
 
 1849 template <
unsigned int N, 
class Real>
 
 1850 template <
class C1, 
class KernelIterator, 
class OutIterator>
 
 1853                                        KernelIterator kernels, KernelIterator kernelsEnd,
 
 1856     typedef typename std::iterator_traits<KernelIterator>::value_type KernelArray;
 
 1857     typedef typename KernelArray::value_type KernelValue;
 
 1858     typedef typename std::iterator_traits<OutIterator>::value_type OutArray;
 
 1859     typedef typename OutArray::value_type OutValue;
 
 1861     vigra_precondition((IsSameType<KernelValue, Complex>::value),
 
 1862          "FFTWConvolvePlan::executeMany(): kernels have unsuitable value_type.");
 
 1863     vigra_precondition((IsSameType<OutValue, Complex>::value),
 
 1864          "FFTWConvolvePlan::executeMany(): outputs have unsuitable value_type.");
 
 1866     Shape paddedShape = checkShapesComplex(in.shape(), kernels, kernelsEnd, outs),
 
 1867           diff = paddedShape - in.shape(),
 
 1869           right = in.shape() + left;
 
 1871     detail::fftEmbedArray(in, fourierArray);
 
 1872     forward_plan.execute(fourierArray, fourierArray);
 
 1874     for(; kernels != kernelsEnd; ++kernels, ++outs)
 
 1876         if(useFourierKernel)
 
 1878             fourierKernel = *kernels;
 
 1883             detail::fftEmbedKernel(*kernels, fourierKernel);
 
 1884             forward_plan.execute(fourierKernel, fourierKernel);
 
 1887         fourierKernel *= fourierArray;
 
 1889         backward_plan.execute(fourierKernel, fourierKernel);
 
 1891         *outs = fourierKernel.subarray(left, right);
 
 1897 template <
unsigned int N, 
class Real>
 
 1898 template <
class KernelIterator, 
class OutIterator>
 
 1899 typename FFTWConvolvePlan<N, Real>::Shape
 
 1900 FFTWConvolvePlan<N, Real>::checkShapes(Shape in,
 
 1901                                        KernelIterator kernels, KernelIterator kernelsEnd,
 
 1904     vigra_precondition(kernels != kernelsEnd,
 
 1905         "FFTWConvolvePlan::checkShapes(): empty kernel sequence.");
 
 1908     for(; kernels != kernelsEnd; ++kernels, ++outs)
 
 1910         vigra_precondition(in == outs->shape(),
 
 1911             "FFTWConvolvePlan::checkShapes(): shape mismatch between input and (one) output.");
 
 1912         kernelMax = max(kernelMax, kernels->shape());
 
 1914     vigra_precondition(
prod(kernelMax) > 0,
 
 1915         "FFTWConvolvePlan::checkShapes(): all kernels have size 0.");
 
 1919 template <
unsigned int N, 
class Real>
 
 1920 template <
class KernelIterator, 
class OutIterator>
 
 1921 typename FFTWConvolvePlan<N, Real>::Shape
 
 1922 FFTWConvolvePlan<N, Real>::checkShapesFourier(Shape in,
 
 1923                                                KernelIterator kernels, KernelIterator kernelsEnd,
 
 1926     vigra_precondition(kernels != kernelsEnd,
 
 1927         "FFTWConvolvePlan::checkShapesFourier(): empty kernel sequence.");
 
 1929     Shape complexShape = kernels->shape(),
 
 1930           paddedShape  = fftwCorrespondingShapeC2R(complexShape);
 
 1932     for(
unsigned int k=0; k<N; ++k)
 
 1933         vigra_precondition(in[k] <= paddedShape[k],
 
 1934              "FFTWConvolvePlan::checkShapesFourier(): kernels too small for given input.");
 
 1936     for(; kernels != kernelsEnd; ++kernels, ++outs)
 
 1938         vigra_precondition(in == outs->shape(),
 
 1939             "FFTWConvolvePlan::checkShapesFourier(): shape mismatch between input and (one) output.");
 
 1940         vigra_precondition(complexShape == kernels->shape(),
 
 1941             "FFTWConvolvePlan::checkShapesFourier(): all kernels must have the same size.");
 
 1943     return complexShape;
 
 1946 template <
unsigned int N, 
class Real>
 
 1947 template <
class KernelIterator, 
class OutIterator>
 
 1948 typename FFTWConvolvePlan<N, Real>::Shape
 
 1949 FFTWConvolvePlan<N, Real>::checkShapesComplex(Shape in,
 
 1950                                                KernelIterator kernels, KernelIterator kernelsEnd,
 
 1953     vigra_precondition(kernels != kernelsEnd,
 
 1954         "FFTWConvolvePlan::checkShapesComplex(): empty kernel sequence.");
 
 1956     Shape kernelShape = kernels->shape();
 
 1957     for(; kernels != kernelsEnd; ++kernels, ++outs)
 
 1959         vigra_precondition(in == outs->shape(),
 
 1960             "FFTWConvolvePlan::checkShapesComplex(): shape mismatch between input and (one) output.");
 
 1961         if(useFourierKernel)
 
 1963             vigra_precondition(kernelShape == kernels->shape(),
 
 1964                 "FFTWConvolvePlan::checkShapesComplex(): Fourier domain kernels must have identical size.");
 
 1968             kernelShape = max(kernelShape, kernels->shape());
 
 1971     vigra_precondition(
prod(kernelShape) > 0,
 
 1972         "FFTWConvolvePlan::checkShapesComplex(): all kernels have size 0.");
 
 1974     if(useFourierKernel)
 
 1976         for(
unsigned int k=0; k<N; ++k)
 
 1977             vigra_precondition(in[k] <= kernelShape[k],
 
 1978                  "FFTWConvolvePlan::checkShapesComplex(): kernels too small for given input.");
 
 1983         return fftwBestPaddedShape(in + kernelShape - Shape(1));
 
 2019 template <
unsigned int N, 
class Real>
 
 2047     template <
class C1, 
class C2, 
class C3>
 
 2051                       unsigned int planner_flags = FFTW_ESTIMATE)
 
 2052     : 
BaseType(in, kernel, out, planner_flags)
 
 2067     template <
class C1, 
class C2, 
class C3>
 
 2069                      bool useFourierKernel = 
false,
 
 2070                      unsigned int planner_flags = FFTW_ESTIMATE)
 
 2071     : 
BaseType(inOut, kernel, false, planner_flags)
 
 2073         ignore_argument(useFourierKernel);
 
 2080     template <
class C1, 
class C2, 
class C3>
 
 2084               unsigned int planner_flags = FFTW_ESTIMATE)
 
 2086         vigra_precondition(in.
shape() == out.
shape(),
 
 2087                            "FFTWCorrelatePlan::init(): input and output must have the same shape.");
 
 2098     template <
class C1, 
class C2, 
class C3>
 
 2103         BaseType::executeImpl(in, kernel, out, 
true);
 
 2113 template <
unsigned int N, 
class Real, 
class C1, 
class C2>
 
 2116                  MultiArrayView<N, FFTWComplex<Real>, C2> out)
 
 2118     FFTWPlan<N, Real>(in, out, FFTW_FORWARD).execute(in, out);
 
 2121 template <
unsigned int N, 
class Real, 
class C1, 
class C2>
 
 2124                         MultiArrayView<N, FFTWComplex<Real>, C2> out)
 
 2126     FFTWPlan<N, Real>(in, out, FFTW_BACKWARD).execute(in, out);
 
 2129 template <
unsigned int N, 
class Real, 
class C1, 
class C2>
 
 2132                  MultiArrayView<N, FFTWComplex<Real>, C2> out)
 
 2134     if(in.shape() == out.shape())
 
 2138         FFTWPlan<N, Real>(out, out, FFTW_FORWARD).execute(out, out);
 
 2142         FFTWPlan<N, Real>(in, out).execute(in, out);
 
 2145         vigra_precondition(
false,
 
 2146             "fourierTransform(): shape mismatch between input and output.");
 
 2149 template <
unsigned int N, 
class Real, 
class C1, 
class C2>
 
 2152                         MultiArrayView<N, Real, C2> out)
 
 2155         "fourierTransformInverse(): shape mismatch between input and output.");
 
 2156     FFTWPlan<N, Real>(in, out).execute(in, out);
 
 2350 template <
unsigned int N, 
class Real, 
class C1, 
class C2, 
class C3>
 
 2353             MultiArrayView<N, Real, C2> kernel,
 
 2354             MultiArrayView<N, Real, C3> out)
 
 2356     FFTWConvolvePlan<N, Real>(in, kernel, out).execute(in, kernel, out);
 
 2359 template <
unsigned int N, 
class Real, 
class C1, 
class C2, 
class C3>
 
 2362             MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
 
 2363             MultiArrayView<N, Real, C3> out)
 
 2365     FFTWConvolvePlan<N, Real>(in, kernel, out).execute(in, kernel, out);
 
 2374 template <
unsigned int N, 
class Real, 
class C1, 
class C2, 
class C3>
 
 2377             MultiArrayView<N, FFTWComplex<Real>, C2> kernel,
 
 2378             MultiArrayView<N, FFTWComplex<Real>, C3> out,
 
 2379             bool fourierDomainKernel)
 
 2381     FFTWConvolvePlan<N, Real>(in, kernel, out, fourierDomainKernel).execute(in, kernel, out);
 
 2390 template <
unsigned int N, 
class Real, 
class C1,
 
 2391           class KernelIterator, 
class OutIterator>
 
 2394                 KernelIterator kernels, KernelIterator kernelsEnd,
 
 2397     FFTWConvolvePlan<N, Real> plan;
 
 2398     plan.initMany(in, kernels, kernelsEnd, outs);
 
 2399     plan.executeMany(in, kernels, kernelsEnd, outs);
 
 2408 template <
unsigned int N, 
class Real, 
class C1,
 
 2409           class KernelIterator, 
class OutIterator>
 
 2412                 KernelIterator kernels, KernelIterator kernelsEnd,
 
 2414                 bool fourierDomainKernel)
 
 2416     FFTWConvolvePlan<N, Real> plan;
 
 2417     plan.initMany(in, kernels, kernelsEnd, outs, fourierDomainKernel);
 
 2418     plan.executeMany(in, kernels, kernelsEnd, outs);
 
 2471 template <
unsigned int N, 
class Real, 
class C1, 
class C2, 
class C3>
 
 2474             MultiArrayView<N, Real, C2> kernel,
 
 2475             MultiArrayView<N, Real, C3> out)
 
 2477     FFTWCorrelatePlan<N, Real>(in, kernel, out).execute(in, kernel, out);
 
 2484 #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:1022
FFTWCorrelatePlan()
Create an empty plan. 
Definition: multi_fft.hxx:2032
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:1457
~FFTWPlan()
Destructor. 
Definition: multi_fft.hxx:960
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:1003
const difference_type & shape() const 
Definition: multi_array.hxx:1648
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:1326
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:1273
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:1311
void convolveFFT(...)
Convolve an array with a kernel by means of the Fourier transform. 
Main MultiArray class containing the memory management. 
Definition: multi_array.hxx:2474
void execute(MultiArrayView< N, Real, C1 > in, MultiArrayView< N, FFTWComplex< Real >, C2 > out) const 
Execute a real-to-complex transform. 
Definition: multi_fft.hxx:1036
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:1360
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:971
bool odd(int t)
Check if an integer is odd. 
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: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:880
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:1341
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:900
Definition: multi_fft.hxx:852
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:1397
void convolveFFTComplexMany(...)
Convolve a complex-valued array with a sequence of kernels by means of the Fourier transform...
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:2081
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:2048
MultiArrayView< N, T, StridedArrayTag > permuteStridesDescending() const 
Definition: multi_array.hxx:2173
doxygen_overloaded_function(template<...> void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays. 
Definition: metaprogramming.hxx:116
bool even(int t)
Check if an integer is even. 
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:1509
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:2099
Definition: multi_fft.hxx:2020
void convolveFFTMany(...)
Convolve a real-valued array with a sequence of kernels by means of the Fourier transform. 
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:1251
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:920
Definition: multi_fft.hxx:1195
void correlateFFT(...)
Correlate an array with a kernel by means of the Fourier transform. 
FFTWPlan(FFTWPlan const &other)
Copy constructor. 
Definition: multi_fft.hxx:930
Base class for, and view to, vigra::MultiArray. 
Definition: multi_array.hxx:704
FFTWPlan()
Create an empty plan. 
Definition: multi_fft.hxx:867
FFTWPlan & operator=(FFTWPlan const &other)
Copy assigment. 
Definition: multi_fft.hxx:943
void convolveFFTComplex(...)
Convolve a complex-valued array by means of the Fourier transform. 
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:1296
T sign(T t)
The sign function. 
Definition: mathutil.hxx:591
void execute(MultiArrayView< N, FFTWComplex< Real >, C1 > in, MultiArrayView< N, Real, C2 > out) const 
Execute a complex-to-real transform. 
Definition: multi_fft.hxx:1050
void swap(MultiArray &other)
Definition: multi_array.hxx:3070
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:987
FFTWConvolvePlan()
Create an empty plan. 
Definition: multi_fft.hxx:1214
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:2068
difference_type strideOrdering() const 
Definition: multi_array.hxx:1617
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:1230