power

power(dx, y1, /, axis=0, **kwargs)[source]

Return the spectral decomposition of a real-valued array along an arbitrary axis.

Parameters

  • dx (float or array-like) – Dimension step size in physical units. Used to scale fx.

  • y (array-like) – The input data.

  • axis (int, optional) – Location of the “time” axis.

  • wintype (str or (str, float), optional) – The window specification, passed to window. The resulting weights are used to window the data before carrying out spectral decompositions. See notes for details.

  • nperseg (int, optional) – The time dimension window or segment length, passed to window. If None, windowing is not carried out. See notes for details.

  • detrend ({‘constant’, ‘linear’}, optional) – Passed as the type argument to scipy.signal.detrend. 'constant' removes the mean and 'linear' removes the linear trend.

Returns

  • fx (array-like) – Frequencies. Units are the inverse of the dx units.

  • P (array-like) – Power spectrum. Units are the input units squared.

Notes

The Fourier coefficients are scaled so that total variance is equal to one half the sum of the right-hand coefficients. This is more natural for the real-valued datasets typically used by physical scientists, and matches the convention from Professor Elizabeth Barnes’s objective analysis course notes. This differs from the numpy convention, which scales the coefficients so that total variance is equal to the sum of squares of all coefficients, or twice the right-hand coefficients.

Windowing is carried out by applying the wintype weights to successive time segments of length nperseg (overlapping by one half the window length), taking spectral decompositions of each weighted segment, then taking the average of the result for all segments. Note that non-boxcar windowing reduces the total power amplitude and results in loss of information. It may often be preferable to follow the example of [RH91] and smooth in frequency space with a Gaussian filter after the decomposition has been carried out.

The below example shows that the extent of power reduction resulting from non-boxcar windowing depends on the character of the signal.

>>> import numpy as np
>>> import climopy as climo
>>> state = np.random.RandomState(51423)
>>> w = climo.window(200, 'hanning')
>>> y1 = np.sin(np.arange(0, 8 * np.pi - 0.01, np.pi / 25)) # basic signal
>>> y2 = state.rand(200) # complex signal
>>> for y in (y1, y2):
...     yvar = y.var()
...     Y = (np.abs(np.fft.fft(y)[1:] / y.size) ** 2).sum()
...     Yw = (np.abs(np.fft.fft(y * w)[1:] / y.size) ** 2).sum()
...     print('Boxcar', Y / yvar)
...     print('Hanning', Yw / yvar)
Boxcar 1.0
Hanning 0.375
Boxcar 0.9999999999999999
Hanning 0.5728391743988162

References

RH91

William J. Randel and Isaac M. Held. Phase Speed Spectra of Transient Eddy Fluxes and Critical Layer Absorption. Journal of the Atmospheric Sciences, 48(5):688–697, March 1991. doi:10.1175/1520-0469(1991)048<0688:PSSOTE>2.0.CO;2.

Examples

>>> import numpy as np
>>> import xarray as xr
>>> import climopy as climo
>>> state = np.random.RandomState(51423)
>>> ureg = climo.ureg
>>> x = xr.DataArray(
...     np.arange(1000, dtype=float) * ureg.day, dims=('time',), name='time'
... )
>>> y = xr.DataArray(
...     state.rand(1000, 50) * ureg.K, dims=('time', 'space'), name='variable'
... )
>>> f, P = climo.power(x, y, axis=0)
>>> f.head()
<xarray.DataArray 'time' (f: 5)>
<Quantity([0.001 0.002 0.003 0.004 0.005], '1 / day')>
Dimensions without coordinates: f
Attributes:
    long_name:  frequency
>>> P.head()
<xarray.DataArray 'variable' (f: 5, space: 5)>
<Quantity([[1.18298459e-04 1.09036599e-04 3.76719588e-04 3.36341980e-04
  1.62219721e-04]
 [1.74596838e-04 2.09816475e-04 1.26945159e-04 3.71464785e-04
  4.18144064e-04]
 [6.38233218e-05 1.52714331e-05 8.25923743e-05 3.31492680e-04
  3.06457271e-05]
 [1.71082223e-04 3.95808831e-05 5.23555579e-04 6.15207584e-05
  5.59359839e-05]
 [1.34026670e-04 3.90031444e-05 2.11700762e-04 2.69345283e-05
  6.86978186e-05]], 'kelvin ** 2')>
Coordinates:
  * f        (f) float64 0.001 0.002 0.003 0.004 0.005
Dimensions without coordinates: space