by Mithrandir

Yet another puzzle solved in Haskell. This time we have a raytracing problem. A simplified version but, nonetheless, hard to do in a functional way. At least in a non-monadic fashion.

Enough talk. Let’s see the problem: suppose we have a box with size n by m. Inside the box we have some light sensors – some cells which are initially black and which will turn to another colour as soon as one light ray comes very close to them. The sensors are placed at all integer coordinates in the box (assuming that the maximum allowed coordinates are n and m). Their sensitivity is chosen such that a cell is active only if the ray passes the point where the sensor is placed. Filling the box with sensors, we would like to obtain a picture of their state after several light sources are placed in the box.

To make the problem harder, the walls of the box are treated with a substance which makes them reflect any incident radiation. Except on some points, where all radiation is absorbed.

Our task was to write a function which would receive as input the dimensions of the box, a list of all absorbent spots and a list of the light sources. Each light was given as a triple (position, direction, RGB colour). To simplify things a little, there are onle some directions allowed: a direction is valid if it makes an angle with the Ox axis whose tangent is rational. That is, each direction can be written as vector with integer coordinates. We call these angles rational for brevity.

We begin by defining our box data structure and some several type synonyms to help us understand the code. As you see, we are using again our old network module. In fact, I intend to expand it in the following days following some usage patterns that I’ve observed while writing code which uses it.

import Network
import System.Environment(getArgs)
data Box = Box {
  net :: Network (Int, Int) (Int, Int, Int)
  , n :: Int
  , m :: Int
type Pos = (Int, Int)
type Dir = (Int, Int)
type Col = (Int, Int, Int)
type Dim = (Int, Int)

Now, we need to build the box. In fact, as can be seen from the data structure definition, our box doesn’t hold very much information: only it’s size and the colour of each sensor are stored.

buildBox :: Int -> Int -> Box
buildBox n m = Box box n m
    box = Network (\_ -> (0, 0, 0)) neigh
    neigh p = filter valid (list p)
    list (x,y) = [(x, y-1), (x, y+1), (x-1, y-1), (x-1, y), (x-1, y+1),
      (x+1, y-1), (x+1, y), (x+1, y+1)]
    valid (x, y) = and [x>=0, x<=n, y>=0, y<=n]

Because the box is limited by some walls we need to make sure that we don’t reach invalid coordinates. To do this, we use the valid predicate in which we have used the and function for brevity and readability.

Now, we can start raytracing. First, we need to extract the next sensor which will turn on for the same ray considering the current one and the ray’s direction. Knowing that the sensors are at integer coordinates and that the rays are guaranteed to make only rational angles the code is simpler than in the general case.

nextPosDir :: (Pos, Dir) -> Dim -> (Pos, Dir)
nextPosDir ((x, y), (dx, dy)) (n, m) = inter (np x dx n, np y dy m)
    inter ((a, b), (c, d)) = ((a, c), (b, d))
    np pos delta limit
      | p < 0 = (-p, -delta)
      | p == limit = (limit-2, -delta)
      | p > limit = (2 * limit - p, -delta)
      | otherwise = (p, delta)
        p = pos + delta

Updating one sensor colour is simple: one simple takes the biggest red, green and blue component from each ray that passes through that sensor.

updateColor :: Pos -> Col -> Col -> Col
updateColor (x, y) (r, g, b) (or, og, ob) = (max r or, max g og, max b ob)

Now comes the most tricky part: simulating one single ray and turning on all the relevant sensors. The until function is useful again.

traceRay :: Box -> [Pos] -> (Pos, Dir, Col) -> Box
traceRay b abs l@(p, d, c) = f $ until test next initial
    f (x, y, z) = x
    (np, nd) = nextPosDir (p, d) (n b, m b)
    initial = (b, np, nd)
    next (box@(Box net n m), pos, dir) = (newbox, newpos, newdir)
        oldval = net `valueAt` pos
        newval = updateColor pos c oldval
        newbox = Box (net `update` (pos, newval)) n m
        (newpos, newdir) = nextPosDir (pos, dir) (n, m)
    test (_, pos, dir) = or [pos `elem` abs, and [d == dir, p == pos]]

After we can trace a ray, is simple to trace the rays generated by all light sources:

rayTrace :: Box -> [Pos] -> [(Pos, Dir, Col)] -> Box
rayTrace b _ [] = b
rayTrace b abs (light:lights) = traceRay (rayTrace b abs lights) abs light

Now, our problem is solved. Only 50-60 lines of code and no State Monad used.

To make matters funnier, let’s write the output to a SVG file instead of printing it to stdout. For this, we need to write to a file a standard header. What all those fields are is not my concern right now, one can simply read the W3C documentation if interested.

filehead :: String
filehead = concat ["<?xml version=\"1.0\" standalone=\"no\"?>\n",
  "<!DOCTYPE svg PUBLIC \"-//W3C//DTD SVG 1.1//EN\" ",
  "<svg xmlns=\"http://www.w3.org/2000/svg\" version=\"1.1\" "]

See the concat above? It is used because writing the entire string on a single line will make that line non-readable. Writing it on multiple lines with ++ or \<newline>\ can be error prone and misleading at best. This function, concat does the needed job: construct a big string from a list of smaller strings.

The above header is not really closed. We need to give some image information before ending it: what size will the image be? To do this, we have to know the layout of the box: it’s size.

sizedef :: Box -> String
sizedef b = concat ["width=\"", show $ n b, "cm\" height=\"", (show $ m b),
  "cm\" viewBox =\"0 0 ", show $ 100 * n b, " ", show $ 100 * m b, "\">"]

Now, we turn our attention to printing only a single cell. We will use a SVG rect element and fill all the required fields.

printCell :: Box -> Pos -> String
printCell box (x, y) = concat ["<rect x=\"", show $ 100 * x, "\" y=\"",
  show $ 100 * y, "\" width=\"100\" height=\"100\" stroke=\"none\" ",
  "fill=\"rgb(", show r, ", ", show g, ", ", show b, ")\"/>"]
    (r, g, b) = net box `valueAt` (x, y)

Printing the entire box is easy.

boxprint :: Box -> String
boxprint b = unlines $ map (printCell b) list
    list = [(x, y) | x <- [0..n b - 1], y <- [0..m b -1]]

Now, we print the end element, construct all the tags and write them to file. See the use of another useful function: unlines which adds a end-line marker between each string before concatenating each string. Useful for writing each rect on its line.

fileend :: String
fileend = "</svg>"
makeTags :: Box -> String
makeTags box = unlines [filehead ++ sizedef box, boxprint box, fileend]
writeBox box = writeFile "result.svg" $ makeTags box

Now, we define our test data

main = writeBox $ rayTrace (buildBox 100 100) [(0, 5), (1, 0), (4, 0),
  (32, 0), (48, 0), (64, 0), (72, 99), (45, 99), (3, 99), (0, 4), (0,
  19), (0, 27), (99, 3), (99, 15), (99, 67)]
  [((42, 42), (1, 1), (100, 2, 4)),
  ((24, 6), (1, 5), (0, 0, 255)),
  ((23, 5), (2, 3), (1, 129, 128)),
  ((68, 29), (-1, -5), (255, 255, 0)),
  ((25, 89), (1, -1), (255, 255, 255)),
  ((1, 97), (1, 1), (255, 0, 255))]

And run it.

$ time ./test

real 0m1.507s
user 0m1.468s
sys 0m0.012s

The results are breathtaking. See the following snippet.


That’s all for this problem, I will come back soon with another article.

About these ads