Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Using scipy.signal.spectral.lombscargle for period discovery

The new Scipy v0.11 offers a package for spectral analysis. Unfortunately the documentation is sparse and there aren't many available examples.

As a baby example, I'm trying to do period discovery of a sine wave. Unfortunately it predicts a period of 1 instead of the expected 2pi. Any ideas?

# imports the numerical array and scientific computing packages
import numpy as np
import scipy as sp
from scipy.signal import spectral

# generates 100 evenly spaced points between 1 and 1000
time = np.linspace(1, 1000, 100)

# computes the sine value of each of those points
mags = np.sin(time)

# scales the sine values so that the mean is 0 and the variance is 1 (the documentation specifies that this must be done)
scaled_mags = (mags-mags.mean())/mags.std()

# generates 1000 frequencies between 0.01 and 1
freqs = np.linspace(0.01, 1, 1000)

# computes the Lomb Scargle Periodogram of the time and scaled magnitudes using each frequency as a guess
periodogram = spectral.lombscargle(time, scaled_mags, freqs)

# returns the inverse of the frequence (i.e. the period) of the largest periodogram value
1/freqs[np.argmax(periodogram)]

This returns 1 instead of the expected period of 2pi ~= 1/0.6366. Any ideas?

like image 850
rhombidodecahedron Avatar asked Nov 12 '12 18:11

rhombidodecahedron


1 Answers

Please notice that the last argument of spectral.lombscargle is the angular frequency according to the docstring:

Parameters
----------
x : array_like
Sample times.
y : array_like
Measurement values.
freqs : array_like
Angular frequencies for output periodogram.
like image 56
btel Avatar answered Sep 23 '22 15:09

btel