cutcutcodec.core.signal.psd.welch¶

cutcutcodec.core.signal.psd.welch(signal_1: Tensor, signal_2: Tensor | None = None, band: Real | None = None)[source]¶

Estimate the power spectral density (PSD) ie intercorrelation with the Welch method.

It is based on cutcutcodec.core.signal.psd.intercorr().

The slices are smoothed using a Kaiser window, calculated with cutcutcodec.core.signal.window.kaiser().

\[\begin{split}\begin{cases} \Delta f = \frac{1 + 2 \beta}{n_w} \\ \end{cases}\end{split}\]

Parameters¶

signal_1torch.Tensor: A real stationary time signal.
signal_2default = signal_1: Another real stationary time signal. If provided, the intercorrelation is calculated. The returned signal is therefore a complex signal. If omitted (default), the autocorrelation is calculated, and the returned signal is therefore real positive.
bandfloat, optional: The normlised frequency resolution (Delta f) in \(\left(r \in \left]0, \frac{1}{2}\right[\right)\), for a sample rate of 1. Higher it is, better is the frequency resolution but the greater is the variance.

Returns¶

psdtorch.Tensor: The correlation of signals, estimated using the Welch method. In the case of autocorrelation, this is an estimation of the power spectral density.

Seealso¶

scipy.

Examples¶

>>> import math
>>> import torch
>>> import matplotlib.pyplot as plt
>>> from cutcutcodec.core.signal.psd import welch

Generate a test signal, a 2 Vrms sine wave at 1234 Hz, corrupted by 0.001 V**2/Hz of white noise sampled at 10 kHz.

>>> fs, N, amp, freq = 10e3, 1e5, 2*math.sqrt(2), 1234.0
>>> noise_power = 0.001 * fs / 2
>>> time = torch.arange(N) / fs
>>> x = amp*torch.sin(2*torch.pi*freq*time)
>>> x += torch.randn(time.shape) * math.sqrt(noise_power)

Compute and plot the power spectral density.

>>> psd = welch(x, band=2/1024)
>>> f = torch.fft.rfftfreq(2*psd.shape[-1]-1, 1/fs)
>>> _ = plt.semilogy(f, psd)
>>> _ = plt.xlabel('frequency [Hz]')
>>> _ = plt.ylabel('PSD [V**2/Hz]')
>>> # plt.show()