c GrADS-enabled example model code for ATM562 (Fortran 77 version)
c http://www.atmos.albany.edu/facstaff/rfovell/ATM562/grads_example_code_augmented.f
c version of 10/24/2018

c this template illustrates a procedure to write data to GrADS files:
c  one GrADS .dat file is created, and TWO GrADS .ctl files
c  the first .ctl file is temporary, to help you visualize model
c  output during a long run, or if a model crash prevents creation
c  of the final .ctl file (which is done last)

      program task3

c ----------------------------------------------------------------
c  your non-executable code, including parameter, dimension, data
c  NOTE my GrADS code expects certain arrays, including a base state
c   mean wind (ub) and a 3D array for north-south velocity (v) that
c   can be utilized for other variables for models lacking Coriolis
c  ADDITIONALLY the leapfrog scheme needs three time levels;
c   I designate these as thp, th, and thm, for perturbation theta
c  SO your code may include the following:
c ----------------------------------------------------------------
c base state vectors
      dimension tb(nz),qb(nz),pb(nz),pib(nz),rhou(nz),rhow(nz),ub(nz)
c 2D prognostic arrays - 3 time levels
      dimension thp(nx,nz),th(nx,nz),thm(nx,nz)   ! pot. temp. pert.
      dimension wp(nx,nz),w(nx,nz),wm(nx,nz)      ! vert. vel.
      dimension up(nx,nz),u(nx,nz),um(nx,nz)      ! zonal horiz. vel.
      dimension vp(nx,nz),v(nx,nz),vm(nx,nz)      ! merid. horiz. vel.
      dimension qvp(nx,nz),qv(nx,nz),qvm(nx,nz)   ! water vapor
      dimension pip(nx,nz),pi(nx,nz),pim(nx,nz)   ! ndim pert. pres.

	[...]

c ----------------------------------------------------------------
c  end of non-executable code here
c ----------------------------------------------------------------
	
c ----------------------------------------------------------------
c  your model setup code is here
c ----------------------------------------------------------------
c  * declare your time step, grid spacings (dt,dx,dz)
c  * declare your model end time, plotting interval (timend,nplt)
c  * set up your 1D base state (tb,pib,rhou,rhow,ub)
c  * set your model 2D initial condition (thermal, etc.)
c examples... you need to change these

      dt=5.0                  ! seconds
      timend=1000             ! seconds
      nplt=60                 ! plot every 60 time steps, or 300 sec

	[...]

c ----------------------------------------------------------------
c  GrADS initialization is here
c ----------------------------------------------------------------
c  * set your casename (will create 'casename'.ctl, 'casename'.dat)
c  * create a temporary control (*.ctl) file
c  * write initial condition to data (*.dat) file
c EXAMPLES... you need to change these

      casename='test' ! declare 'character*80 casename' above
      iplt = 0        ! counter for plotting calls
      byteswap = 0    ! if byteswapping not needed, set = 0
      igradscnt =999  ! grads counter - dummy value to start

c  * initial call to write_to_grads_ctl sets up temporary control file

      call write_to_grads_ctl(casename,nx,nz,dt,dx,dz,
     1 igradscnt,nplt,byteswap) ! create the temporary control file

c  * write out model data at initial time - i.e., call dumpgrads

      igradscnt =0    ! reset grads counter
       call dumpgrads(igradscnt,u,v,w,th,pi,qv,rhou,tb,pib,ub,qb,
     1   nx,nz,cpd)   

c ----------------------------------------------------------------
c  model integration loop STARTS here
c ----------------------------------------------------------------
c  * update plot counters
c  * integrate your equations, take care of boundary conditions
c  * set for next time step
c example...

       iplt = iplt + 1


c ----------------------------------------------------------------
c  decide if it's time to plot...
c ----------------------------------------------------------------
c  * if it's time to plot, call dumpgrads and reset plot counter
c example...

      if(iplt.eq.nplt)then
       call dumpgrads(igradscnt,u,v,w,th,pi,qv,rhou,tb,pib,ub,qb,
     1   nx,nz,cpd)                
       iplt=0
      endif



c ----------------------------------------------------------------
c  model integration loop ENDS here
c ----------------------------------------------------------------
c  * go back to start of model integration loop unless it's timend


c ----------------------------------------------------------------
c  main routine ENDS here
c ----------------------------------------------------------------
c  * before you stop model, create FINAL grads control file
c  * this control file will have the proper number of times listed

      call write_to_grads_ctl(casename,nx,nz,dt,dx,dz,
     1 igradscnt,nplt,byteswap)

      stop
      end

c ----------------------------------------------------------------
c  after your main routine, include GrADS subroutines:
c    subroutine dumpgrads
c    subroutine write_to_grads_dat
c    subroutine write_to_grads_ctl
c  * you can also keep them in separate files and link to them
c  * EXAMPLE:   gfortran task3.f grads_routines_augmented.f
c ----------------------------------------------------------------