haskellbrot

gistfile1.hs

import Data.Ratio
import Data.List
import Graphics.GD

type MandelNum = Double

maxIterations :: Int
maxIterations = 255

main :: IO()
main = makeImage

makeImage :: IO()
makeImage = do
   rendering <- imageRender (2048, 1536) 9 (-0.5, 0.0)
   savePngFile "output.png" rendering

makeASCIIArt :: IO()
makeASCIIArt = putStrLn $ asciiRender (100, 100) 6 (-0.7, 0.1)


imageRender :: Size -> Int -> (MandelNum, MandelNum) -> IO(Image)
imageRender size@(width, height) zoom center = do
   image <- newImage size
   drawList imageData width 0 image
   return image
   where
      imageData = concat $ mandelbrotFastImage size zoom center

drawList :: [Int] -> Int -> Int -> Image -> IO()
drawList []             _     _     image = return ()
drawList (pixel:pixels) width index image = drawPixel >> drawRest
   where
      drawPixel  = setPixel point pixelColor image
      drawRest   = drawList pixels width (index+1) image
      point      = (index `mod` width, index `div` width)
      pixelColor
         | pixel == maxIterations = rgb 0 0 0
         | otherwise              = rgb pixel 0 (maxIterations-pixel)

asciiRender :: (Int, Int) -> Int -> (MandelNum, MandelNum) -> String
asciiRender size zoom center = unlines $ seconds symbols
   where
      symbols   = map (map symbolify) image
      image     = mandelbrotFastImage size zoom center
      
      seconds []         = []
      seconds (a:b:rest) = a:(seconds rest)
      
      symbolify value
         | value > 200 = ' '
         | value > 180 = '@'
         | value > 160 = '%'
         | value > 140 = '&'
         | value > 120 = '$'
         | value > 100 = '='
         | value > 80  = '*'
         | value > 60  = '+'
         | value > 40  = ';'
         | value > 20  = '.'
         | value > 15  = '`'
         | otherwise   = ' '

-- Quad-tree pruning of the pixel space to determine any areas 
-- that are completely surrounded by the same value:
-- these are then assumed to be the same value throughout
mandelbrotFastImage :: (Int, Int) -> Int -> (MandelNum, MandelNum) -> [[Int]]
mandelbrotFastImage dimensions@(imgWidth, imgHeight) zoom offset
   | (min imgWidth imgHeight) < 3 = generateImage start end
   | isSurrounded                 = filledImage
   | otherwise                    = reformedImage
   where
      isSurrounded   = (borderFirst /= 0) && all (==borderFirst) borderRest
      
      resolution     = 1 / (fromIntegral (2^zoom))
      coordTransform = toMandelbrotCoords offset resolution dimensions
      start          = coordTransform (left, top)
      end            = coordTransform (right, bottom)
      halfWidth      = imgWidth  `div` 2
      halfHeight     = imgHeight `div` 2
      
      top            = 0
      bottom         = imgHeight-1
      left           = 0
      right          = imgWidth-1
      
      -- splitting the space into 4 even parts
      topLeft        = generateQuad ((left, top), (halfWidth-1, halfHeight-1))
      topRight       = generateQuad ((halfWidth, top), (right, halfHeight-1))
      bottomLeft     = generateQuad ((left, halfHeight), (halfWidth-1, bottom))
      bottomRight    = generateQuad ((halfWidth, halfHeight), (right, bottom))
      
      topRow         = concat $ generateQuad ((left, top), (right, top))
      bottomRow      = concat $ generateQuad ((left, bottom), (right, bottom))
      leftCol        = concat $ generateQuad ((left, top), (left, bottom))
      rightCol       = concat $ generateQuad ((right, top), (right, bottom))
      
      reformedImage  = (zipWith (++) topLeft topRight) ++ (zipWith (++) bottomLeft bottomRight)
      filledImage    = take imgHeight $ repeat (take imgWidth $ repeat borderFirst)
      
      generateImage  = mandelbrot maxIterations resolution
      
      (borderFirst:borderRest) = leftCol ++ topRow ++ bottomRow ++ rightCol
      
      
      generateQuad ((x0, y0), (x1, y1)) = mandelbrotFastImage (x1-x0+1, y1-y0+1) zoom quadOffset
         where
            quadOffset = coordTransform ((x0 + x1) `div` 2, (y0 + y1) `div` 2)


mandelbrotImage :: (Int, Int) -> Int -> (MandelNum, MandelNum) -> [[Int]]
mandelbrotImage dimensions@(imgWidth, imgHeight) zoom offset = mandelbrot maxIterations resolution start end 
   where
      resolution     = 1 / (fromIntegral (2^zoom))
      coordTransform = toMandelbrotCoords offset resolution dimensions
      start          = coordTransform (0, 0)
      end            = coordTransform (imgWidth-1, imgHeight-1)


toMandelbrotCoords :: (MandelNum, MandelNum) -> MandelNum -> (Int, Int) -> (Int, Int) -> (MandelNum, MandelNum)
toMandelbrotCoords (offsetX, offsetY) resolution (width, height) (col, row) = (x, y)
   where
      x = offsetX + (fromIntegral (col - (width  `div` 2)) * resolution) - resolution/2
      y = offsetY + (fromIntegral (row - (height `div` 2)) * resolution) - resolution/2


mandelbrot :: Int -> MandelNum -> (MandelNum, MandelNum) -> (MandelNum, MandelNum)  -> [[Int]]
mandelbrot limit resolution (xStart,yStart) (xEnd, yEnd) = map (map (score limit)) coords
   where
      columnStart  = 1 + (floor $ (xEnd - xStart) / resolution)
      coords       = splitEvery columnStart [(x, y) | y <- ys, x <- xs]
      ys           = [yStart,yStart+resolution..yEnd]
      xs           = [xStart,xStart+resolution..xEnd]


-- todo: cycle detection?
score :: Int -> (MandelNum, MandelNum) -> Int
score limit (x, y)
   | isKnownEscaped = limit
   | otherwise = length $ take limit $ takeWhile isNotEscaped $ escapeAttempts
   where
      coord = (x, y)
      escapeAttempts = iterate (attemptEscape coord) (0, 0)
      escapeRadiusSquared = 4
      isNotEscaped (x, y) = (x*x + y*y < escapeRadiusSquared)
      attemptEscape (x0, y0) (x, y) = (x*x - y*y + x0, 2*x*y + y0)
      
      q = (x - 1/4)^2 + y^2
      isKnownEscaped = q * (q + (x - 1/4)) < (1/4)*y^2

splitEvery :: Int -> [a] -> [[a]]
splitEvery _ [] = []
splitEvery n list = first : (splitEvery n rest)
   where
      (first, rest) = splitAt n list