Gaussian signal generation Edit: Could the following be the answer?


Generate WGN-like-signal which is centered around a set dBm value.
Treat that signal like it was a frequency domain representation of an unknown X time-domain signal.
Do some math to get he unit back into Volts/Hz from dBm/Hz (using reference impedance).
Run inverse FFT.


Original Question:

I am looking for some advise on how to generate an array of samples pulled from normal (or Gaussian) distribution.

But, and here's the deal.
I want to control its PSD or ASD "directly", assuming I am using a refernce impedence of 50 Ohms.

I would like to know Ahead-Of-Time (AOT), before signal gets generated, what STD Dev should be applied to achieve a set PSD or at least average spectral level in a given bandwidth/channel.

This value should be set by a user.

Code below should prove I know (more or less) what PSD or ASD is, I am just looking for an algorithm advise, maybe someone did something similar?

Here's how I calculate the PSD or ASD for a generated white-noise-like signal in a standard way using numpy.

import numpy as np
import matplotlib.pyplot as plt
import logging

if __name__ == "__main__":
    logger = logging.getLogger("Main")
    logging.basicConfig(format='%(levelname)s:%(message)s', level=logging.DEBUG)
    DRAW_PLOTS = False
    mu, sigma = 0, 0.001  # mean and std deviation
    N = 1000
    timestamp = 1  # 10 Hz sampling rate

    # This is in Volts and no timestamp
    s = np.random.normal(mu, sigma, N)
    # Because signal is in Volts
    # our result of fft is in Volts/Hz or would be if we say our timestamp for
    # signal `s` is eqault to 1 second
    # This noise_spectrum can be called Amplitude Spectral Density (ASD)
    # If we want PSD we have to go from ASD to PSD:
    # Square the voltage and also maybe apply some reference impedance
    # ... or don't if you want to end up with Voltz (^2) / Hz
    noise_spectrum = np.fft.fft(s)
    # The following is in dBV^2/Hz
    noise_spectrum_db = np.multiply(
        20,
        np.log10(noise_spectrum)
    )
    # The following is in dBW/Hz (divided by 50 Ohm)
    noise_spectrum_dbw = np.multiply(
        20,
        np.log10(noise_spectrum/50)
    )
    freq = np.fft.fftfreq(s.size, timestamp)
    logger.info(
        f"Avg channel power in dBW: {np.real(np.mean(noise_spectrum_dbw))}"
    )
    logger.info(
        f"Avg channel power in dBm: {np.real(np.mean(noise_spectrum_dbm))}"
    )
```
Below is a function which I wrote long back, when I needed to generate AWGN time-domain samples given Noise PSD in dBm/Hz.

AWGN_NOISE() : Generates Additive White Gaussian Noise of PSD power in dBm/Hz

AWGN has Gaussian PDF with 0 mean and $\sigma^{2} = N_{o}/2$

$NoisePSD_{dBm/Hz} = 10.log_{10}(\frac{N_o}{2.BW})$, Why?

Because, Output Noise Power(in dBm) : $$10.log_{10}(\frac{N_o}{2}) = NoisePSD_{dBm/Hz} + 10.log10(BW)$$

Hence, 
$$\sigma = \sqrt{\frac{N_o}{2}} = \sqrt{BW.10^{\frac{NoisePSD_{dBm/Hz}}{10}}}$$

So, the python function would be as follows:

def awgn_noise(length, power, Bandwidth):
    sigma = np.sqrt(Bandwidth * 10**(power/10))
    noise = np.random.normal(0, sigma, length) 
    return noise


Hope this answers the question.