subroutine gfland_table(rlon,rlat,fol)
c
c     VERSION 1.0.0
c     Last updated: 10 Feb 2020
c     Modified Feb 2023 (KM) for NCO specifications
c
c     This program calculates the fraction of land (fol)
c
c     Input: 
c       rlon,rlat - (deg W neg, and deg N pos)
c
c     Output: 
c       fol - Fraction of land (0 to 1)
c
      parameter (mx=3601,my=1801)
c
      dimension gftab(mx,my)
      character*200 include_path, fntab
c
      data iread /1/
      data iperiodic /1/
c
      save gftab
      save slat,elat,dlat,nlat,slon,elon,dlon,nlon,iread
c
      common /oper/ ioper
c
c      if     (ioper .eq. 1) then
c         call getenv ( "ATCFPROBLTY", include_path )      ! ATCF
c      elseif (ioper .eq. 2) then
c         call getenv ( "PRBLTYINC", include_path )        ! IBM
c      else
c         include_path = "include"                         ! local
          call getenv ("SHIPS_COEF", include_path )
c      endif
      ind = index( include_path, " " ) - 1

c     Set defaults
      fol      = 9999.!-999.
      ierr     = 1
c
      if (iread .eq. 1) then
c        Read the data tables
c
         lutab = 27
         fntab = include_path(1:ind)//"/gfland_table.dat"
         open(file=TRIM(fntab),unit=lutab,form='formatted',
     +        status='old',err=927)
c
         read(lutab,*,err=927,end=927) slon,elon,dlon,nlon,
     +                         slat,elat,dlat,nlat
c
         do j=1,nlat
            read(lutab,100,err=927,end=927)  
     +         (gftab(i,j),i=1,nlon)
  100       format(10(f9.1))
         enddo
c
c        Add a column at lon=360 matching lon=0 
         if (iperiodic .eq. 1) then
            nlon=nlon+1
            elon=elon+dlon
            do j=1,nlat
               gftab(nlon,j)=gftab(1,j)
            enddo
         endif
c
         iread=0
      endif
c    
c     Make sure rlon is between 0 and 360
      rlontemp = mod(rlon+360., 360.)
      rlon = rlontemp

c     Check input values
      if (rlat .lt. slat .or. rlat .gt. elat) then
         ierr=2
         call gfland_err_handling(1,6,ierr,
     +                            slon,elon,slat,elat,rlon,rlat,fol)
      endif
c
      if (rlon .lt. slon .or. rlon .gt. elon) then
         ierr=3
         call gfland_err_handling(1,6,ierr,
     +                            slon,elon,slat,elat,rlon,rlat,fol)
      endif
c
c
c     Calculate indices of table value closest but to the lower left
c     of the requested lat/lon
      ii = 1 + ifix( (rlon-slon)/dlon )
      jj = 1 + ifix( (rlat-slat)/dlat )
c
      if (ii .eq. nlon) ii = nlon-1
      if (jj .eq. nlat) jj = nlat-1
c
c     Calculate normalized x,y values relative to lower left table point
      xlon = slon + dlon*float(ii-1)
      x    = (rlon-xlon)/dlon
c
      ylat = slat + dlat*float(jj-1)
      y    = (rlat-ylat)/dlat
c
      xy = x*y
c
c     Interpolate table values to find fol and afraction
      w00 = 1.0+xy-(x+y)
      w10 = x-xy
      w01 = y-xy
      w11 = xy
c
      fol = w00*gftab(ii  ,jj  ) +
     +      w10*gftab(ii+1,jj  ) +
     +      w01*gftab(ii  ,jj+1) +
     +      w11*gftab(ii+1,jj+1)
c
      return
c
  927 continue
      call gfland_err_handling(1,6,ierr,
     +                         slon,elon,slat,elat,rlon,rlat,fol)
      return 
c
      end

      subroutine gfland_err_handling(ierrtype,luerr,ierr,
     +                               slon,elon,slat,elat,rlon,rlat,fol)
      ! This routine handles error processing for tcliper.
      ! Types of errors handled (ierrtype=):
      !    -1: error in gdland_table
      ! Input
      !    ierrtype - error type to process
      !    luerr    - unit number to write to (assumes file is already open)
      !    ierr     - error code to report

      IMPLICIT NONE
 
      !++ Passed variables
      integer, intent(in) :: ierrtype, luerr, ierr
      real, intent(in) :: slon,elon,slat,elat,rlon,rlat
      real, intent(inout) :: fol

      if (ierrtype .eq. 1) then
         write(luerr,*) 'Error in subroutine gfland_table, ierr=',ierr
         write(luerr,*) 'slon,elon,slat,elat,rlon,rlat: ',
     +                   slon,elon,slat,elat,rlon,rlat 
         fol=0.0
      else
         write(luerr,*) 'Unrecognized error in gfland_table'
         stop
      endif

      return
      end subroutine gfland_err_handling