Course notes for Sept 25 class, Part I

The Phong Reflectance Algorithm

The Phong algorithm, originally developed by Bui Tuong Phong, was the first really good attempt in computer graphics to create a visual approximation to the way light interacts with surface materials.

In reality, most "smooth" surfaces have a microstructure of some sort, consisting of little pits and bumps -- a sort of microscopic landscape. So incoming light doesn't bounce of into quite the mirror direction, but rather into a distribution of directions centered around the mirror. This is known as "specular reflection".

The phong model consists of three components:

  • Ambient: reflectance off the surface due to ambient light bouncing around the scene,

  • Diffuse: reflectance due to Lambertian diffuse illumination (like the light that penetrates a plastic surface and bounces off dye particles before re-emerging, and

  • Specular: reflectance due specular reflection on the surface.

When these three components are put together, we get the following formula, which sums over all light sources in the scene:

Argb   +   i Lrgbi ( Drgb max(0, NLi)   +   Srgb max(0, -WR)p )
  • Li and Lrgbi are the direction and color, respectively of the ith light source;

  • R is the direction into which light from direction Li would reflect if the surface were a perfect mirror;

  • Argb, Drgb and Srgb are the respective colors of the ambient, diffuse and specular portions of reflectance.

  • p is the specular power.

Intersection of a ray with a plane

The statement that a plane L = (a,b,c,d) contains a point p = (x,y,z,1) can be expressed by the linear equation:

Lp = 0
Expressed in component form, this is:
ax + by + cz + d = 0
All points for which this is true are said to be on the plane.

Similarly, we can express the set of all points p that are either on plane L or are on a particular side of the plane as:

Lp ≤ 0
ax + by + cz + d ≤ 0
The shape that such an equation defines is called a half space.

Given a ray (V+tW), we can compute the intersection of the ray with plane L as follows:

L • (V+tW) = 0

a(vx + t wx) + b(vy + t wy) + c(vz + t wz) + d = 0

avx + bvy + cvz + d     +     t wx + t wy + t wz     =     0

(LV) + t (LW) = 0

t   =   - ( LV ) / ( LW )

Now consider instead the half-space defined by:

Lp ≤ 0
In this case, the solution for t becomes:
t   ≤   - ( LV ) / ( LW )

This can result in a number of distinct geometric situations, depending on the sign of different parts of this equation:

If LW is negative, it means that the ray is entering the half space, and the root t describes where the entering ray crosses the plane. If the computed root t in this case is negative, it means the ray began inside the half-space, so the solution will consist of the entire ray, from t=0 to t=∞.

If LW is positive, it means that the ray is exiting the half space, and the root t describes where the exiting ray crosses the plane. If the computed root t in this case is positive, it means the ray began outside the half-space, so the solution will consist of the empty set.

There are also two degenerate cases when LW is zero, which means that the ray direction W is tangent to the surface of the plane. In such cases, if LV ≤ 0, then the entire ray lies within the half space (that is, the solution consists of the entire ray), whereas if LV > 0, then none of the ray lies within the half space (that is, the solution is the null set).

Intersection of a ray with a cube

If a plane defines a half space, then six planes can be arranged to define an intersection of half spaces, so as to form a unit cube:

x ≤ 1
x ≥ -1
y ≤ 1
y ≥ -1
z ≤ 1
z ≥ -1
Or, expressing each of these six equations in its canonical form (a,b,c,d), to represent ax + by + cx + z ≤ 0:
One wonderful thing about ray tracing is that we can use it to render this cube by just observing a very simple truth: Within the one dimensional world of the ray, this union of three dimensional point sets becomes a union of one dimensional point sets (ie: line segments).

Since ray tracing to each of the six half spaces results in one (possibly infinite) line segment, then within the ray, the volume within the cube is simply the intersection of these six line segments.

If the resulting intersection is empty, this indicates that the ray has missed the cube.

We already know, for a given ray V+tW and plane L = (a,b,c,d), that there are four possible cases:

LW < 0 The ray enters the half space at t = - LV / LW
LW > 0 The ray exits the half space at t = - LV / LW
LW = 0 and LV ≤ 0 The entire ray lies inside the half space
LW = 0 and LV > 0 The entire ray lies outside the half space
To find the intersection of the ray with the cube, we take the maximum tenter of the entering rays, and the minimum texit of the exiting rays, and see whether tenter < texit. If so, then the ray has entered the cube at tenter. If not, then the ray has missed the cube.

We handle the two degenerate cases as follows:

  • Any L for which LW = 0 and LV ≤ 0, we can simply ignore that L, since the ray's intersection with this half-space is the entire ray.

  • If there is any L for which LW = 0 and LV > 0, then the ray has missed the cube, since the ray's intersection with this half-space is the null set.

We are going to want to do lighting calculations between the ray and the cube. For this reason we need to keep track of the surface normal, at the point where the ray enters the cube. To do this, we need to do two things:

  • Keep track of which of the six ray/half-plane intersections resulted in tenter, and

  • For that L = (a,b,c,d), use normalize(a,b,c) as the surface normal for our lighting computations.

Transforming the cube: from matrix to cube-equations transform

Things get a lot more interesting when you can transform an object that you are raytracing. When combined with a transformation matrix, a single unit cube can then be used to represent any box, no matter its position, its orientation, or the lengths of its sides.

A note about how matrices are stored in GLSL:

There are two ways a matrix can be stored in a one dimensional array: in row-major order, or in column-major order. The GLSL shader language (which you are using) uses column-major order, so the 16 sequential arguments to the mat4() function are interpreted as follows:

   m0 m4 m8  m12
   m1 m5 m9  m13
   m2 m6 m10 m14
   m3 m7 m11 m15
The operation   matrix * point   will treat the point as a column vector, and will multiply the matrix by that vector as follows:
   m0 m4 m8  m12     p0
   m1 m5 m9  m13  *  p1
   m2 m6 m10 m14     p2
   m3 m7 m11 m15     p3
You can also multiply the other way:   point * matrix. This will have the interpretation:
                      m0 m4 m8  m12
   (p0 p1 p2 p3)  *   m1 m5 m9  m13
                      m2 m6 m10 m14
                      m3 m7 m11 m15

We are going to cover matrices in more depth in the next lecture, meanwhile, here is a simple function you can use in your shader program to create a matrix that does both translation and scale. If you want to be ambitious, feel free to try to also implement rotation on your own before we cover it next week:

//  Function that generates a matrix to translate and scale a point:

    mat4 tsMatrix(vec3 t, vec3 s) {
        return mat4(s.x,0.,0.,0., 0.,s.y,0.,0., 0.,0.,s.z,0., t.x,t.y,t.z,1.);
If you already have two vec3 variables T and S specifying a particular translation and scale, then you can create a transformation matrix with tsMatrix as follows:
   mat4 matrix = tsMatrix(T, S);
You can also use tsMatrix to create the inverse of this transformation:
   mat4 inverseMatrix = tsMatrix(-T/S, 1.0/S);
As we discussed in class, if matrix M translates each point p as follows:
p   →   Mp
then in order to preserve the property   Lp = 0   we need to transform L as follows:
L   →   LM-1
To ray trace to a cube that has been transformed, you should trace your ray to a transformed version of each of its six defining linear inequalities.

In all other ways, your ray/cube algorithm stays exactly the same.

If L' = (a',b',c',d') is the transformed plane that defines the face where your ray enters the transformed cube, then you can get the normal vector for lighting computation as: normalize(a',b',c').

Here is an example of tsMatrix() used in various ways, which you might find helpful:

     vec4 point = vec4(1.0,0.0,0.0,1.0);
     vec4 plane = vec4(1.0,0.0,0.0,-1.0);

     vec3 T = vec3(0.0, 1.0, 0.0); // translate upward by one unit
     vec3 S = vec3(0.5, 0.5, 0.5); // scale down by one half

     mat4 matrix = tsMatrix(T, S);
     vec4 newPoint = matrix * point; // result will be: (0.5,1.0,0.0,1.0)

     mat4 inverseMatrix = tsMatrix(-T/S, 1.0/S);
     vec4 newPlane = plane * inverseMatrix;

Still to come:

  • Ray tracing to general second order surfaces:
    • Standard forms (sphere, cylinder, cone)
    • Transforming quadratic by a matrix
    • Solving the transformed quadratic

  • Light sources not at infinity

  • Fog and atmosphere

  • Simple boolean shape combinations

First part of homework:

  1. Add boxes to your ray tracing scene.

  2. Show that you can ray trace to transformed boxes.