wwz_psd (pyleoclim.utils.wwz_psd)

pyleoclim.utils.wwz_psd(ys, ts, freq=None, freq_method='log', freq_kwargs=None, tau=None, c=0.001, nproc=8, detrend=False, params=['default', 4, 0, 1], gaussianize=False, standardize=False, Neff=3, anti_alias=False, avgs=2, method='default')[source]

Return the psd of a timeseries using wwz method.

Parameters
  • ys (array) – a time series, NaNs will be deleted automatically

  • ts (array) – the time points, if ys contains any NaNs, some of the time points will be deleted accordingly

  • freq (array) – vector of frequency

  • freq_method (str) –

    Method to generate the frequency vector if not set directly. The following options are avialable:
    • log (default)

    • lomb-scargle

    • welch

    • scale

    • nfft

    See utils.wavelet.make_freq_vector for details

  • freq_kwargs (dict) – Arguments for the method chosen in freq_method. See specific functions in utils.wavelet for details

  • tau (array) – the evenly-spaced time points, namely the time shift for wavelet analysis

  • c (float) – the decay constant, the default value 1e-3 is good for most of the cases

  • nproc (int) – the number of processes for multiprocessing

  • detrend (str) –

    None - the original time series is assumed to have no trend; ‘linear’ - a linear least-squares fit to ys is subtracted; ‘constant’ - the mean of ys is subtracted ‘savitzy-golay’ - ys is filtered using the Savitzky-Golay

    filters and the resulting filtered series is subtracted from y.

  • params (list) – The paramters for the Savitzky-Golay filters. The first parameter corresponds to the window size (default it set to half of the data) while the second parameter correspond to the order of the filter (default is 4). The third parameter is the order of the derivative (the default is zero, which means only smoothing.)

  • gaussianize (bool) – If True, gaussianizes the timeseries

  • standardize (bool) – If True, standardizes the timeseries

  • method (string) – ‘Foster’ - the original WWZ method; ‘Kirchner’ - the method Kirchner adapted from Foster; ‘Kirchner_f2py’ - the method Kirchner adapted from Foster with f2py ‘default’ - the Numba version of the Kirchner algorithm will be called. Defaults to default

  • Neff (int) – effective number of points

  • anti_alias (bool) – If True, uses anti-aliasing

  • avgs (int) – flag for whether spectrum is derived from instantaneous point measurements (avgs<>1) OR from measurements averaged over each sampling interval (avgs==1)

Returns

  • psd (array) – power spectral density

  • freq (array) – vector of frequency

  • psd_ar1_q95 (array) – the 95% quantile of the psds of AR1 processes

  • psd_ar1 (array) – the psds of AR1 processes

See also

periodogram()

Estimate power spectral density using a periodogram

mtm()

Retuns spectral density using a multi-taper method

lomb_scargle()

Return the computed periodogram using lomb-scargle algorithm

welch()

Estimate power spectral density using the Welch method

References

Foster, G. (1996). Wavelets for period analysis of unevenly sampled time series. The Astronomical Journal, 112(4), 1709-1729. Kirchner, J. W. (2005). Aliasin in 1/f(alpha) noise spectra: origins, consequences, and remedies. Physical Review E covering statistical, nonlinear, biological, and soft matter physics, 71, 66110.

Examples

>>> from pyleoclim import utils
>>> import matplotlib.pyplot as plt
>>> import numpy as np
>>> # Create a signal
>>> time = np.arange(2001)
>>> f = 1/50
>>> signal = np.cos(2*np.pi*f*time)
>>> # Spectral Analysis
>>> res = utils.wwz_psd(signal, time)
>>> # plot
>>> fig = plt.loglog(
...           res['freq'],
...           res['psd'])
>>> plt.xlabel('Frequency')
>>> plt.ylabel('PSD')
>>> plt.show()

(Source code)