/*
*+
* Name:
* palPlante
* Purpose:
* Topocentric RA,Dec of a Solar-System object from heliocentric orbital elements
* Language:
* Starlink ANSI C
* Type of Module:
* Library routine
* Invocation:
* void palPlante ( double date, double elong, double phi, int jform,
* double epoch, double orbinc, double anode, double perih,
* double aorq, double e, double aorl, double dm,
* double *ra, double *dec, double *r, int *jstat );
* Description:
* Topocentric apparent RA,Dec of a Solar-System object whose
* heliocentric orbital elements are known.
* Arguments:
* date = double (Given)
* TT MJD of observation (JD-2400000.5)
* elong = double (Given)
* Observer's east longitude (radians)
* phi = double (Given)
* Observer's geodetic latitude (radians)
* jform = int (Given)
* Element set actually returned (1-3; Note 6)
* epoch = double (Given)
* Epoch of elements (TT MJD)
* orbinc = double (Given)
* inclination (radians)
* anode = double (Given)
* longitude of the ascending node (radians)
* perih = double (Given)
* longitude or argument of perihelion (radians)
* aorq = double (Given)
* mean distance or perihelion distance (AU)
* e = double (Given)
* eccentricity
* aorl = double (Given)
* mean anomaly or longitude (radians, JFORM=1,2 only)
* dm = double (Given)
* daily motion (radians, JFORM=1 only)
* ra = double * (Returned)
* Topocentric apparent RA (radians)
* dec = double * (Returned)
* Topocentric apparent Dec (radians)
* r = double * (Returned)
* Distance from observer (AU)
* jstat = int * (Returned)
* status: 0 = OK
* - -1 = illegal jform
* - -2 = illegal e
* - -3 = illegal aorq
* - -4 = illegal dm
* - -5 = numerical error
* Authors:
* PTW: Pat Wallace (STFC)
* TIMJ: Tim Jenness (JAC, Hawaii)
* {enter_new_authors_here}
* Notes:
* - 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).
* - The longitude and latitude allow correction for geocentric
* parallax. This is usually a small effect, but can become
* important for near-Earth asteroids. Geocentric positions can be
* generated by appropriate use of routines palEpv (or palEvp) and
* palUe2pv.
* - The elements are with respect to the J2000 ecliptic and equinox.
* - A choice of three different element-set options is 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 elements and 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)
*
* Unused arguments (DM for JFORM=2, AORL and DM for JFORM=3) are not
* accessed.
* - Each of the three element sets defines an unperturbed heliocentric
* orbit. For a given epoch of observation, the position of the body
* in its orbit can be predicted from these elements, which are
* called "osculating elements", using standard two-body analytical
* solutions. However, due to planetary perturbations, a given set
* of osculating elements remains usable for only as long as the
* unperturbed orbit that it describes is an adequate approximation
* to reality. Attached to such a set of elements is a date called
* the "osculating epoch", at which the elements are, momentarily,
* a perfect representation of the instantaneous position and
* velocity of the body.
*
* Therefore, for any given problem there are up to three different
* epochs in play, and it is vital to distinguish clearly between
* them:
*
* . The epoch of observation: the moment in time for which the
* position of the body is to be predicted.
*
* . The epoch defining the position of the body: the moment in time
* at which, in the absence of purturbations, the specified
* position (mean longitude, mean anomaly, or perihelion) is
* reached.
*
* . The osculating epoch: the moment in time at which the given
* elements are correct.
*
* For the major-planet and minor-planet cases it is usual to make
* the epoch that defines the position of the body the same as the
* epoch of osculation. Thus, only two different epochs are
* involved: the epoch of the elements and the epoch of observation.
*
* For comets, the epoch of perihelion fixes the position in the
* orbit and in general a different epoch of osculation will be
* chosen. Thus, all three types of epoch are involved.
*
* For the present routine:
*
* . The epoch of observation is the argument DATE.
*
* . The epoch defining the position of the body is the argument
* EPOCH.
*
* . The osculating epoch is not used and is assumed to be close
* enough to the epoch of observation to deliver adequate accuracy.
* If not, a preliminary call to sla_PERTEL may be used to update
* the element-set (and its associated osculating epoch) by
* applying planetary perturbations.
* - Two important sources for orbital elements are Horizons, operated
* by the Jet Propulsion Laboratory, Pasadena, and the Minor Planet
* Center, operated by the Center for Astrophysics, Harvard.
*
* The JPL Horizons elements (heliocentric, J2000 ecliptic and
* equinox) correspond to SLALIB arguments as follows.
*
* Major planets:
*
* JFORM = 1
* EPOCH = JDCT-2400000.5
* ORBINC = IN (in radians)
* ANODE = OM (in radians)
* PERIH = OM+W (in radians)
* AORQ = A
* E = EC
* AORL = MA+OM+W (in radians)
* DM = N (in radians)
*
* Epoch of osculation = JDCT-2400000.5
*
* Minor planets:
*
* JFORM = 2
* EPOCH = JDCT-2400000.5
* ORBINC = IN (in radians)
* ANODE = OM (in radians)
* PERIH = W (in radians)
* AORQ = A
* E = EC
* AORL = MA (in radians)
*
* Epoch of osculation = JDCT-2400000.5
*
* Comets:
*
* JFORM = 3
* EPOCH = Tp-2400000.5
* ORBINC = IN (in radians)
* ANODE = OM (in radians)
* PERIH = W (in radians)
* AORQ = QR
* E = EC
*
* Epoch of osculation = JDCT-2400000.5
*
* The MPC elements correspond to SLALIB arguments as follows.
*
* Minor planets:
*
* JFORM = 2
* EPOCH = Epoch-2400000.5
* ORBINC = Incl. (in radians)
* ANODE = Node (in radians)
* PERIH = Perih. (in radians)
* AORQ = a
* E = e
* AORL = M (in radians)
*
* Epoch of osculation = Epoch-2400000.5
*
* Comets:
*
* JFORM = 3
* EPOCH = T-2400000.5
* ORBINC = Incl. (in radians)
* ANODE = Node. (in radians)
* PERIH = Perih. (in radians)
* AORQ = q
* E = e
*
* Epoch of osculation = Epoch-2400000.5
* History:
* 2012-03-12 (TIMJ):
* Initial version direct conversion of SLA/F.
* Adapted with permission from the Fortran SLALIB library.
* {enter_further_changes_here}
* Copyright:
* Copyright (C) 2004 Patrick T. Wallace
* Copyright (C) 2012 Science and Technology Facilities Council.
* All Rights Reserved.
* Licence:
* This program is free software; you can redistribute it and/or
* modify it under the terms of the GNU General Public License as
* published by the Free Software Foundation; either version 3 of
* the License, or (at your option) any later version.
*
* This program is distributed in the hope that it will be
* useful, but WITHOUT ANY WARRANTY; without even the implied
* warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
* PURPOSE. See the GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
* MA 02110-1301, USA.
* Bugs:
* {note_any_bugs_here}
*-
*/
#include "pal.h"
void palPlante ( double date, double elong, double phi, int jform,
double epoch, double orbinc, double anode, double perih,
double aorq, double e, double aorl, double dm,
double *ra, double *dec, double *r, int *jstat ) {
double u[13];
/* Transform conventional elements to universal elements */
palEl2ue( date, jform, epoch, orbinc, anode, perih, aorq, e, aorl,
dm, u, jstat );
/* If succcessful, make the prediction */
if (*jstat == 0) palPlantu( date, elong, phi, u, ra, dec, r, jstat );
}