
>>>>>>>>>>>>>>>>>>> FILE: ../ephrd/eph_const_ini_mod.f90 <<<<<<<<<<<<<<<<<<<

module eph_const_ini_mod

use constants_mod

implicit none
save

contains

subroutine eph_const_ini(system)

!===============================================================================
! eph_const_ini initializes a few variables in the constants_mod based on 
!     values that are extracted from the JPL Ephemerides or are otherwise 
!     unknown until run time. Specifically, the following constants are 
!     initialized with this routine.
!
!        au
!        au_per_day
!        gm_conv
!        v_light 
!        earth_moon_ratio
!        earth_moon_ratio_p1
!        sin_eps 
!        cos_eps
!
! Input:
!    system - either 'J2000' or 'B1950'
!-------------------------------------------------------------------------------
! Written       Nov 13, 2000       Steve Chesley
! Modified      Feb 20, 2002       Steve Chesley
!    - Use sec_per_day & au_per_day.
! Modified      Mar 06, 2002       Steve Chesley
!    - Update gm_conv based on new value of au.
! Modified      Sep 17, 2003       Steve Chesley
!    - Update this documentation to reflect parameters added in 2002.
!
! $Id: eph_const_ini_mod.f90,v 1.6 2004/12/01 20:52:29 chesley Exp $
!$==============================================================================

>>>>>>>>>>>>>>>>>>> FILE: ../ephrd/open_eph_mod.f90 <<<<<<<<<<<<<<<<<<<

module open_eph_mod

use constants_mod
use get_unit_mod
use eph_const_ini_mod
use intpar_mod

implicit none

contains

!===============================================================================
subroutine open_eph( object, jd_min, jd_max, tref, quiet, de_name, pleph_name, &
   sbeph_name, system, jdeqnx, j2000, jg_sep, file_10, req_num_sb, &
   req_sb_name, asteroid, object_name, iau_number, designation)

character (len=*), intent(in)            :: object
real (kind=wp),    intent(out)           :: jd_min, jd_max
real (kind=wp),    intent(in),  optional :: tref
logical,           intent(in),  optional :: quiet
character (len=*), intent(out), optional :: de_name
character (len=*), intent(out), optional :: pleph_name
character (len=*), intent(out), optional :: sbeph_name
character (len=*), intent(out), optional :: system
real (kind=wp),    intent(out), optional :: jdeqnx
logical,           intent(out), optional :: j2000, jg_sep, file_10
integer,           intent(out), optional :: req_num_sb
character (len=*), intent(out), optional :: req_sb_name(:)
logical,           intent(out), optional :: asteroid
integer,           intent(out), optional :: iau_number
character (len=*), intent(out), optional :: object_name, designation

!-------------------------------------------------------------------------------

! open_eph: Open the planetary and small bodies ephemerides as
!     appropriate based on values set in the VARLIST. Print the list of
!     perturbing bodies unless the 'quiet' flag is set. Return the combined
!     ephemerides start and stop times (JED) and optionally return other
!     descriptive parameters related to the ephemerides and the object.

!
! NOTE: The object's ORB file must be loaded into varlist, typically using 
!  rdorbf().
!
! The planetary ephemerides are taken from the first found among the following:
!    1) VARLIST variable PEFILE
!    2) Environment variable $PEFILE
!    3) $CAET_DATA/pleph.nio
!
! The small body ephemerides are taken from the first found among the following:
!    1) VARLIST variable SBFILE
!    2) Environment variable $SBFILE
!    3) $CAET_DATA/sbeph.nio
!
! The VARLIST array SBPERT can be used to limit the number of perturbing
! small bodies.  If SBPERT(1) = 'T' (case insensitive) then all small
! bodies in the file will be used, if SBPERT(1) = 'F' (case insensitive)
! then none will be used.  If only a subset of the small bodies in the
! file should be used then the array SBPERT() should contain the names
! of the perturbing subset (e.g., 'Vesta', 'Halley', 'P/Halley', '2061
! XX').  If the object name matches one of the small body names (case
! insensitive), then that body will be removed from the small_bodies
! array and will not act as a perturber.
!
! The overlapping time span (jd_min -> jd_max) of the planetary and
! small body (if used) files is returned.  If there is no overlap, this
! routine aborts.
!
! The reference system defined by the planetary is used to initialize
! related constants and variables by calling eph_const_ini(system).
!
!-------------------------------------------------------------------------------
! Input:
!   object       Object line from ORB file.
!
! Output:
!   jd_min       Earliest ephemeris time (JED).
!   jd_max       Latest ephemeris time (JED).
!
! Optional Input:
!   tref         Reference time (JED).
!   quiet        Flag set true to suppress 'Perturbing Bodies' log-style output.
!
! Optional Output:
!   dename       The full (7 character) name of the planetary ephemerides 
!                (e.g. 'DE-0405').
!   pleph_name   The (5 character) name of the planetary ephemerides.
!   sbeph_name   The (12 character) name of the small body ephemerides.
!   system       Reference system ('J2000' or 'B1950').
!   j2000        Flag set true if the ephemeris is J2000 (B1950 assumed 
!                otherwise).
!   jg_sep       Flag set true if the user requested separation of 
!                Jupiter / Galilean Sats.
!   file_10      Flag set true if the user has requested "FILE 10" output.
!   req_num_sb   The number of small-bodies requested.
!   req_sb_name  Array of the small-body names requested.
!   asteroid     Flag set true if the object is not a comet.
!   object_name  The object name (as returned by PARSOB).
!   iau_number   The object IAU number (as returned by PARSOB).
!   designation  The object designation (as returned by PARSOB).
!
!-------------------------------------------------------------------------------
!  Written:  2001-Sep-28  by  Alan B. Chamberlin
!
!  Modified: 2001-Nov-09  by  Alan B. Chamberlin
!   - redesigned the interface to better match a similar module written by 
!     Steve Chesley called "open_ephem_mod".  Included most of Steve's 
!     excellent documentation!
!   - added "open_eph_def" subroutine to initialize all internal VARLIST 
!     variables.
!   - added initialization of the Galilean satellite theory.
!   - added many optional output arguments and coverted some args to optional.
!  Modified: 2001-Nov-13  by  Steven R. Chesley
!   - Use arbitrary unit number when opening GSTFILE.
!   - VLDEF the SBPERT array at the same dimension as max_num_sb.
!   - Use errmsg to report warning when the object's name is also on the 
!     perturber list.
!   - Bug fix in 'asteroid' assignment statement.
!  Modified: 2001-Dec-17  by  Steven R. Chesley
!   - Tweak documentation.
!  Modified: 2002-Jan-24  by  Steven R. Chesley
!   - Do not issue warning when SB perturbers are not used in present solution
!     so long as they are used in current solution. Do issue warning whenever 
!     SB perturbers are _not_ requested.
!  Modified: 2002-Jan-25  by  Steven R. Chesley
!   - Do not issue warning when PENAM is empty.
!  Modified: 2002-Mar-07  by  Steven R. Chesley
!   - Changed intpar call to intpar_ini to enable use of upgraded F90 
!     propagator routines.
!  Modified: 2003-Mar-06  by  Steven R. Chesley
!   - Moved call to eph_const_ini to before intpar_ini to allow the
!     constants to be updated for intpar.
!  Modified: 2003-Mar-26  by  Alan B. Chamberlin
!   - added 'dename' and 'jdeqnx' to the list of optional arguments 
!  Modified: 2004-Dec-01  by  Steven R. Chesley
!   - Changed format of documentation to facilitate auto-extraction
!  Modified: 2004-Dec-17  by  Steven R. Chesley
!   - Used quiet optional arg to also control output of warnings
!  Modified: 2007-Aug-6  by  Steven R. Chesley
!   - Increased size of sb_name from 40 to 120 characters
!
!  $Id: open_eph_mod.f90,v 1.9 2007/08/14 22:09:16 chesley Exp $
!-------------------------------------------------------------------------------
! TO DO:
!  * consider adding a call to getsta with optional returned parameters
!$==============================================================================

!^==============================================================================
subroutine open_eph_def

!-------------------------------------------------------------------------------
! open_eph_def:  Define all VARLIST variables unique to open_eph.
!$------------------------------------------------------------------------------

>>>>>>>>>>>>>>>>>>> FILE: ../ephrd/peread.f <<<<<<<<<<<<<<<<<<<

      SUBROUTINE PEOPEN( PEFILE, DENAM, JDREF, PESTRT, PESTOP)
      CHARACTER*(*)      PEFILE, DENAM
      DOUBLE PRECISION   JDREF, PESTRT, PESTOP
C
C     SUBROUTINE PECACH( PCSTRT, PCSTOP)
      DOUBLE PRECISION   PCSTRT, PCSTOP
C
C     SUBROUTINE PSTATH( T, TREF, LIST, PV)
      DOUBLE PRECISION   T, TREF, PV(6,10)
      INTEGER            LIST(10)
C
C     SUBROUTINE PSTATB( T, TREF, LISTB, PVB)
      INTEGER            LISTB(0:10)
      DOUBLE PRECISION   PVB(6,0:10)
C
C     SUBROUTINE PENUT ( T, TREF, REQNUT, NUT)
      INTEGER            REQNUT
      DOUBLE PRECISION   NUT(4)
C
C     SUBROUTINE PECON ( NAM, NV, VAL)
      CHARACTER*(*)      NAM(*)
      INTEGER            NV
      DOUBLE PRECISION   VAL(*)
C
C=======================================================================
C              PEREAD: Planetary Ephemeris Readers Package
C
C-----------------------------------------------------------------------
C  PEOPEN (Planetary Ephemeris OPEN) opens the planetary ephemeris file
C  and returns the following information about the file: the development
C  number, the reference epoch, and the time span of the ephemeris.
C  PEOPEN must be called before any other routine in this package is
C  called.
C
C  Input:
C
C   PEFILE  File name of the planetary ephemeris file.
C
C  Outputs:
C
C   DENAM   Development Ephemeris name for the file (7 characters)
C           (eg, 'DE-0118').
C
C   JDREF   Reference epoch of the equinox, JED.
C
C   PESTRT  Start epoch for the ephemeris, JED.
C
C   PESTOP  Stop epoch for the ephemeris, JED.
C
C-----------------------------------------------------------------------
C  PECACH (Planetary Ephemeris CACHe) defines the limits of the cache
C  for planetary ephemeris data records.  The cache is filled with
C  records for the requested time interval when the planetary ephemeris
C  is opened.  The cache is filled from the start epoch towards the stop
C  epoch, which may be either forward or backward in time; if there is
C  not enough cache space for the requested interval, it will be
C  truncated on the 'stop' end.  If this entry point is not called, the
C  cache is not used.
C
C  Inputs:
C
C   PCSTRT  Start epoch for the cache, JED.
C
C   PCSTOP  Stop epoch for the cache, JED.
C
C-----------------------------------------------------------------------
C  PSTATH (Planetary ephemeris STATe, Heliocentric) reads and interpo-
C  lates the planetary ephemeris for the inertial-frame components of
C  the heliocentric positions and velocities of specified Solar System
C  bodies.  Bodies are specified by setting appropriate elements of the
C  request list LIST to nonzero values.  Elements 1 through 9 of the
C  list refer to the barycenters of the 9 planetary systems (eg, 3
C  indicates the Earth-Moon barycenter); element 10 indicates the Moon.
C  The user can request either position only or position and velocity
C  for any of these bodies by setting the appropriate request list
C  element to 1 or 2, respectively.  The states of all bodies except the
C  Moon are heliocentric; the state of the Moon is geocentric.  If
C  states relative to the Solar System barycenter are required, use
C  routine PSTATB).
C
C  Inputs:
C
C   T              Time (ephemeris days wrt reference time).
C
C   TREF           Reference time (JED).  For best accuracy, this should
C                  have a fractional part of ".0" or ".5".
C
C   LIST(1..10)    List of indicators specifying the interpolations
C                  desired for each of the 10 bodies:
C                  LIST(i)= 0  no interpolation for body i,
C                           1  position only for body i, and
C                           2  position and velocity for body i.
C
C  Output:
C
C   PV(1..6, 1..10)  State components for the specified bodies.  The
C                    first subscript varies over the components,
C                    which are ordered X, Y, Z, VX, VY, VZ, where V
C                    indicates velocity; the second subscript
C                    indicates the body.  The units are km and km/s.
C
C-----------------------------------------------------------------------
C  PSTATB (Planetary ephemeris STATe, Barycentric) performs the same
C  function as PSTATH, but returns states relative to the Solar System
C  barycenter instead of the Sun.  (The state of the Moon, however is
C  still geocentric.) In addition, this routine returns the state of the
C  Sun relative to the Solar System barycenter, if requested via the 0th
C  element of LISTB.  Note that even though the calling sequence for
C  this routine looks the similar to that for PSTATE, LISTB must contain
C  11 elements (0..10), not 10, and PVB should have 11 columns (0..10).
C  Note also that the order of the planets in the request and state
C  vectors has changed in this version: the Sun is body zero instead of
C  body 11 as is previous versions.
C
C  Inputs:
C
C   T              Time (ephemeris days wrt reference time).
C
C   TREF           Reference time (JED).  For best accuracy, this should
C                  have a fractional part of ".0" or ".5".
C
C   LISTB(0..10)   List of indicators specifying the interpolations
C                  desired for each of 11 bodies:
C                  LISTB(i)= 0  no interpolation for body i,
C                            1  position only for body i, and
C                            2  position and velocity for body i.
C
C  Output:
C
C   PVB(1..6, 0..10)  State components for the specified bodies.  The
C                     first subscript varies over the components,
C                     which are ordered X, Y, Z, VX, VY, VZ, where V
C                     indicates velocity; the second subscript
C                     indicates the body.  The units are km and km/s.
C
C-----------------------------------------------------------------------
C  PENUT (Planetary ephemeris Earth NUTation) reads and interpolates the
C  planetary ephemeris for Earth nutations, which consist of two angles,
C  a nutation in longitude, and a nutation in obliquity.  The user can
C  request the nutation angles alone, or the angles and their first
C  derivatives, by setting the input REQNUT to 1 or 2, respectively.
C
C  Inputs:
C
C   T              Time (ephemeris days wrt reference time).
C
C   TREF           Reference time (JED).  For best accuracy, this should
C                  have a fractional part of ".0" or ".5".
C
C   REQNUT    Indicator specifying desired interpolations:
C             REQNUT= 0  for no interpolation
C                     1  for nutation angles only, and
C                     2  for angles and first derivatives.
C
C  Output:
C
C   NUT(1..4)   Nutation output array.  The order is NL, NO, DNL, DNO,
C               where NL and NO are the nutations in longitude and
C               obliquity, respectively, and D indicates the time
C               derivative.  The units are rad and rad/s.
C
C-----------------------------------------------------------------------
C  PECON (Planetary Ephemeris CONstants) reads the values of specified
C  constants from the SSDPS groups on the planetary ephemeris file.
C  You specify the constants you are interested in by placing their
C  names in the array NAM; PECON looks up the corresponding values and
C  returns them in the array VAL.  Some of the useful constants are:
C    AU      Length of the Astronomical Unit, km
C    CLIGHT  Speed of light, km/s
C    GMi     Gravitational constant GM for body i, in AU^3/d^2,
C            where i= 1, 2, B, 4, ..., 9, or S, where B indicates the
C            Earth-Moon barycenter, and S indicates the Sun)
C    EMRAT   Earth-Moon ratio (mass of Earth over mass of Moon)
C    PESTRT  Start time of ephemeris, JED
C    PESTOP  Stop time of ephemeris, JED
C    BETA, GAMMA  Relativity parameters (unitless)
C
C  A single constant can be conveniently retrieved by specifying the
C  name directly in the call, specifying a 1 for the second argument
C  and the name of the variable in the third.  For example, the length
C  of the astronomical unit is read into into variable AU as follows:
C            CALL PECON( 'AU', 1, AU)
C  PECON is designed so that it only actually reads the file the first
C  time it is called: the names and values of all constants are saved in
C  memory and used on subsequent calls.  Thus, it is not inefficient to
C  call PECON many times, once for each constant.
C
C  Inputs:
C
C   NAM(1..NV)  Array of CHARACTER*6 names of the constants whose values
C               are desired.
C
C   NV          Number of names in NAM (ie number of values to be read).
C
C  Output:
C
C   VAL(1..NV)  Values of the specified constants.
C
C$======================================================================

>>>>>>>>>>>>>>>>>>> FILE: ../ephrd/sbread.f <<<<<<<<<<<<<<<<<<<

      SUBROUTINE SBOPEN( SBFILE, SBNAM, JDREF, ESTRT, ESTOP, NSB)
      CHARACTER*(*)      SBFILE, SBNAM
      DOUBLE PRECISION   JDREF, ESTRT, ESTOP
      INTEGER            NSB
C
C     SUBROUTINE SBHPV ( T, TREF, RQLIST, PV)
      DOUBLE PRECISION   T, TREF, PV(6,*)
      INTEGER            RQLIST(*)
C
C     SUBROUTINE SBCON ( NAM, ISB, VAL)
      CHARACTER*(*)      NAM
      INTEGER            ISB
      DOUBLE PRECISION   VAL(*)
C
C     SUBROUTINE SBSIZE( NSB)
C
C     SUBROUTINE SBFIND( BNAM, NN, SBLIST)
      CHARACTER*(*)      BNAM(*)
      INTEGER            NN, SBLIST(*)
C
C     SUBROUTINE SBGET(  SBLIST, NN, BNAM)
C
C     SUBROUTINE SBONAM( ISB, NAM)
C
C=======================================================================
C              SBREAD: Small Body Ephemeris Reader Package
C
C  The original (version 1) small-body ephemeris file format was changed
C  to version 2.  This version of the reader can handle either version.
C  In version 2 files, only the "combined" section of the file is used.
C
C-----------------------------------------------------------------------
C  SBOPEN (Small Body ephemeris OPEN) opens a small body ephemeris
C  file and returns the following information about the file:
C
C  SBOPEN must be called before any other routine in this package is
C  called.
C
C  Input:
C
C   SBFILE  File name of the small body ephemeris file.
C
C  Outputs:
C
C   SBNAM   Ephemeris file ID.
C   JDREF   Reference epoch of the equinox, JED.
C   ESTRT   Start epoch for the ephemeris, JED.
C   ESTOP   Stop epoch for the ephemeris, JED.
C   NSB     Number of small bodies in the file.
C
C-----------------------------------------------------------------------
C  SBHPV (Small Body ephemeris Heliocentric Position and Velocity)
C  returns the heliocentric position and velocity of each requested
C  body in the file.
C
C  Inputs:
C
C   T              Time (ephemeris days wrt reference time).
C   TREF           Reference time (JED).
C   RQLIST(1..NSB) List of indicators specifying the interpolations
C                  desired for each small body:
C                  RQLIST(i)= 0  no interpolation for body i,
C                             1  position only for body i, and
C                             2  position and velocity for body i.
C
C  Output:
C
C   PV(1..6, 1..NSB) State components for the specified bodies.  The
C                    first subscript varies over the components,
C                    which are ordered X, Y, Z, VX, VY, VZ, where V
C                    indicates velocity; the second subscript
C                    indicates the body.  The units are km and km/s.
C
C-----------------------------------------------------------------------
C  SBCON (Small Body ephemeris CONstants) reads values of specified
C  constants from the ID, SMALL_BODY, and CONSTANTS groups of the small
C  body ephemeris NAVIO file.
C
C  Example:  SBCON( 'GM', 2, GM(1) )  will return GM for the 2nd small
C            body in the file in GM(1).
C  Example:  SBCON( 'GM', 0, GM(1) )  will return GM for all small
C            bodies (NB) in the file in array GM(1..NB).
C
C  RESTRICTIONS:  'GM' or 'RD' only!
C
C  Inputs:
C
C   NAM     Name of constant whose value is desired.
C   ISB     Index within file of small body whose constant(s) is
C           requested.  If ISB = 0, values are returned for all small
C           bodies on the file.
C
C  Output:
C
C   VAL(1..NSB)  Value(s) of the specified constant(s).
C
C-----------------------------------------------------------------------
C  SBSIZE (Small Body ephemeris SIZE) returns the number of small
C  bodies in the ephemeris file.
C
C  Output:
C
C   NSB     Number of small bodies in the file.
C
C-----------------------------------------------------------------------
C  SBFIND (Small Body ephemeris FIND) returns the request list SBLIST,
C  appropriate for the Small Body ephemeris file, corresponding to a
C  given array of small body names.  If a name is not found, an error
c  message is output and execution is aborted.
C
C  NOTE: The name is compared with the "OBJ-NAME" item of the NAVIO
C        Small Body Ephemeris File.  The "P/, C/, D/" prefix for comet
C        names is optional.
C
C  Input:
C
C   BNAM(1..NN)  Array of names for each requested small body.
C                Examples:  'Vesta', 'Halley', 'P/Halley'
C   NN           Number of names in BNAM.
C
C  Output:
C
C   SBLIST(1..NSB)  Request array for each small body in the file.
C                   If the small body (i) was requested in BNAM, then
C                   SBLIST(i) is 1, otherwise SBLIST(i) is 0.  SBLIST
C                   can be used as RQLIST in subsequent calls to SBHPV,
C                   for example.
C
C-----------------------------------------------------------------------
C  SBGET (Small Body ephemeris GET) returns all object names specified
C  in the request list SBLIST.  The names returned are those from the
C  "OBJ-NAME" item of the NAVIO Small Body ephemeris file.
C
C  Input:
C
C   SBLIST(1..NSB)  Request array for each small body in the file.
C                   If SBLIST(i) > 0, body (i) is requested.
C
C  Output:
C
C   NN           Number of names in BNAM.
C   BNAM(1..NN)  Array of object names for each requested small body.
C
C-----------------------------------------------------------------------
C  SBONAM (Small Body ephemeris Output NAMe) returns the object name of
C  small body given its index ISB.
C
C  Input:
C
C    ISB   Index of requested small body within the file.
C
C  Output:
C
C    NAM   Name of the requested small body (NAVIO OBJ-NAME item).
C          If ISB is invalid, NAM is blank.
C
C$======================================================================
