Functions
geocentric.c File Reference

Reduce heliocentric position to geocentric right ascension and declination. More...

#include <math.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include "aephem.h"

Functions

void ae_geocentric (double jd_tt, double q[], double e[], double v_e[], double *ra, double *dec, double *dist)
 Given heliocentric coordinates, reduce to geocentric coordinates. More...
 
void ae_geocentric_from_orbit (double jd_tt, const struct ae_orbit_t *o_orb, const struct ae_orbit_t *q_orb, double *ra, double *dec, double *dist)
 Determine geocentric coordinates using orbital elements. More...
 
int ae_geocentric_from_jpl (struct ae_jpl_handle_t *jh, double jd_tt, int obj_num, double *ra, double *dec, double *dist)
 Determine geocentric coordinates using JPL ephemerides. More...
 
void ae_geocentric_from_cat (double jd_tt, double e[], double v_e[], const struct ae_star_t *star, double *ra, double *dec)
 Reduce a star's catalogue coordinates to its apparent place. More...
 
void ae_geocentric_sun_from_orbit (double jd_tt, const struct ae_orbit_t *o_orb, double *ra, double *dec, double *dist)
 Calculate the geocentric coordinates of the sun from orbital elements. More...
 
void ae_geocentric_moon_from_orbit (double jd_tt, const struct ae_orbit_t *o_orb, double *ra, double *dec, double *dist)
 Calculate the geocentric coordinates of the moon from orbital elements. More...
 

Detailed Description

Reduce heliocentric position to geocentric right ascension and declination.

Function Documentation

void ae_geocentric ( double  jd_tt,
double  q[],
double  e[],
double  v_e[],
double *  ra,
double *  dec,
double *  dist 
)

Given heliocentric coordinates, reduce to geocentric coordinates.

This routine reduces heliocentric coordinates to equitorial coordinates. Unlike ae_geocentric_from_orbit() or ae_geocentric_from_jpl(), it does not automatically compute the heliocentric coordinates. This means that light travel-time is not automatically accounted for! The input position of the planet is what is used.

Corrections done during the reduction are: light deflection, annual aberration, precession and nutation.

This routine should be thought of as a macro. One can call the functions that perform the corrections listed above individually. This routine simply wraps them into one call.

Parameters
jd_ttThe Julian date in TT.
qThe heliocentric rectangular coordinates of the object being observed, in AU.
eThe heliocentric rectangular coordinates of the earth, in AU.
v_eThe velocity of the earth, in AU per day.
raFor returning the right ascension, in degrees.
decFor returning the declination, in degrees.
distFor returning the geocentric distance between earth and the object, in AU. Set to NULL if you do not require this value.

References AE_ADD_ABERRATION, ae_annual_aberration(), AE_FROM_J2000, ae_precess(), ae_rect_to_polar(), ae_relativity(), and aes_nutate().

Referenced by ae_geocentric_from_jpl(), ae_geocentric_from_orbit(), ae_geocentric_sun_from_orbit(), and aes_subobs_point().

void ae_geocentric_from_cat ( double  jd_tt,
double  e[],
double  v_e[],
const struct ae_star_t star,
double *  ra,
double *  dec 
)

Reduce a star's catalogue coordinates to its apparent place.

This routine puts the coordinates in J2000.0 (if not already there), and performs precession, proper motion, parallax, light deflection, aberration and nutation.

This routine should be thought of as a macro. One can call the functions that perform the corrections listed above individually. This routine simply wraps them into one call.

Parameters
jd_ttThe Julian date in TT.
eThe heliocentric rectangular coordinates of the earth, in AU. These can be obtained from ae_kepler() or ae_jpl_get_coords().
v_eThe velocity of the earth, in AU per day. This can be obtained from ae_v_orbit() or ae_jpl_get_coords().
starThe catalogue object. If it is in the FK4 B1950.0 coordinates, this routine will convert it to FK5 J2000.0 coordinates.
raFor returning the apparent (geocentric) right ascension, in degrees.
decFor returning the apparent (geocentric) declination, in degrees.

References AE_ADD_ABERRATION, ae_annual_aberration(), AE_B1950, AE_DTR, ae_fk4_to_fk5(), AE_FROM_J2000, ae_precess(), ae_rect_to_polar(), ae_relativity(), AE_STR, AE_TO_J2000, aes_nutate(), ae_star_t::dec, ae_star_t::epoch, ae_star_t::mudec, ae_star_t::mura, ae_star_t::px, ae_star_t::ra, and ae_star_t::v.

int ae_geocentric_from_jpl ( struct ae_jpl_handle_t jh,
double  jd_tt,
int  obj_num,
double *  ra,
double *  dec,
double *  dist 
)

Determine geocentric coordinates using JPL ephemerides.

This routine reduces the heliocentric equitorial coordinates, obtained by calling ae_jpl_get_coords(), to the geocentric right ascension and declination. The position is corrected for light travel time, light deflection, annual aberration, precession and nutation.

It is assumed that the JPL planetary ephemeris file is being used, since it automatically polls the file for Earth's position and velocity. If you are using a different ephemeris file, you cannot use this routine.

This routine should be thought of as a macro. One can call ae_jpl_get_coords() and the functions that perform the corrections listed above individually. This routine simply wraps them into one call.

To obtain geocentric coordinates from an orbital elements struct, use ae_geocentric_from_orbit().

Parameters
jhA handle initialised by ae_jpl_init().
jd_ttThe Julian date in TT.
obj_numThe object for which to get coordinates. If a planetary ephemeris file is being used, one can use the predefined ae_ss_bodies_t constants for clarity.
raFor returning the geocentric right ascension, in degrees.
decFor returning the geocentric declination, in degrees.
distFor returning the geocentric distance between earth and the object, in AU. Set to NULL if you do not require this value.
Returns
0 on success; on failure, a negative error code from ae_retcode_t.

References ae_geocentric(), ae_jpl_get_coords(), ae_light_t_jpl(), and AE_SS_EARTH.

void ae_geocentric_from_orbit ( double  jd_tt,
const struct ae_orbit_t o_orb,
const struct ae_orbit_t q_orb,
double *  ra,
double *  dec,
double *  dist 
)

Determine geocentric coordinates using orbital elements.

This routine reduces the heliocentric equitorial coordinates, obtained by calling ae_kepler(), to the geocentric right ascension and declination. The position is corrected for light travel time, light deflection, annual aberration, precession and nutation.

This routine should be thought of as a macro. One can call ae_kepler() and the functions that perform the corrections listed above individually. This routine simply wraps them into one call.

To obtain geocentric coordinates from a JPL ephemeris file, use ae_geocentric_from_jpl().

Parameters
jd_ttThe Julian date in TT.
o_orbOrbital elements of the observer.
q_orbThe orbital elements of the object.
raFor returning the geocentric right ascension, in degrees.
decFor returning the geocentric declination, in degrees.
distFor returning the geocentric distance between earth and the object, in AU. Set to NULL if you do not require this value.

References ae_geocentric(), ae_kepler(), ae_light_t_orbit(), and ae_v_orbit().

void ae_geocentric_moon_from_orbit ( double  jd_tt,
const struct ae_orbit_t o_orb,
double *  ra,
double *  dec,
double *  dist 
)

Calculate the geocentric coordinates of the moon from orbital elements.

Corrects for precession, nutation, annual aberration. Neglects light travel time.

Note that if using JPL ephemeris, one can use ae_geocentric_from_jpl() to do this calculation since the moon is an object in those ephemerides.

The Astronomical Almanac (Section D, Daily Polynomial Coefficients) seems to omit the annual aberration at one point, even though the reference ephemeris is inertial . . .

Parameters
jd_ttThe Julian date in TT.
o_orbOrbital elements of the observer.
raFor returning the geocentric right ascension, in degrees.
decFor returning the geocentric declination, in degrees.
distFor returning the geocentric distance between earth and the object, in AU. Set to NULL if you do not require this value.

References AE_DTR, ae_epsilon(), AE_FROM_J2000, ae_gmoon(), ae_kepler(), ae_nutate(), ae_nutation_lon_ob(), AE_R_EARTH_AU, ae_rect_to_polar(), and AE_STR.

void ae_geocentric_sun_from_orbit ( double  jd_tt,
const struct ae_orbit_t o_orb,
double *  ra,
double *  dec,
double *  dist 
)

Calculate the geocentric coordinates of the sun from orbital elements.

Corrects for precession, nutation, annual aberration.

Note that if using JPL ephemeris, one can use ae_geocentric_from_jpl() to do this calculation since the sun is an object in those ephemerides.

Parameters
jd_ttThe Julian date in TT.
o_orbOrbital elements of the observer.
raFor returning the geocentric right ascension, in degrees.
decFor returning the geocentric declination, in degrees.
distFor returning the geocentric distance between earth and the object, in AU. Set to NULL if you do not require this value.

References ae_geocentric(), ae_kepler(), and ae_v_orbit().