diff options
Diffstat (limited to 'src/slalib/planel.f')
| -rw-r--r-- | src/slalib/planel.f | 117 |
1 files changed, 117 insertions, 0 deletions
diff --git a/src/slalib/planel.f b/src/slalib/planel.f new file mode 100644 index 0000000..064e51c --- /dev/null +++ b/src/slalib/planel.f @@ -0,0 +1,117 @@ + SUBROUTINE sla_PLANEL (DATE, JFORM, EPOCH, ORBINC, ANODE, PERIH, + : AORQ, E, AORL, DM, PV, JSTAT) +*+ +* - - - - - - - +* P L A N E L +* - - - - - - - +* +* Heliocentric position and velocity of a planet, asteroid or comet, +* starting from orbital elements. +* +* Given: +* DATE d date, Modified Julian Date (JD - 2400000.5) +* JFORM i choice of element set (1-3; Note 3) +* EPOCH d epoch of elements (TT MJD) +* ORBINC d inclination (radians) +* ANODE d longitude of the ascending node (radians) +* PERIH d longitude or argument of perihelion (radians) +* AORQ d mean distance or perihelion distance (AU) +* E d eccentricity +* AORL d mean anomaly or longitude (radians, JFORM=1,2 only) +* DM d daily motion (radians, JFORM=1 only) +* +* Returned: +* PV d(6) heliocentric x,y,z,xdot,ydot,zdot of date, +* J2000 equatorial triad (AU,AU/s) +* JSTAT i status: 0 = OK +* -1 = illegal JFORM +* -2 = illegal E +* -3 = illegal AORQ +* -4 = illegal DM +* -5 = numerical error +* +* Called: sla_EL2UE, sla_UE2PV +* +* Notes +* +* 1 DATE is the instant for which the prediction is required. It is +* in the TT timescale (formerly Ephemeris Time, ET) and is a +* Modified Julian Date (JD-2400000.5). +* +* 2 The elements are with respect to the J2000 ecliptic and equinox. +* +* 3 Three different element-format options are available: +* +* Option JFORM=1, suitable for the major planets: +* +* EPOCH = epoch of elements (TT MJD) +* ORBINC = inclination i (radians) +* ANODE = longitude of the ascending node, big omega (radians) +* PERIH = longitude of perihelion, curly pi (radians) +* AORQ = mean distance, a (AU) +* E = eccentricity, e (range 0 to <1) +* AORL = mean longitude L (radians) +* DM = daily motion (radians) +* +* Option JFORM=2, suitable for minor planets: +* +* EPOCH = epoch of elements (TT MJD) +* ORBINC = inclination i (radians) +* ANODE = longitude of the ascending node, big omega (radians) +* PERIH = argument of perihelion, little omega (radians) +* AORQ = mean distance, a (AU) +* E = eccentricity, e (range 0 to <1) +* AORL = mean anomaly M (radians) +* +* Option JFORM=3, suitable for comets: +* +* EPOCH = epoch of perihelion (TT MJD) +* ORBINC = inclination i (radians) +* ANODE = longitude of the ascending node, big omega (radians) +* PERIH = argument of perihelion, little omega (radians) +* AORQ = perihelion distance, q (AU) +* E = eccentricity, e (range 0 to 10) +* +* 4 Unused elements (DM for JFORM=2, AORL and DM for JFORM=3) are +* not accessed. +* +* 5 The reference frame for the result is with respect to the mean +* equator and equinox of epoch J2000. +* +* 6 The algorithm was originally adapted from the EPHSLA program of +* D.H.P.Jones (private communication, 1996). The method is based +* on Stumpff's Universal Variables. +* +* Reference: Everhart, E. & Pitkin, E.T., Am.J.Phys. 51, 712, 1983. +* +* P.T.Wallace Starlink 18 March 1999 +* +* Copyright (C) 1999 Rutherford Appleton Laboratory +*- + + IMPLICIT NONE + + DOUBLE PRECISION DATE + INTEGER JFORM + DOUBLE PRECISION EPOCH,ORBINC,ANODE,PERIH,AORQ,E,AORL,DM,PV(6) + INTEGER JSTAT + + DOUBLE PRECISION U(13) + INTEGER J + + + +* Validate elements and convert to "universal variables" parameters. + CALL sla_EL2UE(DATE,JFORM, + : EPOCH,ORBINC,ANODE,PERIH,AORQ,E,AORL,DM,U,J) + +* Determine the position and velocity. + IF (J.EQ.0) THEN + CALL sla_UE2PV(DATE,U,PV,J) + IF (J.NE.0) J=-5 + END IF + +* Wrap up. + JSTAT = J + + END |