In this post, we will develop a simple ray tracing model, which can be a very useful tool in many domains. For instance in acoustics, optics, fluids dynamics... Why so many different domains? Well, this is related to the fact that under certain hypotheses (high frequency and paraxial approximations) waves (light, sound, water...) behave like rays.
There's many good introductions to ray tracing, but here we will define things as we go. The way we want to trace rays is very basic. Essentially, we want to be able to:
To do this, we will model the geometry as a sequence of oriented segments, with wave speed properties on each side of the segment. Any type of geometry will thus be represented as a list of segments. Also, we will only work in 2d.
A single segment might be modelled as follows:
class Segment:
"A class representing segments."
def __init__(self, A, B, normal_vector, celerity_front, celerity_back):
"Inits a segment from A to B, with a normal vector and front and back celerities."
self.A = A
self.B = B
self.normal_vector = normal_vector
self.celerity_front = celerity_front
self.celerity_back = celerity_back
A sphere would thus be represented by a sequence of segments. Any geometry could thus just be a list of segments.
So how does one trace a ray? There's three things here:
Let's tackle these. But first, we need to define a ray.
A ray needs an origin and an end.
class Ray:
"A class representing a ray."
def __init__(self, origin, end):
self.origin = origin
self.end = end
Finally, we also need sources of rays, which are defined by a point and a direction, which is an unit vector.
class Source:
"A ray source."
def __init__(self, origin, direction):
self.origin = origin
self.direction = direction
Finally, we also need points:
from collections import namedtuple
Point = namedtuple('Point', 'x, y')
Point(1, 2)
Point(x=1, y=2)
And vectors:
Vector = namedtuple('Vector', 'u_x, u_y')
Vector(3., 2.4)
Vector(u_x=3.0, u_y=2.4)
Once we have a source, we can trace it until a segment with which it intersects. The intersection point allows us to return a ray.
def trace_ray(source, geometric_object):
"Traces a ray from a source to a geometric object."
intersected_segment = min(geometric_object.segments, key=lambda seg: tof(source, seg))
intersection_point = intersection(source, intersected_segment)
return Ray(source.origin, intersection_point)
# TODO: tof function, intersection function
Let's now turn to the intersection function. Our rays are always straight lines. So it is easy to come up with a way for intersecting a line.
def intersection(source, segment):
"Computes the intersection point between a source and a segment."
# they're parallel: no intersection or infinity
if scalar(source.direction, segment.normal_vector) == 0.:
return None
else:
File "<ipython-input-9-30eed575786f>", line 6 else: ^ SyntaxError: unexpected EOF while parsing
A = Point(0, 0)
B = Point(5, 0)
S = Point(1, -4)
import matplotlib.pyplot as plt
%matplotlib inline
plt.style.use('bmh')
def draw_points(points, labels=None):
"Draws a list of points with labels."
if labels is None:
labels = ['Point {}'.format(ind) for ind, val in enumerate(points)]
for p, label in zip(points, labels):
plt.plot(p.x, p.y, 'o', label=label)
draw_points([A, B, S], ['A', 'B', 'S'])
plt.legend()
plt.xlim(-1, 7)
plt.ylim(-5, 1)
(-5, 1)
AB =
File "<ipython-input-14-318c53a6be6f>", line 1 AB = ^ SyntaxError: invalid syntax
Information about line intersections can be found here:
http://geomalgorithms.com/a05-_intersect-1.html
I thing this is the way to go.
Drawing rays that intersect a straight wall.
s = Point(0, 0)
east = Vector(1, 0)
source = Source(s, direction=east, aperture=(-30, 30), n=11)
top = Point(1, 1)
bottom = Point(1, -1)
wall = Segment(top, bottom)
rays = Intersect(source, wall)
plot(rays)
Drawing rays that get focused by spherical shell.
s = Point(0, 0)
east = Vector(1, 0)
source = Source(s, direction=east, aperture=(-30, 30), n=11)
top = Point(1, 1)
bottom = Point(1, -1)
wall = Segment(top, bottom)
rays = Intersect(source, wall)
plot(rays)