Add basic noises

lisainstrument/__init__.py

@@ -5,3 +5,9 @@
 from .meta import __version__
 from .meta import __author__
 from .meta import __email__
+
+from .noises import white
+from .noises import violet
+from .noises import pink
+from .noises import red
+from .noises import infrared

lisainstrument/noises.py

@@ -3,6 +3,65 @@
 """
 Noises module.
 
+Implements basic random noise generators, and use them to implement instrumental noises.
+
 Authors:
     Jean-Baptiste Bayle <j2b.bayle@gmail.com>
 """
+
+import logging
+import numpy
+import pyplnoise
+
+from numpy import pi, sqrt
+from lisaconstants import c
+
+
+def white(fs, size, psd):
+    """Generate a white """
+    logging.debug("Generating white noise (fs=%s Hz, size=%s, psd=%s)", fs, size, psd)
+    if not psd:
+        logging.debug("Vanishing power spectral density, bypassing noise generation")
+        return 0
+    generator = pyplnoise.WhiteNoise(fs, psd / 2)
+    return generator.get_series(size)
+
+
+def violet(fs, size, psd):
+    """Generate a violet noise in f^2."""
+    logging.debug("Generating violet noise (fs=%s Hz, size=%s, psd=%s)", fs, size, psd)
+    if not psd:
+        logging.debug("Vanishing power spectral density, bypassing noise generation")
+        return 0
+    white_noise = white(fs, size, psd)
+    return numpy.gradient(white_noise, 1 / fs) / (2 * pi)
+
+
+def pink(fs, size, psd):
+    """Generate a pink noise in f^(-1)."""
+    logging.debug("Generating pink noise (fs=%s Hz, size=%s, psd=%s)", fs, size, psd)
+    if not psd:
+        logging.debug("Vanishing power spectral density, bypassing noise generation")
+        return 0
+    generator = pyplnoise.PinkNoise(fs, 1 / size, fs / 2)
+    return sqrt(psd / 2) * generator.get_series(size)
+
+
+def red(fs, size, psd):
+    """Generate a red noise (also Brownian or random walk) in f^(-2)."""
+    logging.debug("Generating red noise (fs=%s Hz, size=%s, psd=%s)", fs, size, psd)
+    if not psd:
+        logging.debug("Vanishing power spectral density, bypassing noise generation")
+        return 0
+    generator = pyplnoise.RedNoise(fs, 1 / size)
+    return sqrt(psd / 2) * generator.get_series(size)
+
+
+def infrared(fs, size, psd):
+    """Generate an infrared noise in f^(-4)."""
+    logging.debug("Generating infrared noise (fs=%s Hz, size=%s, psd=%s)", fs, size, psd)
+    if not psd:
+        logging.debug("Vanishing power spectral density, bypassing noise generation")
+        return 0
+    red_noise = red(fs, size, psd)
+    return numpy.cumsum(red_noise) * (2 * pi / fs)