Course notes for Sept 25 class, Part II

Ray tracing to general second order surfaces

First a special case: ray tracing to an unit radius cylindrical tube

The implicit function for a unit diameter infinitely long cylindrical tube along the z axis is:

x2 + y2 - 1 ≤ 0
We can solve the ray equation manually, like we did for the sphere, as follows:

  • Substituting in V+tW for (x,y,z):
    (Vx + t Wx)2 + (Vy + t Wy)2 - 1 ≤ 0
  • Expanding out:
    (Wx2 + Wy2) t2 + (2 Vx Wx + 2 Vy Wy) t + (Vx2 + Vy2 - 1) ≤ 0
  • Which just equals:
    (WW) t2 + 2(VW) t + (VV - 1) ≤ 0
  • So to solve the quadratic equation in t:
    A = WW
    B = 2 (VW)
    C = VV - 1

But this would not give us the ability to transform the cylinder, nor would it generalize to other quadratic shapes, including cones, paraboloids and hyperboloids.

So instead, let's work on the general case of being able to transform any implicit quadratic function by 4×4 matrix, and then ray trace to it.

The general second order implicit surface:

Any quadratic implicit surface can be expressed as follows:

ax2 + by2 + cz2 + dyz + ezx + fxy + gx + hy + iz + j ≤ 0
The ten coefficients of this equation can be stored as a vector:
[a, b, c, d, e, f, g, h, i, j]
The unit sphere, the infinite unit cylinder and the infinite unit cone, respectively, can all be expressed in such a form:
Sphere: x2 + y2 + z2 - 1 ≤ 0 [1, 1, 1, 0, 0, 0, 0, 0, 0, -1]
Cylinder: x2 + y2 - 1 ≤ 0 [1, 1, 0, 0, 0, 0, 0, 0, 0, -1]
Cone: x2 + y2 - z2 ≤ 0 [1, 1, -1, 0, 0, 0, 0, 0, 0, 0]
Raytracing to a general second order implicit surface:

Ray tracing to the general form just requires substituting the ray into this equation to get a quadratic equation in t, exactly as we did above in the case of the cylinder:

  • Substituting in V+tW for (x,y,z):
    a(Vx+tWx)2 + b(Vy+tWy)2 + c(Vz+tWz)2 + d(Vy+tWy)(Vz+tWz) + e(Vz+tWz)(Vx+tWx) + f(Vx+tWx)(Vy+tWy) + g(Vx+tWx) + h(Vy+tWy) + i(Vz+tWz) + j ≤ 0
  • Expanding out:
    a(Vx+tWx)(Vx+tWx) + b(Vy+tWy)(Vy+tWy) + c(Vz+tWz)(Vz+tWz) + d(Vy+tWy)(Vz+tWz) + e(Vz+tWz)(Vx+tWx) + f(Vx+tWx)(Vy+tWy) + g(Vx+tWx) + h(Vy+tWy) + i(Vz+tWz) + j ≤ 0
  • Rearranging terms:
    (aWx2 + bWy2 + cWz2 + dWyWz + eWzWx + fWxWy) t2 + (2aVxWx + 2bVyWy + 2cVzWz + d(VyWz+VzWy) + e(VzWx+VxWz) + f(VxWy+VyWx) + gWx + hWy + iWz) t + (aVx2 + bVy2 + cVz2 + dVyVz + eVzVx + fVxVy + gVx + hVy + iVz + j) ≤ 0
  • So to solve the quadratic equation in t:
    A   = aWx2 + bWy2 + cWz2 + dWyWz + eWzWx + fWxWy
    B   = 2(aVxWx + bVyWy + cVzWz) + d(VyWz+VzWy) + e(VzWx+VxWz) + f(VxWy+VyWx) + gWx + hWy + iWz
    C   = aVx2 + bVy2 + cVz2 + dVyVz + eVzVx + fVxVy + gVx + hVy + iVz + j

Applying a linear transform to a second order surface

Let us first review the simpler case of linear surfaces: Half space P = (a, b, c, d) contains point p = (x,y,z), if and only if ax + by + cz + d ≤ 0. This can also be expressed as the product:

a b c d   •   x
  ≤   0

or, in other words, L pT ≤ 0. To find the plane that contains transformed points MpT, you need to replace L with LM-1.

Similarly, we can express quadratic equation

ax2 + by2 + cz2 + dyz + ezx + fxy + gx + hy + iz + j ≤ 0
as the following double product:

x y z 1   •  
af eg
0b dh
00 ci
00 0j
  •   x
  ≤   0

or in other words, p Q pT = 0, where Q denotes the 10 quadratic coefficients arranged into the above 4×4 matrix.

This means that if you want to find the three-variable quadratic equation that contains transformed points MpT, you need to replace Q by   (M-1)T Q M-1   because that will give you:

(M pT)T     (M-1)T Q M-1     (M pT)     =

pT MT     (M-1)T Q M-1     M pT     =

p Q pT

when you create the product of

(M-1)T P M-1
you can end up with non-zero values in all sixteen slots of your 4×4 matrix. This is because in addition to x*y, x*z, x*1, y*z, y*1, z*1 terms, you also end up with y*x, z*y, 1*z, z*y, 1*z, 1*z terms, as follows:
x*x x*y x*z x*1
y*x y*y y*z y*1
z*x z*y z*z z*1
1*x 1*y 1*z 1*1
Since multiplication of numbers is commutative, you can just add the six numbers in the lower left (shown in magenta above) into their mirror locations on the upper right, as follows:
00 KO+L
00 0P
This will leave you with just the ten unique coefficients of the quadratic in three variables, in the upper right of your transformed matrix.

Computing the surface normal:

Once you have computed your transformed values for [a,b,c,d,e,f,g,h,i,j], you can use high school differential calculus to compute its vector valued derivative:

f(x,y,z) = ax2 + by2 + cz2 + dyz + ezx + fxy + gx + hy + iz + j

f'(x,y,z) = [  2ax + ez + fy + g  ,  2by + dz + fx + h  ,  2cz + dy + ex + i  ]

Just normalize this derivative to get the normal at the surface point:
N = normalize(f'(x,y,z))

Light sources not at infinity

Light sources do not need to be at infinity. A local light source can add interesting chiarioscuro effects to your scene.

Unlike light sources at infinity, a local light source at some position (x,y,z) in your scene will more brightly illuminate surfaces that are near to it, and less brightly illuminate surface that are farther from it.

If a local light source were an infinitely small point, then intensity drop-off from that point would be 1/r2, where r is the distance of the surface from (x,y,z). In reality, light sources, such as desk lamps, tend to have some non-zero extent, and so their intensity drops off with distance more slowly than 1/r2.

A simple computational trick to model this slower dropoff is to define a power p, where 1.0 ≤ p ≤ 2.0, and define the dropoff function as 1/rp.

Another visual trick you might use with local light sources is to give a very high ambient component to the surface material of one of your objects (eg: a sphere), and then place a local light source inside this object. The object will then appear to be a light source that illuminates the other objects in your scene.

Fog and atmosphere

Fog and atmosphere is very easy to implement and can make your ray traced scenes look a lot cooler.

In the presence of fog, objects gradually fade from view, as a function of their distance from your eye. Uniform fog (ie: fog with the same density everywhere) is quite easy to model and render. The key observation is that for any uniform fog there is some distance D that light will travel before half of the light has been scattered away by fog particles. The other half of the light continues along unobstructed.

Now consider what happens when the light which has made it through the fog travels further by another equal distance D. Half of the remaining light will now be scattered away, and half of the remaining light will make it through the second distance. In other words, the light that makes it through a distance 2D unobstructed is (1/2)2. And in general, the amount of light that makes it through a distance kD is (1/2)k.

To calculate uniform fog which scatters half of its light every multiple of D, you need to measure the distance d traveled by the ray (ie: the distance from the ray origin to the surface point that the ray hits). The proportion of this light that makes it through the fog unscattered will be t=(1/2)d/D. The effect of the fog can be shown by mixing t(surfaceColor) + (1-t)(fogColor), where surfaceColor is the color that you would have calculated had there been no fog. Some uniform shade of light gray, such as (0.8,0.8,0.8), generally works well for the fog color.

Simple boolean shape combinations

In class, we briefly went over the general case of boolean modeling, but trying to implement the general case would be too much for one week.

So for now let's just talk about the simplest case (and a very useful one): intersection of convex shapes.

A convex shape is one that is not reentrant. That is, once a ray has entered the shape, the ray can then exit the shape at most once, and then never reenter that shape again. For example, a cube and a cylinder are convex, but a donut is not, since a single ray can enter and then exit the same donut twice.

Unlike the general case, the rule when raytracing for intersection of two convex shapes is very simple:

  • Let A and B be two shapes.
  • Let V+tW be a ray.
  • Let the entry and exit roots of A be tA0 and tA1, respectively.
  • Let the entry and exit roots of B be tB0 and tB1, respectively.
The intersection A∩B is found as follows:
  • The ray will enter A∩B at max(tA0,tB0)
  • The ray will exit A∩ B at min(tA1,tB1)
Note: if max(tA0,tB0) > min(tA1,tB1), this means the ray does not intersect A∩B.


To ray trace to a unit cylinder (a tube with front and back end-caps):

  • Start with the infinite unit tube:   x2+y2-1 ≤ 0
  • Then intersect with half space:   x - 1 ≤ 0
  • Then intersect with half space:   -x - 1 ≤ 0

Homework for October 1, part II:

  • Implement ray tracing to linearly transformed general quadratic surfaces.

  • Use the above, together with boolean intersection, to tray trace to cylinders with end caps.

  • Extra credit: implement fog and local lights. Implement other quadratic shapes. Also feel free to implement anything else you are inspired to.

As usual, try to make a cool animated scene -- something that you yourself find fun or inspiring.