;  Grads control file
; ---------------
; dset snow.grd
; title snow depth  winter
; options sequential big_endian
; undef -9999.00               
; xdef 29 linear 20. 5.
; ydef 8 linear 35.  5.
; zdef 1 linear 1 1
; tdef 60 linear jan1936 1yr
; vars 1
; p 0 0 precipitation 
; endvars
; ---------------

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
  
  ;PLOT  = True 
  ;NC    = True

  diri  = "/import/home/ashendricks/above/data/ndvi/"
  fili  = "polar_82-17_MaxNDVI3g.grd"                      ; DSET
  fName = diri+fili       

  nlat  = 481                               ; YDEF
  mlon  = 4320                              ; XDEF

  year  = 1982                            ; TDEF
  ntim  = 36                              ; time steps
  ;nmos  = 12
                                          ; not required
  time  = new (ntim, float  ,"No_FillValue") ; generate a "time" variable
  ;date  = new (ntim, integer,"No_FillValue") ; generate a "date" variable
  
  n     = -1
  ;do nmo=1,nmos
  ;   YRM= year*10000 + nmo*100
  ;   ndm= days_in_month(year, nmo)  
  do iy=1,ntim  
       n = n+1
       time(n) = year + n     
  ;     date(n) = YRM + ndy                ; YYYYMMDD
  end do
  ;end do
  time!0         = "time"
  time@long_name = "time"
  time@units     = "year"           ; "yyyy.fraction_of_year"
  time&time      = time

  ;date!0         = "time"
  ;date@long_name = "date"
  ;date@units     = "??????"
  ;date&time      = time
                                          ; generate lat/lon
  lon = fspan(-179.9583,(-179.9583+((mlon-1)*0.0833)),mlon)
  lon!0 = "lon"
  lon@long_name = "longitude"
  lon@units     = "degrees_east"
  print (lon)

  lat = fspan(49.9528,(49.9528+((nlat-1)*0.0833)),nlat)
  lat!0 = "lat"
  lat@long_name = "latitude"
  lat@units     = "degrees_north"
                                          ; create an array to contain data
  UNDEF = -999.                           ; UNDEF
  x     = new ( (/ntim,nlat,mlon/), float, UNDEF)
  x!0   = "time"
  x!1   = "lat"
  x!2   = "lon"
  x&time=  time
  x&lat =  lat
  x&lon =  lon

  x@long_name = "ndvi"         ; VARS
  x@units     = "index"                 
                                    
 ;setfileoption("bin","ReadByteOrder","Native")       ; a033 default
 setfileoption("bin","ReadByteOrder","LittleEndian")
 ;setfileoption("bin","ReadByteOrder","BigEndian")

                                          ; read each record: store in x
  do nt=0,ntim-1                          ; the ::-1 changes the latitude order
     x(nt,:,:) = fbinrecread(fName, nt, (/nlat,mlon/), "float")
  end do

  printVarSummary(x)
  print ("min(x)="+min(x))
  print ("max(x)="+max(x))

  if (NC) then
      nline= inttochar(10)

      diro = "/import/home/ashendricks/above/data/ndvi/"
      filo = fili + ".nc"
      system ("/bin/rm -f "+diro+filo)  ; remove any pre-existing file

      ncdf = addfile (diro+filo, "c")
     ;setfileoption(ncdf,"DefineMode",True)               ; a033 [most efficient]
      
      globeAtt         = 1      ; global [file] attributes
      globeAtt@title   = "MaxNDVI82-17"
      globeAtt@source  = fili 

      globeAtt@story   = nline + \
                         "An NCL user sent a Grads file and .ctl."+nline
     ;globeAtt@NCL     = nline + \
     ;                   "/fs/cgd/home0/shea/ncld/ncld2/ucla_bordoni/rdGrads.ncl_GPCP1"+nline
      globeAtt@creation_date= systemfunc ("date" ) 

      fileattdef( ncdf, globeAtt )
                                    
      dimNames = (/"time", "lat", "lon" /)  
      dimSizes = (/  -1  ,  nlat,  mlon /) 
      dimUnlim = (/ True , False, False /)   
      filedimdef(ncdf, dimNames  , dimSizes,  dimUnlim )

      filevardef   (ncdf, "time"  , typeof(time) , "time" )
      filevarattdef(ncdf, "time", time)
     
      filevardef   (ncdf, "lat", "float", "lat")
      filevarattdef(ncdf, "lat", lat)
     
      filevardef   (ncdf, "lon", "float", "lon")
      filevarattdef(ncdf, "lon", lon)
     
      ;filevardef   (ncdf, "date", typeof(date), "time")
      ;filevarattdef(ncdf, "date", date)
      
      filevardef(ncdf, "GPCP"  , typeof(x) , (/"time", "lat" , "lon"/) )
      filevarattdef(ncdf, "GPCP", x)

      ncdf->time   = (/ time /)
      ;ncdf->date   = (/ date /)
      ncdf->lat    = (/ lat /)
      ncdf->lon    = (/ lon /)
      ncdf->ndvi   = (/ x /) 
  end if

end