ifftCore :: Data Index -> Pull1 (Complex Float) -> Pull1 (Complex Float)
ifftCore n = map (/ complex (i2f (2^(n+1))) 0) . composeOn stage (reverse (0...n))
  where
    stage k vec = indexed1 (length vec) ixf
      where
        ixf i = condition (testBit i k) (twid * (b - a)) (a+b)
          where
            a    = vec !! i
            b    = vec !! (i `xor` k2)
            twid = cis (pi * i2f (lsbs k i) / i2f k2)
            k2   = 1 .<<. k

bitRev :: Type a => Data Index -> Pull1 a -> Pull1 a
bitRev n = composeOn riffle (1...n)

riffle :: Syntax a => Data Index -> Pull DIM1 a -> Pull DIM1 a
riffle k = permute (const $ rotBit k)

-- Helper functions
composeOn :: (Syntax a) => (b -> a -> a) -> Pull DIM1 b -> a -> a
composeOn f v i = fromZero $ fold (flip f) i v

rotBit :: Data Index -> Data Index -> Data Index
rotBit 0 _ = P.error "rotBit: k should be at least 1"
rotBit k i = lefts .|. rights
  where
    ir = i .>>. 1
    rights = ir .&. oneBits k
    lefts  = (((ir .>>. k) .<<. 1) .|. (i .&. 1)) .<<. k
ifftCore :: Data Index -> Pull1 (Complex Float) -> Pull1 (Complex Float)
ifftCore n = map (/ complex (i2f (2^(n+1))) 0) . composeOn stage (reverse (0...n))
  where
    stage k vec = indexed1 (length vec) ixf
      where
        ixf i = condition (testBit i k) (twid * (b - a)) (a+b)
          where
            a    = vec !! i
            b    = vec !! (i `xor` k2)
            twid = cis (pi * i2f (lsbs k i) / i2f k2)
            k2   = 1 .<<. k

bitRev :: Type a => Data Index -> Pull1 a -> Pull1 a
bitRev n = composeOn riffle (1...n)

riffle :: Syntax a => Data Index -> Pull DIM1 a -> Pull DIM1 a
riffle k = permute (const $ rotBit k)

-- Helper functions
composeOn :: (Syntax a) => (b -> a -> a) -> Pull DIM1 b -> a -> a
composeOn f v i = fromZero $ fold (flip f) i v

rotBit :: Data Index -> Data Index -> Data Index
rotBit 0 _ = P.error "rotBit: k should be at least 1"
rotBit k i = lefts .|. rights
  where
    ir = i .>>. 1
    rights = ir .&. oneBits k
    lefts  = (((ir .>>. k) .<<. 1) .|. (i .&. 1)) .<<. k
ifftCore n = map (/ complex (i2f (2^(n+1))) 0) . composeOn stage (reverse (0...n))
  where
    stage k vec = indexed1 (length vec) ixf
      where
        ixf i = condition (testBit i k) (twid * (b - a)) (a+b)
          where
            a    = vec !! i
            b    = vec !! (i `xor` k2)
            twid = cis (pi * i2f (lsbs k i) / i2f k2)
            k2   = 1 .<<. k

bitRev :: Type a => Data Index -> Pull1 a -> Pull1 a
bitRev n = composeOn riffle (1...n)

riffle :: Syntax a => Data Index -> Pull DIM1 a -> Pull DIM1 a
riffle k = permute (const $ rotBit k)

-- Helper functions
composeOn :: (Syntax a) => (b -> a -> a) -> Pull DIM1 b -> a -> a
composeOn f v i = fromZero $ fold (flip f) i v

rotBit :: Data Index -> Data Index -> Data Index
rotBit 0 _ = P.error "rotBit: k should be at least 1"
rotBit k i = lefts .|. rights
  where
    ir = i .>>. 1
    rights = ir .&. oneBits k
    lefts  = (((ir .>>. k) .<<. 1) .|. (i .&. 1)) .<<. k