Random Number Generators

Pseudo Random

To generate random numbers, the simple solution would be to use a simple pseudo random number generator, in the code we have a Linear Congruential Generator with the same parameters as glibc. It works alright, though it's expected the quality to be lower and period shorter than it could be. The generator state is maintained per pixel and initialized at the start with random(). Here's the code:

typedef unsigned int RNG;

__device float path_rng(RNG *rng)
    *rng = (1103515245*(*rng) + 12345);
    return (float)*rng/(float)0xFFFFFFFF;

Quasi Random

We can do much better with quasi random number generators, some useful papers on the topic:

We use Sobol sequences, since they fit our needs quite well. We need something that runs well on the CPU and GPU, supports high path depth and can do adaptive sampling. Let's start by trying to get a clear picture of what we need.


Ideally we would like a function that looks like this:

__device float path_rng(int pass, int pixel_x, int pixel_y, int dimension);

This would give us the random number for a given sample without having to maintain any kind of state. The pass here is the rendering pass, which is increased each time a new sample is taken for that pixel. The pixel_x and pixel_y values are self explanatory, but dimension needs a bit more explanation.

We give each random number used in tracing a path a well defined unique dimension. We have two dimensions for the pixel filter, two for depth of field, and then for each bounce, one to pick a BSDF, two to sample the BSDF, ... .

enum PathTraceDimension {
    PRNG_FILTER_U = 0,
    PRNG_FILTER_V = 1,
    PRNG_LENS_U = 2,
    PRNG_LENS_V = 3,
    PRNG_BASE_NUM = 4,

    PRNG_BSDF_U = 0,
    PRNG_BSDF_V = 1,
    PRNG_BSDF = 2,
    PRNG_LIGHT = 3,
    PRNG_LIGHT_U = 4,
    PRNG_LIGHT_V = 5,
    PRNG_LIGHT_F = 6,

From this we can easily compute the dimension each time we need to get a random number from the generator, for example:


Sobol Sequences

Sobol sequences have the property that they can be used for progressive sampling. Unlike the Halton sequence that can also be used for this, they are not correlated in higher dimensions, and so do not need to be scrambled. Any number from any dimension can be queried without per path precomputation.

Each dimension has it's own sequence, and when rendering the i-th pass, we get element i from the sequence. A sequence is defined by a 32x32 binary matrix, and getting the i-th element in the sequence corresponds to multiplying the matrix by i. With binary operations this ends up being pretty quick:

float sobol(unsigned int vectors[][32], unsigned int dimension, unsigned int i)
    unsigned int result = 0;

    for(unsigned int j = 0; i; i >>= 1, j++)
        if(i & 1)
            result ^= vectors[dimension][j];

    return result * (1.0f/(float)0xFFFFFFFF);

These matrices are not as simple to compute, but luckily the hard work has been done for us, and the data to generate them up to dimension 21201 is available online.

We currently use 4 dimensions initially to sample the subpixel location and lens and 8 numbers per bounce, so that limits us to a maximum path depth of 2649.

At the moment we are using the same sequences in each pixel, decorrelated using a Cranly-Patterson rotation (a simple random 2D shift different for each pixel). This may not be optimal, and there's some experiments in the code with full screen nets, but they give too much correlation.