load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"   
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"   

begin

  fname = "template_ps.nc"
  fout  = addfile(fname,"c")

  fileAtt = True
  fileAtt@MissingValue  = -999
  fileAtt@FillValue     = -999
  fileAtt@date  = systemfunc("date")
  fileAtt@group = "UCLA"
  fileAtt@author= "Bjorn Stevens"
  fileAtt@case  = "Base Case"
  fileAtt@title = "GCSS 7 Profile Statistics"

  fileattdef(fout,fileAtt)

  ntimes   = 9
  nwpoints = 300
  ntpoints = 300
  dimNames = (/"time","zt","zw"/)
  dimSizes = (/ntimes, ntpoints, nwpoints/)
  dimUnlim = (/False,False,False/)
  filedimdef(fout,dimNames,dimSizes,dimUnlim)

  zt = fspan(2.5,1497.5,ntpoints)
  zw = fspan(0.0,1495.0,nwpoints)
  time = fspan(0., 14400., ntimes)

  time@_FillValue = -999
  time@units = "[s]"
  time@longname = "End of Averaging Interval"
  time!0 = "time"
  time&time = time

  zw@_FillValue = -999
  zw@units = "[m]"
  zw@longname = "Height of w points"
  zw!0 = "zw"
  zw&zw = zw

  zt@_FillValue = -999
  zt@units = "[m]"
  zt@longname = "Height of temperature points"
  zt!0 = "zt"
  zt&zt = zt

  fout->time = time
  fout->zt = zt
  fout->zw = zw

  template_w = new ( (/ntimes, nwpoints/), "float")
  template_w!0 = "time"
  template_w!1 = "zw"
  template_w&time = time
  template_w&zw = zw
  template_t = new ( (/ntimes, ntpoints/), "float")
  template_t!0 = "time"
  template_t!1 = "zt"
  template_t&time = time
  template_t&zt = zt

  u = template_t
  u@longname = "zonal wind"
  u@units = "[ms:S:-1:N:]"
  fout->u = u

  v = template_t
  v@longname = "meridional wind"
  v@units = "[ms:S:-1:N:]"
  fout->v = v

  thetal = template_t
  thetal@longname = "liquid water potential temperature"
  thetal@units = "[K]"
  fout->thetal = thetal

  rt = template_t
  rt@longname = "total water mixing ratio"
  rt@units = "[gkg:S:-1:N:]"
  fout->rt = rt

  rl = template_t
  rl@longname = "liquid water mixing ratio"
  rl@units = "[gkg:S:-1:N:]"
  fout->rl = rl

  dn0 = template_t
  dn0@longname = "density"
  dn0@units = "[kgm:S:-3:N:]"
  fout->dn0 = dn0

  u_var = template_t
  u_var@longname = "resolved variance of zonal wind"
  u_var@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->u_var = u_var

  v_var = template_t
  v_var@longname = "resolved variance of meridional wind"
  v_var@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->v_var = v_var

  w_var = template_w
  w_var@longname = "resolved variance of vertical wind"
  w_var@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->w_var = w_var

  w_skw = template_w
  w_skw@longname = "resolved third moment of vertical wind"
  w_skw@units = "[m:S:-3:N:s:S:-3:N:]"
  fout->w_skw = w_skw

  sfs_tke = template_w
  sfs_tke@longname = "subfilter variance of kinetic energy"
  sfs_tke@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->sfs_tke = sfs_tke

  theta_var = template_t
  theta_var@longname = "resolved variance of theta_l"
  theta_var@units = "[K:S:2:N:]"
  fout->theta_var = theta_var

  rt_var = template_t
  rt_var@longname = "resolved variance of total water"
  rt_var@units = "[g:S:-2:N:kg:S:-3:N:]"
  fout->rt_var = rt_var

  rl_var = template_t
  rl_var@longname = "resolved variance of liquid water"
  rl_var@units = "[g:S:-2:N:kg:S:-2:N:]"
  fout->rl_var = rl_var

  rad_flx = template_w
  rad_flx@longname = "Radiative flux"
  rad_flx@units = "[Ws:S:-2:N:]"
  fout->rad_flx = rad_flx

  tot_tw = template_w
  tot_tw@longname = "Total theta_l flux"
  tot_tw@units = "[Ws:S:-2:N:]"
  fout->tot_tw = tot_tw

  sfs_tw = template_w
  sfs_tw@longname = "Subfilter theta_l flux"
  sfs_tw@units = "[Ws:S:-2:N:]"
  fout->sfs_tw = sfs_tw

  tot_rw = template_w
  tot_rw@longname = "Total r_t flux"
  tot_rw@units = "[Ws:S:-2:N:]"
  fout->tot_rw = tot_rw

  sfs_rw = template_w
  sfs_rw@longname = "Subfilter r_t flux"
  sfs_rw@units = "[Ws:S:-2:N:]"
  fout->sfs_rw = sfs_rw

  tot_uw = template_w
  tot_uw@longname = "Total u-flux"
  tot_uw@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->tot_uw = tot_uw

  sfs_uw = template_w
  sfs_uw@longname = "Subfilter u-flux"
  sfs_uw@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->sfs_uw = sfs_uw

  tot_vw = template_w
  tot_vw@longname = "Total v-flux"
  tot_vw@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->tot_vw = tot_vw

  sfs_vw = template_w
  sfs_vw@longname = "Subfilter v-flux"
  sfs_vw@units = "[m:S:-2:N:s:S:-2:N:]"
  fout->sfs_vw = sfs_vw

  E = template_w
  E@longname = "Resolved turbulent kinetic energy"
  E@units = "[m:S:2:N:s:S:-2:N:]"
  fout->E = E

  shr_prd = template_w
  shr_prd@longname = "Resolved shear production"
  shr_prd@units = "[m:S:2:N:s:S:-3:N:]"
  fout->shr_prd = shr_prd

  boy_prd = template_w
  boy_prd@longname = "Resolved buoyancy production"
  boy_prd@units = "[m:S:2:N:s:S:-3:N:]"
  fout->boy_prd = boy_prd

  transport = template_w
  transport@longname = "Resolved transport (turbulent plus pressure)"
  transport@units = "[m:S:2:N:s:S:-3:N:]"
  fout->transport = transport

  dissipation = template_w
  dissipation@longname = "Dissipation "
  dissipation@units = "[m:S:2:N:s:S:-3:N:]"
  fout->dissipation = dissipation

  storage = template_w
  storage@longname = "TKE change across averaging interval"
  storage@units = "[m:S:2:N:s:S:-3:N:]"
  fout->storage = storage

  sfs_boy = template_w
  sfs_boy@longname = "Subfilter buoyancy production"
  sfs_boy@units = "[m:S:2:N:s:S:-3:N:]"
  fout->sfs_boy = sfs_boy

  e = template_w
  e@longname = "Subfilter turbulent kinetic energy"
  e@units = "[m:S:2:N:s:S:-3:N:]"
  fout->e = e

end