subroutine xtr_trc(lon,lat,dpth,map_t,imax,kmax)
c
c=======================================================================
c                                                                    ===
c  This subroutine extracts tracers from the model grid and then     ===
c  interpolates them to the sub-domain grid (lon,lat,dpth).          ===
c                                                                    ===
c  On Input:                                                         ===
c                                                                    ===
c     LON     sub-domain T-point longitude (degrees west, real array)===
c     LAT     sub-domain T-point latitude (degrees north, real array)===
c     DPTH    sub-domain depths (centimeter, real array)             ===
c     KMAX    number of points in the z-direction to interpolate     ===
c             (integer)                                              ===
c     IMAX    number of points in the x-direction to interpolate     ===
c             (integer)                                              ===
c                                                                    ===
c  On Output:                                                        ===
c                                                                    ===
c     MAP_T   sub-domain tracers (real array)                        ===
c                                                                    ===
c  Calls:  DEPTHSLAB, LL2XY, OPICK                                   ===
c                                                                    ===
c=======================================================================
c
c-----------------------------------------------------------------------
c  Define global data.
c-----------------------------------------------------------------------
c
#include <cdefs.h>
#include <param.h>
#include <pconst.h>
#include <fullwd.h>
#include <moddat.h>
#include <rhomean.h>
c
c-----------------------------------------------------------------------
c  Define local  data.
c-----------------------------------------------------------------------
c
      integer i,ic,imax,ip,j,jold,k,kkc,kmax,l,m,net,nst
#ifndef barotropic
     *        ,kc
#endif
      integer icell(4)
      FLOAT
     *      x,xdis,y,ydis,z
      FLOAT
     *      depth(imt,km,2),dpth(ximtkm),lat(ximt),lon(ximt),
     *      map_t(ximtkm,nt),tcell(4,nt),trc(imt,km,2,nt),
     *      tzcell(km,4,nt),zcell(km,4)
c
c=======================================================================
c  Begin executable code.
c=======================================================================
c
c-----------------------------------------------------------------------
c  Extract and linearly interpolate tracers to sub-domain grid.
c-----------------------------------------------------------------------
c
      jold=0
      do 110 l=1,imax
c
c   Find indices for horizontal grid cell.
c
        call ll2xy (lon(l),lat(l),coord,imt,jmt,gridx,gridy,rlngd,rlatd,
     *              delx,dely,thetad,x,y)
        i=int(x)
        j=int(y)
        xdis=x-FLoaT(i)
        ydis=y-FLoaT(j)
c
        icell(1)=i
        icell(2)=i+1
        icell(3)=i+1
        icell(4)=i
c
c   Get slab of data containing grid cell.
c
        if(j.ne.jold) then
          do 10 m=1,nt
            nst=1+imtkm*(m-1)
            net=nst+imtkm-1
            call opick(labs(ndisk),nslab,(j-1)*nslab+1,nst,net,
     *                 trc(1,1,1,m))
            call opick(labs(ndisk),nslab,(j  )*nslab+1,nst,net,
     *                 trc(1,1,2,m))
  10      continue
c
c  Get depths at slabs u-points.
c
          call depthslab(j  ,tgrid,depth(1,1,1))
          call depthslab(j+1,tgrid,depth(1,1,2))
          j=jold
        endif
c
c  Pass information from slabs to a vertical column.
c
        do 20 k=1,km
          zcell(k,1)=depth(icell(1),k,1)
          zcell(k,2)=depth(icell(2),k,1)
          zcell(k,3)=depth(icell(3),k,2)
          zcell(k,4)=depth(icell(4),k,2)
        do 20 m=1,nt
          tzcell(k,1,m)=trc(icell(1),k,1,m)
          tzcell(k,2,m)=trc(icell(2),k,1,m)
          tzcell(k,3,m)=trc(icell(3),k,2,m)
          tzcell(k,4,m)=trc(icell(4),k,2,m)
  20    continue
c
c  Interpolation.
c
        do 100 k=1,kmax
c
c  Vertical interpolation.
c
          ip=k+(l-1)*kmax
          z=abs(dpth(ip))*m2cm
          do 80 ic=1,4
            if(z.le.zcell(1,ic)) then
              do 30 m=1,nt
                tcell(ic,m)=tzcell(1,ic,m)
  30          continue
            elseif(z.ge.zcell(km,ic))then
              do 40 m=1,nt
                tcell(ic,m)=tzcell(km,ic,m)
  40          continue
            else
#ifndef barotropic
              do 50 kc=1,km-1
                if((zcell(kc,ic).lt.z).and.
     *             (z.le.zcell(kc+1,ic))) then
                  kkc=kc
                  goto 60
                endif
  50          continue
  60          continue
#endif
              do 70 m=1,nt
                tcell(ic,m)=((z-zcell(kkc,ic))*tzcell(kkc+1,ic,m)+
     *                       (zcell(kkc+1,ic)-z)*tzcell(kkc,ic,m))/
     *                      (zcell(kkc+1,ic)-zcell(kkc,ic))
  70          continue
            endif
  80      continue
c
c  Horizontal interpolation.
c
          do 90 m=1,nt
            map_t(ip,m)=(c1-ydis)*((c1-xdis)*tcell(1,m)+
     *                             xdis*tcell(2,m))+
     *                  ydis*((c1-xdis)*tcell(4,m)+xdis*tcell(3,m))
            if(m.eq.2) map_t(ip,m)=map_t(ip,m)+smean
  90      continue
 100    continue
 110  continue
      return
      end