-- Hoogle documentation, generated by Haddock
-- See Hoogle, http://www.haskell.org/hoogle/


-- | Bindings to the FFTW library.
--   
--   Bindings to the FFTW library.
--   
--   Provides high performance discrete fourier transforms in arbitrary
--   dimensions. Include transforms of complex data, real data, and real to
--   real transforms.
@package fft
@version 0.1.8.6

module Math.FFT.Base

-- | Types of transforms. Used to control <a>dftShape</a>.
data DFT
CC :: DFT
RC :: DFT
CR :: DFT
CRO :: DFT
RR :: DFT

-- | Tuple of transform dimensions and non-transform dimensions of the
--   array.
type TSpec = ([IODim], [IODim])

-- | Real to Real transform kinds.
data Kind
R2HC :: Kind
HC2R :: Kind
DHT :: Kind
REDFT00 :: Kind
REDFT10 :: Kind
REDFT01 :: Kind
REDFT11 :: Kind
RODFT00 :: Kind
RODFT01 :: Kind
RODFT10 :: Kind
RODFT11 :: Kind

-- | Determine which direction of DFT to execute.
data Sign
DFTForward :: Sign
DFTBackward :: Sign

-- | The <a>Flag</a> type is used to influence the kind of plans which are
--   created. To specify multiple flags, use a bitwise <a>.|.</a>.
newtype Flag
Flag :: FFTWFlag -> Flag
[unFlag] :: Flag -> FFTWFlag

-- | Our API is polymorphic over the real data type. FFTW, at least in
--   principle, supports single precision <a>Float</a>, double precision
--   <a>Double</a> and long double <tt>CLDouble</tt> (presumable?).
class (Storable a, RealFloat a) => FFTWReal a
plan_guru_dft :: FFTWReal a => CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr (Complex a) -> Ptr (Complex a) -> FFTWSign -> FFTWFlag -> IO Plan
plan_guru_dft_r2c :: FFTWReal a => CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr a -> Ptr (Complex a) -> FFTWFlag -> IO Plan
plan_guru_dft_c2r :: FFTWReal a => CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr (Complex a) -> Ptr a -> FFTWFlag -> IO Plan
plan_guru_r2r :: FFTWReal a => CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr a -> Ptr a -> Ptr FFTWKind -> FFTWFlag -> IO Plan

-- | This lock must be taken during <i>planning</i> of any transform. The
--   FFTW library is not thread-safe in the planning phase. Thankfully, the
--   lock is not needed during the execute phase.
lock :: MVar ()
withLock :: IO a -> IO a

-- | Default flag. For most transforms, this is equivalent to setting
--   <a>measure</a> and <a>preserveInput</a>. The exceptions are complex to
--   real and half-complex to real transforms.
nullFlag :: Flag

-- | Allows FFTW to overwrite the input array with arbitrary data; this can
--   sometimes allow more efficient algorithms to be employed.
--   
--   Setting this flag implies that two memory allocations will be done,
--   one for work space, and one for the result. When <a>estimate</a> is
--   not set, we will be doing two memory allocations anyway, so we set
--   this flag as well (since we don't retain the work array anyway).
destroyInput :: Flag

-- | <a>preserveInput</a> specifies that an out-of-place transform must not
--   change its input array. This is ordinarily the default, except for
--   complex to real transforms for which <a>destroyInput</a> is the
--   default. In the latter cases, passing <a>preserveInput</a> will
--   attempt to use algorithms that do not destroy the input, at the
--   expense of worse performance; for multi-dimensional complex to real
--   transforms, however, no input-preserving algorithms are implemented so
--   the Haskell bindings will set <a>destroyInput</a> and do a transform
--   with two memory allocations.
preserveInput :: Flag

-- | Instruct FFTW not to generate a plan which uses SIMD instructions,
--   even if the memory you are planning with is aligned. This should only
--   be needed if you are using the guru interface and want to reuse a plan
--   with memory that may be unaligned (i.e. you constructed the
--   <a>CArray</a> with <a>unsafeForeignPtrToCArray</a>).
unaligned :: Flag

-- | The header claims that this flag is documented, but in reality, it is
--   not. I don't know what it does and it is here only for completeness.
conserveMemory :: Flag

-- | <a>estimate</a> specifies that, instead of actual measurements of
--   different algorithms, a simple heuristic is used to pick a (probably
--   sub-optimal) plan quickly. With this flag, the input/output arrays are
--   not overwritten during planning.
--   
--   This is the only planner flag for which a single memory allocation is
--   possible.
estimate :: Flag

-- | <a>measure</a> tells FFTW to find an optimized plan by actually
--   computing several FFTs and measuring their execution time. Depending
--   on your machine, this can take some time (often a few seconds).
--   <a>measure</a> is the default planning option.
measure :: Flag

-- | <a>patient</a> is like <a>measure</a>, but considers a wider range of
--   algorithms and often produces a "more optimal" plan (especially for
--   large transforms), but at the expense of several times longer planning
--   time (especially for large transforms).
patient :: Flag

-- | <a>exhaustive</a> is like <a>patient</a> but considers an even wider
--   range of algorithms, including many that we think are unlikely to be
--   fast, to produce the most optimal plan but with a substantially
--   increased planning time.
exhaustive :: Flag
unSign :: Sign -> FFTWSign
unKind :: Kind -> FFTWKind

-- | Verify that a plan is valid. Throws an exception if not.
check :: Plan -> IO ()

-- | Confirm that the plan is valid, then execute the transform.
execute :: Plan -> IO ()

-- | In-place normalization outside of IO. You must be able to prove that
--   no reference to the original can be retained.
unsafeNormalize :: (Ix i, Shapable i, Fractional e, Storable e) => [Int] -> CArray i e -> CArray i e

-- | Normalized general complex DFT
dftG :: (FFTWReal r, Ix i, Shapable i) => Sign -> Flag -> [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | Normalized general complex to real DFT where the last transformed
--   dimension is logically even.
dftCRG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Normalized general complex to real DFT where the last transformed
--   dimension is logicall odd.
dftCROG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Multi-dimensional forward DFT.
dftN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | Multi-dimensional inverse DFT.
idftN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | Multi-dimensional forward DFT of real data.
dftRCN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i (Complex r)

-- | Multi-dimensional inverse DFT of Hermitian-symmetric data (where only
--   the non-negative frequencies are given).
dftCRN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i r

-- | Multi-dimensional inverse DFT of Hermitian-symmetric data (where only
--   the non-negative frequencies are given) and the last transformed
--   dimension is logically odd.
dftCRON :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i r
fzr :: b -> [a] -> [(a, b)]
drr :: (FFTWReal r, Ix i, Shapable i) => Kind -> [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional real to real transform. The result is not
--   normalized.
dftRRN :: (FFTWReal r, Ix i, Shapable i) => [(Int, Kind)] -> CArray i r -> CArray i r

-- | Multi-dimensional real to half-complex transform. The result is not
--   normalized.
dftRHN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional half-complex to real transform. The result is not
--   normalized.
dftHRN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Discrete Hartley Transform. The result is not
--   normalized.
dhtN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 1 discrete cosine transform.
dct1N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 2 discrete cosine transform. This is commonly
--   known as <i>the</i> DCT.
dct2N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 3 discrete cosine transform. This is commonly
--   known as <i>the</i> inverse DCT. The result is not normalized.
dct3N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 4 discrete cosine transform.
dct4N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 1 discrete sine transform.
dst1N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 2 discrete sine transform.
dst2N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 3 discrete sine transform.
dst3N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 4 discrete sine transform.
dst4N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | 1-dimensional complex DFT.
dft :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i (Complex r)

-- | 1-dimensional complex inverse DFT. Inverse of <a>dft</a>.
idft :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i (Complex r)

-- | 1-dimensional real to complex DFT.
dftRC :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i (Complex r)

-- | 1-dimensional complex to real DFT with logically even dimension.
--   Inverse of <a>dftRC</a>.
dftCR :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i r

-- | 1-dimensional complex to real DFT with logically odd dimension.
--   Inverse of <a>dftRC</a>.
dftCRO :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i r

-- | 1-dimensional real to half-complex DFT.
dftRH :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional half-complex to real DFT. Inverse of <a>dftRH</a> after
--   normalization.
dftHR :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Discrete Hartley Transform. Self-inverse after
--   normalization.
dht :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 1 discrete cosine transform.
dct1 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 2 discrete cosine transform. This is commonly known
--   as <i>the</i> DCT.
dct2 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 3 discrete cosine transform. This is commonly known
--   as <i>the</i> inverse DCT.
dct3 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 4 discrete cosine transform.
dct4 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 1 discrete sine transform.
dst1 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 2 discrete sine transform.
dst2 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 3 discrete sine transform.
dst3 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 4 discrete sine transform.
dst4 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r
has :: Flag -> Flag -> Bool
infix 7 `has`

-- | Try to transform a CArray with only one memory allocation (for the
--   result). If we can find a way to prove that FFTW already has a
--   sufficiently good plan for this transform size and the input will not
--   be overwritten, then we could call have a version of this that does
--   not require <a>estimate</a>. Since this is not currently the case, we
--   require <a>estimate</a> to be set. Note that we do not check for the
--   <a>preserveInput</a> flag here. This is because the default is to
--   preserve input for all but the C-&gt;R and HC-&gt;R transforms.
--   Therefore, this function must not be called for those transforms,
--   unless <a>preserveInput</a> is set.
transformCArray :: (Ix i, Storable a, Storable b) => Flag -> CArray i a -> (i, i) -> (FFTWFlag -> Ptr a -> Ptr b -> IO Plan) -> CArray i b

-- | Transform a CArray with two memory allocations. This is entirely safe
--   with all transforms, but it must allocate a temporary array to do the
--   planning in.
transformCArray' :: (Ix i, Storable a, Storable b) => Flag -> CArray i a -> (i, i) -> (FFTWFlag -> Ptr a -> Ptr b -> IO Plan) -> CArray i b

-- | All the logic for determining shape of resulting array, and how to do
--   the transform.
dftShape :: (Ix i, Shapable i, Storable e) => DFT -> [Int] -> CArray i e -> ((i, i), TSpec)

-- | A simple helper.
withTSpec :: TSpec -> (CInt -> Ptr IODim -> CInt -> Ptr IODim -> IO a) -> IO a

-- | A generally useful list utility
adjust :: (a -> a) -> Int -> [a] -> [a]

-- | Complex to Complex DFT, un-normalized.
dftGU :: (FFTWReal r, Ix i, Shapable i) => Sign -> Flag -> [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | Real to Complex DFT.
dftRCG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i r -> CArray i (Complex r)

-- | Complex to Real DFT. The first argument determines whether the last
--   transformed dimension is logically odd or even. <a>True</a> implies
--   the dimension is odd.
dftCRG_ :: (FFTWReal r, Ix i, Shapable i) => Bool -> Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Complex to Real DFT where last transformed dimension is logically
--   even.
dftCRGU :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Complex to Real DFT where last transformed dimension is logically odd.
dftCROGU :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Real to Real transforms.
dftRRG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [(Int, Kind)] -> CArray i r -> CArray i r

-- | Queries the FFTW cache. The <a>String</a> can be written to a file so
--   the wisdom can be reused on a subsequent run.
exportWisdomString :: IO String

-- | Add wisdom to the FFTW cache. Returns <a>True</a> if it is successful.
importWisdomString :: String -> IO Bool

-- | Tries to import wisdom from a global source, typically
--   <tt><i>etc</i>fftw/wisdom</tt>. Returns <a>True</a> if it was
--   successful.
importWisdomSystem :: IO Bool
type FFTWFlag = CUInt
c_measure :: FFTWFlag
c_destroy_input :: FFTWFlag
c_unaligned :: FFTWFlag
c_conserve_memory :: FFTWFlag
c_exhaustive :: FFTWFlag
c_preserve_input :: FFTWFlag
type FFTWSign = CInt
c_patient :: FFTWFlag
c_forward :: FFTWSign
c_estimate :: FFTWFlag
c_backward :: FFTWSign
type FFTWKind = CInt
c_r2hc :: FFTWKind
c_hc2r :: FFTWKind
c_dht :: FFTWKind
c_redft00 :: FFTWKind
c_redft10 :: FFTWKind
c_redft01 :: FFTWKind
c_redft11 :: FFTWKind
c_rodft00 :: FFTWKind

-- | Corresponds to the <tt>fftw_iodim</tt> structure. It completely
--   describes the layout of each dimension, before and after the
--   transform.
c_rodft10 :: FFTWKind
data IODim
IODim :: Int -> Int -> Int -> IODim

-- | Logical size of dimension
[nIODim] :: IODim -> Int

-- | Stride along dimension in input array
[isIODim] :: IODim -> Int

-- | Stride along dimension in output array
[osIODim] :: IODim -> Int
c_rodft01 :: FFTWKind
c_rodft11 :: FFTWKind

-- | A plan is an opaque foreign object.
type Plan = Ptr FFTWPlan
type FFTWPlan = ()

-- | Plan a complex to complex transform using the guru interface.
cf_plan_guru_dft :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr (Complex Float) -> Ptr (Complex Float) -> FFTWSign -> FFTWFlag -> IO Plan

-- | Plan a real to complex transform using the guru interface.
cf_plan_guru_dft_r2c :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr Float -> Ptr (Complex Float) -> FFTWFlag -> IO Plan

-- | Plan a complex to real transform using the guru interface.
cf_plan_guru_dft_c2r :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr (Complex Float) -> Ptr Float -> FFTWFlag -> IO Plan

-- | Plan a real to real transform using the guru interface.
cf_plan_guru_r2r :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr Float -> Ptr Float -> Ptr FFTWKind -> FFTWFlag -> IO Plan

-- | Plan a complex to complex transform using the guru interface.
c_plan_guru_dft :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr (Complex Double) -> Ptr (Complex Double) -> FFTWSign -> FFTWFlag -> IO Plan

-- | Plan a real to complex transform using the guru interface.
c_plan_guru_dft_r2c :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr Double -> Ptr (Complex Double) -> FFTWFlag -> IO Plan

-- | Plan a complex to real transform using the guru interface.
c_plan_guru_dft_c2r :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr (Complex Double) -> Ptr Double -> FFTWFlag -> IO Plan

-- | Plan a real to real transform using the guru interface.
c_plan_guru_r2r :: CInt -> Ptr IODim -> CInt -> Ptr IODim -> Ptr Double -> Ptr Double -> Ptr FFTWKind -> FFTWFlag -> IO Plan

-- | Simple plan execution
c_execute :: Plan -> IO ()
c_execute_dft :: Plan -> Ptr (Complex Double) -> Ptr (Complex Double) -> IO ()
c_execute_dft_r2c :: Plan -> Ptr Double -> Ptr (Complex Double) -> IO ()
c_execute_dft_c2r :: Plan -> Ptr (Complex Double) -> Ptr Double -> IO ()
c_execute_r2r :: Plan -> Ptr Double -> Ptr Double -> IO ()
c_export_wisdom_string :: IO CString
c_import_wisdom_string :: CString -> IO CInt
c_import_wisdom_system :: IO CInt

-- | Frees memory allocated by <tt>fftw_malloc</tt>. Currently, we only
--   need this to free the wisdom string.
c_free :: Ptr a -> IO ()
instance GHC.Show.Show Math.FFT.Base.DFT
instance GHC.Classes.Eq Math.FFT.Base.DFT
instance GHC.Show.Show Math.FFT.Base.Kind
instance GHC.Classes.Eq Math.FFT.Base.Kind
instance GHC.Show.Show Math.FFT.Base.Sign
instance GHC.Classes.Eq Math.FFT.Base.Sign
instance Data.Bits.Bits Math.FFT.Base.Flag
instance GHC.Num.Num Math.FFT.Base.Flag
instance GHC.Show.Show Math.FFT.Base.Flag
instance GHC.Classes.Eq Math.FFT.Base.Flag
instance Math.FFT.Base.FFTWReal GHC.Types.Float
instance Math.FFT.Base.FFTWReal GHC.Types.Double


-- | This module exposes an interface to FFTW, the Fastest Fourier
--   Transform in the West.
--   
--   These bindings present several levels of interface. All the higher
--   level functions (<a>dft</a>, <a>idft</a>, <a>dftN</a>, ...) are easily
--   derived from the general functions (<a>dftG</a>, <a>dftRCG</a>, ...).
--   Only the general functions let you specify planner flags. The higher
--   levels all set <a>estimate</a> so you should not have to wait through
--   time consuming planning (see below for more).
--   
--   The simplest interface is the one-dimensional transforms. If you
--   supply a multi-dimensional array, these will only transform the first
--   dimension. These functions only take one argument, the array to be
--   transformed.
--   
--   At the next level, we have multi-dimensional transforms where you
--   specify which dimensions to transform in and the array to transform.
--   For instance
--   
--   <pre>
--   b = dftRCN [0,2] a
--   </pre>
--   
--   is the real to complex transform in dimensions 0 and 2 of the array
--   <tt>a</tt> which must be at least rank 3. The array <tt>b</tt> will be
--   complex valued with the same extent as <tt>a</tt> in every dimension
--   except <tt>2</tt>. If <tt>a</tt> had extent <tt>n</tt> in dimension
--   <tt>2</tt> then the <tt>b</tt> will have extent <tt>a <a>div</a> 2 +
--   1</tt> which consists of all non-negative frequency components in this
--   dimension (the negative frequencies are conjugate to the positive
--   frequencies because of symmetry since <tt>a</tt> is real valued).
--   
--   The real to real transforms allow different transform kinds in each
--   transformed dimension. For example,
--   
--   <pre>
--   b = dftRRN [(0,DHT), (1,REDFT10), (2,RODFT11)] a
--   </pre>
--   
--   is a Discrete Hartley Transform in dimension 0, a discrete cosine
--   transform (DCT-2) in dimension 1, and distrete sine transform (DST-4)
--   in dimension 2 where the array <tt>a</tt> must have rank at least 3.
--   
--   The general interface is similar to the multi-dimensional interface,
--   takes as its first argument, a bitwise <tt>.|.</tt> of planning
--   <a>Flag</a>s. (In the complex version, the sign of the transform is
--   first.) For example,
--   
--   <pre>
--   b = dftG DFTBackward (patient .|. destroy_input) [1,2] a
--   </pre>
--   
--   is an inverse DFT in dimensions 1 and 2 of the complex array
--   <tt>a</tt> which has rank at least 3. It will use the patient planner
--   to generate a (near) optimal transform. If you compute the same type
--   of transform again, it should be very fast since the plan is cached.
--   
--   Inverse transforms are typically normalized. The un-normalized inverse
--   transforms are <a>dftGU</a>, <a>dftCRGU</a> and <a>dftCROGU</a>. For
--   example
--   
--   <pre>
--   b = dftCROGU measure [0,1] a
--   </pre>
--   
--   is an un-normalized inverse DFT in dimensions 0 and 1 of the complex
--   array <tt>a</tt> (representing the non-negative frequencies, where the
--   negative frequencies are conjugate) which has rank at least 2. Here,
--   dimension 1 is logically odd so if <tt>a</tt> has extent <tt>n</tt> in
--   dimension 1, then <tt>b</tt> will have extent <tt>(n - 1) * 2 + 1</tt>
--   in dimension 1. It is more common that the logical dimension is even,
--   in which case we would use <a>dftCRGU</a> in which case <tt>b</tt>
--   would have extent <tt>(n - 1) * 2</tt> in dimension <tt>1</tt>.
--   
--   The FFTW library separates transforms into two steps. First you
--   compute a plan for a given transform, then you execute it. Often the
--   planning stage is quite time-consuming, but subsequent transforms of
--   the same size and type will be extremely fast. The planning phase
--   actually computes transforms, so it overwrites its input array. For
--   many C codes, it is reasonable to re-use the same arrays to compute a
--   given transform on different data. This is not a very useful paradigm
--   from Haskell. Fortunately, FFTW caches its plans so if try to generate
--   a new plan for a transform size which has already been planned, the
--   planner will return immediately. Unfortunately, it is not possible to
--   consult the cache, so if a plan is cached, we may use more memory than
--   is strictly necessary since we must allocate a work array which we
--   expect to be overwritten during planning. FFTW can export its cached
--   plans to a string. This is known as wisdom. For high performance work,
--   it is a good idea to compute plans of the sizes you are interested in,
--   using aggressive options (i.e. <a>patient</a>), use
--   <a>exportWisdomString</a> to get a string representing these plans,
--   and write this to a file. Then for production runs, you can read this
--   in, then add it to the cache with <a>importWisdomString</a>. Now you
--   can use the <a>estimate</a> planner so the Haskell bindings know that
--   FFTW will not overwrite the input array, and you will still get a high
--   quality transform (because it has wisdom).
module Math.FFT

-- | Determine which direction of DFT to execute.
data Sign
DFTForward :: Sign
DFTBackward :: Sign

-- | Real to Real transform kinds.
data Kind
R2HC :: Kind
HC2R :: Kind
DHT :: Kind
REDFT00 :: Kind
REDFT10 :: Kind
REDFT01 :: Kind
REDFT11 :: Kind
RODFT00 :: Kind
RODFT01 :: Kind
RODFT10 :: Kind
RODFT11 :: Kind

-- | Allows FFTW to overwrite the input array with arbitrary data; this can
--   sometimes allow more efficient algorithms to be employed.
--   
--   Setting this flag implies that two memory allocations will be done,
--   one for work space, and one for the result. When <a>estimate</a> is
--   not set, we will be doing two memory allocations anyway, so we set
--   this flag as well (since we don't retain the work array anyway).
destroyInput :: Flag

-- | <a>preserveInput</a> specifies that an out-of-place transform must not
--   change its input array. This is ordinarily the default, except for
--   complex to real transforms for which <a>destroyInput</a> is the
--   default. In the latter cases, passing <a>preserveInput</a> will
--   attempt to use algorithms that do not destroy the input, at the
--   expense of worse performance; for multi-dimensional complex to real
--   transforms, however, no input-preserving algorithms are implemented so
--   the Haskell bindings will set <a>destroyInput</a> and do a transform
--   with two memory allocations.
preserveInput :: Flag

-- | <a>estimate</a> specifies that, instead of actual measurements of
--   different algorithms, a simple heuristic is used to pick a (probably
--   sub-optimal) plan quickly. With this flag, the input/output arrays are
--   not overwritten during planning.
--   
--   This is the only planner flag for which a single memory allocation is
--   possible.
estimate :: Flag

-- | <a>measure</a> tells FFTW to find an optimized plan by actually
--   computing several FFTs and measuring their execution time. Depending
--   on your machine, this can take some time (often a few seconds).
--   <a>measure</a> is the default planning option.
measure :: Flag

-- | <a>patient</a> is like <a>measure</a>, but considers a wider range of
--   algorithms and often produces a "more optimal" plan (especially for
--   large transforms), but at the expense of several times longer planning
--   time (especially for large transforms).
patient :: Flag

-- | <a>exhaustive</a> is like <a>patient</a> but considers an even wider
--   range of algorithms, including many that we think are unlikely to be
--   fast, to produce the most optimal plan but with a substantially
--   increased planning time.
exhaustive :: Flag

-- | 1-dimensional complex DFT.
dft :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i (Complex r)

-- | 1-dimensional complex inverse DFT. Inverse of <a>dft</a>.
idft :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i (Complex r)

-- | Multi-dimensional forward DFT.
dftN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | Multi-dimensional inverse DFT.
idftN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | Normalized general complex DFT
dftG :: (FFTWReal r, Ix i, Shapable i) => Sign -> Flag -> [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | Complex to Complex DFT, un-normalized.
dftGU :: (FFTWReal r, Ix i, Shapable i) => Sign -> Flag -> [Int] -> CArray i (Complex r) -> CArray i (Complex r)

-- | 1-dimensional real to complex DFT.
dftRC :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i (Complex r)

-- | 1-dimensional complex to real DFT with logically even dimension.
--   Inverse of <a>dftRC</a>.
dftCR :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i r

-- | 1-dimensional complex to real DFT with logically odd dimension.
--   Inverse of <a>dftRC</a>.
dftCRO :: (FFTWReal r, Ix i, Shapable i) => CArray i (Complex r) -> CArray i r

-- | Multi-dimensional forward DFT of real data.
dftRCN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i (Complex r)

-- | Multi-dimensional inverse DFT of Hermitian-symmetric data (where only
--   the non-negative frequencies are given).
dftCRN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i r

-- | Multi-dimensional inverse DFT of Hermitian-symmetric data (where only
--   the non-negative frequencies are given) and the last transformed
--   dimension is logically odd.
dftCRON :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i (Complex r) -> CArray i r

-- | Real to Complex DFT.
dftRCG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i r -> CArray i (Complex r)

-- | Normalized general complex to real DFT where the last transformed
--   dimension is logically even.
dftCRG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Normalized general complex to real DFT where the last transformed
--   dimension is logicall odd.
dftCROG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Complex to Real DFT where last transformed dimension is logically
--   even.
dftCRGU :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | Complex to Real DFT where last transformed dimension is logically odd.
dftCROGU :: (FFTWReal r, Ix i, Shapable i) => Flag -> [Int] -> CArray i (Complex r) -> CArray i r

-- | 1-dimensional real to half-complex DFT.
dftRH :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional half-complex to real DFT. Inverse of <a>dftRH</a> after
--   normalization.
dftHR :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Discrete Hartley Transform. Self-inverse after
--   normalization.
dht :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 1 discrete cosine transform.
dct1 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 2 discrete cosine transform. This is commonly known
--   as <i>the</i> DCT.
dct2 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 3 discrete cosine transform. This is commonly known
--   as <i>the</i> inverse DCT.
dct3 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 4 discrete cosine transform.
dct4 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 1 discrete sine transform.
dst1 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 2 discrete sine transform.
dst2 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 3 discrete sine transform.
dst3 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | 1-dimensional Type 4 discrete sine transform.
dst4 :: (FFTWReal r, Ix i, Shapable i) => CArray i r -> CArray i r

-- | Multi-dimensional real to half-complex transform. The result is not
--   normalized.
dftRHN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional half-complex to real transform. The result is not
--   normalized.
dftHRN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Discrete Hartley Transform. The result is not
--   normalized.
dhtN :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 1 discrete cosine transform.
dct1N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 2 discrete cosine transform. This is commonly
--   known as <i>the</i> DCT.
dct2N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 3 discrete cosine transform. This is commonly
--   known as <i>the</i> inverse DCT. The result is not normalized.
dct3N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 4 discrete cosine transform.
dct4N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 1 discrete sine transform.
dst1N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 2 discrete sine transform.
dst2N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 3 discrete sine transform.
dst3N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional Type 4 discrete sine transform.
dst4N :: (FFTWReal r, Ix i, Shapable i) => [Int] -> CArray i r -> CArray i r

-- | Multi-dimensional real to real transform. The result is not
--   normalized.
dftRRN :: (FFTWReal r, Ix i, Shapable i) => [(Int, Kind)] -> CArray i r -> CArray i r

-- | Real to Real transforms.
dftRRG :: (FFTWReal r, Ix i, Shapable i) => Flag -> [(Int, Kind)] -> CArray i r -> CArray i r

-- | Add wisdom to the FFTW cache. Returns <a>True</a> if it is successful.
importWisdomString :: String -> IO Bool

-- | Tries to import wisdom from a global source, typically
--   <tt><i>etc</i>fftw/wisdom</tt>. Returns <a>True</a> if it was
--   successful.
importWisdomSystem :: IO Bool

-- | Queries the FFTW cache. The <a>String</a> can be written to a file so
--   the wisdom can be reused on a subsequent run.
exportWisdomString :: IO String
