Manifest.hs

{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE TypeOperators #-}
module Data.Array.Accelerate.FFTW.Manifest (
   Transform,
   dft,
   idft,
   dft2d,
   idft2d,
   ) where

import qualified Data.Array.Accelerate.FFTW.Private as Batch
import qualified Data.Array.Accelerate as A
import Data.Array.Accelerate.FFTW.Private (Transform)
import Data.Array.Accelerate ((:.))

import Data.Array.Accelerate.Data.Complex (Complex)
import Foreign.Storable.Complex ()

import qualified Numeric.Netlib.Class as Class


{- $setup
>>> import qualified Data.Array.Accelerate.Utility.Lift.Exp as Exp
>>> import qualified Data.Array.Accelerate.Interpreter as AI
>>> import qualified Data.Array.Accelerate as A
>>> import Data.Array.Accelerate.Utility.Lift.Exp (expr)
>>> import Data.Array.Accelerate.Data.Complex (Complex((:+)))
>>> import Data.Array.Accelerate ((:.)((:.)))
>>>
>>> roundComplex :: A.Exp (Complex Float) -> A.Exp (Complex Int)
>>> roundComplex =
>>>    Exp.modify (expr:+expr) $ \(r:+i) -> A.round r :+ (A.round i :: A.Exp Int)
>>>
>>> roundComplexArray ::
>>>    (A.Shape sh) => A.Array sh (Complex Float) -> A.Array sh (Complex Int)
>>> roundComplexArray = AI.run1 (A.map roundComplex)
-}


{- |
>>> roundComplexArray $ dft (A.fromList (A.Z:.3:.3) [1000,0,0, 0,1000,0, 0,0,1000])
Matrix (Z :. 3 :. 3)
  [ 1000 :+ 0,        1000 :+ 0,        1000 :+ 0,
    1000 :+ 0, (-500) :+ (-866),    (-500) :+ 866,
    1000 :+ 0,    (-500) :+ 866, (-500) :+ (-866)]
-}
dft ::
   (A.Shape sh, Class.Real a) =>
   Transform (sh:.Int) (Complex a)
dft = Batch.runComplexTransform $ Class.switchReal Batch.dft Batch.dft

{- |
>>> roundComplexArray $ idft $ dft (A.fromList (A.Z:.3:.3) [1000,0,0, 0,1000,0, 0,0,1000])
Matrix (Z :. 3 :. 3)
  [ 3000 :+ 0,    0 :+ 0,    0 :+ 0,
       0 :+ 0, 3000 :+ 0,    0 :+ 0,
       0 :+ 0,    0 :+ 0, 3000 :+ 0]
-}
idft ::
   (A.Shape sh, Class.Real a) =>
   Transform (sh:.Int) (Complex a)
idft = Batch.runComplexTransform $ Class.switchReal Batch.idft Batch.idft


{- |
>>> roundComplexArray $ dft2d (A.fromList (A.Z:.3:.3) [1000,0,0, 0,0,0, 0,0,0])
Matrix (Z :. 3 :. 3)
  [ 1000 :+ 0, 1000 :+ 0, 1000 :+ 0,
    1000 :+ 0, 1000 :+ 0, 1000 :+ 0,
    1000 :+ 0, 1000 :+ 0, 1000 :+ 0]

>>> roundComplexArray $ dft2d (A.fromList (A.Z:.3:.3) [0,1000,0, 0,0,0, 0,0,0])
Matrix (Z :. 3 :. 3)
  [ 1000 :+ 0, (-500) :+ (-866), (-500) :+ 866,
    1000 :+ 0, (-500) :+ (-866), (-500) :+ 866,
    1000 :+ 0, (-500) :+ (-866), (-500) :+ 866]

>>> roundComplexArray $ dft2d (A.fromList (A.Z:.3:.3) [1000,0,0, 0,1000,0, 0,0,1000])
Matrix (Z :. 3 :. 3)
  [ 3000 :+ 0,    0 :+ 0,    0 :+ 0,
       0 :+ 0,    0 :+ 0, 3000 :+ 0,
       0 :+ 0, 3000 :+ 0,    0 :+ 0]
-}
dft2d ::
   (A.Shape sh, Class.Real a) =>
   Transform (sh:.Int:.Int) (Complex a)
dft2d = Batch.runComplexTransform $ Class.switchReal Batch.dft2d Batch.dft2d

{- |
>>> roundComplexArray $ idft2d $ dft2d (A.fromList (A.Z:.3:.3) [1000,0,0, 0,1000,0, 0,0,1000])
Matrix (Z :. 3 :. 3)
  [ 9000 :+ 0,    0 :+ 0,    0 :+ 0,
       0 :+ 0, 9000 :+ 0,    0 :+ 0,
       0 :+ 0,    0 :+ 0, 9000 :+ 0]
-}
idft2d ::
   (A.Shape sh, Class.Real a) =>
   Transform (sh:.Int:.Int) (Complex a)
idft2d = Batch.runComplexTransform $ Class.switchReal Batch.idft2d Batch.idft2d