; ==============================================================
; eof_2.ncl
; =============================================================
; NCL V6.4.0 has new functions eofunc_n_Wrap and 
; eofunc_ts_n_Wrap that allow you to calculate the EOFs without 
; first having to first reorder the data. See eof_2_640.ncl.
; =============================================================
; These files are loaded by default in NCL V6.2.0 and newer
 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
 load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"

begin
  neof   = 10
;*******************************************
;  Read data
;*******************************************
  diri = "./"
  fili = "gpi.invert.box.ano.nc"

  f    = addfile (diri+fili, "r")
  x    = f->gpi            ; (time,lat,lon) ==> (29,60,180)
  printVarSummary(x) 

;*******************************************
;  EOF 
;*******************************************
  X = x(lat|:,lon|:,time|:)              ; Space x Time

  optEof = True
  eof    = eofunc_Wrap( X, neof, optEof)
  eof_ts = eofunc_ts_Wrap( X, eof, False)

;  eof    = eofunc( X, neof, optEof)
;  eof_ts = eofunc_ts( X, eof, False)
  eof!0 = "time"
  printVarSummary(eof)
  printVarSummary(eof_ts)

;    fbinrecwrite ("eof.dat", -1, eof)

; attempt to write to netcdf
  system("/bin/rm -f eof.nc")   ; remove any pre-existing file
  ncdf     = addfile("eof.nc" ,"c")  ; open output netCDF file
; make time and UNLIMITED dimension
; recommended  for most applications
       filedimdef(ncdf,"time",-1,True)

  ncdf->eof = eof
;end write netcdf

 
;   do n=0,neval-1
;      evts(n,:) = evts(n,:) + evts@ts_mean(n)  ; add time series mean
;   end do                 
                     
;   xEof  = eof2data(ev,evts)       ; reconstruct original array
                                         ; 
;   x2    = eof2data(eof(3,:,:),eof_ts(3,:,:)) ; construct array using
                                           ; first 2 EOFs
;  printVarSummary(x2)
;    fbinrecwrite ("eof_res.dat", -1, x2)
;    fbinrecwrite ("eof.dat", -1, eof) 
;    asciiwrite ("eof.txt", eof)
;   fbinrecwrite ("eof_ts.dat", -xx, eof_ts)
    asciiwrite ("pc1.txt", eof_ts(0,:))
   asciiwrite ("pc2.txt", eof_ts(1,:))
    asciiwrite ("pc3.txt", eof_ts(2,:))
    asciiwrite ("pc4.txt", eof_ts(3,:))
    asciiwrite ("pc5.txt", eof_ts(4,:))

end