program fields
c
c  Test particle fields from the ITGE program.
c
      include 'itg.par'
c
      character*80 str
      real potntial(lz,60,nez)
c
      include 'post.cmn'
      namelist/wdat/ld,md,nd,dt0,nstp,nprnt,nread,
     * lin,x0,y0,z0,pmag,shr,etai,dens,tiovte,rmu1,rmu2,
     * rkap,tkap,imp1,imp2,imkap,etaipar,vy,vz,
     * igradon,iphi00,iflr,mlow,mhi,ihdf,
     * ninterv,movieon,mfr,nfr,nperpmom,qkap,nparmom, 
     * iodd,beta,qkappar,ntrace,ifilter,iexp,etaik,ikx
c
      open(unit=5,file='y11res',status='old',form='formatted')
      open(unit=20,file='itg.fields',status='old')
      open(unit=16,file='itg.in',status='old')

c input files:
c  5 y11res	  a copy of itg.res, the results file from an ITG run.
c 20 itg.fields   phi(x,y)
c  6 itg.in       the original input file for the run
c
      read(16,wdat)
      read(5,1)ld,ldb,nmin,nmax,nd,lin,ne,nfreq,ntotal,md
      read(5,1)nparmom,nperpmom,iodd,iflr,iphi00,ifilter,iexp,igradon
      read(5,1)(mmin(n),n=1,nd)
      read(5,1)(mmax(n),n=1,nd)
      read(5,2)tim,shr,etai,rmu1,rmu2,rkap,tiovte,vy,vz
      read(5,2) dt0,dt,x0,y0,z0
c
 1     format(1x,3i10)
 2     format(1x,3e22.14)
c
      pi=3.141592653589793
      call opngks

      ipage=1
      call set(0.0,1.0,0.0,1.0,0.0,1.0,0.0,7.5,1)
      call plchhq(.5,7.2,'Gyrofluid ITG Simulation',25.,0.,0.)

      if (lin.eq.0) then
      call plchlq(.1,6.0,'E x B nonlinearities included. ',12.,0.,-1.)
      else
      call plchlq(.1,6.0,'This is a linear run. ',12.,0.,-1.)
      endif
      if (igradon.eq.1) then
      call plchlq(.1,5.8,'Background pressure gradient maintained. '
     * ,12.,0.,-1.)
      else
      call plchlq(.1,5.8,'Background pressure allowed to relax. '
     * ,12.,0.,-1.)
      endif
      write(str,111) 'x0 = ',x0
 111  format(a5,f6.2)
      call plchlq(.1,5.6,str,12.,0.,-1.)
      write(str,111) 'y0 = ',y0
      call plchlq(.1,5.4,str,12.,0.,-1.)
      write(str,111) 'z0 = ',z0      
      call plchlq(.1,5.2,str,12.,0.,-1.)
      write(str,102) ld
 102  format('Grid points in the radial direction = ',i3)
      call plchlq(.1,5.0,str,12.,0.,-1.)
      call plchlq(.1,4.8,'Poloidal modes:  ',12.,0.,-1.)
      write(str,103) 'Mmin = ',mmin(nd)
 103  format(a7,i2)
      call plchlq(.2,4.6,str,12.,0.,-1.)
      write(str,103) 'Mmax = ',mmax(nd)
      call plchlq(.2,4.4,str,12.,0.,-1.)
      call plchlq(.1,4.2,'Toroidal modes:  ',12.,0.,-1.)
      write(str,103) 'Nmin = ',nmin
      call plchlq(.2,4.0,str,12.,0.,-1.)
      write(str,103) 'Nmax = ',nmax
      call plchlq(.2,3.8,str,12.,0.,-1.)
      write(str,104) dt0
 104  format('Time step = ',f5.2)
      call plchlq(.1,3.4,str,12.,0.,-1.)
      write(str,105) tim
 105  format('Time at end of run = ',f7.2)
      call plchlq(.1,3.2,str,12.,0.,-1.)
      write(str,106) shr
 106  format('Ln / Ls = ',f6.3)
      call plchlq(.1,2.8,str,12.,0.,-1.)
      write(str,107) etai
 107  format('ETAi = ',f4.1)
      call plchlq(.1,2.6,str,12.,0.,-1.)
      write(str,108) tiovte
 108  format('Ti/Te = ',f4.1)
      call plchlq(.1,2.4,str,12.,0.,-1.)

      if (rmu1.eq.rmu2) then
      write(str,109) rmu1
 109  format('Perpendicular viscosities = ',f5.2)
      call plchlq(.1,2.2,str,12.,0.,-1.)
      else
      write(str,110) rmu1
      call plchlq (.1,2.2,str,12.,0.,-1.)
      write(str,116) rmu2
      call plchlq (.1,2.0,str,12.,0.,-1.)
      write(str,117) rkap
      call plchlq (.1,1.8,str,12.,0.,-1.)
 110  format('Perpendicular viscosity in the density
     * equation = ',f5.2)
 116  format('Perpendicular viscosity in the momentum
     * equation = ',f5.2)
 117  format('Perpendicular viscosity in the temperature
     * equation = ',f5.2)
      endif
      write(str,115) nparmom
      call plchlq (.1,2.0,str,12.,0.,-1.)
 115  format('Number of parallel moments: ',i2)
      write(str,112) nperpmom
      call plchlq (.1,1.8,str,12.,0.,-1.)
 112  format('Number of perpendicular moments: ',i2)
      call finish(ipage)
c
c second page of inputs
c
      call set(0.0,1.0,0.0,1.0,0.0,1.0,0.0,7.5,1)
      call plchhq(.5,7.2,'Input Continued',25.,0.,0.)

      if (vy.eq.0. .and. vz.eq.0.) then
         call plchlq(.1,6.0,'Sheared velocity flow not included',
     *    12.,0.,-1.)
      else
      call plchlq(.1,6.0,'Sheared velocity flow included',12.,0.,-1.)
      write(str,1013) 'vy = ',vy
      call plchlq(.1,5.8,str,12.,0.,-1.)
      write(str,1013) 'vz = ',vz
      call plchlq(.1,5.6,str,12.,0.,-1.)
 1013 format(a5,f12.6)
      endif

      write(str,103) 'iodd = ',iodd
      call plchlq(.1,5.4,str,12.,0.,-1.)
      write(str,103) 'iflr = ',iflr
      call plchlq(.1,5.2,str,12.,0.,-1.)
      write(str,203) 'iphi00 = ',iphi00
 203  format(a9,i2)
      call plchlq(.1,5.0,str,12.,0.,-1.)
      write(str,204) 'ifilter = ',ifilter
      call plchlq(.1,4.8,str,12.,0.,-1.)
 204  format(a10,i2)
      if(nspecies.ne.1) then
      call plchlq(.1,4.6,'Non-adiabatic electrons included',
     * 12.,0.,-1.)
      write(str,205) 'Mass ratio = ',rmime
      call plchlq(.1,4.4,str,12.,0.,-1.)
 205  format(a13,f6.1)
      write(str,206) 'ETAe = ',etae
      call plchlq(.1,4.2,str,12.,0.,-1.)
 206  format(a7,f4.1)
      endif
      start=movieon*dt0
      write(str,301) 'Movie start time = ',start
      call plchlq(.1,4.0,str,12.,0.,-1.)
 301  format(a19,f8.2)
      step=dt0*ninterv
      write(str,302) 'Movie time step = ',step
      call plchlq(.1,3.8,str,12.,0.,-1.)
 302  format(a18,f5.2)
      call finish(ipage)
c
c "potntial" is the electrostatic potential in one z-plane, at one time.
c 
      nsteps=(nstp-movieon)/ninterv
      read(20,2)(((potntial(l,m,n),l=1,ld),m=1,60),n=1,nsteps)
c
c Now you have the potential vs. x, y and t!
c
      call clsgks
      call exit(0)
      end