c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
      function iyddif(iyd0,iyd1)
      dimension ndmon(12)
c                 J  F  M  A  M  J  J  A  S  O  N  D
      data ndmon/31,28,31,30,31,30,31,31,30,31,30,31/
c
c Return number of days from iyd0 to iyd1:
c
      if (iyd0.eq.iyd1) then
        iyddif = 0
        return
      endif
      if (iyd0.gt.iyd1) then
        write(6,"('iyddif: iyd0 > iyd1 (iyd0,1=',2i6,')')")
     +    iyd0,iyd1
        iyddif = -1
        return 
      endif
      if (ichiyd(iyd0).ne.0) then
        write(6,"('>>> iyddif: bad iyd0=',i6)") iyd0
        iyddif = -1
        return
      endif
      if (ichiyd(iyd1).ne.0) then
        write(6,"('>>> iyddif: bad iyd1=',i6)") iyd1
        iyddif = -1
        return
      endif
      iyr0 = iyd0/1000
      iyr1 = iyd1/1000
      iyddif = 0
      if (iyr1.gt.iyr0) then
        do i=iyr0,iyr1-1
          if (mod(i,4).ne.0) then
            iyddif = iyddif+365
          else
            iyddif = iyddif+366
          endif
        enddo
      endif
      ida0 = iyd0-iyr0*1000
      ida1 = iyd1-iyr1*1000
      iyddif = iyddif+(ida1-ida0)
      return
      end
c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
      function ichiyd(iyd)
      dimension ndmon(12)
c                 J  F  M  A  M  J  J  A  S  O  N  D
      data ndmon/31,28,31,30,31,30,31,31,30,31,30,31/
c
c Check validity of yyddd iyd:
c
      ichiyd = 0
      iyr = iyd/1000
      if (iyr.lt.1.or.iyr.gt.99) then
        write(6,"('>>> ichiyd: bad year=',i4)") iyr
        ichiyd = 1
        return
      endif
      iday = iyd-iyr*1000
      leap = 0
      ndmon(2) = 28
      if (mod(iyr,4).eq.0) then
        leap = 1 
        ndmon(2) = 29
      endif
      if (iday.lt.1.or.(leap.le.0.and.iday.gt.365).or.
     +                 (leap.gt.0.and.iday.gt.366)) then
        write(6,"('>>> ichiyd: bad iday=',i5,' (iyr=',i2,')')") 
     +    iday,iyr
        ichiyd = 1
        return
      endif
      return
      end
c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
      function inciyd(iyd,ndays)
c
c Increment iyd (yyddd) by ndays (ndays may be +/-):
c
      if (ndays.eq.0) then
        inciyd = iyd
        return
      endif
      if (ichiyd(iyd).ne.0) then
        write(6,"('>>> inciyd: bad iyd=',i6)") iyd
        return
      endif
      leap = 0
      iyr = iyd/1000
      if (mod(iyr,4).eq.0) leap = 1
      iday = iyd-iyr*1000
      inc = 1
      if (ndays.lt.0) inc = -1
      do i=1,abs(ndays)
        iday = iday+inc
        if (leap.le.0.and.iday.gt.365) then
          iday = 1
          iyr = iyr+1
        endif
        if (leap.gt.0.and.iday.gt.366) then
          iday = 1
          iyr = iyr+1
        endif
        if (iday.lt.1) then
          iyr = iyr-1
          if (mod(iyr,4).eq.0) then
            iday = 366
          else
            iday = 365
          endif
        endif
      enddo
      inciyd = iyr*1000+iday
      return
      end
c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
      subroutine isec2hms(isec,ih,im,is)
c
c Given integer seconds isec (not including days), return 
c integer hour, minute, seconds in ih,im,is: 
c
      ih = isec/3600
      nsec = mod(isec,3600)
      im = nsec/60
      is = mod(nsec,60)
      return
      end
c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
      function lenstr(str)
      character*(*) str
c
c Return index to last non-blank char in str
c
      length = len(str)
      do i=length,1,-1
        if (str(i:i).ne.' ') then
          lenstr = i
          return
        endif
      enddo
      lenstr = 0
      return
      end