# Re: Power Spectrum

From: Dennis Shea (shea AT cgd.ucar.edu)
Date: Fri May 06 2005 - 19:04:19 MDT

> Could anyone tell me how the power spectrum is calculated in the
>function specx_anal€¯?
>

```---
NCL usees FFTPACK to compute fourier info. The documentation
for the forward fft is at:
ftp://ftp.ucar.edu/dsl/lib/fftpack/ezfftf.f
---
[1] detrend the series [optional]
[2] taper the series   [optional]
[3] calculate the variance of the detrended/tapered series
[4] forward fft on series
[5] square the coef [periodogram  ~2 dof]
[6] smooth the periodogram estimates
[7] normalize [6] so that the area under the curve is
equal to the variance calculated in [3]

area_under_curve =  SUM { S(f)*df*frac }

where frac=1.0 except at the beginning
and end values where frac=0.5
====================
>       If my understanding based on the manual is right, the power spectrum
>       My concern is if  spcx=|F(w)|  or  spcx= |F(w)|^2.  or
>spcx=(|F(w)|^2) / (N^2) ? Here F(w) is the Fourier Transform of x(t), | |
means absolute value, N is the length of the F(w).
>
Neither ... The returned apc estimates have been normailized
to get the area under the curve.

FYI: The documentation for the forward fft is at:
ftp://ftp.ucar.edu/dsl/lib/fftpack/ezfftf.f
---
You could do the above in NCL in a step-by-step process:
========================================================
In the last few days, I have been twice asked how
NCL's specx_anal "works". It is pretty basic and
for the most part can be well simulated by use of several
of NCL's built-in functions.
[a] input a series x[*] with no missing values.

[b] remove the mean an, optionally, detrend the series.

http://ngwww.ucar.edu/ngdoc/ng/ref/ncl/functions/dtrend.html

[c] taper the result of [b]

http://ngwww.ucar.edu/ngdoc/ng/ref/ncl/functions/taper.html

[d] Calculate the variance of [c]

[d] perform an FFT [ezfftf]  on [c]

http://ngwww.ucar.edu/ngdoc/ng/ref/ncl/functions/fft.html

[e] square the real (rc) and imaginary (ic) coef from [d]. This gives
the periodogram.

periodogram = rc^2 + ic^2    ; array syntax

The periodgram has a bandwidth of 1/N (=df [delta_frequency])
but no statistical significance.

[f] perform a running average over the periodogram values using
either runave or wrunave. This creates a spectrum.

http://ngwww.ucar.edu/ngdoc/ng/ref/ncl/functions/runave.html
http://ngwww.ucar.edu/ngdoc/ng/ref/ncl/functions/wrunave.html

A runave of say 5 pts creates a bandwidth of 5*df [approx]
and approx 10 degrees of freedom.

[g] normalize the area under the curve [f] to the input variance.
Remember the first and last bandwidths are only 0.5*df.

The numbers will not match 'spec_anal' exactly but should be very close.

++++++++++++++++++++++++++++++++++++++++++++++++++
good luck
D
_______________________________________________
ncl-talk mailing list
ncl-talk@ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
```

This archive was generated by hypermail 2b29 : Sun May 08 2005 - 22:05:31 MDT