How to deal with refraction when the rays start inside of a nested object

I am building a simple raytracer for educational purposes and want to add refraction to objects. Using Snells Law, I am able to create a new ray recursively at the intersection points. The raytracer supports currently only spheres and I use a scene where I have multiple spheres nested inside of each other with different refraction indices.

If I start a ray from outside of the spheres, everything seems simple. You start with the refraction index of the scene, and as soon as you hit the first sphere, refract the ray using the refraction index before and the refraction index of the material of the sphere until you hit the next sphere and so on. Using the normals of the intersection I can determine whether I enter or leave the sphere.

However, I don't understand how I should handle sphere leaves and what to do if the ray doesn't start in the outer part of the scene.

• Can I just take a stack of the refraction indices and go one layer up as soon as I leave a sphere?
• How can I determine with what refraction index I have to start if I start inside of the spheres?

Example

You hav开发者_运维知识库e three spheres, with refraction indices 0.9, 1.1 and 0.8 from outer to inner. Air index is 1.0

• Your camera is outside of the sphere and points at the center of the sphere:

  • start index is 1.0, you first hit the outer sphere with index 0.9 and refract from 1.0 to 0.9 and save that your ray is now in 0.9 material
  • you hit the middle sphere and notice the material constant of 1.1, since you have saved the 0.9, you know that you have to refract from 0.9 to 1.1 and save the 1.1 in addition to the 0.9
  • you hit the inner sphere and refract from 1.1 to 0.8 and you have save until now 0.9, 1.1 and 0.8
  • you hit the inner sphere again (this time you exit it, so you check your saved values and know that you have to switch back to 1.1)
  • ... until you are outside

• Problem now, when the camera is inside the sphere. You won't know to what refraction index you have to switch.

I have a similar ray tracer (written in Python) and stumbled over the same problem: in order to correctly work out the physics one must know the refractive index at each side of the intersection boundary. This took quite a while to solve elegantly, but in the end I went with this solution/design:

Design

1) Scene -- I have a master scene object (basically an array of all objects in the scene), you will probably have something similar. It stores geometrical objects.

Methods:

• intersection_points(ray) - returns a list of all intersection points, sorted by distance from the ray.
• intersection_objects(ray) - returns a list of all intersection objects, sorted by distance from the ray.
• containing_object(ray) - returns the object that contains the ray.
• objects() - returns a list of all the objects in arbitrary order.

Note: The scene adds an extra object to the list: the Scene_Boundary. This is a giant box (or Sphere) that encapsulate the whole scene i.e. EVERYTHING is inside this boundary.

2) Objects -- Make the geometric objects (e.g. your sphere) implement these methods.

Methods:

• contains(ray) - returns True is the ray origin is inside the object, False if on the surface and False if outside
• ray_is_on_surface(ray) - returns True is the ray is on the surface only, otherwise False.
• intersection_points(ray) - returns the intersection point(s) that the ray makes with the object
• surface_normal(ray) - returns the surface normal vector of the surface which the ray struck (this will help with Fresnel reflection and refraction)

For optical calculations the objects must also have a refractive index.

Instance variables:

• refractive_index

Boundary Problem

The problem we want to solve: what is the refractive index inside (n1) and outside (n2) of the boundary? To do this we follow this procedure:

1) Trace the ray through the whole scene:

sphere # origin = (0,0,0), radius = 1
ray  # origin = (0,0,0), direction = (0,0,1) Note: the ray is inside the sphere
scene.add_object(sphere)
ipoints = scene.intersection_points(ray) #  [ (0,0,1), (0,0,10) ]
iobjects = scene.intersection_objects(ray) # [ Sphere, Scene_Boundary]

Remember these are sorted by distance from the ray origin. The last item in ipoints and iobjects is the intersection the ray makes with the scene boundary. We will use this later!

2) n1 is found simply by finding the containing object, e.g.:

obj1 = scene.containing_object(ray) # Scene_Boundary
n1 = obj1.refractive_index() # n1 = 1. Scene_Boundary always has refractive index of Air

3) n2 is found by looking one object ahead in the iobject list, e.g. in pseudocode:

index = iobjects.index_of_object(obj1)
obj2 = iobjects[index+1]
n2 = obj2.refractive_index() # n2 = 1.5 e.g. Glass

4) Get the surface normal for later use:

normal = obj1.surface_normal(ray)

You have all the information you need to calculate the correct reflection and refraction. This is general enough to work even if the ray is outside the object, but occasionally I did need to implement some logical filtering to make the algorithm more robust, but that's basically it!

Reflection and Refraction

You can reflect a vector just by knowing the surface normal. In Python using numpy I do it like this,

def reflect_vector(normal, vector):
   d = numpy.dot(normal, vector)
   return vector - 2 * d * normal

Refraction (as discussed) needs n1 and n2 values:

def fresnel_refraction(normal, vector, n1, n2):
    n = n1/n2
    dot = np.dot(norm(vector), norm(normal))
    c = np.sqrt(1 - n**2 * (1 - dot**2))
    sign = 1
    if dot < 0.0:
        sign = -1
    refraction = n * vector + sign*(c - sign*n*dot) * normal
    return norm(refraction)

Finally you will need to calculate the reflection coefficient for the ray where angle is the angle between the ray direction and the surface normal (this assumed your ray is 'unpolarised'). Compare this with a random number between 0 and 1 to decide if reflection occurs.

def fresnel_reflection(angle, n1, n2):
    assert 0.0 <= angle <= 0.5*np.pi, "The incident angle must be between 0 and 90 degrees to calculate Fresnel reflection."
    # Catch TIR case
    if n2 < n1:
        if angle > np.arcsin(n2/n1):
            return 1.0

    Rs1 = n1 * np.cos(angle) - n2 * np.sqrt(1 - (n1/n2 * np.sin(angle))**2)
    Rs2 = n1 * np.cos(angle) + n2 * np.sqrt(1 - (n1/n2 * np.sin(angle))**2)
    Rs = (Rs1/Rs2)**2
    Rp1 = n1 * np.sqrt(1 - (n1/n2 * np.sin(angle))**2) - n2 * np.cos(angle)
    Rp2 = n1 * np.sqrt(1 - (n1/n2 * np.sin(angle))**2) + n2 * np.cos(angle)
    Rp = (Rp1/Rp2)**2
    return 0.5 * (Rs + Rp)

Final comments

This all comes from my Python optical ray tracing project which is not yet released(!), but you can check here for some detail: http://daniel.farrell.name/freebies/pvtrace. I like Python! There are a number of Python ray tracing projects listed here, http://groups.google.com/group/python-ray-tracing-community/web/list-of-python-statistical-ray-tracers . Finally, be careful with fractional refractive indices in your example, the equation will breakdown.

Update

Screenshot of this implemented in my ray tracer, available at http://github.com/danieljfarrell/pvtrace

Posting this from a physics point of view and not a raytracing implementation point of view :P.

Snell's law states that the ratio of the sines of the incident angle and the refractive angle are equal to the inverse of the ratio of the refractive index of the two mediums on either side of the boundary.

Thus, when you have a ray approaching a new material and you want to know the angle in the new material, you need to know the angle that the ray is hitting the new material, the index of refraction of the new material, and the index of refraction of the material the ray is currently in.

As you say refraction works fine moving into the sphere, you must know the index of refraction of each sphere and your scene already.

I'd say creating a stack of refractive indices would be a good way to deal with going into a bunch of nested materials, as you're going to have to touch all the refractive indexes that you push onto the stack again as you move out of the nested set of spheres.

As to determining what refractive index you have to start with as you leave the spheres, you're always saying sin(theta1)/sin(theta2) = [refractive index of 2]/[refractive index of 1]. Thus, you need the refractive index of the material that you're currently in and the index of the material that you're going to be moving towards.

Apologies if I misunderstood your question, but I hope that helps!