Showing with 81 additions and 54 deletions.
  1. +79 −52 pyfftw.pyx
  2. +2 −2 setup.py
131 changes: 79 additions & 52 deletions pyfftw.pyx
View file
Expand Up @@ -19,10 +19,12 @@
import numpy as np
cimport numpy as np
from libc.stdlib cimport malloc, free
from libc cimport limits

from pyfftw cimport _fftw_iodim, fftw_iodim, fftwf_plan
cimport pyfftw


directions = {'FFTW_FORWARD': FFTW_FORWARD,
'FFTW_BACKWARD': FFTW_BACKWARD}

Expand Down Expand Up @@ -118,37 +120,49 @@ cdef void _fftwl_destroy_plan(void *_plan):

fftwl_destroy_plan(<fftwl_plan>_plan)

# fftw_schemes is a dictionary with a mapping from keys to tuples
# that take the format:
# (planner index, executor index, destroyer index,
# input dtype, output dtype, array validator function)
# fftw_schemes is a dictionary with a mapping from a keys,
# which are a tuple of the string representation of numpy
# dtypes to a scheme name.
#
# scheme_functions is a dictionary of functions, either
# an index to the array of functions in the case of
# 'planner', 'executor' and 'destroyer' or a callable
# in the case of 'validator'
#
# The indices refer to the relevant functions for each scheme,
# The array indices refer to the relevant functions for each scheme,
# the tables to which are defined below.
#
# array validator function is a callable that has the following
# signature:
# The 'validator' function is a callable for validating the arrays
# that has the following signature:
# bool callable(in_array, out_array)
# and checks that the arrays a valid pair. If it is set to None,
# and checks that the arrays are a valid pair. If it is set to None,
# then the default check is applied, which confirms that the arrays
# have the same shape.
fftw_schemes = {
'complex128': (0, 0, 0, 'complex128', 'complex128', None),
'complex64': (1, 1, 1, 'complex64', 'complex64', None),
'complex_longdouble': (2, 2, 2, 'clongdouble', 'clongdouble', None)}
(np.dtype('complex128'), np.dtype('complex128')): 'complex128',
(np.dtype('complex64'), np.dtype('complex64')): 'complex64',
(np.dtype('clongdouble'), np.dtype('clongdouble')): 'complex_longdouble'}

scheme_functions = {
'complex128': {'planner': 0, 'executor':0,
'destroyer':0, 'validator':None},
'complex64': {'planner':1, 'executor':1,
'destroyer':1, 'validator':None},
'complex_longdouble': {'planner':2, 'executor':2,
'destroyer':2, 'validator':None}}

# Planner table (with the same number of elements as there are schemes).
cdef fftw_generic_plan_guru planners[3]

cdef fftw_generic_plan_guru * _build_planner_list():

planners[fftw_schemes['complex128'][0]] = \
planners[scheme_functions['complex128']['planner']] = \
<fftw_generic_plan_guru>&_fftw_plan_guru_dft

planners[fftw_schemes['complex64'][0]] = \
planners[scheme_functions['complex64']['planner']] = \
<fftw_generic_plan_guru>&_fftwf_plan_guru_dft

planners[fftw_schemes['complex_longdouble'][0]] = \
planners[scheme_functions['complex_longdouble']['planner']] = \
<fftw_generic_plan_guru>&_fftwl_plan_guru_dft

return planners
Expand All @@ -158,13 +172,13 @@ cdef fftw_generic_execute executors[3]

cdef fftw_generic_execute * _build_executor_list():

executors[fftw_schemes['complex128'][1]] = \
executors[scheme_functions['complex128']['executor']] = \
<fftw_generic_execute>&_fftw_execute_dft

executors[fftw_schemes['complex64'][1]] = \
executors[scheme_functions['complex64']['executor']] = \
<fftw_generic_execute>&_fftwf_execute_dft

executors[fftw_schemes['complex_longdouble'][1]] = \
executors[scheme_functions['complex_longdouble']['executor']] = \
<fftw_generic_execute>&_fftwl_execute_dft

return executors
Expand All @@ -174,13 +188,13 @@ cdef fftw_generic_destroy_plan destroyers[3]

cdef fftw_generic_destroy_plan * _build_destroyer_list():

destroyers[fftw_schemes['complex128'][2]] = \
destroyers[scheme_functions['complex128']['destroyer']] = \
<fftw_generic_destroy_plan>&_fftw_destroy_plan

destroyers[fftw_schemes['complex64'][2]] = \
destroyers[scheme_functions['complex64']['destroyer']] = \
<fftw_generic_destroy_plan>&_fftwf_destroy_plan

destroyers[fftw_schemes['complex_longdouble'][2]] = \
destroyers[scheme_functions['complex_longdouble']['destroyer']] = \
<fftw_generic_destroy_plan>&_fftwl_destroy_plan

return destroyers
Expand Down Expand Up @@ -243,25 +257,22 @@ cdef class ComplexFFTW:
self.__dims = NULL
self.__howmany_dims = NULL

cdef bytes scheme = str('none')

for each_scheme in fftw_schemes:
if input_array.dtype == fftw_schemes[each_scheme][3] and \
output_array.dtype == fftw_schemes[each_scheme][4]:
scheme = each_scheme
break

if (scheme == 'none' ):
try:
scheme = fftw_schemes[
(input_array.dtype, output_array.dtype)]
except KeyError:
raise ValueError('The output array and input array dtypes '+\
'do not correspond to a valid fftw scheme.')

if fftw_schemes[each_scheme][5] == None:

functions = scheme_functions[scheme]

if functions['validator'] == None:
if not (output_array.shape == input_array.shape):
raise ValueError('The output array should be the same '+\
'shape as the input array for the given array '+\
'dtypes.')
else:
if not fftw_schemes[each_scheme][5](input_array, output_array):
if not functions['validator'](input_array, output_array):
raise ValueError('The input array and output array are '+\
'invalid complementary shapes for their dtypes.')

Expand All @@ -287,9 +298,9 @@ cdef class ComplexFFTW:
cdef fftw_generic_destroy_plan *destroyers = _build_destroyer_list()
cdef fftw_generic_execute *executors = _build_executor_list()

self.__fftw_planner = planners[fftw_schemes[scheme][0]]
self.__fftw_execute = executors[fftw_schemes[scheme][1]]
self.__fftw_destroy = destroyers[fftw_schemes[scheme][2]]
self.__fftw_planner = planners[functions['planner']]
self.__fftw_execute = executors[functions['executor']]
self.__fftw_destroy = destroyers[functions['destroyer']]

self.__flags = 0
for each_flag in flags:
Expand Down Expand Up @@ -330,8 +341,27 @@ cdef class ComplexFFTW:
self.__output_strides = \
np.array(output_array.strides)/output_array.itemsize

# Fill in the stride and shape information
# Make sure that the arrays are not too big for fftw
cdef int i
for i in range(0, len(self.__input_shape)):
if self.__input_shape[i] >= <Py_ssize_t> limits.INT_MAX:
raise ValueError('Dimensions of the input array must be ' +\
'less than ', str(limits.INT_MAX))

if self.__input_strides[i] >= <Py_ssize_t> limits.INT_MAX:
raise ValueError('Strides of the input array must be ' +\
'less than ', str(limits.INT_MAX))

for i in range(0, len(self.__output_shape)):
if self.__output_shape[i] >= <Py_ssize_t> limits.INT_MAX:
raise ValueError('Dimensions of the output array must be ' +\
'less than ', str(limits.INT_MAX))

if self.__output_strides[i] >= <Py_ssize_t> limits.INT_MAX:
raise ValueError('Strides of the output array must be ' +\
'less than ', str(limits.INT_MAX))

# Fill in the stride and shape information
for i in range(0, self.__rank):
self.__dims[i]._n = self.__input_shape[_axes][i]
self.__dims[i]._is = self.__input_strides[_axes][i]
Expand All @@ -346,8 +376,8 @@ cdef class ComplexFFTW:
self.__plan = self.__fftw_planner(
self.__rank, <fftw_iodim *>self.__dims,
self.__howmany_rank, <fftw_iodim *>self.__howmany_dims,
<void *>self.__input_array.data,
<void *>self.__output_array.data,
<void *>np.PyArray_DATA(self.__input_array),
<void *>np.PyArray_DATA(self.__output_array),
self.__direction, self.__flags)

if self.__plan == NULL:
Expand All @@ -363,15 +393,15 @@ cdef class ComplexFFTW:
The currently supported schemes are as follows:
+-----------------+------------------+----------------+
| ``input_array`` | ``output_array`` | Relative shape |
+=================+==================+================+
| ``complex64`` | ``complex64`` | Equal size |
+-----------------+------------------+----------------+
| ``complex128`` | ``complex128`` | Equal size |
+-----------------+------------------+----------------+
| ``clongdouble`` | ``clongdouble`` | Equal size |
+-----------------+------------------+----------------+
+-----------------+------------------+---------------------------------------------+
| ``input_array`` | ``output_array`` | Requirements on array shapes |
+=================+==================+=============================================+
| ``complex64`` | ``complex64`` | ``input_array.shape == output_array.shape`` |
+-----------------+------------------+---------------------------------------------+
| ``complex128`` | ``complex128`` | ``input_array.shape == output_array.shape`` |
+-----------------+------------------+---------------------------------------------+
| ``clongdouble`` | ``clongdouble`` | ``input_array.shape == output_array.shape`` |
+-----------------+------------------+---------------------------------------------+
``clongdouble`` typically maps directly to ``complex256``
or ``complex192``, dependent on platform.
Expand All @@ -380,9 +410,6 @@ cdef class ComplexFFTW:
can be planned for any set of strides on either the input
or the output.
Equal size means the ``shape`` attribute of the arrays is
the same.
``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 an
Expand Down Expand Up @@ -517,8 +544,8 @@ cdef class ComplexFFTW:
Execute the planned operation.
'''
self.__fftw_execute(self.__plan,
<void *>self.__input_array.data,
<void *>self.__output_array.data)
<void *>np.PyArray_DATA(self.__input_array),
<void *>np.PyArray_DATA(self.__output_array))


cpdef n_byte_align_empty(shape, n, dtype='float64', order='C'):
Expand Down
4 changes: 2 additions & 2 deletions setup.py
View file
Expand Up @@ -44,10 +44,10 @@

setup_args = {
'name': 'pyFFTW',
'version': '0.5.0',
'version': '0.5.1',
'author': 'Henry Gomersall',
'author_email': 'heng@kedevelopments.co.uk',
'description': 'A pythonic wrapper around FFTW, presenting a unified interface for all the supported transforms.',
'description': 'A pythonic wrapper around FFTW, the FFT library, presenting a unified interface for all the supported transforms.',
'url': 'http://hgomersall.github.com/pyFFTW/',
'long_description': long_description,
'classifiers': [
Expand Down