Jon Harrop <use
...@jdh30.plus.com> writes:
>> Where is the code actually used for benchmarking, BTW?
> Here's the OCaml:
[... big snip ... ]
Ok, I took the trouble and translated your Ocaml code to SML for use
under SML/NJ. (MLton needs a minor adjustment.) The code is only a
few lines longer than the original and should come in under 100loc
when ignoring comments and blank lines.
On an old 400 MHz Celeron, the program runs almost 30% faster than the
original SML code that you posted earlier. I'd be interested to hear
how it fares on the Athlon...
Matthias
------------------
(* This is a straight translation of Jon Harrop's OCaml ray tracer code
* for the "Great Computer Language Shootout", with very minor modifications.
* Transliteration for SML/NJ done by Matthias Blume. *)
(* The following comment was in the original Ocaml code. It does not
* refer to differences between the Ocaml and SML/NJ versions: *)
(* This implementation differs from the original in several ways:
Uses an implicit scene, generated as a ray is traced rather than being
precalculated and stored explicitly in a tree.
Specialized shadow-ray intersection functions. *)
structure Ray:sig val main:string*string list->OS.Process.status end = struct
val delta = Math.sqrt 2.22044604925031e~16
(* 3D vector and associated functions *)
type vec = real * real * real
infix 7 *| fun s *| (x, y, z) = (s*x, s*y, s*z) : vec
infix 6 |+| fun (x, y, z) |+| (x', y', z') = (x+x', y+y', z+z') : vec
infix 6 |-| fun (x, y, z) |-| (x', y', z') = (x-x', y-y', z-z') : vec
infix 7 |*| fun (x, y, z) |*| (x', y', z') = x*x'+y*y'+z*z' : real
fun unitise r = (1.0 / Math.sqrt (r |*| r)) *| r
fun ~|(x,y,z) = (~x,~y,~z) : vec
(* Calculate the parametric intersection of the given ray with the given
sphere. *)
fun ray_sphere (orig, dir, center, radius) =
let val v = center |-| orig
val b = v |*| dir
val disc = b * b - v |*| v + radius * radius
in if disc < 0.0 then Real.posInf
else let val disc = Math.sqrt disc
val t2 = b + disc
in if t2 < 0.0 then Real.posInf
else let val t1 = b - disc in if t1 > 0.0 then t1 else t2 end
end
end
(* Calculate whether or not the given ray intersects the given sphere. *)
fun ray_sphere' (orig, dir, center, radius) =
let val v = center |-| orig
val b = v |*| dir
val disc = b * b - v |*| v + radius * radius
in if disc < 0.0 then false else b + Math.sqrt disc >= 0.0 end
(* Ratio of the radii of one level of spheres to the next. *)
val s = 6.0 / Math.sqrt 12.0
(* Find the first intersection point of the given ray with the scene. *)
fun intersect (level, orig, dir) =
let fun of_scene (center, radius, lambda, normal, level) =
if level = 1 then
let val lambda' = ray_sphere (orig, dir, center, radius)
in if lambda' >= lambda then (lambda, normal)
else (lambda', unitise (orig |+| lambda' *| dir |-| center))
end
else if ray_sphere (orig, dir, center, 3.0 * radius) >= lambda
then (lambda, normal)
else let val accu = of_scene (center, radius, lambda, normal, 1)
val r = 0.5 * radius val l = level - 1 val r' = s * r
fun h (dx, dz, (lam, norm)) =
of_scene (center |+| (dx, r', dz), r, lam, norm, l)
val mr' = ~r'
in h (r', mr', h (r', r', h (mr', r', h (mr', mr', accu))))
end
in of_scene ((0.0, ~1.0, 0.0), 1.0, Real.posInf, (0.0, 0.0, 0.0), level) end
(* Find if the given ray intersects the scene. *)
fun intersect' (level, orig, dir) =
let fun of_scene (center, radius, level) =
if level = 1 then ray_sphere' (orig, dir, center, radius) else
(* Exploit short-circuit evaluation of boolean comparisons to
* terminate this function early. *)
ray_sphere' (orig, dir, center, (3.0 * radius)) andalso
(of_scene (center, radius, 1) orelse
let val r = 0.5 * radius val l = level-1 val r' = s * r
in of_scene (center |+| (~r', r', ~r'), r, l) orelse
of_scene (center |+| (r', r', ~r'), r, l) orelse
of_scene (center |+| (~r', r', r'), r, l) orelse
of_scene (center |+| (~r', r', r'), r, l)
end)
in of_scene ((0.0, ~1.0, 0.0), 1.0, level) end
(* Trace a single ray by casting it into the scene and, if it intersects
anything, casting a second ray toward the light to determine occlusion.
*)
fun ray_trace (l, light, orig, dir) =
let val (lam, n) = intersect (l, orig, dir)
in if lam >= Real.posInf then 0.0
else let val g = ~ (n |*| light)
(* If we are on the shadowed side of a sphere then don't bother casting
a
shadow ray as we know it will intersect the same sphere. *)
in if g<=0.0 orelse intersect' (l, orig|+|lam*|dir|+|delta*|n, ~|light)
then 0.0 else g
end
end
(* Ray trace the scene at the given resolution. *)
fun doit n (* Resolution *) =
let val n2 = n div 2 val rn = real n val ns = Int.toString n
val light = unitise (~1.0, ~3.0, 2.0) (* Light direction *)
val level = 6 val ss = 4 (* ss: oversampling *)
val rss' = 1.0 / real ss val rss'sq255 = rss' * rss' * 255.0
fun yl y = if y < 0 then OS.Process.success else xl (0, y)
and xl (x, y) = if x >= n then yl (y-1) else dxl (0.0, 0, x, y)
and dxl (g, dx, x, y) =
if dx < ss then dyl (g, 0, dx, x, y)
else (print (String.str (chr (Real.trunc (0.5 + g * rss'sq255))));
xl (x+1, y))
and dyl (g, dy, dx, x, y) =
if dy >= ss then dxl (g, dx+1, x, y)
else let val orig = (0.0, 0.0, ~4.0)
fun cvt (w, dw) = real (w-n2) + real dw * rss'
val dir = unitise (cvt (x, dx), cvt (y, dy), rn)
in dyl (g+ray_trace (level, light, orig, dir), dy+1, dx, x, y)
end
in print (concat ["P5\n", ns, " ", ns, "\n255\n"]); yl (n-1) end
fun main (_, arg :: _) = doit (getOpt (Int.fromString arg, 256))
| main _ = doit 256
end