fftw_3.html

FFTW, a collection of fast C routines to compute the Discrete Fourier Transform in one or more dime

<P>
Specifically, for a real transform of dimensions
n<sub>1</sub> x n<sub>2</sub> x ... x n<sub>d</sub>,
the complex data is an
n<sub>1</sub> x n<sub>2</sub> x ... x (n<sub>d</sub>/2+1)
array of <CODE>fftw_complex</CODE> values in row-major order (with the
division rounded down).  That is, we only store the lower half (plus one
element) of the last dimension of the data from the ordinary complex
transform.  (We could have instead taken half of any other dimension,
but implementation turns out to be simpler if the last, contiguous,
dimension is used.)


<P>
<A NAME="IDX178"></A>
<A NAME="IDX179"></A>
Since the complex data is slightly larger than the real data, some
complications arise for in-place transforms.  In this case, the final
dimension of the real data must be padded with extra values to
accommodate the size of the complex data--two extra if the last
dimension is even and one if it is odd.  That is, the last dimension of
the real data must physically contain
2 * (n<sub>d</sub>/2+1)
<CODE>fftw_real</CODE> values (exactly enough to hold the complex data).
This physical array size does not, however, change the <EM>logical</EM>
array size--only
n<sub>d</sub>
values are actually stored in the last dimension, and
n<sub>d</sub>
is the last dimension passed to <CODE>rfftwnd_create_plan</CODE>.




<H3><A NAME="SEC38">Strides in In-place RFFTWND</A></H3>

<P>
<A NAME="IDX180"></A>
<A NAME="IDX181"></A>
The fact that the input and output datatypes are different for rfftwnd
complicates the meaning of the <CODE>stride</CODE> and <CODE>dist</CODE> parameters
of in-place transforms--are they in units of <CODE>fftw_real</CODE> or
<CODE>fftw_complex</CODE> elements?  When reading the input, they are
interpreted in units of the datatype of the input data.  When writing
the output, the <CODE>istride</CODE> and <CODE>idist</CODE> are translated to the
output datatype's "units" in one of two ways, corresponding to the two
most common situations in which <CODE>stride</CODE> and <CODE>dist</CODE> parameters
are useful.  Below, we refer to these "translated" parameters as
<CODE>ostride_t</CODE> and <CODE>odist_t</CODE>.  (Note that these are computed
internally by rfftwnd; the actual <CODE>ostride</CODE> and <CODE>odist</CODE>
parameters are ignored for in-place transforms.)


<P>
First, there is the case where you are transforming a number of
contiguous arrays located one after another in memory.  In this
situation, <CODE>istride</CODE> is <CODE>1</CODE> and <CODE>idist</CODE> is the product of
the physical dimensions of the array.  <CODE>ostride_t</CODE> and
<CODE>odist_t</CODE> are then chosen so that the output arrays are contiguous
and lie on top of the input arrays.  <CODE>ostride_t</CODE> is therefore
<CODE>1</CODE>.  For a real-to-complex transform, <CODE>odist_t</CODE> is
<CODE>idist/2</CODE>; for a complex-to-real transform, <CODE>odist_t</CODE> is
<CODE>idist*2</CODE>.


<P>
The second case is when you have an array in which each element has
<CODE>nc</CODE> components (e.g. a structure with <CODE>nc</CODE> numeric fields),
and you want to transform all of the components at once.  Here,
<CODE>istride</CODE> is <CODE>nc</CODE> and <CODE>idist</CODE> is <CODE>1</CODE>.  For this
case, it is natural to want the output to also have <CODE>nc</CODE>
consecutive components, now of the output data type; this is exactly
what rfftwnd does.  Specifically, it uses an <CODE>ostride_t</CODE> equal to
<CODE>istride</CODE>, and an <CODE>odist_t</CODE> of <CODE>1</CODE>.  (Astute readers will
realize that some extra buffer space is required in order to perform
such a transform; this is handled automatically by rfftwnd.)


<P>
The general rule is as follows.  <CODE>ostride_t</CODE> equals <CODE>istride</CODE>.
If <CODE>idist</CODE> is <CODE>1</CODE> and <CODE>idist</CODE> is less than
<CODE>istride</CODE>, then <CODE>odist_t</CODE> is <CODE>1</CODE>.  Otherwise, for a
real-to-complex transform <CODE>odist_t</CODE> is <CODE>idist/2</CODE> and for a
complex-to-real transform <CODE>odist_t</CODE> is <CODE>idist*2</CODE>.




<H3><A NAME="SEC39">Destroying a Multi-dimensional Plan</A></H3>


<PRE>
#include &#60;rfftw.h&#62;

void rfftwnd_destroy_plan(rfftwnd_plan plan);
</PRE>

<P>
<A NAME="IDX182"></A>


<P>
The function <CODE>rfftwnd_destroy_plan</CODE> frees the plan <CODE>plan</CODE>
and releases all the memory associated with it.  After destruction,
a plan is no longer valid.




<H3><A NAME="SEC40">What RFFTWND Really Computes</A></H3>
<P>
<A NAME="IDX183"></A>


<P>
The conventions that we follow for the real multi-dimensional transform
are analogous to those for the complex multi-dimensional transform. In
particular, the forward transform has a negative sign in the exponent
and neither the forward nor the backward transforms will perform any
normalization.  Computing the backward transform of the forward
transform will multiply the array by the product of its dimensions (that
is, the logical dimensions of the real data).  The forward transform is
real-to-complex and the backward transform is complex-to-real.


<P>
<A NAME="IDX184"></A>
<A NAME="IDX185"></A>
The Gods forbade using HTML to display mathematical formulas.  Please
see the TeX or Postscript version of this manual for the proper
definition of the n-dimensional real Fourier transform that RFFTW
uses.  For completeness, we include a bitmap of the TeX output below:
<P><center><IMG SRC="equation-4.gif" ALIGN="top"></center>




<H2><A NAME="SEC41">Wisdom Reference</A></H2>

<P>
<A NAME="IDX186"></A>


<H3><A NAME="SEC42">Exporting Wisdom</A></H3>


<PRE>
#include &#60;fftw.h&#62;

void fftw_export_wisdom(void (*emitter)(char c, void *), void *data);
void fftw_export_wisdom_to_file(FILE *output_file);
char *fftw_export_wisdom_to_string(void);
</PRE>

<P>
<A NAME="IDX187"></A>
<A NAME="IDX188"></A>
<A NAME="IDX189"></A>


<P>
These functions allow you to export all currently accumulated
<CODE>wisdom</CODE> in a form from which it can be later imported and
restored, even during a separate run of the program. (See Section <A HREF="fftw_2.html#SEC13">Words of Wisdom</A>.)  The current store of <CODE>wisdom</CODE> is not
affected by calling any of these routines.


<P>
<CODE>fftw_export_wisdom</CODE> exports the <CODE>wisdom</CODE> to any output
medium, as specified by the callback function
<CODE>emitter</CODE>. <CODE>emitter</CODE> is a <CODE>putc</CODE>-like function that
writes the character <CODE>c</CODE> to some output; its second parameter is
the <CODE>data</CODE> pointer passed to <CODE>fftw_export_wisdom</CODE>.  For
convenience, the following two "wrapper" routines are provided:


<P>
<CODE>fftw_export_wisdom_to_file</CODE> writes the <CODE>wisdom</CODE> to the
current position in <CODE>output_file</CODE>, which should be open with write
permission.  Upon exit, the file remains open and is positioned at the
end of the <CODE>wisdom</CODE> data.


<P>
<CODE>fftw_export_wisdom_to_string</CODE> returns a pointer to a
<CODE>NULL</CODE>-terminated string holding the <CODE>wisdom</CODE> data. This
string is dynamically allocated, and it is the responsibility of the
caller to deallocate it with <CODE>fftw_free</CODE> when it is no longer
needed.


<P>
All of these routines export the wisdom in the same format, which we
will not document here except to say that it is LISP-like ASCII text
that is insensitive to white space.




<H3><A NAME="SEC43">Importing Wisdom</A></H3>


<PRE>
#include &#60;fftw.h&#62;

fftw_status fftw_import_wisdom(int (*get_input)(void *), void *data);
fftw_status fftw_import_wisdom_from_file(FILE *input_file);
fftw_status fftw_import_wisdom_from_string(const char *input_string);
</PRE>

<P>
<A NAME="IDX190"></A>
<A NAME="IDX191"></A>
<A NAME="IDX192"></A>


<P>
These functions import <CODE>wisdom</CODE> into a program from data stored by
the <CODE>fftw_export_wisdom</CODE> functions above. (See Section <A HREF="fftw_2.html#SEC13">Words of Wisdom</A>.)
The imported <CODE>wisdom</CODE> supplements rather than replaces any
<CODE>wisdom</CODE> already accumulated by the running program (except when
there is conflicting <CODE>wisdom</CODE>, in which case the existing wisdom is
replaced).


<P>
<CODE>fftw_import_wisdom</CODE> imports <CODE>wisdom</CODE> from any input medium,
as specified by the callback function <CODE>get_input</CODE>. <CODE>get_input</CODE>
is a <CODE>getc</CODE>-like function that returns the next character in the
input; its parameter is the <CODE>data</CODE> pointer passed to
<CODE>fftw_import_wisdom</CODE>. If the end of the input data is reached
(which should never happen for valid data), it may return either
<CODE>NULL</CODE> (ASCII 0) or <CODE>EOF</CODE> (as defined in <CODE>&#60;stdio.h&#62;</CODE>).
For convenience, the following two "wrapper" routines are provided:


<P>
<CODE>fftw_import_wisdom_from_file</CODE> reads <CODE>wisdom</CODE> from the
current position in <CODE>input_file</CODE>, which should be open with read
permission.  Upon exit, the file remains open and is positioned at the
end of the <CODE>wisdom</CODE> data.


<P>
<CODE>fftw_import_wisdom_from_string</CODE> reads <CODE>wisdom</CODE> from the
<CODE>NULL</CODE>-terminated string <CODE>input_string</CODE>.


<P>
The return value of these routines is <CODE>FFTW_SUCCESS</CODE> if the wisdom
was read successfully, and <CODE>FFTW_FAILURE</CODE> otherwise. Note that, in
all of these functions, any data in the input stream past the end of the
<CODE>wisdom</CODE> data is simply ignored (it is not even read if the
<CODE>wisdom</CODE> data is well-formed).




<H3><A NAME="SEC44">Forgetting Wisdom</A></H3>


<PRE>
#include &#60;fftw.h&#62;

void fftw_forget_wisdom(void);
</PRE>

<P>
<A NAME="IDX193"></A>


<P>
Calling <CODE>fftw_forget_wisdom</CODE> causes all accumulated <CODE>wisdom</CODE>
to be discarded and its associated memory to be freed. (New
<CODE>wisdom</CODE> can still be gathered subsequently, however.)




<H2><A NAME="SEC45">Memory Allocator Reference</A></H2>


<PRE>
#include &#60;fftw.h&#62;

void *(*fftw_malloc_hook) (size_t n);
void (*fftw_free_hook) (void *p);
</PRE>

<P>
<A NAME="IDX194"></A>
<A NAME="IDX195"></A>
<A NAME="IDX196"></A>
<A NAME="IDX197"></A>


<P>
Whenever it has to allocate and release memory, FFTW ordinarily calls
<CODE>malloc</CODE> and <CODE>free</CODE>.  
If <CODE>malloc</CODE> fails, FFTW prints an error message and exits.  This
behavior may be undesirable in some applications. Also, special
memory-handling functions may be necessary in certain
environments. Consequently, FFTW provides means by which you can install
your own memory allocator and take whatever error-correcting action you
find appropriate.  The variables <CODE>fftw_malloc_hook</CODE> and
<CODE>fftw_free_hook</CODE> are pointers to functions, and they are normally
<CODE>NULL</CODE>.  If you set those variables to point to other functions,
then FFTW will use your routines instead of <CODE>malloc</CODE> and
<CODE>free</CODE>.  <CODE>fftw_malloc_hook</CODE> must point to a <CODE>malloc</CODE>-like
function, and <CODE>fftw_free_hook</CODE> must point to a <CODE>free</CODE>-like
function.




<H2><A NAME="SEC46">Thread safety</A></H2>

<P>
<A NAME="IDX198"></A>
<A NAME="IDX199"></A>
Users writing multi-threaded programs must concern themselves with the
<EM>thread safety</EM> of the libraries they use--that is, whether it is
safe to call routines in parallel from multiple threads.  FFTW can be
used in such an environment, but some care must be taken because certain
parts of FFTW use private global variables to share data between calls.
In particular, the plan-creation functions share trigonometric tables
and accumulated <CODE>wisdom</CODE>.  (Users should note that these comments
only apply to programs using shared-memory threads.  Parallelism using
MPI or forked processes involves a separate address-space and global
variables for each process, and is not susceptible to problems of this
sort.)


<P>
The central restriction of FFTW is that it is not safe to create
multiple plans in parallel.  You must either create all of your plans
from a single thread, or instead use a semaphore, mutex, or other
mechanism to ensure that different threads don't attempt to create plans
at the same time.  The same restriction also holds for destruction of
plans and importing/forgetting <CODE>wisdom</CODE>.  Once created, a plan may
safely be used in any thread.


<P>
The actual transform routines in FFTW (<CODE>fftw_one</CODE>, etcetera) are
re-entrant and thread-safe, so it is fine to call them simultaneously
from multiple threads.  Another question arises, however--is it safe to
use the <EM>same plan</EM> for multiple transforms in parallel?  (It would
be unsafe if, for example, the plan were modified in some way by the
transform.)  We address this question by defining an additional planner
flag, <CODE>FFTW_THREADSAFE</CODE>.
<A NAME="IDX200"></A>
When included in the flags for any of the plan-creation routines,
<CODE>FFTW_THREADSAFE</CODE> guarantees that the resulting plan will be
read-only and safe to use in parallel by multiple threads.


<P><HR><P>
Go to the <A HREF="fftw_1.html">first</A>, <A HREF="fftw_2.html">previous</A>, <A HREF="fftw_4.html">next</A>, <A HREF="fftw_10.html">last</A> section, <A HREF="fftw_toc.html">table of contents</A>.
</BODY>
</HTML>