Introducing our performance example - Ray tracing

We will be using ray tracing for most of our performance journey with Python

The purpose of this newsletter is discussing Python performance issues. In order to discuss the many perspectives related to performance we are going introduce an application that will serve as our example for many articles to come.

Our example will be a ray tracer. A ray tracer is a graphics rendering technique to generate photo-realistic images. It is widely used in modern games and in Hollywood film CGI.

Ray tracing is a complex technology, but Peter Shirley has written a set of books, that introduce the basic concepts in an amazingly simple way. We will use his perspective as a simple starting point to ray tracing that will allow us to illustrate Python’s performance issues.

Remember that our guiding light is writing performant Python code, not the specific problem - ray tracing. Sometimes will will gloss over ray tracing related issues. If you are interested in that topic, I strongly recommend you to have a look at Peter Shirley’s books, and the accompanying website.

Over time we will be discussing many different implementations of ray tracing in Python… and in C as sometimes we really need a low-level language to implement the most performance critical parts of our code. C is the most common language to do this.

The first implementation that we will present here is based on NumPy. But don’t look at it as a great approach (it’s not). It’s just a starting point that will help us on our journey to write for efficient Python code.

You can find this code in Our repository. The code of this article is actually in a single small Python file (which also includes the code for the next article).

The main ray tracer function

Here is the function that runs the ray tracer for a complete image:

def run_raytrace():
    nx = 2000
    ny = 1000
    array = np.empty((nx, ny, 3), dtype=np.uint8)
    lower_left_corner = np.array([-2, -1, -1], dtype=np.float64)
    horizontal = np.array([4, 0, 0], dtype=np.float64)
    vertical = np.array([0, 2, 0], dtype=np.float64)
    origin = np.array([0, 0, 0], dtype=np.float64)
    for j in range(ny):
        for i in range(nx):
            u = i / nx
            v = j / ny
            ray = np.vstack((
                origin,
                lower_left_corner + horizontal * u + vertical * v
            ))
            col = color_sphere(ray)
            ir = int(255.99 * col[0])
            ig = int(255.99 * col[1])
            ib = int(255.99 * col[2])
            array[i, ny - j - 1] =  (ir, ig, ib)
    return array

We have an image of size 2000x1000. Our image - hosted on a NumPy array - will be stored in the array variable. Each point in the image has 3 components - RGB - hence the 3 on the definition of array.

We then have a couple of for loops (a bad idea that we will discuss in later installments), to compute all pixels in our image.

Our ray will have an origin of [0, 0, 0] and then will have a direction depending on the position on our image - which is coded with a 3D vector. Hence a ray is a matrix of 2 (origin, direction) vectors (each with 3 components). For details about the underlying design check Peter Shirley’s books.

We then call a function, color_sphere, to return the RGB components for a certain pixel. The function - which also defines the objects in the scene - can be found below.

By the way, our first image is a red sphere on a blue horizon gradient, like this:

Determining the color of each pixel

To determine the color of each pixel we will use this coloring function:

def color_sphere(ray):
    if hit_sphere(np.array([0, 0, -1]), 0.5, ray):
        return np.array([1, 0, 0], dtype=np.float64)
    unit_direction = ray[1] / np.linalg.norm(ray[1])
    t = 0.5 * (unit_direction[1] + 1.0)
    return np.array([1.0, 1.0, 1.0]) * (1 - t) + np.array([0.5, 0.7, 1.0]) * t

The first two lines check if we are hitting our sphere on our ray. The last three lines compute a background blue gradient in case we are not hitting the sphere. Again, as ray tracing is not the issue of our discussion, you are referred to Peter Shirley’s books if you want more details.

The function to detect if we are hitting a sphere is:

def hit_sphere(center, radius, ray):
    oc = ray[0] - center
    a = np.dot(ray[1], ray[1])
    b = 2.0 * np.dot(oc, ray[1])
    c = np.dot(oc, oc) - radius * radius
    discriminant = b * b - 4 * a * c
    return discriminant > 0

That’s it, that is really all our code to generate the image above.

Housekeeping

In order to call this code we will be using typer, a library to generate command line interfaces based on typing information. You might be thinking: but the code above has no type annotations, then what is the point of using a library based on typing? Great question, something we will delve into at a later stage.

We will also use Pillow to generate images from arrays.

So, to run the code above we will have the following to wrap:

import numba
import numpy as np
import typer
from PIL import Image

app = typer.Typer()
@app.command()

def run_non_numba():
    array = run_raytrace()
    img = Image.fromarray(np.transpose(array, (1, 0, 2)))
    img.save("out.png")

app()

We will be using Poetry to maintain our code.

If you get the git repository you can run this code with:

poetry run python -m numba_raytracer.basic run-non-numba

Next steps

Here we introduced a bare bones Ray tracer based on NumPy that we will use to explore performance issues in the Python ecosystem

Our next post will discuss a Numba version of this code - which will be substantially faster than the version presented here.