Easy distribution ray tracer

I am having my graphics class at Westminster College (not a new job-- I am teaching a course for them as an adjunct-- so go buy my apps!) turn their Whitted-style ray tracer into a distribution ray tracer.   I think how we are doing it is the easiest way possible, but I still managed to butcher the topic in lecture so I am posting this mainly for my students.   But if anybody tries it or sees something easier please weigh in.   Note these are not "perfect" distributions, but neither is uniform or phong really :)

Suppose we have a basic ray tracer where we send a shadow and reflection ray.

The basic code (where rgb colors geometric vectors are vec3) will look like:
vec3 rayColor(ray r, int depth) {
    vec3 color(0)
    if (depth > maxDepth) return color
    if (hitsScene(r)) {  // somehow this passes back state like hitpoint p
         if (Rs > 0) { // Rs is rgb specular reflectance
                ray s(p, reflect(r, N))
                color += rayColor(s, depth+1)
         }
         if (Rd > 0) { //Rd is rgb diffuse reflectance
                 ray shadow(p, L)   // L is the direction to the one light
                 if (not hitsScene(shadow)) {
                       color += Rd*lightColor*max(0, dot(N,L) )
                 }
          }
      }
      else 
             color = background(r)
      return color
}
               
           

The three spheres we add to get fuzzy effect.

Now let's assume we want to fuzz things up to get soft shadows and glossy reflections, and that we want diffuse rays for bounce lighting (even one level of depth would help).   We just need to pick random points in three sphere, one for each effect.   Let's assume we have a function rs() that returns a random point in a unit sphere centered at the origin (see last blog post).   All we need to do is randomly perturb each shadow and specular reflection ray, and we can generate another diffuse ray that will get bounce light.  This is:

vec3 rayColor(ray r, int depth) {
    vec3 color(0)
    if (depth > maxDepth) return color
    if (hitsScene(r)) {  // somehow this passes back state like hitpoint p
         if (Rs > 0) { // Rs is rgb specular reflectance
                // radius_specular is a material parameter where 0 = perfect
                ray s(p, reflect(r, N)+ radius_specular*rs())
                color += rayColor(s, depth+1)
         }
         if (Rd > 0) { //Rd is rgb diffuse reflectance 
                ray shadow(p, L + radius_shadow*rs())   // L is the direction to the one light
                if (not hitsScene(shadow)) {
                       color += Rd*lightColor*max(0, dot(N,L) )
                }
                ray diffuse(p, N + rs()) // assumes N is unit length
                color += Rd*rayColor(diffuse, depth+1)
          }
      }
      else 
             color = background(r)
      return color
}