From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>

Date: Fri Dec 02 2011 - 08:33:34 MST

Received on Fri Dec 2 08:33:42 2011

Date: Fri Dec 02 2011 - 08:33:34 MST

ezfftf is a direct call to the FFTPACK routine at:

http://www.netlib.org/fftpack/ezfftf.f

or, with better documentation

ftp://ftp.ucar.edu/dsl/lib/fftpack/ezfftf.f

===

As with any FFT, the issue is with how the coefficients are normalized.

For ezfftf:

var_zfft2 = 0.5*sum(z(0:nfreq-2)) + z(nfreq-1)

This is a *biased* estimate of the variance. [ie: 1/N and not 1/(N-1) ]

===

I was a bit surprised. I thought the variance would be

var_zfft2 = 0.5*sum(z)

but the core documentation indicates differently.

Attached is a test script.

On 12/2/11 2:01 AM, Wolfgang Langhans wrote:

*> Hi,
*

*>
*

*> I am a bit confused about the outcome of the ezfftf Fourier Forward
*

*> analysis. Should not the sum over all squared norms equal the total
*

*> biased variance of the data series? How can I understand the different
*

*> results obtained in the little example below? How comes that the sum
*

*> amounts almost twice the actual variance?
*

*>
*

*> Thanks for your help in understanding this!
*

*> Wolfgang
*

*>
*

*> ;======================================================
*

*> low = -5.0
*

*> high = 5.0
*

*> con = (high - low) / 32766.0 ; 32766.0 forces a 0.0 to 1.0 range
*

*> Z = new(100,float)
*

*> do i=0,dimsizes(Z)-1
*

*> Z(i) = low + con * rand()
*

*> end do
*

*>
*

*> zmean = avg(Z)
*

*> zvarbiased = 1./dimsizes(Z) * sum( (Z-zmean)^2 ) ;total biased
*

*> variance of original data
*

*>
*

*> zfft = ezfftf(Z)
*

*> z = zfft(0,:)^2 + zfft(1,:)^2
*

*> nfreq = dimsizes(z)
*

*> totvar = sum(z(0:nfreq-1)) ; sum over all squared norms of the
*

*> complex Fourier transforms
*

*>
*

*> print("Variance: "+ zvarbiased + " Sum over spectrum: "+totvar )
*

*> ;======================================================
*

*> _______________________________________________
*

*> ncl-talk mailing list
*

*> List instructions, subscriber options, unsubscribe:
*

*> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
*

_______________________________________________

ncl-talk mailing list

List instructions, subscriber options, unsubscribe:

http://mailman.ucar.edu/mailman/listinfo/ncl-talk

- text/plain attachment: tst_ezfftf.ncl

*
This archive was generated by hypermail 2.1.8
: Fri Dec 02 2011 - 16:10:18 MST
*