subroutine yd2jd( yd, jd )
!
!     Given Gregorian calendar date of the form yyyyddd.dd, returns
!     the Julian Day Number (and fraction).  (Year "1 B.C." is
!     Gregorian year 0, year "2 B.C. is Gregorian year -1, etc.)
!
!     B. Knapp, 95.09.25
!
!     Caution:  This routine should not be used with historical dates
!     of historical events prior to the adoption of the Gregorian
!     calendar (1582 to 1918, depending on country); instead, for
!     such dates a *Julian* calendar-date-to-JD # service should be
!     used.
!
C
C     RCS DATA
C     
C     $Header$
C     
C     $Log$
C
C
      implicit none
!
!     Input:
      real*8 yd
!
!     Output:
      real*8 jd
!
!     Local:
      real*8 ydabs,d
      integer*4 y
!
!
!     Extract year and day number
      ydabs = abs( yd )
      y = int( ydabs/1000.0d0 )
      d = mod( ydabs, 1000.0d0 )
!
!     Add 10000 years (to handle negative years)
      if ( yd .ge. 0 ) then
         y = 10000+y
      else
         y = 10000-y
      endif
!
!     Compute the Julian day number by adding the appropriate number
!     of days for the elapsed years and leap-year days to the JD #
!     for day 1.0, year 0, and then removing the 10000-year offset.
      jd = 1721059.5d0 + 365.0d0*y +
     &     dble( (y-1)/4 - (y-1)/100 + (y-1)/400 ) - 3652425.d0 + d
      return
      end