c
c------------------------------------------------------------------
c Begin file /home/sting/foster/plt/clrcon.f
c------------------------------------------------------------------
c
      subroutine clrcon(zz,idimx,nx,ny,xc,xd,yc,yd,cint,cmin,cmax,
     +                  iaxes,lnumx,numx,mnrx,labx,fmtx,nfmtx,
     +                  lnumy,numy,mnry,laby,fmty,nfmty,
     +                  lbar,barl,barr,barb,bart,ip)
c
c Use conpack to make color fill contour plot of zz(idimx,ny):
c
      include 'tgcmparam.h'
      include 'color.h'
      include 'cpwrk.h'
      include 'lblbar.h'
c
      dimension zz(idimx,ny)
      dimension numx(lnumx),numy(lnumy)
      character*16 labx,laby
      character*8 fmtx,fmty
      data nclrs/14/ iclu/0/, iclci/1/
      external drawcl,colram
      character*2 mapproj
c
      call cpsetr('XC1 - X COORD AT INDEX 1',xc)
      call cpsetr('XCM - X COORD AT INDEX M',xd)
      call cpsetr('YC1 - Y COORD AT INDEX 1',yc)
      call cpsetr('YCN - Y COORD AT INDEX N',yd)
c
c With CLS=14 (and CIS at its default of 0), conpack should
c   pick 'nice' contour levels such that there will be at least
c   14 levels
c
      icls = 13
 100  continue
      call cpsetr('ORV',1.e12)

      call cpgeti('MAP',imap)
      call cpgeti('SET',iset)
      call mapgtc('PR',mapproj)
      call mapgti('EL',iel)
c
c Hardwire contour interval and/or levels if desired:
      if (cint.eq.0.) then
        call cpseti('CLS -- contour level selector',icls)
      else
        call cpseti('CLS -- contour level selector',1)
        call cpsetr('CMN -- contour minimum',cmin)
        call cpsetr('CMX -- contour maximum',cmax)
        call cpsetr('CIS -- contour interval specifier',cint)
      endif
c
c Initialize conpack and pick contour levels:
c (ncl = number of contour levels)
c
      call cprect(zz,idimx,nx,ny,rwrk,lrwrk,iwrk,liwrk)
      call cppkcl(zz,rwrk,iwrk)
      call cpgeti('NCL -- number of contour levels',ncl)
c
c NCL is number of contour levels and NCLRS is number of colors
c (excluding black and white). We need NCL+1 colors for NCL levels.
c Therefore, the max number of contour levels our color table can
c handle is NCLRS-1. If NCL exceeds this, we reset the conpack level
c selector ICLS, and try it again (unless we are hardwiring interval):
c
      if (ncl.gt.nclrs-1.and.cint.eq.0.) then
c       write(6,"('>>> clrcon: ncl > nclrs-1: ncl=',i3,' nclrs=',
c    +    i3,' will try again ')") ncl,nclrs
        icls = icls-1
        goto 100
      elseif (ncl.gt.nclrs-1) then
        write(6,"('>>> clrcon warning: ncl > nclrs-1: ncl=',i3,
     +    ' nclrs=',i3)") ncl,nclrs
      endif 
c
c Loop through contour levels:
c
      do 200 i=1,ncl
        call cpseti('PAI -- parameter array index',i)
        call cpseti('CLC -- contour line color index',iclci)
        call cpseti('LLC -- line label color index',iwhite)
        call cpgetr('CLV -- contour level',clev)
        call cpseti('CLU -- contour level use flag',iclu)
        call cpgeti('AIA - area identifier above',iiaia)
        call cpgeti('AIB - area identifier below',iiaib)
        lfin(i) = ixclr(iiaib)
        lfin(i+1) = ixclr(iiaia)
        write(labclv(i),labform(ip)) clev    
 200  continue
c
c Initialize area map:
      call arinam(iama,liama)
c
c Add contour lines to area map:
      call cpclam(zz,rwrk,iwrk,iama)
c
c Scan and color the area map:
      call arscam(iama,xcra,ycra,lcra,iaia,igia,10,colram)
c
c Draw contour lines:
      call cpcldr(zz,rwrk,iwrk)
c
c Draw info labels:
      call cplbdr(zz,rwrk,iwrk)
c
c Label axes:
c
      if (iaxes.gt.0)
     +  call labaxes(lnumx,numx,mnrx,labx,fmtx,nfmtx,
     +               lnumy,numy,mnry,laby,fmty,nfmty)
c
c Color label bar:
c lbar = 0 -> do not draw color bar
c lbar = 1 -> draw vertical color bar
c lbar = 2 -> draw horizontal color bar
c
      if (lbar.gt.0) then
c
c Make a color scale bar for color filled global contours:
c (Bar is positioned horizontally below the plot)
c 
        if (lbar.eq.1) then ! vertical
          ihov = 1
          wsfb = .25
          hsfb = 1.
        else                ! horizontal
          ihov = 0 
          wsfb = 1.
c         hsfb = .25 
          hsfb = .50 
        endif
        call cpgeti('NCL - Number of contour levels',ncl)
        nbox = ncl+1
        iftp = 1
        nlbs = ncl
        lbab = 1
        call lblbar(ihov,barl,barr,barb,bart,nbox,wsfb,hsfb,lfin,
     +    iftp,labclv,nlbs,lbab)
      endif
c
      return
      end