Hi,

I thought that was transparent !!

However there are a couple of differences from Hart's basic method:

1. It allows for RIFS - that's why the code to find the next transform is looking up in an array.

2. It allows for depth restrictions on transforms - that's why there's a look up in an array to find if a given transform is allowed at a given depth.

Here's some pseudo-code, I think I've got it correct. Hope it helps - the only bit that I'm hoping I've got correct is the loop to go to the next transform at the current depth and if we've done them all then to go up a depth.

Here's the pseudo-code for a single viewing ray:

calculate the base point for the ray
calculate the normalised direction vector for the ray
get the front clip (distance along ray from base point)
get the back clip offset (ditto)
set the "found" distance to a big value

if we've hit the bounding sphere

  store base point for depth 0
  store direction vectors for depth 0
  set current transform to first transform
  (normally number of transforms - 1)
  set scale at depth 0 to 1
  set depth to 0

  while depth>=0

    store current transform as transform for this depth
    current scale = scale at this depth*scale for current transform
    increment depth
    transform base point and direction vector storing values as "current" and for the new depth
    set flag to false (=="goto next transform at current depth"
    store current scale as scale for new depth

    if we're within the bailout distance

      if nearest intersection with sphere is <= distance to back clip and farthest intersection with sphere is >= front clip

        if reached max depth or (depth>min depth and reached minscale)

          if new nearest distance<found distance

            set found distance and back clip to new distance
            store any colouring information

          else

            current transform = first transform (normally number of transforms - 1)
            set flag to true (=="goto next depth" i.e. NOT "goto next transform at current depth")

          endif

        endif (reached max depth etc)

      endif (intersection with sphere within scan)

    endif (within bailout distance)

    if flag false (== "goto next transform at current depth")

      depth = depth - 1
      current transform = transform[depth]
      while flag==false

        current transform = current transform - 1
        if current transform<0

          depth = depth - 1
          if depth<0

            f = true

          else

            current transform = transform[depth]

          endif

        else

          f = true

        endif

      endwhile

      if depth>=0

        set current base point and direction vectors to values from depth

      endif

    endif

  endwhile (depth>=0)

  if found distance<big value (set at start)

    point = base point + direction vector * found

  endif

x2=1.0/(dx*dx+dy*dy+dz*dz) = 1.0

ah sorry, all clear now (i am used to code with const-declarations, which helps a lot in matters like these)! however, i think you can still avoid those operations if you know:

1. that (dx,dy,dz) is normalised in the beginning
2. 1 / |M| for all the xform matrices M

then you could compute (dx,dy,dz) = M(dx,dy,dz) * Mscale, where Mscale is inverse of detM, to simplify some of the later computations.