Re: Re: function in NCL for 10-day high pass filter

From: Dennis Shea (shea AT
Date: Thu Oct 27 2005 - 15:05:01 MDT

> As to my knowledge, there are *at least" two ways in NCL you can use for
band-filtering analysis (also appropriate for low-frequencey filter and
high-frequncey filter).
> 1) Fourier harmonics anal: using *ezfftf* and *ezfftb* together, which is
most tranditional and simpliest way to extract the signals you need from a time
> 2) Lancos filtering:, which
also works well.

Yes, this is correct. You can do filtering via FFTs.

It is "easy" but there are some caveats.

*Conceptually*, you do the the following
  [1] Perform a forward ftt on a series, "x"
  [2] Set the fourier coefficients of unwanted
      frequencies to 0.0
  [3] Reconstruct, the series via ezfftb.
      x = random_uniform(-1, 1, 24) ; generate a series
      fc = ezfftf( x ) ; (2,12)
      fc(:,0:4) = 0.0 ; set low freq to 0.0
      fc(:,8:11) = 0.0 ; set high freq to 0.0
      xBand = ezfftb(fc) ; reconstruct band passed series
*Reality* is a bit different. :-(

  [1] The time series are not infinitely long or cyclic
  [2] Anytime you truncate the coefficients you get Gibbs Phenomena
  [3] Not all the numbers you get back are interpretable.
  (1) The FFT assumes periodic data [based on sines and cosines].
      Since time series are not cyclic, the series should be *tapered*
      before doing the FFT to minimize "leakage"
  (2) Gibbs phenomena will affect results. This may/may-not be important.
      It depends on the series.
  (3) When weights are applied to the time series [eg, wgt_runave],
      then if it is, say, an 11-point running average, the first
      and last 5 points are missing. When you use the FFT approach,
      numbers are returned at *every* location. Gee, isn't that
      great!!! The problem is that a certain number of points
      at the beginning and end of the reconstructed series
      are bogus.
      Let the "x" be anomailies
      x = random_uniform(-1, 1, 24) ; generate a series
      xTaper = taper(x, 0.1, 0)
      fc = ezfftf( xTaper ) ; (2,12)
      fc(:,0:4) = 0.0 ; set low freq to 0.0
      fc(:,8:11) = 0.0 ; set high freq to 0.0
      xBand = ezfftb(fc) ; reconstruct band passed series
      Ignore the first and last "n" points. You have to
      know what "n" is. One way, create a filtered series
      that was generated via the wgt_runave. Then do
      the same via an fft and see where they 'match up'.
      I recommend setting the bogus points to _FillValue
      or only explicitly the "good" points.
Good luck

ncl-talk mailing list

This archive was generated by hypermail 2b29 : Thu Oct 27 2005 - 18:24:17 MDT