Fourier Transforms (baseband_tasks.fourier
)¶
Introduction¶
The Fourier transform module contains classes that wrap various fast Fourier
transform (FFT) packages, in particular numpy.fft
and pyfftw.FFTW
. The
purpose of the module is to give the packages a common interface, and to allow
individual transforms to be defined once, then re-used multiple times. This is
especially useful for FFTW, which achieves its fast transforms through prior
planning.
The module currently does not support Hermitian Fourier transforms - frequency-domain values are always treated as complex.
Using the Fourier Module¶
To make FFTs, easiest is to use the fft_maker
factory. In our examples, we will use it with the numpy fft back-end:
>>> from baseband_tasks.fourier import fft_maker
>>> fft_maker.set('numpy')
<ScienceState fft_maker: NumpyFFTMaker()>
The set()
method allows one to
choose any of the FFT maker classes -
e.g. NumpyFFTMaker
or
PyfftwFFTMaker
. Package-level options, such
as the flags to FFTW
, can be passed as **kwargs
.
To create a transform, we pass the time-dimension data array shape and dtype, transform direction (‘forward’ or ‘backward’), transform axis (if the data is multi-dimensional), normalization convention and sample rate:
>>> import numpy as np
>>> import astropy.units as u
>>> fft = fft_maker((1000,), 'float64', direction='forward', ortho=True,
... sample_rate=1.*u.kHz)
Here, we have chosen orthogonal normalization, which normalizes both the
frequency and time-domain outputs by \(1 / n^{1/2}\), where \(n\) is
the length of the time-domain array. We can now perform the transform by
calling fft
:
>>> y = np.sin(2. * np.pi * np.arange(1000))
>>> Y = fft(y)
Y
is the Fourier transform of y
. To obtain the inverse, we use the
inverse
method in fft
:
>>> ifft = fft.inverse()
>>> y_copy = y.copy()
>>> yn = ifft(Y)
>>> np.allclose(y, y_copy)
True
Note that we compare to a copy of the input; if possible for a given
Fourier implementation (e.g., in pyfftw` but not in `numpy`), the ``inverse
implementation reuses input and output arrays of the forward transform to
save memory, so at the end on would have yn is y
.
To show information about the transform, we can simply print the instance:
>>> fft
<NumpyFFT direction=forward,
axis=0, ortho=True, sample_rate=1.0 kHz
Time domain: shape=(1000,), dtype=float64
Frequency domain: shape=(501,), dtype=complex128>
To obtain the sample frequencies of Y
, we can use the
frequency
property:
>>> fft.frequency[:10]
<Quantity [0. , 0.001, 0.002, 0.003, 0.004, 0.005, 0.006, 0.007, 0.008,
0.009] kHz>
For multi-dimensional arrays, the sample frequencies are for the transformed axis.
Reference/API¶
baseband_tasks.fourier Package¶
Fourier transform module.
This module provides a standard interface to various FFT packages.
Routine listings¶
fft_maker
: primary interface for creating FFT instances.
Implementation Notes¶
For each packages, there is a corresponding *FFTMaker
class, which
holds default information needed for creating an FFT instance. For
instance, for PyfftwFFTMaker
, this holds flags
, threads
, etc.
These *FFTMaker
instances in turn can be used to create *FFT
instances which are set up to do the FFT on data with a given shape, in
a given direction, etc.
Classes¶
FFT factory class utilizing |
|
|
FFT factory class utilizing the |
|
Create an FFT, with a settable default engine. |
Class Inheritance Diagram¶
baseband_tasks.fourier.base Module¶
Base classes and tools for the fourier module.
Implementation Notes¶
The base classes provide common code for adding new FFT engines.
In particular, the FFTMakerBase
class is subclassed to create
*FFTMaker
classes (where *
stands for a package such as pyfftw
or numpy
). These classes are initialized with default information
needed for creating an FFT instance. For instance, for
PyfftwFFTMaker
, this holds flags
, threads
, etc. Via the
FFTMakerMeta
meta class, all such maker classes are registered in the
FFT_MAKER_CLASSES
dict, keyed by a lower-case version of the name
(with fftmaker
removed).
These *FFTMaker
instances in turn can be used, via their
__call__
method, to create *FFT
subclass instances, setting
relevant attributes such as time_shape
, time_dtype
, axis
,
etc., and instantiating them for a given direction.
The *FFT
classes themselves are most easily based on FFTBase
,
which defines properties for accessing the various attributes, a default
__init__
method that allows on to create an instance for a given
direction of the FFT, and a FFTBase.__call__()
method for actually
doing the FFT on data. The FFTBase
class also defines convenience
properties, such as FFTBase.frequency
to get the sample frequencies,
and FFTBase.inverse()
for getting the FFT in the inverse
direction.
Selection of a default FFT package is done via fft_maker
, which stores
a default *FFTMaker
instance. It is based on
astropy.utils.state.ScienceState
, but adds a fft_maker.system_default
factory that one can access by setting the state to None
, as well as
__new__
method which allows one to use the default *FFTMaker
instance to create an *FFT
instance.
Classes¶
Base class for all FFT factories. |
|
|
Framework for single pre-defined FFT and its associated metadata. |
|
Create an FFT, with a settable default engine. |
|
Registry of FFT maker classes. |
Variables¶
Dict for storing FFT maker classes, indexed by their name or prefix. |
Class Inheritance Diagram¶
baseband_tasks.fourier.numpy Module¶
FFT maker and class using the numpy.fft
routines.
Classes¶
|
Single pre-defined FFT based on |
FFT factory class utilizing |
Class Inheritance Diagram¶
baseband_tasks.fourier.pyfftw Module¶
FFT maker and class using pyfftw
routines.
Implementation Notes¶
The code for PyfftwFFTBase
is relatively complex to ensure that the
input and output arrays are re-used, even between the forward and
backward transforms (if created using PyfftwFFTBase.inverse()
)
Classes¶
|
Single pre-defined FFT based on |
|
FFT factory class utilizing the |