pyfftw
- The core¶
The core of pyfftw
consists of the FFTW
class,
wisdom functions and a couple of
utility functions for dealing with aligned
arrays.
This module represents the full interface to the underlying FFTW
library. However, users may find it easier to
use the helper routines provided in pyfftw.builders
.
FFTW Class¶
-
class
pyfftw.
FFTW
(input_array, output_array, axes=(-1, ), direction='FFTW_FORWARD', flags=('FFTW_MEASURE', ), threads=1, planning_timelimit=None)¶ FFTW is a class for computing the complex N-Dimensional DFT or inverse DFT of an array using the FFTW library. The interface is designed to be somewhat pythonic, with the correct transform being inferred from the dtypes of the passed arrays.
On instantiation, the dtypes and relative shapes of the input array and output arrays are compared to the set of valid (and implemented) FFTW schemes. If a match is found, the plan that corresponds to that scheme is created, operating on the arrays that are passed in. If no scheme can be created, then
ValueError
is raised.The actual FFT or iFFT is performed by calling the
execute()
method.The arrays can be updated by calling the
update_arrays()
method.The created instance of the class is itself callable, and can perform the execution of the FFT, both with or without array updates, returning the result of the FFT. Unlike calling the
execute()
method, calling the class instance will also optionally normalise the output as necessary. Additionally, calling with an input array update will also coerce that array to be the correct dtype.See the documentation on the
__call__()
method for more information.Arguments:
input_array
andoutput_array
should be numpy arrays. The contents of these arrays will be destroyed by the planning process during initialisation. Information on supported dtypes for the arrays is given below.axes
describes along which axes the DFT should be taken. This should be a valid list of axes. Repeated axes are only transformed once. Invalid axes will raise anIndexError
exception. This argument is equivalent to the same argument innumpy.fft.fftn()
, except for the fact that the behaviour of repeated axes is different (numpy.fft
will happily take the fft of the same axis if it is repeated in theaxes
argument). Rudimentary testing has suggested this is down to the underlying FFTW library and so unlikely to be fixed in these wrappers.direction
should be a string and one of'FFTW_FORWARD'
or'FFTW_BACKWARD'
, which dictate whether to take the DFT (forwards) or the inverse DFT (backwards) respectively (specifically, it dictates the sign of the exponent in the DFT formulation).Note that only the Complex schemes allow a free choice for
direction
. The direction must agree with the the table below if a Real scheme is used, otherwise aValueError
is raised.
flags
is a list of strings and is a subset of the flags that FFTW allows for the planners:'FFTW_ESTIMATE'
,'FFTW_MEASURE'
,'FFTW_PATIENT'
and'FFTW_EXHAUSTIVE'
are supported. These describe the increasing amount of effort spent during the planning stage to create the fastest possible transform. Usually'FFTW_MEASURE'
is a good compromise. If no flag is passed, the default'FFTW_MEASURE'
is used.'FFTW_UNALIGNED'
is supported. This tells FFTW not to assume anything about the alignment of the data and disabling any SIMD capability (see below).'FFTW_DESTROY_INPUT'
is supported. This tells FFTW that the input array can be destroyed during the transform, sometimes allowing a faster algorithm to be used. The default behaviour is, if possible, to preserve the input. In the case of the 1D Backwards Real transform, this may result in a performance hit. In the case of a backwards real transform for greater than one dimension, it is not possible to preserve the input, making this flag implicit in that case. A little more on this is given below.'FFTW_WISDOM_ONLY'
is supported. This tells FFTW to raise an error if no plan for this transform and data type is already in the wisdom. It thus provides a method to determine whether planning would require additional effor or the cached wisdom can be used. This flag should be combined with the various planning-effort flags ('FFTW_ESTIMATE'
,'FFTW_MEASURE'
, etc.); if so, then an error will be raised if wisdom derived from that level of planning effort (or higher) is not present. If no planning-effort flag is used, the default of'FFTW_ESTIMATE'
is assumed. Note that wisdom is specific to all the parameters, including the data alignment. That is, if wisdom was generated with input/output arrays with one specific alignment, using'FFTW_WISDOM_ONLY'
to create a plan for arrays with any different alignment will cause the'FFTW_WISDOM_ONLY'
planning to fail. Thus it is important to specifically control the data alignment to make the best use of'FFTW_WISDOM_ONLY'
.
The FFTW planner flags documentation has more information about the various flags and their impact. Note that only the flags documented here are supported.
threads
tells the wrapper how many threads to use when invoking FFTW, with a default of 1.planning_timelimit
is a floating point number that indicates to the underlying FFTW planner the maximum number of seconds it should spend planning the FFT. This is a rough estimate and corresponds to calling offftw_set_timelimit()
(or an equivalent dependent on type) in the underlying FFTW library. IfNone
is set, the planner will run indefinitely until all the planning modes allowed by the flags have been tried. See the FFTW planner flags page for more information on this.
Schemes
The currently supported schemes are as follows:
Type input_array.dtype
output_array.dtype
Direction Complex complex64
complex64
Both Complex complex128
complex128
Both Complex clongdouble
clongdouble
Both Real float32
complex64
Forwards Real float64
complex128
Forwards Real longdouble
clongdouble
Forwards Real1 complex64
float32
Backwards Real1 complex128
float64
Backwards Real1 clongdouble
longdouble
Backwards 1 Note that the Backwards Real transform for the case in which the dimensionality of the transform is greater than 1 will destroy the input array. This is inherent to FFTW and the only general work-around for this is to copy the array prior to performing the transform. In the case where the dimensionality of the transform is 1, the default is to preserve the input array. This is different from the default in the underlying library, and some speed gain may be achieved by allowing the input array to be destroyed by passing the
'FFTW_DESTROY_INPUT'
flag.clongdouble
typically maps directly tocomplex256
orcomplex192
, andlongdouble
tofloat128
orfloat96
, dependent on platform.The relative shapes of the arrays should be as follows:
- For a Complex transform,
output_array.shape == input_array.shape
- For a Real transform in the Forwards direction, both the following
should be true:
output_array.shape[axes][-1] == input_array.shape[axes][-1]//2 + 1
- All the other axes should be equal in length.
- For a Real transform in the Backwards direction, both the following
should be true:
input_array.shape[axes][-1] == output_array.shape[axes][-1]//2 + 1
- All the other axes should be equal in length.
In the above expressions for the Real transform, the
axes
arguments denotes the unique set of axes on which we are taking the FFT, in the order passed. It is the last of these axes that is subject to the special case shown.The shapes for the real transforms corresponds to those stipulated by the FFTW library. Further information can be found in the FFTW documentation on the real DFT.
The actual arrangement in memory is arbitrary and the scheme can be planned for any set of strides on either the input or the output. The user should not have to worry about this and any valid numpy array should work just fine.
What is calculated is exactly what FFTW calculates. Notably, this is an unnormalized transform so should be scaled as necessary (fft followed by ifft will scale the input by N, the product of the dimensions along which the DFT is taken). For further information, see the FFTW documentation.
The FFTW library benefits greatly from the beginning of each DFT axes being aligned on the correct byte boundary, enabling SIMD instructions. By default, if the data begins on such a boundary, then FFTW will be allowed to try and enable SIMD instructions. This means that all future changes to the data arrays will be checked for similar alignment. SIMD instructions can be explicitly disabled by setting the FFTW_UNALIGNED flags, to allow for updates with unaligned data.
byte_align()
andempty_aligned()
are two methods included with this module for producing aligned arrays.The optimum alignment for the running platform is provided by
pyfftw.simd_alignment
, though a different alignment may still result in some performance improvement. For example, if the processor supports AVX (requiring 32-byte alignment) as well as SSE (requiring 16-byte alignment), then if the array is 16-byte aligned, SSE will still be used.It’s worth noting that just being aligned may not be sufficient to create the fastest possible transform. For example, if the array is not contiguous (i.e. certain axes are displaced in memory), it may be faster to plan a transform for a contiguous array, and then rely on the array being copied in before the transform (which
pyfftw.FFTW
will handle for you when accessed through__call__()
).-
N
¶ The product of the lengths of the DFT over all DFT axes. 1/N is the normalisation constant. For any input array A, and for any set of axes, 1/N * ifft(fft(A)) = A
-
simd_aligned
¶ Return whether or not this FFTW object requires simd aligned input and output data.
-
input_alignment
¶ Returns the byte alignment of the input arrays for which the
FFTW
object was created.Input array updates with arrays that are not aligned on this byte boundary will result in a ValueError being raised, or a copy being made if the
__call__()
interface is used.
-
output_alignment
¶ Returns the byte alignment of the output arrays for which the
FFTW
object was created.Output array updates with arrays that are not aligned on this byte boundary will result in a ValueError being raised.
-
flags
¶ Return which flags were used to construct the FFTW object.
This includes flags that were added during initialisation.
-
input_array
¶ Return the input array that is associated with the FFTW instance.
-
output_array
¶ Return the output array that is associated with the FFTW instance.
-
input_shape
¶ Return the shape of the input array for which the FFT is planned.
-
output_shape
¶ Return the shape of the output array for which the FFT is planned.
-
input_strides
¶ Return the strides of the input array for which the FFT is planned.
-
output_strides
¶ Return the strides of the output array for which the FFT is planned.
-
input_dtype
¶ Return the dtype of the input array for which the FFT is planned.
-
output_dtype
¶ Return the shape of the output array for which the FFT is planned.
-
direction
¶ Return the planned FFT direction. Either ‘FFTW_FORWARD’ or ‘FFTW_BACKWARD’.
-
axes
¶ Return the axes for the planned FFT in canonical form. That is, as a tuple of positive integers. The order in which they were passed is maintained.
-
ortho
¶ If
ortho=True
both the forward and inverse transforms are scaled by 1/sqrt(N).
-
normalise_idft
¶ If
normalise_idft=True
, the inverse transform is scaled by 1/N.
-
__call__
()¶ - __call__(input_array=None, output_array=None, normalise_idft=True,
- ortho=False)
Calling the class instance (optionally) updates the arrays, then calls
execute()
, before optionally normalising the output and returning the output array.It has some built-in helpers to make life simpler for the calling functions (as distinct from manually updating the arrays and calling
execute()
).If
normalise_idft
isTrue
(the default), then the output from an inverse DFT (i.e. when the direction flag is'FFTW_BACKWARD'
) is scaled by 1/N, where N is the product of the lengths of input array on which the FFT is taken. If the direction is'FFTW_FORWARD'
, this flag makes no difference to the output array.If
ortho
isTrue
, then the output of both forward and inverse DFT operations is scaled by 1/sqrt(N), where N is the product of the lengths of input array on which the FFT is taken. This ensures that the DFT is a unitary operation, meaning that it satisfies Parseval’s theorem (the sum of the squared values of the transform output is equal to the sum of the squared values of the input). In other words, the energy of the signal is preserved.If either
normalise_idft
orortho
areTrue
, then ifft(fft(A)) = A.When
input_array
is something other than None, then the passed in array is coerced to be the same dtype as the input array used when the class was instantiated, the byte-alignment of the passed in array is made consistent with the expected byte-alignment and the striding is made consistent with the expected striding. All this may, but not necessarily, require a copy to be made.As noted in the scheme table, if the FFTW instance describes a backwards real transform of more than one dimension, the contents of the input array will be destroyed. It is up to the calling function to make a copy if it is necessary to maintain the input array.
output_array
is always used as-is if possible. If the dtype, the alignment or the striding is incorrect for the FFTW object, then aValueError
is raised.The coerced input array and the output array (as appropriate) are then passed as arguments to
update_arrays()
, after whichexecute()
is called, and then normalisation is applied to the output array if that is desired.Note that it is possible to pass some data structure that can be converted to an array, such as a list, so long as it fits the data requirements of the class instance, such as array shape.
Other than the dtype and the alignment of the passed in arrays, the rest of the requirements on the arrays mandated by
update_arrays()
are enforced.A
None
argument to either keyword means that that array is not updated.The result of the FFT is returned. This is the same array that is used internally and will be overwritten again on subsequent calls. If you need the data to persist longer than a subsequent call, you should copy the returned array.
-
update_arrays
(new_input_array, new_output_array)¶ Update the arrays upon which the DFT is taken.
The new arrays should be of the same dtypes as the originals, the same shapes as the originals and should have the same strides between axes. If the original data was aligned so as to allow SIMD instructions (e.g. by being aligned on a 16-byte boundary), then the new array must also be aligned so as to allow SIMD instructions (assuming, of course, that the
FFTW_UNALIGNED
flag was not enabled).The byte alignment requirement extends to requiring natural alignment in the non-SIMD cases as well, but this is much less stringent as it simply means avoiding arrays shifted by, say, a single byte (which invariably takes some effort to achieve!).
If all these conditions are not met, a
ValueError
will be raised and the data will not be updated (though the object will still be in a sane state).
-
execute
()¶ Execute the planned operation, taking the correct kind of FFT of the input array (i.e.
FFTW.input_array
), and putting the result in the output array (i.e.FFTW.output_array
).
-
get_input_array
()¶ Return the input array that is associated with the FFTW instance.
Deprecated since 0.10. Consider using the
FFTW.input_array
property instead.
-
get_output_array
()¶ Return the output array that is associated with the FFTW instance.
Deprecated since 0.10. Consider using the
FFTW.output_array
property instead.
Wisdom Functions¶
Functions for dealing with FFTW’s ability to export and restore plans, referred to as wisdom. For further information, refer to the FFTW wisdom documentation.
-
pyfftw.
export_wisdom
()¶ Return the FFTW wisdom as a tuple of strings.
The first string in the tuple is the string for the double precision wisdom. The second string in the tuple is the string for the single precision wisdom. The third string in the tuple is the string for the long double precision wisdom.
The tuple that is returned from this function can be used as the argument to
import_wisdom()
.
-
pyfftw.
import_wisdom
(wisdom)¶ Function that imports wisdom from the passed tuple of strings.
The first string in the tuple is the string for the double precision wisdom. The second string in the tuple is the string for the single precision wisdom. The third string in the tuple is the string for the long double precision wisdom.
The tuple that is returned from
export_wisdom()
can be used as the argument to this function.This function returns a tuple of boolean values indicating the success of loading each of the wisdom types (double, float and long double, in that order).
-
pyfftw.
forget_wisdom
()¶ Forget all the accumulated wisdom.
Utility Functions¶
-
pyfftw.
simd_alignment
¶ An integer giving the optimum SIMD alignment in bytes, found by inspecting the CPU (e.g. if AVX is supported, its value will be 32).
This can be used as
n
in the arguments forbyte_align()
,empty_aligned()
,zeros_aligned()
, andones_aligned()
to create optimally aligned arrays for the running platform.
-
pyfftw.
byte_align
(array, n=None, dtype=None)¶ Function that takes a numpy array and checks it is aligned on an n-byte boundary, where
n
is an optional parameter. Ifn
is not provided then this function will inspect the CPU to determine alignment. If the array is aligned then it is returned without further ado. If it is not aligned then a new array is created and the data copied in, but aligned on the n-byte boundary.dtype
is an optional argument that forces the resultant array to be of that dtype.
-
pyfftw.
empty_aligned
(shape, dtype='float64', order='C', n=None)¶ Function that returns an empty numpy array that is n-byte aligned, where
n
is determined by inspecting the CPU if it is not provided.The alignment is given by the final optional argument,
n
. Ifn
is not provided then this function will inspect the CPU to determine alignment. The rest of the arguments are as pernumpy.empty()
.
-
pyfftw.
zeros_aligned
(shape, dtype='float64', order='C', n=None)¶ Function that returns a numpy array of zeros that is n-byte aligned, where
n
is determined by inspecting the CPU if it is not provided.The alignment is given by the final optional argument,
n
. Ifn
is not provided then this function will inspect the CPU to determine alignment. The rest of the arguments are as pernumpy.zeros()
.
-
pyfftw.
ones_aligned
(shape, dtype='float64', order='C', n=None)¶ Function that returns a numpy array of ones that is n-byte aligned, where
n
is determined by inspecting the CPU if it is not provided.The alignment is given by the final optional argument,
n
. Ifn
is not provided then this function will inspect the CPU to determine alignment. The rest of the arguments are as pernumpy.ones()
.
-
pyfftw.
is_byte_aligned
()¶ is_n_byte_aligned(array, n=None)
Function that takes a numpy array and checks it is aligned on an n-byte boundary, where
n
is an optional parameter, returningTrue
if it is, andFalse
if it is not. Ifn
is not provided then this function will inspect the CPU to determine alignment.
-
pyfftw.
n_byte_align
(array, n, dtype=None)¶ This function is deprecated:
byte_align
should be used instead.Function that takes a numpy array and checks it is aligned on an n-byte boundary, where
n
is an optional parameter. Ifn
is not provided then this function will inspect the CPU to determine alignment. If the array is aligned then it is returned without further ado. If it is not aligned then a new array is created and the data copied in, but aligned on the n-byte boundary.dtype
is an optional argument that forces the resultant array to be of that dtype.
-
pyfftw.
n_byte_align_empty
(shape, n, dtype='float64', order='C')¶ This function is deprecated:
empty_aligned
should be used instead.Function that returns an empty numpy array that is n-byte aligned.
The alignment is given by the first optional argument,
n
. Ifn
is not provided then this function will inspect the CPU to determine alignment. The rest of the arguments are as pernumpy.empty()
.
-
pyfftw.
is_n_byte_aligned
(array, n)¶ This function is deprecated:
is_byte_aligned
should be used instead.Function that takes a numpy array and checks it is aligned on an n-byte boundary, where
n
is a passed parameter, returningTrue
if it is, andFalse
if it is not.
-
pyfftw.
next_fast_len
(target)¶ Find the next fast transform length for FFTW.
FFTW has efficient functions for transforms of length 2**a * 3**b * 5**c * 7**d * 11**e * 13**f, where e + f is either 0 or 1.
Parameters: target (int) – Length to start searching from. Must be a positive integer. Returns: out – The first fast length greater than or equal to target. Return type: int Examples
On a particular machine, an FFT of prime length takes 2.1 ms:
>>> from pyfftw.interfaces import scipy_fftpack >>> min_len = 10007 # prime length is worst case for speed >>> a = numpy.random.randn(min_len) >>> b = scipy_fftpack.fft(a)
Zero-padding to the next fast length reduces computation time to 406 us, a speedup of ~5 times:
>>> next_fast_len(min_len) 10080 >>> b = scipy_fftpack.fft(a, 10080)
Rounding up to the next power of 2 is not optimal, taking 598 us to compute, 1.5 times as long as the size selected by next_fast_len.
>>> b = fftpack.fft(a, 16384)
Similar speedups will occur for pre-planned FFTs as generated via pyfftw.builders.