u:
p:

# Search

python/noise.html 2020-02-05

# Noise

## White noise

Making noise in python is very simple. To make 1000 gaussian white noise samples do:

```#!/usr/bin/env python3
import numpy as np
n = 1000
np.random.normal(0.0,  1.0, n)
```

These samples are independant and have a gaussian distribution with mean = 0.0 and SD = 1.0:

Note that for noise to be white there is absolutely no requirement to have a gaussian distribution. The following yields (non gaussian) white noise:

```np.random.randint(0, 2, n)
```

## Pink noise

Pink noise requires making a spectrum with power spectral density (PSD) S = 1/f and reverse fourier transforming it. Its quite simple in Python:

```#!/usr/bin/python3
import numpy as np

def spectrum_noise(spectrum_func, samples=1024, rate=44100):
"""
make noise with a certain spectral density
"""
freqs = np.fft.rfftfreq(samples, 1.0/rate)            # real-fft frequencies (not the negative ones)
spectrum = np.zeros_like(freqs, dtype='complex')      # make complex numbers for spectrum
spectrum[1:] = spectrum_func(freqs[1:])               # get spectrum amplitude for all frequencies except f=0
phases = np.random.uniform(0, 2*np.pi, len(freqs)-1)  # random phases for all frequencies except f=0
spectrum[1:] *= np.exp(1j*phases)                     # apply random phases
noise = np.fft.irfft(spectrum)                        # return the reverse fourier transform
noise = np.pad(noise, (0, samples - len(noise)), 'constant') # add zero for odd number of input samples

return noise

def pink_spectrum(f, f_min=0, f_max=np.inf, att=np.log10(2.0)*10):
"""
Define a pink (1/f) spectrum
f     = array of frequencies
f_min = minimum frequency for band pass
f_max = maximum frequency for band pass
att   = attenuation per factor two in frequency in decibel.
Default is such that a factor two in frequency increase gives a factor two in power attenuation.
"""
# numbers in the equation below explained:
#  0.5: take the square root of the power spectrum so that we get an amplitude (field) spectrum
# 10.0: convert attenuation from decibel to bel
#  2.0: frequency factor for which the attenuation is given (octave)
s = f**-( 0.5 * (att/10.0) / np.log10(2.0) )  # apply attenuation
s[np.logical_or(f < f_min, f > f_max)] = 0    # apply band pass
return s
```

The code should be pretty self explenatory. The default attenuation gives standard 1/f pink noise. This is sometimes called 3 dB rolloff.

### Usage

To make a wav-file with this noise in the 5 - 10 kHz band scaled to the maximum volume ofg wav-files, use:

```from scipy.io import wavfile
x = spectrum_noise(lambda x:pink_spectrum(x, 5000, 10000), 4096)
x = np.int16(x/ max(abs(x)) * (2**15 - 1))
wavfile.write("pink.wav", 44100, x)
```

### Modulation

For 44100 samples of high pass noise above 140 Hz with 6 dB rolloff and 2 Hz modulation with random modulation phase, use:

```from scipy.io import wavfile
samples = 44100 # number of samples
rate = 44100    # rate (Hz)
f_mod = 2       # modulation frequency (Hz)
m = 1           # modulation depth
x = spectrum_noise(lambda x:pink_spectrum(x, 140, np.inf, 6.0), samples, rate)
x *= 1.0 + m * np.sin(2*np.pi * ( f_mod * (np.arange(0, samples)/rate + np.random.rand()) )) # modulation
factor = (2**13-1) / np.std(x)  # scale rms to 2^15 / 4 \approx -12 dBFS RMS
#factor = (2**15-1) / max(abs(x)) # scale peak to 2^15 = 0 dBFS peak
x = np.int16(x * factor)  # scale and convert to 16 bit integer

wavfile.write("pink.wav", rate, x)

```

### Finding the level

We now have a file nicely normalised to -12dBFS RMS. But how does this compare to the sound that this noise has to be added to? Plenty of programs that can find an audio file's RMS:

```normalize-audio -n sound.wav
ffmpeg -i sound.wav -filter:a volumedetect -f null /dev/null
ffmpeg-normalize -n sound.wav -f --print-stats -nt rms
```

Note that there are other way of determining the level of a file: peak and ebu.

### Normalisation (scaling)

```ffmpeg-normalize -nt ebu -t -12  -ofmt wav -f pink.wav -o pink2.wav
```

### Mixing

Now that we know how loud the noise has to become, we can add it to the sound:

```sox -m pink.wav sound.wav pink3.wav
sox pink.wav pink2.wav --plot gnuplot -V2 loudness  | gnuplot -p
```

### Misc commands

```cvlc --play-and-exit pink.wav # play with cvlc
python3 -c "from scipy.io import wavfile; print(wavfile.read('pink.wav').shape)" # show number of samples
```