Fractal Rendering [Color Version]

fractal.hs

import Data.Array.Repa ( Array, DIM2, DIM3, Z(..), (:.)(..) )
import qualified Data.Array.Repa          as R
import qualified Data.Array.Repa.IO.DevIL as D
import Data.Word  ( Word8   )
import Data.Fixed ( divMod' )
import Debug.Trace

import qualified Data.Map as M
import qualified Data.Vector.Unboxed as UV

traceS a = trace $ show a

type R     = Float
type R2    = (R, R)
type Angle = R

pixels :: Int
pixels = 640
pixels_half = div pixels 2
scale :: R
scale = 160

a = -1.4
b = 1.6
c = 1.0
d = 0.7
type R3 = (R,R,R)
start :: R3 
start = (0.156452879309654, 0.466939155012369, 0.144256137311459)
clifford_iter :: R -> R3 -> R3 
clifford_iter dt (x,y,z) = (dx, dy, dz)
  where dx = dt * (sin (a*y) - c * cos (a*x))
        dy = dt * (sin (b*x) - d * cos (b*y))
        dz = 0
clifford_num :: R -> R3 -> Int -> R3
clifford_num _ v 0 = v
clifford_num dt v n = clifford_num dt (clifford_iter dt v) (n-1)
clifford_expand :: R -> R3 -> Int -> [R3]
clifford_expand _ v 0 = [v]
clifford_expand dt v n = v : clifford_expand dt (clifford_iter dt v) (n-1)
project :: R3 -> Int
project (x,y,z) = 
  if vy >= pixels || vy < 0 then error $ "vy out of range (" ++ show x ++ " ," ++ show y ++ " ," ++ show z ++ ")" else
  vx * pixels + vy
  where vx = truncate (x*scale) + pixels_half -- + truncate (z * 0.1 * scale)
        vy = truncate (y*scale) + pixels_half -- + truncate (z * 0.1 * scale)
updateFromPositions :: UV.Vector R3 -> [R3] -> UV.Vector R3
updateFromPositions arr ps =
  let alist = M.toList $ foldr (\p map->M.insertWith threePlus (project p) (mapToColor map) map) (M.empty) ps in
  arr UV.// alist where
    development map = fromIntegral $ M.size map
    mapToColor map = (0.02, 0.02 + 0.0000001 * development map, max 0 (0.03 - 0.0000001 * development map))
    threePlus (a1,b1,c1) (a2,b2,c2) = (a1+a2, b1+b2, c1+c2)

vectorToImage :: UV.Vector R3 -> Array DIM3 Word8
vectorToImage vec = R.traverse arr8 (:. 4) chans where
    firstOfThree (a,b,c)  = a
    secondOfThree (a,b,c) = b
    thirdOfThree (a,b,c)  = c
    arr  = R.fromVector (Z :. pixels :. pixels) vec
    rToWord8 = (floor . (*255) . min 1 . max 0)
    arr8 = R.map (\(x,y,z)->(rToWord8 x,rToWord8 y, rToWord8 z)) arr
    chans _ (Z :. _ :. _ :. 3) = 255  -- alpha channel
    chans a (Z :. x :. y :. 0) = firstOfThree  $ a (Z :. x :. y)
    chans a (Z :. x :. y :. 1) = secondOfThree $ a (Z :. x :. y)
    chans a (Z :. x :. y :. 2) = thirdOfThree $ a (Z :. x :. y)

main = do
  let trajectory = clifford_expand 1.0 start 1000000
      zeros = cycle [0]
      whiteVector = UV.fromList $ (take (pixels*pixels) $ zip3 zeros zeros zeros)
      vect = updateFromPositions whiteVector trajectory
      image = vectorToImage vect
  D.runIL $ D.writeImage "outfrac.png" image

TODO : make parameters (a,b,c,d, start) into configurable environment