Brownian Noise

When we attempt to model a physical phenomenon for entertainment (or simulation) purposes, it’s less a question of perfectly fitting the real-world process that creates it, and more a question of getting close enough for the communication we intend. This is certainly true for Brownian motion.

For those who aren’t familiar, Brownian motion is the wobbly oscillation we tend to see in falling snow, or dust motes floating through the air. They repel each other slightly, and carry a form of influence throughout the entire medium; but otherwise act independently. Other great examples include cloud diffusion, or sediment flow.

Obviously this is all over the place in modern games, but if we were to simulate it exactly, we would have serious issues. We would need to keep track of the properties of a crippling number of particles, along with their influence on neighboring particles. This is possible in theory, but not readily in practice.

So instead, we take on what’s called a Lorentzian model (as opposed to Newtonian, which would be keeping track of each individual particle), and just consider the net flow!

I’m going to go over some modern C code for what’s called Fractional (or Fractal) Brownian Motion, fBm for short, here. Then, I’m going to explain it, bit by bit. The only thing I’m going to leave in the abstract is the noise function, which can be any random number generator of your choice; an LCG, /dev/random, whatever you would like.

Every modern language I know of seems to have at least one LCG random function built into its libraries, so this shouldn’t be a problem for you; just understand that noise(p) takes the point p as a parameter and uses it to generate a difficult-to-predict value.

(If you must generate one, and all you have is an LCG expecting a scalar number, I suggest simply multiplying each dimension by some factor, running it through the same LCG with its value as the seed, and popping the result into a return vector. It doesn’t have to be perfect, just convincing.)

We’ll need a handful of parameters to get this to work in a customizable fashion. Aside from the point we’re trying to generate this noise at, we’ll also want the fractal increment, the lacunarity, and the number of noise octaves. That probably threw many readers for a loop, as it isn’t at all obvious what those mean; so let’s discuss them in turn.

Octaves

Fractal Brownian motion is technically not one layer of noise, but several, summed, and usually normalized (which just means spanned from zero to one). Each component noise is a wave of a specific frequency.

If all of the waves were of identical frequency, we would lose the character of the fBm, and just have what equated to a scaled noise pattern. So, with each component wave, we need to increase the frequency by a certain amount.

“Octave” is a term most are more familiar with musically—the number after a letter note, like C6, or D7, is its relative octave on a piano. The catch is, a difference of one octave is a doubling of frequency—as middle C is 262 Hz, the C one octave above it (C5) is 524 Hz. The same principle applies for our octaves in noise.

(In fact, you actually can tie this to something like an LV2 and create some weird, by which I mean awesome, sound effects; but that’s much harder to demonstrate on a web page, so we’ll be sticking with graphics for now!)

The more octaves you add to your noise function, the more detailed and winding the result is. But, this brings us to our next point. We could just assume that each layer of noise is an octave higher in its sampling, but that does feel a little limiting. What if it was half an octave? Or two octaves? This has a noticeable affect on the form of the noise, so we’ll want it to be adjustable. This is the…

Lacunarity

So, if we wanted to go purely by octave, we would simply use a frequency equal to 2^x, where x is some integer. We might have a few constants and factors in there to muddy things up to our liking, but this would be the core of it. While there are plenty of arguments for simplicity’s sake for sticking with 2^x styled values, when we aren’t necessarily going for simplicity, we aren’t beholden to it at all.

I’m going to mathematically refer to that two number as L. It’s our lacunarity, and it doesn’t have to be two. Lacunarity is an interesting term for it, as it’s actually derived from the Latin lacuna which means “lake”. As lacunarity increases, the size of the gaps in our fractal increase in turn, so it’s less mathematically defined as ontogenetically defined—that is, how it looks, not what it is.

(In fact, lacunarity actually has a very context-dependent definition; in the study of fractals, fractals are often classed by a lacunarity property which can be geometrically measured, with something like the box counting algorithm. This is not necessarily always related to our use of it here.)

The last element to look at, given our understanding of the need for multiple octaves, lacunarity, and a noise generator, is what’s going on in x. There’s nothing saying that it has to be an integer, or even a typical flat value. This ends up being critical for fBm, too. We refer to x with the…

Fractal Increment

Or, for short, H. Generally, in our octave multiplier L^x, x=-H × octave number. So, where i is an integer in the range of zero to the number of octaves, each noise value is multiplied by an amplifier of L^-Hi.

(I believe the choice for H for the fractal increment relates to its similarity to the fractal Hausdorff dimension, but can’t be certain of this at the moment. We’ll use it out of tradition.)

The rationale behind the negative is to avoid the typically inconvenient issue of a positive power, as these don’t converge to zero at infinity. I won’t lie, sometimes that’s actually what we want; but usually we don’t and it leads to very strange behaviors. Feel free to set your increment to a negative number, giving us a positive exponent, if you want to see the result!

This gives us quite a few elements which we can slide over for the sake of controlling our Brownian noise. Let’s look at our final math for it.

Final Mathematics

So, we’re iterating up to k octaves. For each octave, we will be adjusting our sample point, conveniently by multiplying it by the lacunarity. We’ll acquire noise from that sample point, and then multiply it by a factor of L^-iH, where i is our current octave number. That is, we’ll be doing this:

Fractal Brownian Motion, in a Nutshell

Don’t worry if that looks scary, we’ll go over it in basic procedural code. It isn’t so bad.

We’ll start with something easy to understand, but not quite perfect, and work from there.

#include <math.h>

double fBm( Vector point, double H, double lacunarity, double octaves )
{
    double value = 0.0;
    
    for (int i = 0; i < octaves; i++) {
        value += Noise( point ) * pow( lacunarity, -H*i );
        point *= lacunarity;
    }
    
    return value;
}

(Again, remember that I’ve not defined what Noise is here; that’s just a function that takes a vector value and returns a random float. There are oodles to choose from, so pick what’s most convenient to you, in your environment. In fact, it could arguably be any function period, but the squirrellier the results, the better.)

Let’s discuss it in detail. We take a basic Vector structure called point, and iterate over octaves with it. For each iteration, we add noise, seeded with our current point, to our final value. But, first we multiply it by our lacunarity to the power of negative increment times the octave number.

That amplitude multiplier of L^-iH ensures that each octave is smaller than the previous octaves in magnitude. Remember that if this were not negative, each value might be successively larger. There are uses for that, but they aren’t typical and aren’t usually useful.

Then, to take care of the variation of our sample point with each octave, we multiply it by lacunarity, so each call to noise will return a different value.

(If you would like to test the behavior of this function and don’t have a Vector class or struct handy, feel free to just use a double for it and watch the behavior of the resulting values over an iteration. You’ll just be seeing it in one dimension, but it will all still be there.)

We aren’t quite done, though. What if we have a non-integer number of octaves? Currently, the behavior of the function is floored to the integer immediately beneath the specified number of octaves. Thankfully, this is relatively simple to fix.

double fBm( Vector point, double H, double lacunarity, double octaves )
{
    double value = 0.0;
    
    int i;
    for (i=0; i<octaves; i++) {
        value += Noise( point ) * pow( lacunarity, -H*i );
        point *= lacunarity;
    }
    
    double remainder = octaves - (int)octaves;
    value += remainder * Noise( point ) * pow( lacunarity, -H*i );
    return value;
}

Let’s discuss my boldfaced additions here. I’ve declared i in the next-higher scope, because we’re going to need to refer back to it.

We get the fractional number of octaves by subtracting the integer value, which is traditionally the floor value, of octaves from it. (Some languages have a frac function that takes care of this, but this is usually what it amounts to internally.)

We then add this remainder value, multiplied by noise on our current reference point (and I remind you, this is different from the last one, as it’s already been multiplied by lacunarity), multiplied by L^-Hi, to our final value.

Technically, this is it—you’ve got unbounded fractal Brownian motion. Dependent on your noise function, this could be a set of continuous vectors, a cloud of densities, or whatever you would like; the fBm part is primarily in the octaves method. However, if you’re going for something like densities, you may want to bound this value to [0, 1], so you can properly display it on a monitor.

How you choose to bound it likely depends on what you’re doing with it, but by dividing by the sum of maximum possible values, we can ensure that it’s always between zero and one.

double fBm( Vector point, double H, double lacunarity, double octaves )
{
    double value = 0.0;
    double sum = 0.0;
    
    int i;
    for (i=0; i<octaves; i++) {
        value += Noise( point ) * pow( lacunarity, -H*i );
        sum += pow(lacunarity, -H*i);
        point *= lacunarity;
    }
    
    double remainder = octaves - (int)octaves;
    value += remainder * Noise( point ) * pow( lacunarity, -H*i );
    sum += remainder * pow(lacunarity, -H*i);
    
    value /= sum;
    
    return value;
}

Of course, you may just want it centered around the median or mean, which would involve keeping track of maximums and minimums for each value, but this is easily done with another variables and a comparison operator. It really comes down to the end user!

I would like to thank the F. Kenton Musgrave for introducing me to this function on a technical level. Most of the changes to it are largely cosmetic, and the chapter is a great read in any case. Musgrave, F. Kenton, (2002) Texturing and Modeling — A Procedural Approach, Third Edition, An Introduction to Fractals, pg. 429-445

The Results (in OSL)

I’ve rewritten this function in OSL (Open Shading Language) and implemented it in Cycles, to demonstrate its results. OSL conveniently includes a noise function which accepts any form of vector or color as a parameter. The final material is thus:

#include <stdosl.h>

shader fBM(vector uv = 0, float octaves = 1, float lacunarity = 2, float H = 1, output float hue = 0)
{
    float value, remainder;
    int i;
    
    value = 0.0;
    vector p = uv;
    
    for(i = 0; i < octaves; i++) {
        //*additive* cascade
        value += noise(p) * pow(lacunarity, -H * i);
        p *= lacunarity;
    }
    
    remainder = octaves - (int)octaves;
    if(remainder) {
        value += remainder * noise(p) * pow(lacunarity, -H * i);
    }
    
    hue = value;
}

It’s likely that, if you’re familiar with C, you understand most of this too. I’ll go over the minor differences—each parameter to the shader has to have a default value, which is listed after the parameter with an equals sign and then its respective value. The output parameter just lists an output field from our shader, which is later plugged into an HSV node to produce a color.

I haven’t bothered to normalize it in script, as hue is circular anyway and will normalize itself. (Processor cycles are precious, after all!)

If we set UV to the object coordinate, octaves to nine, lacunarity to a basic two, and H to 0.5, we get this:

fBm demo’d over a hue space

Drop our lacunarity to one, and we have this:

A simple drop of lacunarity can dramatically change the effect!

If we drop our increment, H, to 0.1, we get a much noisier image.

Dropping our increment increases the frequency of change dramatically.

To demonstrate why our increment is multiplied by a negative, let’s drop it to a negative value itself to cancel it out. If we drop from 0.5 to -0.5, we end up with unbounded values, which while it sounds cool ultimately just gives us grey noise of infinite frequency.

Believe it or not, this is extremely fine color noise from a negative increment.

I promise that there are uses for this kind of behavior, but as you might guess, there are precious few.

Lastly, dropping the number of octaves simply removes fine detail from the fBm, and gives us a much smoother cloud that follows the same pattern.

Fewer octaves intuitively results in more regularity.

Hopefully, you have an idea how fractal Brownian motion actually works now. You can alter its input parameters to skew it in different ways, or apply it to any number of values. If we wanted a normalized gradient field, as an example, we could simply pipe our “hue” into the polar coordinates of the vector. It’s applicable in any number of dimensions—which is part of why I’m not picking what noise is for you—and can provide a convincing emulation of any number of real-world fields.