DiscreteFourierTransformInverse

class odl.trafos.fourier.DiscreteFourierTransformInverse(*args, **kwargs)[source]

Bases: odl.trafos.fourier.DiscreteFourierTransformBase

Plain backward DFT, only evaluating the trigonometric sum.

This operator calculates the inverse DFT:

f[k] = 1/prod(N) * sum_j( f_hat[j] * exp(+- 1j*2*pi * j*k/N) )

without any further shifting or scaling compensation. See the Numpy FFT documentation, the pyfftw API documentation or What FFTW really computes for further information.

See also

DiscreteFourierTransform
FourierTransformInverse
numpy.fft.ifftn

n-dimensional inverse FFT routine

numpy.fft.irfftn

n-dimensional half-complex inverse FFT

odl.trafos.backends.pyfftw_bindings.pyfftw_call

apply an FFTW transform

References

Attributes
adjoint

Adjoint transform, equal to the inverse.

axes

Axes along the FT is calculated by this operator.

domain

Set of objects on which this operator can be evaluated.

halfcomplex

Return True if the last transform axis is halved.

impl

Backend for the FFT implementation.

inverse

Inverse Fourier transform.

is_functional

True if this operator’s range is a Field.

is_linear

True if this operator is linear.

range

Set in which the result of an evaluation of this operator lies.

sign

Sign of the complex exponent in the transform.

Methods

_call(self, x, out, \*\*kwargs)

Implement self(x, out[, **kwargs]).

clear_fftw_plan(self)

Delete the FFTW plan of this transform.

derivative(self, point)

Return the operator derivative at point.

init_fftw_plan(self[, planning_effort])

Initialize the FFTW plan for this transform for later use.

norm(self[, estimate])

Return the operator norm of this operator.

__init__(self, range, domain=None, axes=None, sign='+', halfcomplex=False, impl=None)[source]

Initialize a new instance.

Parameters
rangeDiscretizedSpace

Range of the inverse Fourier transform. If its DiscretizedSpace.exponent is equal to 2.0, this operator has an adjoint which is equal to the inverse.

domainDiscretizedSpace, optional

Domain of the inverse Fourier transform. If not given, the domain is determined from range and the other parameters as a uniform_discr with unit cell size and exponent p / (p - 1) (read as ‘inf’ for p=1 and 1 for p=’inf’).

axessequence of ints, optional

Dimensions in which a transform is to be calculated. None means all axes.

sign{‘-‘, ‘+’}, optional

Sign of the complex exponent.

halfcomplexbool, optional

If True, interpret the last axis in axes as the negative frequency part of the transform of a real signal and calculate a “half-complex-to-real” inverse FFT. In this case, the domain has by default the shape floor(N[i]/2) + 1 in this axis i. Otherwise, domain and range have the same shape. If range is a complex space, this option has no effect.

impl{‘numpy’, ‘pyfftw’}, optional

Backend for the FFT implementation. The ‘pyfftw’ backend is faster but requires the pyfftw package. None selects the fastest available backend.

Examples

Complex-to-complex (default) transforms have the same grids in domain and range:

>>> range = odl.uniform_discr([0, 0], [2, 4], (2, 4),
...                           nodes_on_bdry=True)
>>> ifft = DiscreteFourierTransformInverse(range)
>>> ifft.domain.shape
(2, 4)
>>> ifft.range.shape
(2, 4)

Complex-to-real transforms have a domain grid with shape n // 2 + 1 in the last tranform axis:

>>> range = odl.uniform_discr([0, 0, 0], [2, 3, 4], (2, 3, 4),
...                           nodes_on_bdry=True)
>>> axes = (0, 1)
>>> ifft = DiscreteFourierTransformInverse(
...     range, halfcomplex=True, axes=axes)
>>> ifft.domain.shape   # shortened in the second axis
(2, 2, 4)
>>> ifft.range.shape
(2, 3, 4)