GNU Info

Info Node: (fftw.info)rfftwnd

(fftw.info)rfftwnd


Next: Array Dimensions for Real Multi-dimensional Transforms Prev: rfftwnd_create_plan Up: Real Multi-dimensional Transforms Reference
Enter node , (file) or (file)node

Computing the Real Multi-dimensional Transform
----------------------------------------------

     #include <rfftw.h>
     
     void rfftwnd_real_to_complex(rfftwnd_plan plan, int howmany,
                                  fftw_real *in, int istride, int idist,
                                  fftw_complex *out, int ostride, int odist);
     void rfftwnd_complex_to_real(rfftwnd_plan plan, int howmany,
                                  fftw_complex *in, int istride, int idist,
                                  fftw_real *out, int ostride, int odist);
     
     void rfftwnd_one_real_to_complex(rfftwnd_plan p, fftw_real *in,
                                      fftw_complex *out);
     void rfftwnd_one_complex_to_real(rfftwnd_plan p, fftw_complex *in,
                                      fftw_real *out);

   These functions compute the real multi-dimensional Fourier Transform,
using a plan created by `rfftwnd_create_plan' (Note: Plan Creation for
Real Multi-dimensional Transforms.). (Note that
the plan determines the rank and dimensions of the array to be
transformed.)  The ``rfftwnd_one_'' functions provide a simplified
interface for the common case of single input array of stride 1.
Unlike other transform routines in FFTW, we here use separate functions
for the two directions of the transform in order to correctly express
the datatypes of the parameters.

   *Important:* When invoked for an out-of-place,
`FFTW_COMPLEX_TO_REAL' transform with `rank > 1', the input array is
overwritten with scratch values by these routines.  The input array is
not modified for `FFTW_REAL_TO_COMPLEX' transforms or for
`FFTW_COMPLEX_TO_REAL' with `rank == 1'.

Arguments
.........

   * `plan' is the plan created by `rfftwnd_create_plan'.  (*note Plan
     Creation for Real Multi-dimensional Transforms:
     rfftwnd_create_plan.). In the case of two and three-dimensional
     transforms, it could also have been created by
     `rfftw2d_create_plan' or `rfftw3d_create_plan', respectively.

     `FFTW_REAL_TO_COMPLEX' plans must be used with the
     ``real_to_complex'' functions, and `FFTW_COMPLEX_TO_REAL' plans
     must be used with the ``complex_to_real'' functions.  It is an
     error to mismatch the plan direction and the transform function.

   * `howmany' is the number of transforms to be computed.

   * `in', `istride' and `idist' describe the input array(s).  There
     are `howmany' input arrays; the first one is pointed to by `in',
     the second one is pointed to by `in + idist', and so on, up to `in
     + (howmany - 1) * idist'.  Each input array is stored in row-major
     format (Note: Multi-dimensional Array Format.), and is not
     necessarily contiguous in memory.  Specifically, `in[0]' is the
     first element of the first array, `in[istride]' is the second
     element of the first array, and so on.  In general, the `i'-th
     element of the `j'-th input array will be in position `in[i *
     istride + j * idist]'. Note that, here, `i' refers to an index into
     the row-major format for the multi-dimensional array, rather than
     an index in any particular dimension.

     The dimensions of the arrays are different for real and complex
     data, and are discussed in more detail below (Note: Array
     Dimensions for Real Multi-dimensional Transforms.).

        - *In-place transforms*: For plans created with the
          `FFTW_IN_PLACE' option, the transform is computed
          in-place--the output is returned in the `in' array.  The
          meaning of the `stride' and `dist' parameters in this case is
          subtle and is discussed below (Note: Strides in In-place
          RFFTWND.).

   * `out', `ostride' and `odist' describe the output array(s).  The
     format is the same as that for the input array.  See below for a
     discussion of the dimensions of the output array for real and
     complex data.

        - *In-place transforms*: These parameters are ignored for plans
          created with the `FFTW_IN_PLACE' option.

   The function `rfftwnd_one' transforms a single, contiguous input
array to a contiguous output array.  By definition, the call
     rfftwnd_one_...(plan, in, out)
   is equivalent to
     rfftwnd_...(plan, 1, in, 1, 1, out, 1, 1)


automatically generated by info2www version 1.2.2.9