I'm having difficulty figuring out how to convert from PSD (V^2/Hz) to dBV (1 Volt RMS reference). The goal is to have a 1 Volt RMS input show as 0 dBV.
Essentially I calculated the PSD of a signal with white noise using the Welch method, and need to plot dBV vs. Frequency.
Any help would greatly be appreciated!
First and foremost realize that "1V rms" corresponds to the total power (square root of it) in the signal while what we see in the FFT result is the power in each bin (and each bin acts like a narrow band filter), so given that it wouldn't make sense for a power spectral density given in power per unit frequency (such as dB/Hz) to scale the vertical axis by the total power in the waveform. If the noise is white (spread evenly across frequency), then the power will be distributed equally to each bin, and as we increase the number of bins, the amount in each bin will reduce accordingly (this is the "FFT noise floor"). The typical approach for a power spectral density is to normalize the results once properly compensated for the window (see detail below) by the resolution bandwidth, and then have the vertical units in dB/Hz as a spectral density. For resolution bandwidth of white noise signals I use the equivalent noise bandwidth (ENBW) of the window which formula is given at the bottom of this post.
Instead of using the Welch method consider computing the FFT directly on a properly windowed time domain signal. This method is simpler and suitable for many applications where the spectral averaging benefit of the Welch method is not needed (and the spectrum itself with higher resolution overall can itself still be averaged).
With this approach the following will properly scale all signals for levels of individual bins or total power computation of either signal bins or all noise bins.
First: normalize the DFT by the total number of samples $N$. This way each DFT bin directly represents the magnitude of its complex phasor in time, for example the two bins for a cosine will each be one half consistent with Euler's formula:
$$cos(\omega n) = 0.5 e^{j\omega n} + 0.5 e^{-j \omega n}$$
If the time domain signal is real and you want the DFT to represent only the positive half spectrum (Single-sided spectrum) with the rms magnitude of the time domain signal (0.707 for the case above), the result for the positive bins only is used and increased by $\sqrt{2}$, properly representing the root-mean-square of the positive bin with its complex conjugate.
Next: compensate for windowing losses.
To scale the vertical axis for individual discrete tones you can simply divide they DFT by the coherent gain of the window (detailed below). With this, the bins will be properly scaled within the scalloping loss of the window. This is simple and will visually properly place the peaks at the expected location, but will be in error based on scalloping loss (which can be minimized for better windows), and further would require a different compensation for noise samples or waveforms that have a spectral density that spans multiple bins. In many cases such as for visual appearance and scaling this would be adequate and sufficient.
Noise is adjusted by the non-coherent gain, and is given in dB/Hz as a noise density. The non-coherent gain will properly factor the resolution bandwidth of the window (a rectangular window has a bandwidth of $f_s/N$ in Hz where $f_s$ is the sampling rate in Hz, while any other window will have a wider bandwidth and thus let more noise into each bin.) The equivalent noise bandwidth for the window (ENBW below) is the equivalent brickwall filter that would let the same total power assuming white noise through to each DFT bin. Once we divide each FFT sample (bin) of white noise by the non-coherent gain, then each bin will properly represent the total power over 1 bin of bandwidth ($f_s/N$ in Hz where $f_s$ is the sampling rate in Hz). Thus if we want the plot to show the power in 1 Hz of bandwidth, we need to scale the plot by the bandwidth of 1 bin (reduce by $10Log_{10}(f_s/N)$ to then get a vertical axis properly in units of dB/Hz).
For accurate power measurement of signal and noise from a DFT spectrum, an approach computing the root sum square (complex conjugate square!) of all significant bins associated with a tone is recommended (as windowing will spread the tone across multiple bins), and then this root sum of squares result is then divided by the non-coherent gain of the window, which will provide a very accurate result of the rms level with no scalloping loss effects whatsoever. Using the non-coherent gain properly accounts for the resolution bandwidth of the window, which for all windows other than a rectangular window spans multiple bins. And conveniently all samples are treated identically in processing to determine power (noise and signal) with the same adjustment by non-coherent gain.
To provide the results in dB in either case, since these are magnitude (rms) quantities the calculation would be done using $20\log10()$, except for ENBW which if given in dB would $10\log10()$. (As Dan Freitas has pointed out to me, ENBW is typically not thought of in dB quantities, but it can be convenient when working with data already in dB).
Summary of Window Compensation Factors
Note that the following apply to the non-zero window specifically and do not include additional zero padding that may or may not also be used. For the case of multiplying the waveform by a window (such as Kaiser, Hamming, Hann etc) and zero padding) the window is to match the length of the data prior to zero padding. Zero padding if used will then interpolate more samples in the resulting DFT without changing the resulting magnitudes or the resolution bandwidth.
Coherent Gain (Window Average): $$G_c = \frac{\sum_{n=0}^{N-1} w[n]}{N}$$
Noncoherent Gain (RMS of Window): $$G_c = \sqrt{\frac{\sum_{n=0}^{N-1} w[n]^2}{N}}$$
Scalloping Loss is the worst case variation between a tone at bin center and mid-way between bins, and is computed as follows:
The coherent gain for a tone mid-way between bins is given as:
$$G_{hc} = \frac{|\sum_{n=0}^{N-1} w[n]e^{j \pi n/N}|}{N}$$
And the scalloping loss is the ratio of the coherent gain mid-way between the bins to the coherent gain at bin center (or difference in dB):
$$G_s = \frac{G_{hc}}{G_c}$$
Processing Gain:
$$G_p = \frac{G_c}{G_{nc}}$$
Equivalent Noise Bandwidth:
$$\text{ENBW} = \frac{1}{G_p^2}$$
For further smoothing or averaging of the spectrum, the equivalent of video bandwidth on a spectrum analyzer, I recommend using zero phase post processing filters on the DFT output such as filtfilt which is available in MATLAB, Octave and Python scipy.signal. This will allow for filtering the frequency domain signal without shifting the frequencies.
