; http://www.wavemetrics.com/products/igorpro/dataanalysis/signalprocessing/powerspectra.htm

     x    = (/ 1002, 1017, 1018, 1020, 1018, 1027, \
               1028, 1030, 1012, 1012, 982, 1012, \
               1001, 996, 995, 1011, 1027, 1025, \
               1030, 1016, 996, 1006, 1002, 982 /)*1.0       ; even
             ;;1030, 1016, 996, 1006, 1002, 982, 999 /)*1.0  ; odd

     N    = dimsizes(x)
     df   = 1.0/N
     xAvg = avg(x)
     xVar = variance(x)*(N-1.)/N    ; biased estimate

     x    = x - xAvg

     print ("N="+N+" df="+df+" avg(x)="+xAvg+"  xVar="+xVar)
     print("---")

;******************************************
; cfftf
;******************************************
     cf   = cfftf (x, 0.0, 0) ; imaginary part set to 0.0
     printVarSummary(cf)
     print("---")

     cf  = cf*(1.0/N)         ; normalization for 2-sided fft
     print("cfftf: "+sprintf("%9.5f",cf@frq)+" "+sprintf("%9.3f", cf(0,:))+" "+sprintf("%9.3f", cf(1,:)) )
     print("---")

     cf2   = cf(0,:)^2 + cf(1,:)^2     ; sum after normalization
     cf2sum= sum(cf2)
     print("---")
     print("cf2sum="+cf2sum+"   ratio="+(cf2sum/xVar))
     print("---")

;******************************************
; ezfftf
;******************************************

     cfez = ezfftf (x)
     printVarSummary(cfez)
     print("---")
     print("ezfftf: "+sprintf("%9.3f", cfez(0,:))+" "+sprintf("%9.3f", cfez(1,:)) )
     print("---")

     cfez2   = cfez(0,:)^2 + cfez(1,:)^2
     cfez2sum= sum(cfez2)
     print("cfez2sum="+cfez2sum)
     print("---")

     if (N%2 .eq. 0) then
         N2 = N/2
         cfez_var = 0.5*sum(cfez2(0:N2-2)) + cfez2(N2-1)
     else
         cfez_var = 0.5*sum(cfez2)
     end if
     print("cfez_var="+cfez_var+"   ratio="+(cfez_var/xVar))
     print("---")