ml4gw.spectral

Several implementation details are derived from the scipy csd and welch implementations. For more info, see

https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.welch.html

and

https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.csd.html

Functions

fast_spectral_density(x, nperseg, nstride, ...)

Compute the power spectral density of a multichannel timeseries or a batch of multichannel timeseries, or the cross power spectral density of two such timeseries.

median(x, axis)

Implements a median calculation that matches numpy's behavior for an even number of elements and includes the same bias correction used by scipy's implementation.

normalize_by_psd(X, psd, sample_rate, pad)

spectral_density(x, nperseg, nstride, ...[, ...])

Compute the power spectral density of a multichannel timeseries or a batch of multichannel timeseries.

truncate_inverse_power_spectrum(psd, ...[, ...])

Truncate the length of the time domain response of a whitening filter built using the specified psd so that it has maximum length fduration seconds.

whiten(X, psd, fduration, sample_rate[, ...])

Whiten a batch of timeseries using the specified background one-sided power spectral densities (PSDs), modified to have the desired time domain response length fduration and possibly to highpass/lowpass filter.

ml4gw.spectral.fast_spectral_density(x, nperseg, nstride, window, scale, average='median', y=None)

Compute the power spectral density of a multichannel timeseries or a batch of multichannel timeseries, or the cross power spectral density of two such timeseries. This implementation is non-exact for the two lowest frequency bins, since it implements centering using the mean for the entire timeseries rather than on a per- window basis. The benefit of this is a faster implementation, which might be beneficial for cases where the lowest frequency components are going to be discarded anyway.

Parameters:
  • x (Union[Float[Tensor, 'time'], Float[Tensor, 'channel time'], Float[Tensor, 'batch channel time']]) -- The timeseries tensor whose power spectral density to compute, or for cross spectral density the timeseries whose fft will be conjugated. Can have shape (batch_size, num_channels, length * sample_rate), (num_channels, length * sample_rate), or (length * sample_rate).

  • nperseg (int) -- Number of samples included in each FFT window

  • nstride (int) -- Stride between FFT windows

  • window (Float[Tensor, '{nperseg//2+1}']) -- Window array to multiply by each FFT window before FFT computation. Should have length nperseg // 2 + 1.

  • scale (float) -- Scale factor to multiply the FFT'd data by, related to desired units for output tensor (e.g. letting this equal 1 / (sample_rate * (window**2).sum()) will give output units of density, :math``text{Hz}^-1``.

  • average (str) -- How to aggregate the contributions of each FFT window to the spectral density. Allowed options are 'mean' and 'median'.

  • y (Union[Float[Tensor, 'time'], Float[Tensor, 'channel time'], Float[Tensor, 'batch channel time'], None]) -- Timeseries tensor to compute cross spectral density with x. If left as None, x's power spectral density will be returned. Otherwise, if x is 1D, y must also be 1D. If x is 2D, the assumption is that this represents a single multi-channel timeseries, and y must be either 2D or 1D. In the former case, the cross-spectral densities of each channel will be computed individually, so y must have the same shape as x. Otherwise, this will compute the CSD of each of x's channels with y. If x is 3D, this will be assumed to be a batch of multi-channel timeseries. In this case, y can either be 3D, in which case each channel of each batch element will have its CSD calculated or 2D, which has two different options. If y's 0th dimension matches x's 0th dimension, it will be assumed that y represents a batch of 1D timeseries, and for each batch element this timeseries will have its CSD with each channel of the corresponding batch element of x calculated. Otherwise, it sill be assumed that y represents a single multi-channel timeseries, in which case each channel of y will have its CSD calculated with the corresponding channel in x across _all_ of x's batch elements.

Return type:

Union[Float[Tensor, 'frequency'], Float[Tensor, 'channel frequency'], Float[Tensor, 'batch channel frequency']]

Returns:

Tensor of power spectral densities of x or its cross spectral density with the timeseries in y.

ml4gw.spectral.median(x, axis)

Implements a median calculation that matches numpy's behavior for an even number of elements and includes the same bias correction used by scipy's implementation.

Return type:

']

Parameters:
  • x (Float[Tensor, '... size'])

  • axis (int)

ml4gw.spectral.normalize_by_psd(X, psd, sample_rate, pad)
Parameters:
  • X (Float[Tensor, 'batch num_ifos time'])

  • psd (Float[Tensor, 'num_ifos frequency'])

  • sample_rate (float)

  • pad (int)

ml4gw.spectral.spectral_density(x, nperseg, nstride, window, scale, average='median')

Compute the power spectral density of a multichannel timeseries or a batch of multichannel timeseries. This implementation is exact for all frequency bins, but slower than the fast implementation.

Parameters:
  • x (Union[Float[Tensor, 'time'], Float[Tensor, 'channel time'], Float[Tensor, 'batch channel time']]) -- The timeseries tensor whose power spectral density to compute, or for cross spectral density the timeseries whose fft will be conjugated. Can have shape (batch_size, num_channels, length * sample_rate), (num_channels, length * sample_rate), or (length * sample_rate).

  • nperseg (int) -- Number of samples included in each FFT window

  • nstride (int) -- Stride between FFT windows

  • window (Float[Tensor, '{nperseg//2+1}']) -- Window array to multiply by each FFT window before FFT computation. Should have length nperseg // 2 + 1.

  • scale (float) -- Scale factor to multiply the FFT'd data by, related to desired units for output tensor (e.g. letting this equal 1 / (sample_rate * (window**2).sum()) will give output units of density, \(\text{Hz}^-1\).

  • average (str) -- How to aggregate the contributions of each FFT window to the spectral density. Allowed options are 'mean' and 'median'.

Return type:

Union[Float[Tensor, 'frequency'], Float[Tensor, 'channel frequency'], Float[Tensor, 'batch channel frequency']]

ml4gw.spectral.truncate_inverse_power_spectrum(psd, fduration, sample_rate, highpass=None, lowpass=None)

Truncate the length of the time domain response of a whitening filter built using the specified psd so that it has maximum length fduration seconds. This is meant to mitigate the impact of sharp features in the background PSD causing time domain responses longer than the segments to which the whitening filter will be applied.

Implementation details adapted from here.

Parameters:
  • psd (Float[Tensor, 'num_ifos frequency']) -- The one-sided power spectral density used to construct a whitening filter.

  • fduration (Union[Float[Tensor, 'time'], float]) -- Desired length in seconds of the time domain response of a whitening filter built using this PSD, or a window of this length to taper the edges of the time domain response of the filter. If passed as a float, a Hann window of this length will be used.

  • sample_rate (float) -- Rate at which the time domain data to which the whitening filter will be applied has been sampled.

  • highpass (Optional[float]) -- If specified, will zero out the frequency response of all frequencies below this value in Hz. If left as None, no highpass filtering will be applied.

  • lowpass (Optional[float]) -- If specified, will zero out the frequency response of all frequencies above this value in Hz. If left as None, no lowpass filtering will be applied.

Return type:

Float[Tensor, 'num_ifos frequency']

Returns:

The PSD with its time domain response truncated

to fduration and any filtered frequencies tapered.

ml4gw.spectral.whiten(X, psd, fduration, sample_rate, highpass=None, lowpass=None)

Whiten a batch of timeseries using the specified background one-sided power spectral densities (PSDs), modified to have the desired time domain response length fduration and possibly to highpass/lowpass filter.

Parameters:
  • X (Float[Tensor, 'batch num_ifos time']) -- batch of multichannel timeseries to whiten

  • psd (Float[Tensor, 'num_ifos frequency']) -- PSDs use to whiten the data. The frequency response of the whitening filter will be roughly the inverse of the square root of this PSD, ensuring that data from the same distribution will have approximately uniform power after whitening. If 2D, each batch element in X will be whitened using the same PSDs. If 3D, each batch element will be whitened by the PSDs contained along the 0th dimenion of psd, and so the first two dimensions of X and psd should match.

  • fduration (Union[Float[Tensor, 'time'], float]) -- Desired length in seconds of the time domain response of a whitening filter built using this PSD, or a window of this length to taper the edges of the time domain response of the filter. If passed as a float, a Hann window of this length will be used. Moreover, half of this length will be removed from each edge of the whitened timeseries to account for filter settle-in time.

  • sample_rate (float) -- Rate at which the data in X has been sampled

  • highpass (Optional[float]) -- The frequency in Hz at which to highpass filter the data, setting the frequency response in the whitening filter to 0. If left as None, no highpass filtering will be applied.

  • lowpass (Optional[float]) -- The frequency in Hz at which to lowpass filter the data, setting the frequency response in the whitening filter to 0. If left as None, no lowpass filtering will be applied.

Return type:

Float[Tensor, 'batch num_ifos time']

Returns:

Batch of whitened multichannel timeseries with

fduration / 2 seconds trimmed from each side.