integer k1,k2,nymd1,nhms1,nymd2,nhms2,i

      call TIME_BOUND  ( nymd,nhms, nymd1,nhms1, nymd2,nhms2, k1,k2 )

      zoch(i) = vgz0(k2,ityp(i))*facp + vgz0(k1,ityp(i))*facm

      rdcs(i) = vgrd(k2,ityp(i))*facp + vgrd(k1,ityp(i))*facm

      rootl   = vgrt(k2,ityp(i))*facp + vgrt(k1,ityp(i))*facm

      d       = vgdd(k2,ityp(i))*facp + vgdd(k1,ityp(i))*facm

      do 3000 k2=k1+1,np+1

           call TABLUP(k1,k2,m,n,np,nx,nh,nt,sh2o,swpre,swtem,

           call TABLUP(k1,k2,m,n,np,nx,nh,nt,sh2o,swpre,swtem,

           call TABLUP(k1,k2,m,n,np,nx,nh,nt,sh2o,swpre,swtem,

        call H2OKDIS(ib,m,n,np,k2-1,fkw,gkw,ne,h2oexp,conexp,

          call TABLUP(k1,k2,m,n,np,nx,nc,nt,sco3,scopre,scotem,

          call CO2KDIS(m,n,np,k2-1,co2exp,tco2,trant)

          call TABLUP(k1,k2,m,n,np,nx,no,nt,sco3,scopre,scotem,

          call N2OKDIS(ib,m,n,np,k2-1,n2oexp,tn2o,trant)

          call CH4KDIS(ib,m,n,np,k2-1,ch4exp,tch4,trant)

          call COMKDIS(ib,m,n,np,k2-1,comexp,tcom,trant)

          call CFCKDIS(m,n,np,k2-1,f11exp,tf11,trant)

          call CFCKDIS(m,n,np,k2-1,f12exp,tf12,trant)

          call CFCKDIS(m,n,np,k2-1,f22exp,tf22,trant)

          call B10KDIS(m,n,np,k2-1,h2oexp,conexp,co2exp,n2oexp

         ff=0.5+(0.3739+(0.0076+0.1185*asyal(i,j,k2-1,ib))

     *      *asyal(i,j,k2-1,ib))*asyal(i,j,k2-1,ib)

         taux=taual(i,j,k2-1,ib)*(1.-ssaal(i,j,k2-1,ib)*ff)

        if( k2 .le. ict ) then

            clrhi(i,j)=min(clr(i,j,k2-1),clrhi(i,j))

        elseif( k2 .gt. ict .and. k2 .le. icb ) then

            clrmd(i,j)=min(clr(i,j,k2-1),clrmd(i,j))

        elseif( k2 .gt. icb ) then

            clrlw(i,j)=min(clr(i,j,k2-1),clrlw(i,j))

       if (k2 .eq. k1+1) then

          flcd(i,j,k2)=flcd(i,j,k2)+blayer(i,j,k1)

          flxd(i,j,k2)=flxd(i,j,k2)+blayer(i,j,k1)

         xx=trant(i,j)*dlayer(i,j,k2)

         flcd(i,j,k2) =flcd(i,j,k2)+xx

         flxd(i,j,k2) =flxd(i,j,k2)+xx*fclr(i,j)

     *                   *dlayer(i,j,k2)

           rflx(i,j,k2) = rflx(i,j,k2)+trant(i,j)*fclr(i,j)

     *                   *dlayer(i,j,k2)

           rflc(i,j,k2) = rflc(i,j,k2)+trant(i,j)

      subroutine TABLUP(k1,k2,m,n,np,nx,nh,nt,sabs,spre,stem,w1,p1,

      integer k1,k2,m,n,np,nx,nh,nt,i,j

        x1=sabs(i,j,k2)-sabs(i,j,k1)

        x2=(spre(i,j,k2)-spre(i,j,k1))/x1

        x3=(stem(i,j,k2)-stem(i,j,k1))/x1

      integer i,j,k,ip,iw,it,ib,ik,iq,isb,k1,k2

      INTEGER km, k1, k2

        k2 = myNr

        k2 = kLoc

        DO k=k1,k2