
      program tigcmloc
c
c ti(e)gcm location processor:
c
      include 'tigcmloc.h'
      dimension iasf(13)
      data iasf/13*1/
      data luhist/8/, lumag/99/, iframe/0/, zmflag/999./, r12flag/999./
      data gmflag/999./, luascii/90/
      pointer (pfields,fields(1))
      pointer (putkmx,utkmx(1)), (phtutkmx,htutkmx(1)),
     +        (putjmx,utjmx(1))
      logical exists
      data labshort
     + /'TN      ','UN      ','VN      ','O2      ','O1      ',  
     +  'N4S     ','NO      ','O+      ','N2D     ','TI      ',  
     +  'TE      ','NE      ','O2+     ','W       ','Z       ',  
     +  'POT     ','N2      ','UI      ','VI      ','WI      ',
     +  'FOF2    ','HMF2    ','RHO     ','O/(O2+N2','O/N2    '/
c
      data logterp
c       tn  un  vn  o2  o1  n4s  no  o+  n2d  ti  te  ne  o2+  w  z
     + / 0,  0,  0,  1,  1,   1,  1,  1,   1,  0,  0,  1,  1,  0, 0,
c       pot  n2  ui  vi  wi fof2 hmf2 rho o/(o2+n2)  o/n2
     +    0,  1,  0,  0,  0,   0,   0,  1,       0,     0/
c
c ixgcm = initialize field indices
c
      do i=1,nftot
        ixgcm(i) = i
      enddo
c
c Initialization: set up tgcm grid, get user input, set
c                 some graphics default parameters:
c
      call setgrid(gcmlat,glat1,dlat,jmx, gcmlon,glon1,dlon,imx,
     +             gcmzp,zp1,dzp,kmx)
      call getinp
      if (iwrascii.gt.0) then
        inquire(file='tigcmloc.dat',exist=exists)
        if (exists) istat = ishell('rm tigcmloc.dat')
      endif
c
c Set up gks, conpack, and autograph:
c
      call opngks
      call gsclip(0)
      call gsasf(iasf)
      call gsfais(1)
      call cpsetup(cpspval)
      call cpseti('LLB',1)
      if (icolor.gt.0) then
        call cpseti('LBC',1)
        call cpseti('LLB',3)
      endif
      call setag
c
c Allocate space for needed fields and plot arrays:
c
c fields(imx,kmx,jmx,nfget) (will be deallocated after ut loop)
c utkmx(ntms,kmx,nfplt,nloc) (plot array for locations)
c htutkmx(ntms,kmx,nloc) = heights for extra right hand side axis (pltloc)
c utjmx(ntms,jmx,nfplt,nzphtlon) (plot array for ut vs lat)
c
      call hpalloc(pfields,imx*kmx*jmx*nfget,ier,0)
      if (ier.ne.0) then
        write(6,"('>>> tigcmloc: NOT ENOUGH MEMORY: increase ',
     +    '#QSUB -lM')")
        stop 'memory'
      endif
      if (nloc.gt.0) then
        if (nfplt.gt.0) call hpalloc(putkmx,ntms*kmx*nfplt*nloc,ier,1)
        call hpalloc(phtutkmx,ntms*kmx*nloc,ier,1)
      endif
      if (nzphtlon.gt.0.and.nfplt.gt.0) 
     +  call hpalloc(putjmx,ntms*jmx*nfplt*nzphtlon,ier,1)
c
c Ut loop:
c
      istart = 1
      iprhdr = 1
      do it=1,ntms
        call getgcmf(histvol,nhvols,istart,mtimes(1,it),fields,nfget,
     +    ifget,imx,kmx,jmx,nfhist,iden,luhist,iprhdr,ivol,istiegcm,ier)   
        if (ier.ne.0) then
          write(6,"('tigcmloc: error ',i3,' from getgcmf: it=',i3,
     +      ' mtimes=',i2,':',i2,':',i2)") ier,it,(mtimes(i,it),i=1,3)
          stop 'getgcmf'
        endif
        if (.not.istiegcm.and.(ifplt(ixpot).gt.0.or.ifplt(ixui).gt.0.or.
     +      ifplt(ixvi).gt.0.or.ifplt(ixwi).gt.0))
     +    call dynoff(ifplt(ixpot),ifplt(ixui),ifplt(ixvi),ifplt(ixwi),
     +      histvol(ivol))
        if (it.eq.1) call setflab(istiegcm,iden,flab,nftot)
        call mkfields(fields,utkmx,htutkmx,utjmx,it)
c
c Make line plots at current time (field on x-axis, zp/ht on y-axis):
c
        if (nutline.gt.0) then
          nzprange = izprange(2)-izprange(1)
          do i=1,nutline
            if (mutline(1,i).eq.mtimes(1,it).and.
     +          mutline(2,i).eq.mtimes(2,it).and.
     +          mutline(3,i).eq.mtimes(3,it)) then
              if (nzprange.gt.0) call lineut(fields,utkmx,htutkmx,it,0)
              if (nhtscale.gt.0) call lineut(fields,utkmx,htutkmx,it,1)
            endif
          enddo
        endif
c
c Dirk's doppler:
c     SUBROUTINE my_DOPTUV(LUOUT,UTC,IT,IX,JX,LX,fields)
c
        if (idop.gt.0.and.nloc.gt.0) then
          ludop = 0
          do l=1,nloc
            call getlocix(gloc(1,l),gloc(2,l),ixlat,ixlon)
            call my_doptuv(ludop,ut(it),it,ixlat,ixlon,l,fields)
          enddo
        endif
c
c End ut loop:
        istart = 0
        iprhdr = 0
      enddo
      call hpdeallc(pfields,ier,1)
c
c Contour ut on x-axis, zp/ht on y-axis:
c
      if (nloc.gt.0) then
        nzprange = izprange(2)-izprange(1)
        do l=1,nloc      
          write(6,"(' ')")
          if (iutvert.gt.0) then
            if (nzprange.gt.0) call pltloc(utkmx,htutkmx,0,l)
            if (nhtscale.gt.0) call pltloc(utkmx,htutkmx,1,l)
          endif
c
c Make line plots with ut on x-axis at selected zp and/or ht
c
          call setag
          call agseti('AXIS/LEFT/CONTROL.',1)
          call agseti('AXIS/BOTTOM/CONTROL.',1)
          if (nzphtlin.gt.0) call linezpht(utkmx,htutkmx,l)
        enddo 
      endif
c
c Contour ut on x-axis, lat on y-axis:
c
      if (nzphtlon.gt.0) then
        do l=1,nzphtlon
          write(6,"(' ')")
          call pltutlat(utjmx,l)
        enddo
      endif
c
      if (nloc.gt.0.and.nfplt.gt.0) call hpdeallc(putkmx,ier,1)
c
c Plot Dirk's doppler:
c
      if (idop.gt.0.and.nloc.gt.0) then
        do l=1,nloc
          call pltdop(l,htutkmx)
        enddo
      endif
c
      call clsgks
      if (iwrascii.eq.2.and.lenstr(dirascii).gt.0)
     +  call rcpfile(luascii,'tigcmloc.dat',dirascii)
      if (lenstr(dircgm).ne.0) call rcpfile(0,'gmeta',dircgm)
      stop 'success!'
      end
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
      subroutine dynoff(i1,i2,i3,i4,vol)
      character*(*) vol
      write(6,"(' ')")
      write(6,"('!!!!!!!!!!!!!!!!!!!! NOTICE !!!!!!!!!!!!!!!!!!!')")
      write(6,"(' ')")
      write(6,"('HISTORY VOLUME ',a,' IS NOT A TIEGCM',
     +  ' VOLUME ',/,'I AM TURNING OFF POTENTIAL AND/OR ION ',
     +  'DRIFT PLOTS')") vol(1:lenstr(vol))
      i1 = 0
      i2 = 0
      i3 = 0
      i4 = 0
      write(6,"(' ')")
      write(6,"('!!!!!!!!!!!!!!!!!!!! NOTICE !!!!!!!!!!!!!!!!!!!')")
      write(6,"(' ')")
      return
      end
c
