diff options
Diffstat (limited to 'src/slalib/plante.f')
| -rw-r--r-- | src/slalib/plante.f | 157 |
1 files changed, 157 insertions, 0 deletions
diff --git a/src/slalib/plante.f b/src/slalib/plante.f new file mode 100644 index 0000000..e8f89d9 --- /dev/null +++ b/src/slalib/plante.f @@ -0,0 +1,157 @@ + SUBROUTINE sla_PLANTE (DATE, ELONG, PHI, JFORM, EPOCH, + : ORBINC, ANODE, PERIH, AORQ, E, + : AORL, DM, RA, DEC, R, JSTAT) +*+ +* - - - - - - - +* P L A N T E +* - - - - - - - +* +* Topocentric apparent RA,Dec of a Solar-System object whose +* heliocentric orbital elements are known. +* +* Given: +* DATE d MJD of observation (JD - 2400000.5) +* ELONG d observer's east longitude (radians) +* PHI d observer's geodetic latitude (radians) +* JFORM i choice of element set (1-3; Note 4) +* 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: +* RA,DEC d RA, Dec (topocentric apparent, radians) +* R d distance from observer (AU) +* JSTAT i status: 0 = OK +* -1 = illegal JFORM +* -2 = illegal E +* -3 = illegal AORQ +* -4 = illegal DM +* -5 = numerical error +* +* 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 longitude and latitude allow correction for geocentric +* parallax. This is usually a small effect, but can become +* important for Earth-crossing asteroids. Geocentric positions +* can be generated by appropriate use of routines sla_EVP and +* sla_PLANEL. +* +* 3 The elements are with respect to the J2000 ecliptic and equinox. +* +* 4 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 +* 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 +* 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 +* +* 5 Unused elements (DM for JFORM=2, AORL and DM for JFORM=3) are +* not accessed. +* +* Called: sla_GMST, sla_DT, sla_EPJ, sla_PVOBS, sla_PRENUT, +* sla_PLANEL, sla_DMXV, sla_DCC2S, sla_DRANRM +* +* P.T.Wallace Starlink 17 March 1999 +* +* Copyright (C) 1999 Rutherford Appleton Laboratory +*- + + IMPLICIT NONE + + DOUBLE PRECISION DATE,ELONG,PHI + INTEGER JFORM + DOUBLE PRECISION EPOCH,ORBINC,ANODE,PERIH,AORQ,E, + : AORL,DM,RA,DEC,R + INTEGER JSTAT + +* Light time for unit distance (sec) + DOUBLE PRECISION TAU + PARAMETER (TAU=499.004782D0) + + INTEGER I + DOUBLE PRECISION DVB(3),DPB(3),VSG(6),VSP(6),V(6),RMAT(3,3), + : VGP(6),STL,VGO(6),DX,DY,DZ,D,TL + DOUBLE PRECISION sla_GMST,sla_DT,sla_EPJ,sla_DRANRM + + + +* Sun to geocentre (J2000). + CALL sla_EVP(DATE,2000D0,DVB,DPB,VSG(4),VSG) + +* Sun to planet (J2000). + CALL sla_PLANEL(DATE,JFORM,EPOCH,ORBINC,ANODE,PERIH,AORQ, + : E,AORL,DM,VSP,JSTAT) + +* Geocentre to planet (J2000). + DO I=1,6 + V(I)=VSP(I)-VSG(I) + END DO + +* Precession and nutation to date. + CALL sla_PRENUT(2000D0,DATE,RMAT) + CALL sla_DMXV(RMAT,V,VGP) + CALL sla_DMXV(RMAT,V(4),VGP(4)) + +* Geocentre to observer (date). + STL=sla_GMST(DATE-sla_DT(sla_EPJ(DATE))/86400D0)+ELONG + CALL sla_PVOBS(PHI,0D0,STL,VGO) + +* Observer to planet (date). + DO I=1,6 + V(I)=VGP(I)-VGO(I) + END DO + +* Geometric distance (AU). + DX=V(1) + DY=V(2) + DZ=V(3) + D=SQRT(DX*DX+DY*DY+DZ*DZ) + +* Light time (sec). + TL=TAU*D + +* Correct position for planetary aberration + DO I=1,3 + V(I)=V(I)-TL*V(I+3) + END DO + +* To RA,Dec. + CALL sla_DCC2S(V,RA,DEC) + RA=sla_DRANRM(RA) + R=D + + END |