SuperNOVAS v1.5
The NOVAS C library, made better
Loading...
Searching...
No Matches
Non-equatorial coordinates

Enumerations

enum  novas_nutation_direction { NUTATE_TRUE_TO_MEAN = -1 , NUTATE_MEAN_TO_TRUE }
 Direction constant for nutation(), between mean and true equatorial coordinates. More...
 
enum  novas_wobble_direction { WOBBLE_ITRS_TO_TIRS = 0 , WOBBLE_TIRS_TO_ITRS , WOBBLE_ITRS_TO_PEF , WOBBLE_PEF_TO_ITRS }
 Direction constants for polar wobble corrections via the wobble() function. More...
 

Functions

int ecl2equ (double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, double elon, double elat, double *restrict ra, double *restrict dec)
 Convert ecliptic longitude and latitude to right ascension and declination.
 
short ecl2equ_vec (double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, const double *in, double *out)
 Converts an ecliptic position vector to an equatorial position vector.
 
short equ2ecl (double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, double ra, double dec, double *restrict elon, double *restrict elat)
 Convert right ascension and declination to ecliptic longitude and latitude.
 
short equ2ecl_vec (double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, const double *in, double *out)
 Converts an equatorial position vector to an ecliptic position vector.
 
int equ2gal (double ra, double dec, double *restrict glon, double *restrict glat)
 Converts ICRS right ascension and declination to galactic longitude and latitude.
 
int gal2equ (double glon, double glat, double *restrict ra, double *restrict dec)
 Converts galactic longitude and latitude to ICRS right ascension and declination.
 
int hor_to_itrs (const on_surface *restrict location, double az, double za, double *restrict itrs)
 Converts astrometric (unrefracted) azimuth and zenith angles at the specified observer location to a unit position vector in the Earth-fixed ITRS frame.
 
int itrs_to_hor (const on_surface *restrict location, const double *restrict itrs, double *restrict az, double *restrict za)
 Converts a position vector in the Earth-fixed ITRS frame to astrometric (unrefracted) azimuth and zenith angles at the specified observer location.
 
int novas_e2h_offset (double dra, double ddec, double pa, double *restrict daz, double *restrict del)
 Converts coordinate offsets, from the local equatorial system to local horizontal offsets.
 
double novas_epa (double ha, double dec, double lat)
 Returns the Parallactic Angle (PA) calculated for an RA/Dec location of the sky at a given sidereal time.
 
int novas_h2e_offset (double daz, double del, double pa, double *restrict dra, double *restrict ddec)
 Converts coordinate offsets, from the local horizontal system to local equatorial offsets.
 
int novas_hor_to_app (const novas_frame *restrict frame, double az, double el, RefractionModel ref_model, enum novas_reference_system sys, double *restrict ra, double *restrict dec)
 Converts an observed azimuth and elevation coordinate to right ascension (R.A.) and declination coordinates expressed in the coordinate system of choice.
 
double novas_hpa (double az, double el, double lat)
 Returns the Parallactic Angle (PA) calculated for a horizontal Az/El location of the sky.
 
int novas_los_to_xyz (const double *los, double lon, double lat, double *xyz)
 Converts a 3D line-of-sight vector (δφ, δθ δr) to a rectangular equatorial (δx, δy, δz) vector.
 
int novas_uvw_to_xyz (const double *uvw, double ha, double dec, double *xyz)
 Converts equatorial u,v,w projected (absolute or relative) coordinates to rectangular telescope x,y,z coordinates (in ITRS) to for a specified line of sight.
 
int novas_xyz_to_los (const double *xyz, double lon, double lat, double *los)
 Converts a 3D rectangular equatorial (δx, δy, δz) vector to a polar (δφ, δθ δr) vector along a line-of-sight.
 
int novas_xyz_to_uvw (const double *xyz, double ha, double dec, double *uvw)
 Converts rectangular telescope x,y,z (absolute or relative) coordinates (in ITRS) to equatorial u,v,w projected coordinates for a specified line of sight.
 
int wobble (double jd_tt, enum novas_wobble_direction direction, double xp, double yp, const double *in, double *out)
 Corrects a vector in the ITRS (rotating Earth-fixed system) for polar motion, and also corrects the longitude origin (by a tiny amount) to the Terrestrial Intermediate Origin (TIO).
 

Detailed Description

Enumeration Type Documentation

◆ novas_nutation_direction

Direction constant for nutation(), between mean and true equatorial coordinates.

See also
nutation()
Enumerator
NUTATE_TRUE_TO_MEAN 

Change from true equator to mean equator (i.e.

undo wobble corrections). You may use any non-zero value as well.

NUTATE_MEAN_TO_TRUE 

Change from mean equator to true equator (i.e. apply wobble corrections)

◆ novas_wobble_direction

Direction constants for polar wobble corrections via the wobble() function.

See also
novas_app_to_hor(), novas_hor_to_app(), wobble(), NOVAS_WOBBLE_DIRECTIONS
Enumerator
WOBBLE_ITRS_TO_TIRS 

use for wobble() to change from ITRS (Earth-fixed) to TIRS (pseudo Earth-fixed).

It includes TIO longitude correction.

Since
1.4
WOBBLE_TIRS_TO_ITRS 

use for wobble() to change from TIRS (pseudo Earth-fixed) to ITRS (Earth-fixed).

It includes TIO longitude correction.

Since
1.4
WOBBLE_ITRS_TO_PEF 

use for wobble() to change from ITRS (Earth-fixed) Pseudo Earth Fixed (PEF).

It does not include TIO longitude correction. Otherwise, it's the same as WOBBLE_ITRS_TO_TIRS

WOBBLE_PEF_TO_ITRS 

use for wobble() to change from Pseudo Earth Fixed (PEF) to ITRS (Earth-fixed).

It does not include TIO longitude correction. Otherwise, it's the same as WOBBLE_TIRS_TO_ITRS

Function Documentation

◆ ecl2equ()

int ecl2equ ( double jd_tt,
enum novas_equator_type coord_sys,
enum novas_accuracy accuracy,
double elon,
double elat,
double *restrict ra,
double *restrict dec )

Convert ecliptic longitude and latitude to right ascension and declination.

To convert GCRS ecliptic coordinates (mean ecliptic and equinox of J2000.0), set 'coord_sys' to NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set to anything, since J2000.0 is assumed. Otherwise, all input coordinates are dynamical at'jd_tt'.

Parameters
jd_tt[day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2])
coord_sysThe astrometric reference system of the coordinates. If 'coord_sys' is NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to J2000 ecliptic coordinates.
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
elon[deg] Ecliptic longitude in degrees, referred to specified ecliptic and equinox of date.
elat[deg] Ecliptic latitude in degrees, referred to specified ecliptic and equinox of date.
[out]ra[h] Right ascension in hours, referred to specified equator and equinox of date.
[out]dec[deg] Declination in degrees, referred to specified equator and equinox of date.
Returns
0 if successful, or else 1 if the value of 'coord_sys' is invalid.
Since
1.0
Author
Attila Kovacs
See also
ecl2equ_vec(), equ2ecl()

References ecl2equ_vec().

◆ ecl2equ_vec()

short ecl2equ_vec ( double jd_tt,
enum novas_equator_type coord_sys,
enum novas_accuracy accuracy,
const double * in,
double * out )

Converts an ecliptic position vector to an equatorial position vector.

To convert ecliptic coordinates (mean ecliptic and equinox of J2000.0) to GCRS RA and dec to, set 'coord_sys' to NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set to anything, since J2000.0 is assumed. Otherwise, all input coordinates are dynamical at 'jd_tt'.

Parameters
jd_tt[day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2])
coord_sysThe astrometric reference system type of the coordinates
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
inPosition vector, referred to specified ecliptic and equinox of date.
[out]outPosition vector, referred to specified equator and equinox of date. It can be the same vector as the input.
Returns
0 if successful, -1 if either vector argument is NULL or the accuracy is invalid, or else 1 if the value of 'coord_sys' is invalid.
See also
ecl2equ(), equ2ecl_vec()

References e_tilt(), frame_tie(), J2000_TO_ICRS, mean_obliq(), NOVAS_FULL_ACCURACY, NOVAS_GCRS_EQUATOR, NOVAS_MEAN_EQUATOR, NOVAS_REDUCED_ACCURACY, and NOVAS_TRUE_EQUATOR.

◆ equ2ecl()

short equ2ecl ( double jd_tt,
enum novas_equator_type coord_sys,
enum novas_accuracy accuracy,
double ra,
double dec,
double *restrict elon,
double *restrict elat )

Convert right ascension and declination to ecliptic longitude and latitude.

To convert GCRS RA and dec to ecliptic coordinates (mean ecliptic and equinox of J2000.0), set 'coord_sys' to NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set to anything, since J2000.0 is assumed. Otherwise, all input coordinates are dynamical at 'jd_tt'.

Parameters
jd_tt[day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2])
coord_sysThe astrometric reference system of the coordinates. If 'coord_sys' is NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to J2000 ecliptic coordinates.
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
ra[h] Right ascension in hours, referred to specified equator and equinox of date.
dec[deg] Declination in degrees, referred to specified equator and equinox of date.
[out]elon[deg] Ecliptic longitude in degrees, referred to specified ecliptic and equinox of date.
[out]elat[deg] Ecliptic latitude in degrees, referred to specified ecliptic and equinox of date.
Returns
0 if successful, or else 1 if the value of 'coord_sys' is invalid.
See also
equ2ecl_vec(), ecl2equ()

References equ2ecl_vec().

◆ equ2ecl_vec()

short equ2ecl_vec ( double jd_tt,
enum novas_equator_type coord_sys,
enum novas_accuracy accuracy,
const double * in,
double * out )

Converts an equatorial position vector to an ecliptic position vector.

To convert ICRS RA and dec to ecliptic coordinates (mean ecliptic and equinox of J2000.0), set 'coord_sys' to NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set to anything, since J2000.0 is assumed. Otherwise, all input coordinates are dynamical at 'jd_tt'.

Parameters
jd_tt[day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2])
coord_sysThe astrometric reference system type of the coordinates.
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
inPosition vector, referred to specified equator and equinox of date.
[out]outPosition vector, referred to specified ecliptic and equinox of date. It can be the same vector as the input. If 'coord_sys' is NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to J2000 ecliptic coordinates.
Returns
0 if successful, -1 if either vector argument is NULL or the accuracy is invalid, or else 1 if the value of 'coord_sys' is invalid.
See also
equ2ecl(), ecl2equ_vec()

References e_tilt(), frame_tie(), ICRS_TO_J2000, mean_obliq(), NOVAS_FULL_ACCURACY, NOVAS_GCRS_EQUATOR, NOVAS_MEAN_EQUATOR, NOVAS_REDUCED_ACCURACY, and NOVAS_TRUE_EQUATOR.

◆ equ2gal()

int equ2gal ( double ra,
double dec,
double *restrict glon,
double *restrict glat )

Converts ICRS right ascension and declination to galactic longitude and latitude.

REFERENCES:

  1. Hipparcos and Tycho Catalogues, Vol. 1, Section 1.5.3.
Parameters
ra[h] ICRS right ascension in hours.
dec[deg] ICRS declination in degrees.
[out]glon[deg] Galactic longitude in degrees.
[out]glat[deg] Galactic latitude in degrees.
Returns
0 if successful, or -1 if either of the output pointer arguments are NULL.
See also
gal2equ()

◆ gal2equ()

int gal2equ ( double glon,
double glat,
double *restrict ra,
double *restrict dec )

Converts galactic longitude and latitude to ICRS right ascension and declination.

REFERENCES:

  1. Hipparcos and Tycho Catalogues, Vol. 1, Section 1.5.3.
Parameters
glon[deg] Galactic longitude in degrees.
glat[deg] Galactic latitude in degrees.
[out]ra[h] ICRS right ascension in hours.
[out]dec[deg] ICRS declination in degrees.
Returns
0 if successful, or -1 if either of the output pointer arguments are NULL.
Since
1.0
Author
Attila Kovacs
See also
equ2gal()

◆ hor_to_itrs()

int hor_to_itrs ( const on_surface *restrict location,
double az,
double za,
double *restrict itrs )

Converts astrometric (unrefracted) azimuth and zenith angles at the specified observer location to a unit position vector in the Earth-fixed ITRS frame.

Parameters
locationGeodetic (ITRF / GRS80) observer location on Earth
az[deg] astrometric (unrefracted) azimuth angle at observer location [0:360]. It may be NULL if not required.
za[deg] astrometric (unrefracted) zenith angle at observer location [0:180]. It may be NULL if not required.
[out]itrsUnit 3-vector direction in Earth-fixed ITRS frame
Returns
0 if successful, or else -1 if the location or the input vector is NULL.
Since
1.0
Author
Attila Kovacs
See also
itrs_to_hor(), itrs_to_cirs(), itrs_to_tod(), refract()

◆ itrs_to_hor()

int itrs_to_hor ( const on_surface *restrict location,
const double *restrict itrs,
double *restrict az,
double *restrict za )

Converts a position vector in the Earth-fixed ITRS frame to astrometric (unrefracted) azimuth and zenith angles at the specified observer location.

Parameters
locationGeodetic (ITRF / GRS80) observer location on Earth.
itrs3-vector position in Earth-fixed ITRS frame
[out]az[deg] astrometric (unrefracted) azimuth angle at observer location [0:360]. It may be NULL if not required.
[out]za[deg] astrometric (unrefracted) zenith angle at observer location [0:180]. It may be NULL if not required.
Returns
0 if successful, or else -1 if the location or the input vector is NULL.
Since
1.0
Author
Attila Kovacs
See also
hor_to_itrs(), cirs_to_itrs(), tod_to_itrs(), refract_astro()

◆ novas_e2h_offset()

int novas_e2h_offset ( double dra,
double ddec,
double pa,
double *restrict daz,
double *restrict del )

Converts coordinate offsets, from the local equatorial system to local horizontal offsets.

Converting between local flat projections and spherical coordinates usually requires a WCS projection.

REFERENCES:

  1. Calabretta, M.R., & Greisen, E.W., (2002), Astronomy & Astrophysics, 395, 1077-1122.
Parameters
dra[arcsec] Projected ffset position in the apparent true-of-date R.A. direction. E.g. The projected offset between two RA coordinates at a same reference declination, is δRA = (RA2 - RA1) * cos(Dec0).
ddec[arcsec] Projected offset position in the apparent true-of-date declination direction.
pa[deg] Parallactic Angle
[out]daz[arcsec] Output offset position in the local azimuth direction. It can be a pointer to one of the input coordinates, or NULL if not required.
[out]del[arcsec] Output offset position in the local elevation direction. It can be a pointer to one of the input coordinates, or NULL if not required.
Returns
0
Since
1.3
Author
Attila Kovacs
See also
novas_h2e_offset(), novas_epa()

References novas_h2e_offset().

◆ novas_epa()

double novas_epa ( double ha,
double dec,
double lat )

Returns the Parallactic Angle (PA) calculated for an RA/Dec location of the sky at a given sidereal time.

The PA is the angle between the local horizontal coordinate directions and the local true-of-date equatorial coordinate directions, at the given location and time. The polar wobble is not included in the calculation.

The Parallactic Angle is sometimes referrred to as the Vertical Position Angle (VPA). Both define the same quantity.

Parameters
ha[h] Hour angle (LST - RA) i.e., the difference between the Local (apparent) Sidereal Time and the apparent (true-of-date) Right Ascension of observed source.
dec[deg] Apparent (true-of-date) declination of observed source
lat[deg] Geodetic latitude of observer
Returns
[deg] Parallactic Angle (PA). I.e., the clockwise position angle of the elevation direction w.r.t. the declination axis in the equatorial system. Same as the clockwise position angle of the declination direction w.r.t. the elevation axis, in the horizontal system.
Since
1.3
Author
Attila Kovacs
See also
novas_hpa(), novas_e2h_offset()

◆ novas_h2e_offset()

int novas_h2e_offset ( double daz,
double del,
double pa,
double *restrict dra,
double *restrict ddec )

Converts coordinate offsets, from the local horizontal system to local equatorial offsets.

Converting between local flat projections and spherical coordinates usually requires a WCS projection.

REFERENCES:

  1. Calabretta, M.R., & Greisen, E.W., (2002), Astronomy & Astrophysics, 395, 1077-1122.
Parameters
daz[arcsec] Projected offset position in the azimuth direction. The projected offset between two azimuth positions at the same reference elevation is δAz = (Az2 - Az1) * cos(El0).
del[arcsec] projected offset position in the elevation direction
pa[deg] Parallactic Angle
[out]dra[arcsec] Output offset position in the local true-of-date R.A. direction. It can be a pointer to one of the input coordinates, or NULL if not required.
[out]ddec[arcsec] Output offset position in the local true-of-date declination direction. It can be a pointer to one of the input coordinates, or NULL if not required.
Returns
0
Since
1.3
Author
Attila Kovacs
See also
novas_e2h_offset(), novas_hpa()

◆ novas_hor_to_app()

int novas_hor_to_app ( const novas_frame *restrict frame,
double az,
double el,
RefractionModel ref_model,
enum novas_reference_system sys,
double *restrict ra,
double *restrict dec )

Converts an observed azimuth and elevation coordinate to right ascension (R.A.) and declination coordinates expressed in the coordinate system of choice.

The observer must be located on the surface of Earth, or else the call will return with an error. The caller may optionally supply a refraction model of choice to calculate an appropriate elevation angle that includes a refraction correction for Earth's atmosphere. If no such model is provided, the provided elevation value will be assumed to be an astrometric elevation without a refraction correction.

Parameters
frameObserver frame, defining the time and place of observation (on Earth).
az[deg] Observed azimuth angle. It may be NULL if not required.
el[deg] Observed elevation angle. It may be NULL if not required.
ref_modelAn appropriate refraction model, or NULL to assume unrefracted elevation. Depending on the refraction model, you might want to make sure that the weather parameters were set when the observing frame was defined.
sysAstronomical coordinate system in which the output is R.A. and declination values are to be calculated.
[out]ra[h] Calculated apparent right ascension (R.A.) coordinate
[out]dec[deg] Calculated apparent declination coordinate
Returns
0 if successful, or else an error from itrs_to_tod() or itrs_to_cirs(), or -1 (errno will indicate the type of error).
Since
1.1
Author
Attila Kovacs
See also
novas_app_to_hor()
novas_app_to_geom(), novas_standard_refraction(), novas_optical_refraction(), novas_radio_refraction(), novas_wave_refraction()

References novas_timespec::fjd_tt, hor_to_itrs(), novas_timespec::ijd_tt, NOVAS_AIRBORNE_OBSERVER, NOVAS_CIRS, NOVAS_GCRS, NOVAS_ICRS, NOVAS_ITRS, NOVAS_J2000, NOVAS_MOD, NOVAS_OBSERVER_ON_EARTH, NOVAS_REFRACT_OBSERVED, NOVAS_TIRS, NOVAS_TOD, spin(), vector2radec(), wobble(), WOBBLE_ITRS_TO_PEF, and WOBBLE_ITRS_TO_TIRS.

◆ novas_hpa()

double novas_hpa ( double az,
double el,
double lat )

Returns the Parallactic Angle (PA) calculated for a horizontal Az/El location of the sky.

The PA is the angle between the local horizontal coordinate directions and the local true-of-date equatorial coordinate directions at the given location. The polar wobble is not included in the calculation.

The Parallactic Angle is sometimes referrred to as the Vertical Position Angle (VPA). Both define the same quantity.

Parameters
az[deg] Azimuth angle
el[deg] Elevation angle
lat[deg] Geodetic latitude of observer
Returns
[deg] Parallactic Angle (PA). I.e., the clockwise position angle of the declination direction w.r.t. the elevation axis in the horizontal system. Same as the the clockwise position angle of the elevation direction w.r.t. the declination axis in the equatorial system.
Since
1.3
Author
Attila Kovacs
See also
novas_epa(), novas_h2e_offset()

◆ novas_los_to_xyz()

int novas_los_to_xyz ( const double * los,
double lon,
double lat,
double * xyz )

Converts a 3D line-of-sight vector (δφ, δθ δr) to a rectangular equatorial (δx, δy, δz) vector.

Parameters
los[arb.u.] Line-of-sight 3-vector (δφ, δθ δr).
lon[deg] Line-of-sight longitude.
lat[deg] Line-of-sight latitude.
[out]xyz[arb.u.] Output rectangular equatorial 3-vector (δx, δy, δz), in the same units as the input. It may be the same vector as the input.
Returns
0 if successful, or else -1 if either vector argument is NULL (errno will be set to EINVAL).
Since
1.3
Author
Attila Kovacs
See also
novas_xyz_to_los(), novas_uvw_to_xyz()

◆ novas_uvw_to_xyz()

int novas_uvw_to_xyz ( const double * uvw,
double ha,
double dec,
double * xyz )

Converts equatorial u,v,w projected (absolute or relative) coordinates to rectangular telescope x,y,z coordinates (in ITRS) to for a specified line of sight.

u,v,w are Cartesian coordinates (u,v) along the local equatorial R.A. and declination directions as seen from a direction on the sky (w). As such, they are effectively ITRS-based line-of-sight (LOS) coordinates.

x,y,z are Cartesian coordinates w.r.t the Greenwich meridian in the ITRS frame. The directions are x: long=0, lat=0; y: long=90, lat=0; z: lat=90.

Parameters
xyz[arb.u.] Absolute or relative u,v,w coordinates (double[3]).
ha[h] Hourangle (LST - RA) i.e., the difference between the Local (apparent) Sidereal Time and the apparent (true-of-date) Right Ascension of observed source.
dec[deg] Apparent (true-of-date) declination of source
[out]uvw[arb.u.] Converted x,y,z coordinates (double[3]) in the same unit as uvw. It may be the same vector as the input.
Returns
0 if successful, or else -1 if either vector argument is NULL (errno will be set to EINVAL)
Since
1.3
Author
Attila Kovacs
See also
novas_xyz_to_uvw()

References novas_los_to_xyz().

◆ novas_xyz_to_los()

int novas_xyz_to_los ( const double * xyz,
double lon,
double lat,
double * los )

Converts a 3D rectangular equatorial (δx, δy, δz) vector to a polar (δφ, δθ δr) vector along a line-of-sight.

Parameters
xyz[arb.u.] Rectangular equatorial 3-vector (δx, δy, δz).
lon[deg] Line-of-sight longitude.
lat[deg] Line-of-sight latitude.
[out]los[arb.u.] Output line-of-sight 3-vector (δφ, δθ δr), in the same units as the input. It may be the same vector as the input.
Returns
0 if successful, or else -1 if either vector argument is NULL (errno will be set to EINVAL).
Since
1.3
Author
Attila Kovacs
See also
novas_los_to_xyz(), novas_xyz_to_uvw()

◆ novas_xyz_to_uvw()

int novas_xyz_to_uvw ( const double * xyz,
double ha,
double dec,
double * uvw )

Converts rectangular telescope x,y,z (absolute or relative) coordinates (in ITRS) to equatorial u,v,w projected coordinates for a specified line of sight.

x,y,z are Cartesian coordinates w.r.t the Greenwich meridian, in the ITRS frame. The directions are x: long=0, lat=0; y: long=90, lat=0; z: lat=90.

u,v,w are Cartesian coordinates (u,v) along the local equatorial R.A. and declination directions as seen from a direction on the sky (w). As such, they are effectively ITRS-based line-of-sight (LOS) coordinates.

Parameters
xyz[arb.u.] Absolute or relative x,y,z coordinates (double[3]).
ha[h] Hourangle (LST - RA) i.e., the difference between the Local (apparent) Sidereal Time and the apparent (true-of-date) Right Ascension of observed source.
dec[deg] Apparent (true-of-date) declination of source
[out]uvw[arb.u.] Converted u,v,w coordinates (double[3]) in same units as xyz. It may be the same vector as the input.
Returns
0 if successful, or else -1 if either vector argument is NULL (errno will be set to EINVAL)
Since
1.3
Author
Attila Kovacs
See also
novas_uvw_to_xyz()

References novas_xyz_to_los().

◆ wobble()

int wobble ( double jd_tt,
enum novas_wobble_direction direction,
double xp,
double yp,
const double * in,
double * out )

Corrects a vector in the ITRS (rotating Earth-fixed system) for polar motion, and also corrects the longitude origin (by a tiny amount) to the Terrestrial Intermediate Origin (TIO).

The ITRS vector is thereby transformed to the Terrestrial Intermediate Reference System (TIRS), or equivalently the Pseudo Earth Fixed (PEF), based on the true (rotational) equator and TIO; or vice versa. Because the true equator is the plane orthogonal to the direction of the Celestial Intermediate Pole (CIP), the components of the output vector are referred to z and x axes toward the CIP and TIO, respectively.

NOTES:

  1. You may obtain Earth orientation parameters from the IERS site and Bulletins. For sub-mas accuracy you should augment the (interpolated) published values for diurnal and semi-diurnal librations and the effect of ocean tides, e.g. by using novas_diurnal_eop() or novas_diurnal_eop_at_time().
  2. Generally, this function should not be called if global pole offsets were set via cel_pole() and then used via place() or one of its variants to calculate Earth orientation corrected (TOD or CIRS) apparent coordinates. In such cases, calling wobble() would apply duplicate corrections. It is generally best to forgo using cel_pole() going forward, and instead apply Earth orinetation corrections with wobble() only when converting vectors between the Earth-fixed ITRS and TIRS frames.
  3. The Earth orientation parameters xp, yp should be provided in the same ITRF realization as the observer location for an Earth-based observer. You can use novas_itrf_transform_eop() to convert the EOP values as necessary.

REFERENCES:

  1. Kaplan, G. H. et. al. (1989). Astron. Journ. 97, 1197-1210.
  2. Lambert & Bizouard (2002), Astronomy and Astrophysics 394, 317-321.
  3. IERS Conventions 2010, Chapter 5, https://iers-conventions.obspm.fr/content/chapter5/icc5.pdf
  4. IERS Conventions 2010, Chapter 8, https://iers-conventions.obspm.fr/content/chapter8/icc8.pdf
Parameters
jd_tt[day] Terrestrial Time (TT) based Julian date.
directionWOBBLE_ITRS_TO_TIRS (0) or WOBBLE_TIRS_TO_ITRS (1) to include corrections for the TIO's longitude (IAU 2006 method); or else WOBBLE_ITRS_TO_PEF (2) or WOBBLE_PEF_TO_ITRS (3) to correct for dx, dy but not for the TIO's longitude (old, pre IAU 2006 method). Negative values default to WOBBLE_TIRS_TO_ITRS.
xp[arcsec] Conventionally-defined X coordinate of Celestial Intermediate Pole with respect to ITRS pole. As measured or else the interpolated published value from IERS, possibly augmented for diurnal variations caused by librations ocean tides if precision below the milliarcsecond level is required.
yp[arcsec] Conventionally-defined Y coordinate of Celestial Intermediate Pole with respect to ITRS pole, in arcseconds. As measured or else the interpolated published value from IERS, possibly augmented for diurnal variations caused by librations ocean tides if precision below the milliarcsecond level is required.
inInput position vector, geocentric equatorial rectangular coordinates, in the original system defined by 'direction'
[out]outOutput Position vector, geocentric equatorial rectangular coordinates, in the final system defined by 'direction'. It can be the same vector as the input.
Returns
0 if successful, or -1 if the direction is invalid output vector argument is NULL.
See also
novas_diurnal_eop(), novas_diurnal_eop_at_time(), novas_itrf_transfor_eop()

References NOVAS_WOBBLE_DIRECTIONS, WOBBLE_ITRS_TO_PEF, WOBBLE_ITRS_TO_TIRS, and WOBBLE_TIRS_TO_ITRS.