Calculates the power spectral density.

PSD(series, len)

series |
- |
Any series or multi-column table. |

len |
- |
Optional. An integer, the FFT length. Defaults to the length of the input series. If |

A frequency domain series.

W1: gsin(100,1.0,0.2)*5;setvunits("V")

W2: psd(W1)

W3: psd(W1, 1024)

max(W2) occurs at 0.2 Hz. with amplitude 1.0*100*(5^2)/2 = 1250.

mean(W1*W1) == 12.50

sum(W1*W1) == 1250

sum(W2) == 1250

area(W2) == 8.33

area(W3) == 12.50

The true PSD for the series is continuous. Because the PSD function returns samples of the continuous power spectral density, we must use more samples for the area to approximate the mean squared power as demonstrated by W3.

PSD (Power Spectral Density) is now a built-in function rather than a macro. The PSD is normalized so the area of the PSD is ideally equal to the mean of the input series squared, i.e.:

area(PSD(s)) == mean(s*s) ideally

The PSD is calculated by the FFT and has the following form:

psd(s) = 2*deltax(s)*(mag(fft(s))^2)/length(s)

with a total of int(len/2) + 1 frequency values from 0 to Fs/2 Hz., where Fs is the sampling rate of the data (i.e. rate(s)). The first value (DC component) and the last value (at Fs/2, the Nyquist frequency) are not scaled by 2 to preserve Parseval’s theorem.

For Fs equal to 1.0, by Parseval’s theorem, the sum of the PSD terms equals the sum of the series squared, i.e.:

sum(psd(s)) == sum(s*s)

A sinewave of amplitude A, frequency F, sample interval T (where T = 1/Fs), and length L, yields a PSD with an amplitude of T L A^{2}/2 at frequency F. If the input series is in Volts, the resulting PSD has units of V^{ 2}/Hz.

If len is larger than the length of series, the series is zero padded to length len before calculating the PSD. If len is less than the series length, the series is truncated to length len. If not specified, len defaults to the length of series.

The PSD is often calculated by averaging overlapped segments of the input series. The Advanced DSP Module includes functions to perform such a calculation.

See POWSPEC to compute the power spectrum.

See SPECTRUM to compute a magnitude normalized FFT.

See SPECGRAM to compute a joint time-frequency distribution.

See WINFUNC for a list of windowing functions and amplitude correction schemes.

See INVPSD to construct a time series from a PSD.

See DADiSP/FFTXL to optimize the underlying FFT computation.

Oppenheim & Shafer

Digital Signal Processing

Prentice Hall, 1975

pp. 548-556

Oppenheim & Shafer

Discrete-Time Signal Processing

Prentice Hall, 1989

pp. 730-742

Digital Signal Processing Committee

Programs for Digital Signal Processing

IEEE Press, 1979

Section 2.1-1 - 2.1-10

S. Lawrence Marple, Jr.

Digital Spectral Analysis with Applications

Prentice Hall, 1987

pp. 152-158