COLORED_NOISE
1/F^alpha Power Law Noise Generation

COLORED_NOISE is a FORTRAN90 library which generates sequences that simulate 1/f^alpha power law noise. This includes white noise (alpha = 0), pink noise (alpha = 1) and brown noise or Brownian motion (alpha = 2), but also values of alpha between 0 and 2.

The original code listing by Kasdin referenced a number of functions from the Numerical Recipes library (FOUR1, FREE_VECTOR, GASDEV, RAN1, REALFT, VECTOR). Numerical Recipes is a proprietary library whose components cannot be freely distributed. Moreover, the referenced functions have some peculiarities of implementation (the FFT is not normalized; the FFT data must be a power of 2 in order; the code uses single precision real arithmetic).

In the code presented here, references to Numerical Recipes functions have been replaced by references to, and the source code of, nonproprietary code. In particular, the Fourier transform is implemented by a "slow Fourier transform" method, and by simple uniform and normal random number generators. You are welcome to make an efficient code by replacing these routines. The purpose of this posting is primarily to demonstrate the method.

Licensing:

The computer code and data files described and made available on this web page are distributed under the GNU LGPL license.

Languages:

COLORED_NOISE is available in a C version and a C++ version and a FORTRAN90 version and a MATLAB version and a Python version.

Related Data and Programs:

CORRELATION, a FORTRAN90 library which contains examples of statistical correlation functions.

NORMAL, a FORTRAN90 library which computes elements of a sequence of pseudorandom normally distributed values.

PINK_NOISE, a FORTRAN90 library which computes a "pink noise" signal obeying a 1/f power law.

SDE, a FORTRAN90 library which illustrates the properties of stochastic differential equations, and common algorithms for their analysis, by Desmond Higham;

SFTPACK, a FORTRAN90 library which implements the "slow" Fourier transform, intended as a teaching tool and comparison with the fast Fourier transform.

STOCHASTIC_RK, a FORTRAN90 library which applies a Runge-Kutta scheme to a stochastic differential equation.

UNIFORM, a FORTRAN90 library which computes elements of a uniform pseudorandom sequence.

Reference:

1. Martin Gardner,
   White and brown music, fractal curves and one-over-f fluctuations,
   Scientific American,
   Volume 238, Number 4, April 1978, pages 16-32.
2. Jeremy Kasdin,
   Discrete Simulation of Colored Noise and Stochastic Processes and 1/f^a Power Law Noise Generation,
   Proceedings of the IEEE,
   Volume 83, Number 5, 1995, pages 802-827.
3. Edoardo Milotti,
   1/f noise: a pedagogical review,
   arXiv:physics/0204033.
4. Sophocles Orfanidis,
   Introduction to Signal Processing,
   Prentice-Hall, 1995,
   ISBN: 0-13-209172-0,
   LC: TK5102.5.O246.
5. William Press,
   Flicker Noises in Astronomy and Elsewhere,
   Comments on Astrophysics,
   Volume 7, Number 4, 1978, pages 103-119.
6. Miroslav Stoyanov, Max Gunzburger, John Burkardt,
   Pink Noise, 1/f^alpha Noise, and Their Effect on Solutions of Differential Equations,
   International Journal for Uncertainty Quantification,
   Volume 1, Number 3, pages 257-278, 2011.

Source Code:

• colored_noise.f90, the source code.

Examples and Tests:

• colored_noise_test.f90 a sample calling program.
• colored_noise_test.txt, the output file.

ALPHA is a sequence of 9 sets of data, for ALPHA = 0.00 to ALPHA = 2.00, each time computing a sequence of 128 values, and using the same sequence of random numbers. This makes the effect of ALPHA very clear.

• make_plots.m a MATLAB program which reads the "txt" files and makes plots.
• alpha_0.00.txt values for ALPHA = 0.00.
• alpha_0.00.png, a plot.
• alpha_0.25.txt values for ALPHA = 0.25.
• alpha_0.25.png, a plot.
• alpha_0.50.txt values for ALPHA = 0.50.
• alpha_0.50.png, a plot.
• alpha_0.75.txt values for ALPHA = 0.75.
• alpha_0.75.png, a plot.
• alpha_1.00.txt values for ALPHA = 1.00.
• alpha_1.00.png, a plot.
• alpha_1.25.txt values for ALPHA = 1.25.
• alpha_1.25.png, a plot.
• alpha_1.50.txt values for ALPHA = 1.50.
• alpha_1.50.png, a plot.
• alpha_1.75.txt values for ALPHA = 1.75.
• alpha_1.75.png, a plot.
• alpha_2.00.txt values for ALPHA = 2.00.
• alpha_2.00.png, a plot.

List of Routines:

• F_ALPHA generates a 1/F^ALPHA noise sequence.
• R8_NORMAL_01 returns a unit pseudonormal R8.
• R8_UNIFORM_01 returns a unit pseudorandom R8.
• R8VEC_PRINT prints an R8VEC.
• R8VEC_PRINT_PART prints "part" of an R8VEC.
• R8VEC_SFTB computes a "slow" backward Fourier transform of real data.
• R8VEC_SFTF computes a "slow" forward Fourier transform of real data.
• TIMESTAMP prints the current YMDHMS date as a time stamp.

You can go up one level to the FORTRAN90 source codes.