Understanding Voss-McCartney pink noise generation algorithm

Question

I'm implementing the Voss-McCartney pink noise generation algorithm.

If you follow the link above, you can read:

from James McCartney 2 Sep 1999 21:00:30 -0600:

The top end of the spectrum wasn't as good. The cascade of sin(x)/x shapes that I predicted in my other post was quite obvious. Ripple was only about 2dB up to Fs/8 and 4dB up to Fs/5. The response was about 5dB down at Fs/4 (one of the sin(x)/x nulls), and there was a deep null at Fs/2. (These figures are a bit rough. More averaging would have helped.)

You can improve the top octave somewhat by adding a white noise generator at the same amplitude as the others. Which fills in the diagram as follows:

xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
x x x x x x x x x x x x x x x x
 x   x   x   x   x   x   x   x
   x       x       x       x
       x               x
               x

It'll still be bumpy up there, but the nulls won't be as deep.

If I understand it well, this algorithm generates pink noise by adding random (white?) noise sources at different frequencies¹

However, I don't fully understand the explanation given in the quote above for the extra white noise generator on the "top row". Can someone clarify how/why it improves the algorithm? Does that make it a good algorithm for pink noise generation for audio applications? Especially, shouldn't I discard the first samples until all the "rows" were mixed into the signal (in the ASCII art quoted above, that would mean discarding 15 first samples)?

---

_{¹ I'm not sure of the wording here. Do not hesitate to correct me if I'm wrong}

Answer 1

First, let’s separate the Voss part from the McCartney. The first generates a 1/f distribution of random numbers, with power inversely proportional to frequency—pink noise. McCartney proposed a change that gives a more flat computational load.

A brief overview of Voss:

Start with random number generators, one for each bit in a binary counter. Let’s consider a simple case of three:

000
001
010
011
100
101
110
111
(and continuing with rollover back to 000)

The counter increments by 1 for each new output sample. There are three random number generators, one for each column. When the digit in a column changes, its corresponding generator produces a new random number. Otherwise, the generator maintains its previous value.

As the final step, for each cycle, all generators are summed to produce the output of the process.

The computational load is uneven. Every other cycle has only a single generator change, while every fourth time all three change (011 -> 100, and 111 -> 000). Practical pink source need more generators, so the worst case gets worse.

McCartney’s method reworks the process for an even load.

Additional thoughts:

I believe that Voss’ goal was not in generating an audio-rate pink noise signal, but a pink sequence that mimics the qualities of slower natural changes, including changes in loudness and pitch of speech and music.

I learned of the algorithm in Friendly Stories About Computers / Synthesizers, by John Simonton Jr., in which he cites the Gardner article and Voss, and lays out the algorithm that he attributes to Voss: Five, “four-side dice”—using a 5-bit counter, each bit-change of which generates a random integer value 0-3. The output is a value 0-15, but having the quality that changes will typically be small, but less often the changes can be large—at times rolling all dice.

For instance, you might use this as an index into a musical scale. In that way, most often the next note is near the last one played, scale-wise, but at times there is the change of bigger leaps, as more of the dice roll. I used this algorithm about 30 years ago as an example included with my HyperMIDI product (MIDI scripting for HyperCard). For this type of use, it seems that McCartney’s load balancing defeats the purpose.

I have not looked at the quality versus computation for an audio noise source, but the method of using a 3 dB/octave filter on white noise works very well, which requires just a few poles and zero.

Answer 2

So let's look at what the author of the article you linked to says further down; Output samples are on the top row, and are the sum of all the other rows at that time.

Output  /---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\/---\
        \___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/\___/
Row -1  /---/---/---/---/---/---/---/---/---/---/---/---/---/---/---/---/---/---

        ___/___/___/___/___/___/___/___/___/___/___/___/___/___/___/___/___/___/
Row 0   /--------/--------/--------/--------/--------/--------/--------/--------/--------

        ________/________/________/________/________/________/________/________/________/
Row 1   --------------/------------------/------------------/------------------/--------------
        ______________/__________________/__________________/__________________/______________
Row 2   ------------------------/--------------------------------------/------------------------
        ________________________/______________________________________/________________________
Row 3   --------------------------------------------/--------------------------------------------
        ____________________________________________/____________________________________________
Row 4   ------------------------------------------------------------------------------------/----
        ____________________________________________________________________________________/____

This means that the above diagram has multiple different white sequences, which they only occasionally change – let's formalize that. Start with only the two top rows:

• Row -1 is simply white noise
• Row 0 is white noise, interpolated by a factor of 2 with a 2-sample-boxcar-filter / sample-and-hold. That gives that noise an (aliased) sinc shape, which is essentially a low-pass shape

Row 1…N do the same, with the sincs becomming narrower by factors of 2.

Thinking about the discrete PSD of this:

• Row -1 has a constant discrete PSD
• Row 0 adds sinc(2f)²-shaped power to that
• Row 1 adds sinc(4f)²-shaped power to that
• and so on

All in all, I don't have a proof that this becomes perfectly pink at hand, it probably doesn't within finite observation, but it's kind of intuitive to think that close to 0 Hz, all the main lobes of these sinc²s add up, and with every doubling of frequency, you get closer to the zeros of more sinc²s.

The proposed algorithm really doesn't seem so elegant – generating good (discrete) white (pseuderandom) noise is actually surprisingly hard for longer observational windows (which is what you need if you want to assess the quality of something), and hence, having a pseudorandom generator¹ run at asymptotically twice the sampling rate seems more effort then letting it run at the sampling rate and then using an appropriate low-pass filter that approximates the desired spectral shape (in this case, $\lvert H(f)\rvert \propto \frac1f$); at least on modern CPUs, which have excellent SIMD instructions (i.e. highly optimized for running filters, not so much for running pseudo-random noise generators), the difference between holding and adding up many noise values and doing a FIR is that the FIR requires multiplication of held values with constants (the filter taps) – and since that can typically done in a fused multiply-accumulate operation.

Now, on an ASIC or FPGA, things might look different; if the amplitude distribution of the noise doesn't matter (i.e. there's no need to add up anything but uniformly drawn, uncorrelated samples), then you can actually save on complexity by doing the "simpler" thing, i.e. logical operations needed to generate e.g. XOROSHIRO128** would very likely be clocked much higher than the multipliers needed for a nice FIR filter.

---

¹you don't need multiple generators – you just ask that one white one more often; white samples are uncorrelated in every subsampling!