This week we're going to start implementing the actual raytracer. The goal for today's session is to bring us to the point where we can render spheres after the phong reflection model.
Some Raytracing Basics
Alright. This is how the most naive raytracer works:
All you need is an eye-point in space, objects (in our case a sphere), an image plane and one or more light sources.
Originating from the eye point we're going to shoot a ray through each of the pixels in the image and see if it hits the sphere - or to put it in math-lingo:
Given a ray consisting of an eye point e and a direction vector d, we're looking for a scalar value t so that the distance of the point e+(d*t) to the center of the sphere c is equal to its radius r.
All the points on the ray are described by the line equation
. The points on the sphere's surface are described by the sphere equation
. To get the points that the surface of the sphere and the ray have in common, we insert the first into the second equation:
. For simplified reading we declare
and resolve the equation to t:
. Note that the square root is not defined for negative values. If the part under the square root is less than zero, the ray doesn't intersect the sphere.
So what we need first is code that lets us perform operations on vectors.
vectr.clj
(ns vectr)
(defstruct vectr :x :y :z)
(defn make-vectr
"Creates a new vector out of its x, y and z coordinate"
[x y z]
(struct vectr (float x) (float y) (float z)) )
(defn vectr-add
"Adds two vectors"
[a b]
(make-vectr (+ (:x a) (:x b))
(+ (:y a) (:y b))
(+ (:z a) (:z b))))
(defn vectr-subtract
"Subtracts the second vector off the first"
[a b]
(make-vectr (- (:x a) (:x b))
(- (:y a) (:y b))
(- (:z a) (:z b))))
(defn vectr-dot
"Computes the dot product of two vectors"
[a b]
(+ (* (:x a) (:x b))
(* (:y a) (:y b))
(* (:z a) (:z b))))
(defn vectr-scale
"Scales the a vector"
[vectr scalar]
(make-vectr (* scalar (:x vectr))
(* scalar (:y vectr))
(* scalar (:z vectr))))
(defn vectr-length-square
"Computes the square of the vectors length"
[vectr]
(+ (* (:x vectr) (:x vectr))
(* (:y vectr) (:y vectr))
(* (:z vectr) (:z vectr))))
(defn vectr-length
"Computes the vectors length"
[vectr]
(Math/sqrt (vectr-length-square vectr)))
(defn vectr-normalize
"Scales the vector, so that its length is 1"
[vectr]
(vectr-scale vectr (/ 1 (vectr-length vectr))))
(defn vectr-cross
"Computes the cross product of two vectors"
[a b]
(make-vectr (- (* (:y a) (:z b))
(* (:z a) (:y b)))
(- (* (:z a) (:x b))
(* (:x a) (:z b)))
(- (* (:x a) (:y b))
(* (:y a) (:x b)))))
Now we can construct a ray ...
ray.clj
(ns ray
(:use vectr))
(defstruct ray :o :d)
(defn make-ray
[from to]
(struct ray from (vectr-normalize (vectr-subtract to from))))
(defn ray-point-at
[ray t]
(vectr-add (:o ray) (vectr-scale (:d ray) t)))
.. and a sphere.
sphere.clj
(ns sphere
(:use vectr))
(defstruct sphere :c :r)
(defn make-sphere
"Creates a sphere out of its center and radius"
[center radius]
(struct sphere center (float radius)))
(defn sphere-intersect
"Tests whether a ray intersects a sphere.
Returns the intersection points or nil"
[sphere ray]
(let [v (vectr-subtract (:o ray) (:c sphere))
v-dot-d (vectr-dot v (:d ray))
q (- (vectr-dot v v)
(Math/pow (:r sphere) 2.0))
sub-sqrt (- (* v-dot-d v-dot-d) q)]
(if (> sub-sqrt 0.0)
(let [sqrt (Math/sqrt sub-sqrt)]
[(- (* v-dot-d -1) sqrt)
(+ (* v-dot-d -1) sqrt)])
nil )))
We can now use these structures to replace our stub code in cray.clj to perform a real intersection test between a ray and a sphere.
cray.clj
(let [img (make-image 400 400)
eye (make-vectr 200 200 -500 )
sphere (make-sphere (make-vectr 200 200 300) 100)]
(image-every-pixel
img
(fn [img x y w h]
(let [pixel (make-vectr x (- h 1 y) 0 )
ray (make-ray eye pixel)]
(image-set-pixel!
img
x y
(if (not (nil? (sphere-intersect sphere ray)))
blue-color
black-color)))))
Executing this, we again get a black image with a blue circle in the middle - but this time we did it the proper way :).
cray $ clj cray.clj out.png
The Phong Reflection Model
Now, a blue circle isn't that exciting. What we really want is a nicely shaded sphere. In order to do this we need to introduce a shading model. One of the most popular ones is the Phong Reflection Model, named after Bui Tuong Phong. Although not being physically correct, it provides an acceptable tradeoff between accuratesse and algorithmic complexity.
It divides an intersection point's color into three separate parts:
-
Ambient
This portion of the color is constant all over the object, independent of the intersection point's surface normale. The blue circle that we rendered up to now is basically correct phong shading with 100% ambient portion, no diffuse and no specular part.
-
Diffuse
This part depends on the surface normale of the intersection point between the ray and the object. The object's color is scaled by the dot product of the reflection vector and the vector that points from the intersection point towards a light source. If these are orthogonal, the dot product is 0 and the diffuse part is black. If the vectors have the same direction, their dot product is 1 and the diffuse part is the object's color itself.
-
Specular
The specular component is reponsible for the little highlights on shiny surfaces, such as glass or plastic. If the reflection vector points almost directly to a light source, a portion of the light source's color is added.
Jump to the illustration
The weight of these components can be played with to define an objects material properties. A low ambient value with medium diffuse and high specular values looks like plastic. A higher ambient value with a high diffuse and zero specular value resembles a pebble.
What computational tools do we need to proceed?
We need to be able to add and scale colors.
color.clj
(defn color-clamp
"Makes sure that the components of the color are in legal range"
[clr]
(make-color
(max (min (:r clr) 1.0) 0.0)
(max (min (:g clr) 1.0) 0.0)
(max (min (:b clr) 1.0) 0.0 )))
(defn color-add
"Adds some colors"
[& args]
(color-clamp
(make-color
(reduce (fn [v obj] (+ (:r obj) v ) ) 0 args)
(reduce (fn [v obj] (+ (:g obj) v ) ) 0 args)
(reduce (fn [v obj] (+ (:b obj) v ) ) 0 args) ) )
)
(defn color-scale
"Dims a color by a scalar value. color*0==black, color*1==color"
[clr in-coeff]
(let [coeff (float in-coeff)]
(make-color (* coeff (:r clr))
(* coeff (:g clr))
(* coeff (:b clr)) ) ) )
A structure to hold the properties of a light source (position and color)
light.clj
(ns light)
(defstruct light :pos :col)
(defn make-light
"Creates a light source out of its position and color"
[pos col]
(struct light pos col))
Then we need a way to encapsulate material properties and some predefinitions for the comfort of use.
material.clj
(ns material
(:use color))
(defstruct material :col :amb :diff :spec :phong)
(defn make-material
[col amb diff spec phong]
(struct material col amb diff spec phong) )
(defn material-plastic
[col]
(make-material col 0.2 0.8 0.8 64) )
(def material-red-plastic (material-plastic red-color))
(def material-green-plastic (material-plastic green-color))
(def material-blue-plastic (material-plastic blue-color))
Now we can take one of our materials and apply it to our sphere.
sphere.clj
(defstruct sphere :c :r :material)
(defn make-sphere
"Creates a sphere out of its center, radius and material"
[center radius material]
(struct sphere center (float radius) material))
Given an intersection point, we need to compute the surface normale. This is easy for a sphere. We just subtract the center of the sphere off the intersection point and normalize the resulting vector.
(defn sphere-surface-normal
"Computes the surface-normal for a given intersection point"
[sphere point]
(do
(vectr-scale (vectr-subtract point (:c sphere)) (/ 1 (:r sphere)))))
Last but not least we need to put our puzzle together, compute the reflection vector r out or the ray's direction v and the surface normal n at the intersection point.
Once we have this, we can finally add shading to our sphere.
cray.clj
(defn phong-compose
[point ray normal material light]
(let [reflection (compute-reflection-ray ray point normal)
light-vectr (vectr-normalize (vectr-subtract point
(:pos light)))
half-vectr (vectr-normalize (vectr-add light-vectr (:d ray)))
diffuse (max 0.0 (vectr-dot light-vectr
(vectr-scale normal -1) ) )
specular (Math/pow (max 0.0
(vectr-dot half-vectr
(vectr-scale normal -1)))
(:phong material))
color (:col material)]
(color-add (color-scale color (:amb material))
(color-scale color (* diffuse (:diff material)))
(color-scale (:col light) (* specular
(:spec material))))))
(defn compute-color
[sphere light ray]
(let [intersects (sphere-intersect sphere ray)
material (:material sphere)]
(if (nil? intersects)
black-color
(let [point (ray-point-at ray (first intersects))
normal (sphere-surface-normal sphere point)]
(phong-compose point ray normal material light)))))
Finally we replace the simplistic intersection test in our main loop with compute-color
cray.clj
(let [pixel (make-vectr x (- h 1 y) 0 )
ray (make-ray eye pixel)]
(image-set-pixel!
img
x y
(compute-color sphere light ray))))))
Lets run the raytracer again and have a look at the result (Click to enlarge).
This looks more like a sphere, doesn't it? By a slight modification of the phong-compose-function, we can make the raytracer render only one of the components of the phong model for better understanding:
+
+
=
Jump to the explanation
Okay, pfeww. That's it for today. Next time we clean up the code a little. Right now, our computation model is pretty limited:
- Our objects, lights, materials, etc. are hardcoded in the main loop. A dynamic world-object is missing
- The compute-color-function knows only about about spheres. But all it should know about is generic objects. We need some polymorphism here.
If you wanna hack around in the sources, rather than copy and paste the code from the site you can clone my cray-repository at github