annotate src/fftw-3.3.5/doc/FAQ/fftw-faq.html/section3.html @ 127:7867fa7e1b6b

Current fftw source
author Chris Cannam <cannam@all-day-breakfast.com>
date Tue, 18 Oct 2016 13:40:26 +0100
parents
children
rev   line source
cannam@127 1 <!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 3.2//EN">
cannam@127 2 <html>
cannam@127 3 <head><title>
cannam@127 4 FFTW FAQ - Section 3
cannam@127 5 </title>
cannam@127 6 <link rev="made" href="mailto:fftw@fftw.org">
cannam@127 7 <link rel="Contents" href="index.html">
cannam@127 8 <link rel="Start" href="index.html">
cannam@127 9 <link rel="Next" href="section4.html"><link rel="Previous" href="section2.html"><link rel="Bookmark" title="FFTW FAQ" href="index.html">
cannam@127 10 </head><body text="#000000" bgcolor="#FFFFFF"><h1>
cannam@127 11 FFTW FAQ - Section 3 <br>
cannam@127 12 Using FFTW
cannam@127 13 </h1>
cannam@127 14
cannam@127 15 <ul>
cannam@127 16 <li><a href="#fftw2to3" rel=subdocument>Q3.1. Why not support the FFTW 2 interface in FFTW
cannam@127 17 3?</a>
cannam@127 18 <li><a href="#planperarray" rel=subdocument>Q3.2. Why do FFTW 3 plans encapsulate the input/output arrays and not just
cannam@127 19 the algorithm?</a>
cannam@127 20 <li><a href="#slow" rel=subdocument>Q3.3. FFTW seems really slow.</a>
cannam@127 21 <li><a href="#slows" rel=subdocument>Q3.4. FFTW slows down after repeated calls.</a>
cannam@127 22 <li><a href="#segfault" rel=subdocument>Q3.5. An FFTW routine is crashing when I call it.</a>
cannam@127 23 <li><a href="#fortran64" rel=subdocument>Q3.6. My Fortran program crashes when calling FFTW.</a>
cannam@127 24 <li><a href="#conventions" rel=subdocument>Q3.7. FFTW gives results different from my old
cannam@127 25 FFT.</a>
cannam@127 26 <li><a href="#nondeterministic" rel=subdocument>Q3.8. FFTW gives different results between runs</a>
cannam@127 27 <li><a href="#savePlans" rel=subdocument>Q3.9. Can I save FFTW's plans?</a>
cannam@127 28 <li><a href="#whyscaled" rel=subdocument>Q3.10. Why does your inverse transform return a scaled
cannam@127 29 result?</a>
cannam@127 30 <li><a href="#centerorigin" rel=subdocument>Q3.11. How can I make FFTW put the origin (zero frequency) at the center of
cannam@127 31 its output?</a>
cannam@127 32 <li><a href="#imageaudio" rel=subdocument>Q3.12. How do I FFT an image/audio file in <i>foobar</i> format?</a>
cannam@127 33 <li><a href="#linkfails" rel=subdocument>Q3.13. My program does not link (on Unix).</a>
cannam@127 34 <li><a href="#linkheader" rel=subdocument>Q3.14. I included your header, but linking still
cannam@127 35 fails.</a>
cannam@127 36 <li><a href="#nostack" rel=subdocument>Q3.15. My program crashes, complaining about stack
cannam@127 37 space.</a>
cannam@127 38 <li><a href="#leaks" rel=subdocument>Q3.16. FFTW seems to have a memory leak.</a>
cannam@127 39 <li><a href="#allzero" rel=subdocument>Q3.17. The output of FFTW's transform is all zeros.</a>
cannam@127 40 <li><a href="#vbetalia" rel=subdocument>Q3.18. How do I call FFTW from the Microsoft language du
cannam@127 41 jour?</a>
cannam@127 42 <li><a href="#pruned" rel=subdocument>Q3.19. Can I compute only a subset of the DFT outputs?</a>
cannam@127 43 <li><a href="#transpose" rel=subdocument>Q3.20. Can I use FFTW's routines for in-place and out-of-place matrix
cannam@127 44 transposition?</a>
cannam@127 45 </ul><hr>
cannam@127 46
cannam@127 47 <h2><A name="fftw2to3">
cannam@127 48 Question 3.1. Why not support the FFTW 2 interface in FFTW
cannam@127 49 3?
cannam@127 50 </A></h2>
cannam@127 51
cannam@127 52 FFTW 3 has semantics incompatible with earlier versions: its plans can
cannam@127 53 only be used for a given stride, multiplicity, and other
cannam@127 54 characteristics of the input and output arrays; these stronger
cannam@127 55 semantics are necessary for performance reasons. Thus, it is
cannam@127 56 impossible to efficiently emulate the older interface (whose plans can
cannam@127 57 be used for any transform of the same size). We believe that it
cannam@127 58 should be possible to upgrade most programs without any difficulty,
cannam@127 59 however.
cannam@127 60 <h2><A name="planperarray">
cannam@127 61 Question 3.2. Why do FFTW 3 plans encapsulate the input/output arrays
cannam@127 62 and not just the algorithm?
cannam@127 63 </A></h2>
cannam@127 64
cannam@127 65 There are several reasons:
cannam@127 66 <ul>
cannam@127 67 <li>It was important for performance reasons that the plan be specific to
cannam@127 68 array characteristics like the stride (and alignment, for SIMD), and
cannam@127 69 requiring that the user maintain these invariants is error prone.
cannam@127 70
cannam@127 71 <li>In most high-performance applications, as far as we can tell, you are
cannam@127 72 usually transforming the same array over and over, so FFTW's semantics
cannam@127 73 should not be a burden.
cannam@127 74 <li>If you need to transform another array of the same size, creating a
cannam@127 75 new plan once the first exists is a cheap operation.
cannam@127 76
cannam@127 77 <li>If you need to transform many arrays of the same size at once, you
cannam@127 78 should really use the <code>plan_many</code> routines in FFTW's &quot;advanced&quot;
cannam@127 79 interface.
cannam@127 80 <li>If the abovementioned array characteristics are the same, you are
cannam@127 81 willing to pay close attention to the documentation, and you really
cannam@127 82 need to, we provide a &quot;new-array execution&quot; interface to
cannam@127 83 apply a plan to a new array.
cannam@127 84 </ul>
cannam@127 85
cannam@127 86 <h2><A name="slow">
cannam@127 87 Question 3.3. FFTW seems really slow.
cannam@127 88 </A></h2>
cannam@127 89
cannam@127 90 You are probably recreating the plan before every transform, rather
cannam@127 91 than creating it once and reusing it for all transforms of the same
cannam@127 92 size. FFTW is designed to be used in the following way:
cannam@127 93
cannam@127 94 <ul>
cannam@127 95 <li>First, you create a plan. This will take several seconds.
cannam@127 96
cannam@127 97 <li>Then, you reuse the plan many times to perform FFTs. These are fast.
cannam@127 98
cannam@127 99 </ul>
cannam@127 100 If you don't need to compute many transforms and the time for the
cannam@127 101 planner is significant, you have two options. First, you can use the
cannam@127 102 <code>FFTW_ESTIMATE</code> option in the planner, which uses heuristics
cannam@127 103 instead of runtime measurements and produces a good plan in a short
cannam@127 104 time. Second, you can use the wisdom feature to precompute the plan;
cannam@127 105 see <A href="#savePlans">Q3.9 `Can I save FFTW's plans?'</A>
cannam@127 106 <h2><A name="slows">
cannam@127 107 Question 3.4. FFTW slows down after repeated
cannam@127 108 calls.
cannam@127 109 </A></h2>
cannam@127 110
cannam@127 111 Probably, NaNs or similar are creeping into your data, and the
cannam@127 112 slowdown is due to the resulting floating-point exceptions. For
cannam@127 113 example, be aware that repeatedly FFTing the same array is a diverging
cannam@127 114 process (because FFTW computes the unnormalized transform).
cannam@127 115
cannam@127 116 <h2><A name="segfault">
cannam@127 117 Question 3.5. An FFTW routine is crashing when I call
cannam@127 118 it.
cannam@127 119 </A></h2>
cannam@127 120
cannam@127 121 Did the FFTW test programs pass (<code>make check</code>, or <code>cd tests; make bigcheck</code> if you want to be paranoid)? If so, you almost
cannam@127 122 certainly have a bug in your own code. For example, you could be
cannam@127 123 passing invalid arguments (such as wrongly-sized arrays) to FFTW, or
cannam@127 124 you could simply have memory corruption elsewhere in your program that
cannam@127 125 causes random crashes later on. Please don't complain to us unless
cannam@127 126 you can come up with a minimal self-contained program (preferably
cannam@127 127 under 30 lines) that illustrates the problem.
cannam@127 128
cannam@127 129 <h2><A name="fortran64">
cannam@127 130 Question 3.6. My Fortran program crashes when calling
cannam@127 131 FFTW.
cannam@127 132 </A></h2>
cannam@127 133
cannam@127 134 As described in the manual, on 64-bit machines you must store the
cannam@127 135 plans in variables large enough to hold a pointer, for example
cannam@127 136 <code>integer*8</code>. We recommend using <code>integer*8</code> on 32-bit machines as well, to simplify porting.
cannam@127 137
cannam@127 138 <h2><A name="conventions">
cannam@127 139 Question 3.7. FFTW gives results different from my old
cannam@127 140 FFT.
cannam@127 141 </A></h2>
cannam@127 142
cannam@127 143 People follow many different conventions for the DFT, and you should
cannam@127 144 be sure to know the ones that we use (described in the FFTW manual).
cannam@127 145 In particular, you should be aware that the
cannam@127 146 <code>FFTW_FORWARD</code>/<code>FFTW_BACKWARD</code> directions correspond to signs of -1/+1 in the exponent of the DFT definition.
cannam@127 147 (<i>Numerical Recipes</i> uses the opposite convention.)
cannam@127 148 <p>
cannam@127 149 You should also know that we compute an unnormalized transform. In
cannam@127 150 contrast, Matlab is an example of program that computes a normalized
cannam@127 151 transform. See <A href="#whyscaled">Q3.10 `Why does your inverse transform return a scaled
cannam@127 152 result?'</A>.
cannam@127 153 <p>
cannam@127 154 Finally, note that floating-point arithmetic is not exact, so
cannam@127 155 different FFT algorithms will give slightly different results (on the
cannam@127 156 order of the numerical accuracy; typically a fractional difference of
cannam@127 157 1e-15 or so in double precision).
cannam@127 158 <h2><A name="nondeterministic">
cannam@127 159 Question 3.8. FFTW gives different results between
cannam@127 160 runs
cannam@127 161 </A></h2>
cannam@127 162
cannam@127 163 If you use <code>FFTW_MEASURE</code> or <code>FFTW_PATIENT</code> mode, then the algorithm FFTW employs is not deterministic: it depends on
cannam@127 164 runtime performance measurements. This will cause the results to vary
cannam@127 165 slightly from run to run. However, the differences should be slight,
cannam@127 166 on the order of the floating-point precision, and therefore should
cannam@127 167 have no practical impact on most applications.
cannam@127 168
cannam@127 169 <p>
cannam@127 170 If you use saved plans (wisdom) or <code>FFTW_ESTIMATE</code> mode, however, then the algorithm is deterministic and the results should be
cannam@127 171 identical between runs.
cannam@127 172 <h2><A name="savePlans">
cannam@127 173 Question 3.9. Can I save FFTW's plans?
cannam@127 174 </A></h2>
cannam@127 175
cannam@127 176 Yes. Starting with version 1.2, FFTW provides the
cannam@127 177 <code>wisdom</code> mechanism for saving plans; see the FFTW manual.
cannam@127 178
cannam@127 179 <h2><A name="whyscaled">
cannam@127 180 Question 3.10. Why does your inverse transform return a scaled
cannam@127 181 result?
cannam@127 182 </A></h2>
cannam@127 183
cannam@127 184 Computing the forward transform followed by the backward transform (or
cannam@127 185 vice versa) yields the original array scaled by the size of the array.
cannam@127 186 (For multi-dimensional transforms, the size of the array is the
cannam@127 187 product of the dimensions.) We could, instead, have chosen a
cannam@127 188 normalization that would have returned the unscaled array. Or, to
cannam@127 189 accomodate the many conventions in this matter, the transform routines
cannam@127 190 could have accepted a &quot;scale factor&quot; parameter. We did not
cannam@127 191 do this, however, for two reasons. First, we didn't want to sacrifice
cannam@127 192 performance in the common case where the scale factor is 1. Second, in
cannam@127 193 real applications the FFT is followed or preceded by some computation
cannam@127 194 on the data, into which the scale factor can typically be absorbed at
cannam@127 195 little or no cost.
cannam@127 196 <h2><A name="centerorigin">
cannam@127 197 Question 3.11. How can I make FFTW put the origin (zero frequency) at
cannam@127 198 the center of its output?
cannam@127 199 </A></h2>
cannam@127 200
cannam@127 201 For human viewing of a spectrum, it is often convenient to put the
cannam@127 202 origin in frequency space at the center of the output array, rather
cannam@127 203 than in the zero-th element (the default in FFTW). If all of the
cannam@127 204 dimensions of your array are even, you can accomplish this by simply
cannam@127 205 multiplying each element of the input array by (-1)^(i + j + ...),
cannam@127 206 where i, j, etcetera are the indices of the element. (This trick is a
cannam@127 207 general property of the DFT, and is not specific to FFTW.)
cannam@127 208
cannam@127 209 <h2><A name="imageaudio">
cannam@127 210 Question 3.12. How do I FFT an image/audio file in
cannam@127 211 <i>foobar</i> format?
cannam@127 212 </A></h2>
cannam@127 213
cannam@127 214 FFTW performs an FFT on an array of floating-point values. You can
cannam@127 215 certainly use it to compute the transform of an image or audio stream,
cannam@127 216 but you are responsible for figuring out your data format and
cannam@127 217 converting it to the form FFTW requires.
cannam@127 218
cannam@127 219 <h2><A name="linkfails">
cannam@127 220 Question 3.13. My program does not link (on
cannam@127 221 Unix).
cannam@127 222 </A></h2>
cannam@127 223
cannam@127 224 The libraries must be listed in the correct order
cannam@127 225 (<code>-lfftw3 -lm</code> for FFTW 3.x) and <i>after</i> your program sources/objects. (The general rule is that if <i>A</i> uses <i>B</i>, then <i>A</i> must be listed before <i>B</i> in the link command.).
cannam@127 226 <h2><A name="linkheader">
cannam@127 227 Question 3.14. I included your header, but linking still
cannam@127 228 fails.
cannam@127 229 </A></h2>
cannam@127 230
cannam@127 231 You're a C++ programmer, aren't you? You have to compile the FFTW
cannam@127 232 library and link it into your program, not just
cannam@127 233 <code>#include &lt;fftw3.h&gt;</code>. (Yes, this is really a FAQ.)
cannam@127 234 <h2><A name="nostack">
cannam@127 235 Question 3.15. My program crashes, complaining about stack
cannam@127 236 space.
cannam@127 237 </A></h2>
cannam@127 238
cannam@127 239 You cannot declare large arrays with automatic storage (e.g. via
cannam@127 240 <code>fftw_complex array[N]</code>); you should use <code>fftw_malloc</code> (or equivalent) to allocate the arrays you want
cannam@127 241 to transform if they are larger than a few hundred elements.
cannam@127 242
cannam@127 243 <h2><A name="leaks">
cannam@127 244 Question 3.16. FFTW seems to have a memory
cannam@127 245 leak.
cannam@127 246 </A></h2>
cannam@127 247
cannam@127 248 After you create a plan, FFTW caches the information required to
cannam@127 249 quickly recreate the plan. (See <A href="#savePlans">Q3.9 `Can I save FFTW's plans?'</A>) It also maintains a small amount of other persistent memory. You can deallocate all of
cannam@127 250 FFTW's internally allocated memory, if you wish, by calling
cannam@127 251 <code>fftw_cleanup()</code>, as documented in the manual.
cannam@127 252 <h2><A name="allzero">
cannam@127 253 Question 3.17. The output of FFTW's transform is all
cannam@127 254 zeros.
cannam@127 255 </A></h2>
cannam@127 256
cannam@127 257 You should initialize your input array <i>after</i> creating the plan, unless you use <code>FFTW_ESTIMATE</code>: planning with <code>FFTW_MEASURE</code> or <code>FFTW_PATIENT</code> overwrites the input/output arrays, as described in the manual.
cannam@127 258
cannam@127 259 <h2><A name="vbetalia">
cannam@127 260 Question 3.18. How do I call FFTW from the Microsoft language du
cannam@127 261 jour?
cannam@127 262 </A></h2>
cannam@127 263
cannam@127 264 Please <i>do not</i> ask us Windows-specific questions. We do not
cannam@127 265 use Windows. We know nothing about Visual Basic, Visual C++, or .NET.
cannam@127 266 Please find the appropriate Usenet discussion group and ask your
cannam@127 267 question there. See also <A href="section2.html#runOnWindows">Q2.2 `Does FFTW run on Windows?'</A>.
cannam@127 268 <h2><A name="pruned">
cannam@127 269 Question 3.19. Can I compute only a subset of the DFT
cannam@127 270 outputs?
cannam@127 271 </A></h2>
cannam@127 272
cannam@127 273 In general, no, an FFT intrinsically computes all outputs from all
cannam@127 274 inputs. In principle, there is something called a
cannam@127 275 <i>pruned FFT</i> that can do what you want, but to compute K outputs out of N the
cannam@127 276 complexity is in general O(N log K) instead of O(N log N), thus saving
cannam@127 277 only a small additive factor in the log. (The same argument holds if
cannam@127 278 you instead have only K nonzero inputs.)
cannam@127 279
cannam@127 280 <p>
cannam@127 281 There are some specific cases in which you can get the O(N log K)
cannam@127 282 performance benefits easily, however, by combining a few ordinary
cannam@127 283 FFTs. In particular, the case where you want the first K outputs,
cannam@127 284 where K divides N, can be handled by performing N/K transforms of size
cannam@127 285 K and then summing the outputs multiplied by appropriate phase
cannam@127 286 factors. For more details, see <A href="http://www.fftw.org/pruned.html">pruned FFTs with FFTW</A>.
cannam@127 287 <p>
cannam@127 288 There are also some algorithms that compute pruned transforms
cannam@127 289 <i>approximately</i>, but they are beyond the scope of this FAQ.
cannam@127 290
cannam@127 291 <h2><A name="transpose">
cannam@127 292 Question 3.20. Can I use FFTW's routines for in-place and
cannam@127 293 out-of-place matrix transposition?
cannam@127 294 </A></h2>
cannam@127 295
cannam@127 296 You can use the FFTW guru interface to create a rank-0 transform of
cannam@127 297 vector rank 2 where the vector strides are transposed. (A rank-0
cannam@127 298 transform is equivalent to a 1D transform of size 1, which. just
cannam@127 299 copies the input into the output.) Specifying the same location for
cannam@127 300 the input and output makes the transpose in-place.
cannam@127 301
cannam@127 302 <p>
cannam@127 303 For double-valued data stored in row-major format, plan creation looks
cannam@127 304 like this: <pre>
cannam@127 305 fftw_plan plan_transpose(int rows, int cols, double *in, double *out)
cannam@127 306 {
cannam@127 307 const unsigned flags = FFTW_ESTIMATE; /* other flags are possible */
cannam@127 308 fftw_iodim howmany_dims[2];
cannam@127 309
cannam@127 310 howmany_dims[0].n = rows;
cannam@127 311 howmany_dims[0].is = cols;
cannam@127 312 howmany_dims[0].os = 1;
cannam@127 313
cannam@127 314 howmany_dims[1].n = cols;
cannam@127 315 howmany_dims[1].is = 1;
cannam@127 316 howmany_dims[1].os = rows;
cannam@127 317
cannam@127 318 return fftw_plan_guru_r2r(/*rank=*/ 0, /*dims=*/ NULL,
cannam@127 319 /*howmany_rank=*/ 2, howmany_dims,
cannam@127 320 in, out, /*kind=*/ NULL, flags);
cannam@127 321 }
cannam@127 322 </pre>
cannam@127 323 (This entry was written by Rhys Ulerich.)
cannam@127 324 <hr>
cannam@127 325 Next: <a href="section4.html" rel=precedes>Internals of FFTW</a>.<br>
cannam@127 326 Back: <a href="section2.html" rev=precedes>Installing FFTW</a>.<br>
cannam@127 327 <a href="index.html" rev=subdocument>Return to contents</a>.<p>
cannam@127 328 <address>
cannam@127 329 <A href="http://www.fftw.org">Matteo Frigo and Steven G. Johnson</A> / <A href="mailto:fftw@fftw.org">fftw@fftw.org</A>
cannam@127 330 - 30 July 2016
cannam@127 331 </address><br>
cannam@127 332 Extracted from FFTW Frequently Asked Questions with Answers,
cannam@127 333 Copyright &copy; 2016 Matteo Frigo and Massachusetts Institute of Technology.
cannam@127 334 </body></html>