A periodic signal can be analyzed in two ways. The most common way to identify periodic signals is to convert a time series into the Fourier domain. The time series will be some value as a function of time, and the Fourier domain will generally be represented by power as a function of frequency. The definition of the Fourier transform that I will use is

F(\nu) = \int_{-\infty}^{+\infty} f(t) e^{-2\pi i\nu\t},

(I apologize, but LaTeX still seems to be disabled in the tiki-wiki package), where f(t) is the time series, \nu is the frequency, and F(\nu) is the Fourier transform. The Fourier transform (as well as concepts such as sampling and windowing) are described well in Numerical Recipes in C (chapter 12 in the 2nd editionhttp://www.nrbook.com/nr3/).

If a periodic signal is known to exist in time-series data, then one can analyze the phase and amplitude of the periodic signal using direct means, such as folding the data about a known period and measuring the phase and amplitude through maximum likelihood or least-squares modeling.

Identifying a Periodic Signal

There are many techniques for identifying a periodic signal. Which technique to use will depend upon how the data was recorded, and what the data represents. However, all techniques require an understanding of some Practical Concepts for Fourier Timing.

Regularly-sampled, densely populated data with identical uncertainties

The best technique to search regularly-sampled data for periodicities is the Fast Fourier Transform (FFT). This is explained in Numerical Recipes, and implemented very well in the FFTW library. The statistical distribution of noise in the power spectrum from a Fast Fourier transform is a chi-squared distribution with two degrees of freedom; the statistics of detecting a signal within noise is described in Vaughan et al. (1994). The best overall reference for using Fourier Transforms is van der Klis (1989), which is unfortunately not electronically available, but is definitely worth going to a library for (and photocopying).

Very sparsely sampled data in the form of discrete events

In X-ray and gamma-ray astronomy, the event rates from a source may be extremely low, with only a few counts arriving during a week-long satellite observation. Nonetheless, if the flux is highly-modulated, and there are at least 50 or so counts in an observing campaign, one may still be able to identify a period. In this case, the Rayleigh (or Z-squared) statistic is an excellent algorithm. It operates on the arrival times of events, and computes a direct Fourier transform. White noise powers are distributed as a chi-squared distribution, as with the FFT. However, it is extremely inefficient if you have more than 10,000 data points (use an FFT for that).

Irregularly sampled, densely populated data with equal uncertainties

If your data is irregularly sampled but dense, such as data from the RXTE All-Sky Monitor, use the Lomb-Scargle periodogram to search for periods. Note that with irregularly-sampled data, the window function will often generate spurious peaks in the power spectrum. Therefore, it is often necessary to estimate the significance of signals using Monte Carlo simulations.

Data of varying accuracy

If the data is not of uniform quality, such that you can assume a uniform variance for all the measurements, then the only way to carry out a periodicity search is to model the data with sinusoids of varying frequencies. This can be done relatively quickly and robustly by folding the data about a fixed period and using a linear least squares minimization.

Non-sinusoidal periods

If the periodic signal is not expected to be sinusoidal, then the Fourier power spectrum will contain multiple harmonics. One could approach the data with the techniques above, by searching an FFT by summing the trial periods and their harmonics, or by modifying Rayleigh statistic to sum multiple harmonics. One could also apply any of a number of techniques in which a time series is folded around trial periods, and the variance of the mean in each bin is used as a measure of whether there is a periodic signal. The idea there is that if there is no signal at the trial period, the mean will be the same in each phase bin, whereas a periodic signal would create a different mean in each bin as you traced out the pulsation's phase. One implementation of this is the phase-dispersion minimization technique of Stellingwerf (1978), although be warned that the statistical distribution he assumes for identifying signals in noise is incorrect.

Studying signals with known periods

If you have identified a period in a dataset, you probably want more information about the signal, such as an accurate estimate of the period and phase to construct and ephemeris, or characterizations of any variations in the period (e.g., caused by the Doppler shift of a pulsar in a binary, or instabilities in the mechanism producing the period).

You can refine a frequency measurement by over-sampling the signal in the Fourier domain. The most common way of doing this is "zero padding" the data, that is, adding extra time bins containing the mean value of the data to the end of the time series, and then taking the Fourier Transform. The resolution of the resulting data is 1/T, where T is the length of the padded data set. However, this does not yield a phase, information about frequency derivatives, or uncertainties on the frequency measurement.

The best way to study a periodic signal is a technique called phase connection (Manchester, R. N., & Taylor, J. H. 1977, Pulsars, San Francisco: W. H. Freeman and Co.). The phase of a periodic signal is just

\phi(t) = \nu t

If the frequency is not constant, however, the phase can be written more generally as

\phi(t) = \phi(0) + \int_0^t \nu(t) dt

The time dependence on the frequency might be introduced by Doppler shifts caused by binary motion, or by a drift in frequency with time. In either case, one can expand the time-dependent frequency around a central frequency,

\nu(t) = \phi(0) + \nu_0 t + \int_0^t \nu^\prime(t) dt

where \nu^\prime is a correction to the frequency. Phase connection works by measuring the phase of the data assuming that \nu_0 is the correct frequency, and the modeling the phase residuals to determine the parameters of a functional form representing \nu^\prime. For instance, if one is looking to measure a frequency drift, one might set the integral equal to

\int_0^t \nu^\prime(t) dt = \delta\phi + \delta\nu t + \frac{1}{2} \dot{\nu} t^2.

Phase connection is implemented by the following steps:

  1. Choose a trial period based on a Fourier power spectrum or other measurement
  2. Split the data into small intervals, and fold each interval about the trial period
  3. Measure the phase of the folded signal in each interval, by fitting it with a sinusoid or by cross-correlating it with a template pulse
  4. Estimate the uncertainty on each phase measurement
  5. If the phase changes with time, this indicates that the trial period is incorrect. By fitting a polynomial (or other appropriate function) to the residual phases, one can obtain a new phase model.

The resulting model can provide both accurate estimates of a signal's phase, frequency, and frequency derivatives (or orbital parameters), and uncertainties on those estimates.

Page last modified on Saturday 30 of January, 2010 16:26:30 EST