c
c------------------------------------------------------------------
c Begin file /home/sting/foster/tigcmloc/pltutlat.f
c------------------------------------------------------------------
c
      subroutine pltutlat(utf,l)
c
c Contour ut on x-axis, lat on y-axis at selected zp or ht:
c
      include 'tigcmloc.h'
      character*56 toplab, fieldlab
      character*80 rec80
      dimension viewport(4),plt(ntms,jmx),
     +  utf(ntms,jmx,nfplt,nzphtlon)
      data viewport /.15,.89,.26,.91/
      data ncalls/0/
      save ncalls
c
      ncalls = ncalls+1
      xmid = 0.5*(viewport(1)+viewport(2))
      rily = -.17
      igmflag = ifix(gmflag)
      if (izphtlon(2,l).ne.igmflag) rily = -.34
      call cpseti('SET',0)
      call cpseti('MAP',0)
      call cpsetr('XC1',utinc(1))
      call cpsetr('XCM',utinc(ntms))
      call cpsetr('ILX',xmid)
      call cpsetr('ILY',rily)
      call cpseti('ILP',0)
      call cpsetr('ILS',.016)
      call cpsetr('YC1',gcmlat(1))
      call cpsetr('YCN',gcmlat(jmx))
      call set(viewport(1),viewport(2),viewport(3),viewport(4),
     +         utinc(1),utinc(ntms),-90.,90.,1)
c
c Field loop: 
c (fof2 and hmf2 are height independent -- they may be plotted
c  only once by this routine, and not at all if plotting at
c  interpolated heights (izphtlon(1,l)=0))
c
      do 200 ip=1,nftot
        if (ifplt(ip).le.0) goto 200
        if ((ip.eq.ixfof2.or.ip.eq.ixhmf2).and.izphtlon(1,l).le.0)
     +    goto 200
c       if (ncalls.gt.1.and.(ip.eq.ixfof2.or.ip.eq.ixhmf2)) goto 200
        plt(:,:) = utf(:,:,ifplt(ip),l)
        if (logutlat.gt.0.and.logterp(ip).gt.0) 
     +    call log10f(plt,ntms*jmx,1.e-20,cpspval)
c
c Contour:
c
        if (icolor.le.0) then
          call contour(plt,ntms,ntms,jmx,cint(ip),cmin(ip),cmax(ip))
        else
          call conclr(plt,ntms,ntms,jmx,cint(ip),cmin(ip),cmax(ip))
        endif
c
c Add axes labels:
c
c       call labrect(utinc,ntms,gcmlat,jmx,'UT (HRS)',
c    +    'LATITUDE (DEG)',0.)
        isltax = 0
        if (izphtlon(2,l).ne.igmflag) isltax = 1
        call labutxy(mtimes,ntms,gcmlat,jmx,'LATITUDE (DEG)',
     +    0.,isltax,zphtlon(2,l))    
c
c Field label at top:
c
        call clearstr(toplab)
        if (logutlat.gt.0.and.logterp(ip).gt.0) then
          write(toplab,"('LOG10 ',a)") flab(ip)
        else
          write(toplab,"(a)") flab(ip)
        endif
        call wrlab(toplab(1:lenstr(toplab)),xmid,viewport(4)+0.07,0.)
        fieldlab = toplab
c
c "location" lab:
c
        call clearstr(toplab)
        xoff = .15
c
c At selected zp and longitude:
        if (izphtlon(1,l).gt.0.and.izphtlon(2,l).ne.igmflag) then
          if (ip.ne.ixfof2.and.ip.ne.ixhmf2) then
            write(toplab,"('ZP=',f5.1,' LON=',f7.2)") 
     +        (zphtlon(i,l),i=1,2)
          else
            write(toplab,"(' LON=',f7.2)") zphtlon(2,l)
          endif
          xoff = .25
c
c At selected zp and zonal means:
        elseif (izphtlon(1,l).gt.0.and.izphtlon(2,l).eq.igmflag) then
          if (ip.ne.ixfof2.and.ip.ne.ixhmf2) then
            write(toplab,"('ZP=',f5.1,' (ZONAL MEANS)')") zphtlon(1,l)
          else
            write(toplab,"('(ZONAL MEANS)')")
          endif
c
c At selected height and longitude:
        elseif (izphtlon(1,l).le.0.and.izphtlon(2,l).ne.igmflag) then
          write(toplab,"('HT=',f6.2,' LON=',f7.2)") 
     +      (zphtlon(i,l),i=1,2)
          xoff = .25
c
c At selected height and zonal means (interpolate first, then take means)
        elseif (izphtlon(1,l).le.0.and.izphtlon(2,l).eq.igmflag) then
          write(toplab,"('HT=',f6.2,' (ZONAL MEANS)')") zphtlon(1,l)
        endif
        call wrlab(toplab(1:lenstr(toplab)),xmid,viewport(4)+0.03,0.)
        call clearstr(rec80)
        write(rec80,"(a)") toplab(1:lenstr(toplab))
c
c History volume label at bottom:
c
        call clearstr(toplab)
        if (nhvols.eq.1) then
          write(toplab,"('HISTORY=',a)") histvol(ivol)
        else
          write(toplab,"('FIRST HISTORY VOL=',a)") histvol(ivol)
        endif
        call wrlab(toplab(1:lenstr(toplab)),xmid,viewport(3)-xoff,
     +    .010)
c
c Wrap it up:
        call frame
        iframe = iframe+1
        if (iwrascii.gt.0) then
          call wrascii(iwrascii,luascii,plt,ntms,jmx,utinc,gcmlat,
     +      'UT (HRS)','LATITUDE (DEG)',histvol(ivol),fieldlab,
     +      iframe,rec80,'tigcmloc',dirascii)
        endif
        if (izphtlon(1,l).gt.0.and.izphtlon(2,l).ne.igmflag) then
          if (ip.ne.ixfof2.and.ip.ne.ixhmf2) then
            write(6,"('utlat frame ',i3,' field ',a,' zp=',f5.1,
     +        ' lon=',f8.2)") iframe,labshort(ip),zphtlon(1,l),
     +        zphtlon(2,l)
          else
            write(6,"('utlat frame ',i3,' field ',a,9x,
     +        ' lon=',f8.2)") iframe,labshort(ip),zphtlon(2,l)
          endif
        elseif (izphtlon(1,l).gt.0.and.izphtlon(2,l).eq.igmflag) then
          if (ip.ne.ixfof2.and.ip.ne.ixhmf2) then
            write(6,"('utlat frame ',i3,' field ',a,' zp=',f5.1,
     +        ' (zonal means)')") iframe,labshort(ip),zphtlon(1,l)
          else
            write(6,"('utlat frame ',i3,' field ',a,9x,
     +        ' (zonal means)')") iframe,labshort(ip)
          endif
        elseif (izphtlon(1,l).le.0.and.izphtlon(2,l).ne.igmflag) then
          write(6,"('utlat frame ',i3,' field ',a,' ht=',f5.1,
     +      ' lon=',f8.2)") iframe,labshort(ip),zphtlon(1,l),
     +      zphtlon(2,l)
        elseif (izphtlon(1,l).le.0.and.izphtlon(2,l).eq.igmflag) then
          write(6,"('utlat frame ',i3,' field ',a,' ht=',f5.1,
     +      ' (zonal means)')") iframe,labshort(ip),zphtlon(1,l)
        endif
c
c End field loop
 200  continue
c
      return
      end