c  This program reads t and s fields and outputs level fields

      parameter (imax=162)
      parameter (jmax=110)

      parameter (maxlev=401)
      dimension tkp(4),skp(4)

      dimension tkpl(maxlev),skpl(maxlev)

C  Open input file
      open(10,file='ts_layer.asc',status='old',readonly)

C  Open output file
      open(11,file='ts_level.dat',status='new',form='unformatted')

      pi=4.*abs(atan(1.))
      rearth=6370.e5
      dypdeg=pi/180.*rearth
      ymin=-29.*dypdeg
      ymax=+25.*dypdeg
      xmin=0.*dypdeg
      xmax=+80.*dypdeg
      dx=(xmax-xmin)/(imax-2)
      dy=(ymax-ymin)/(jmax-2)

5     continue

      read(10,*,end=900)h1b,h2b,h3b,h4b
      read(10,*)djd,tkp,skp

c  Interpolate velocity level fields from convergences
      do m=1,maxlev
c      Find the desired layer
       dep=(m-1)
       if(dep.le.h1b) then
         lay=1
       elseif(dep.le.(h1b+h2b)) then
         lay=2
       elseif(dep.le.(h1b+h2b+h3b)) then
         lay=3
       elseif(dep.le.(h1b+h2b+h3b+h4b)) then
         lay=4
       else
         lay=4
       endif

       tkpl(m)=tkp(lay)
       skpl(m)=skp(lay)

      enddo

      write(11)tkpl
      write(11)skpl

      goto 5
900   close(10)
      close(11)

      end