![]() |
SuperNOVAS v1.5
The NOVAS C library, made better
|
Various transformations between different coordinate systems used in astronomy, such as equatorial, ecliptic, Galactic, or local horizontal coordinate systems. 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 | equ2hor (double jd_ut1, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const on_surface *restrict location, double ra, double dec, enum novas_refraction_model ref_option, double *restrict zd, double *restrict az, double *restrict rar, double *restrict decr) |
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. | |
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. | |
Various transformations between different coordinate systems used in astronomy, such as equatorial, ecliptic, Galactic, or local horizontal coordinate systems.
This module also provides functions otherwise relating to the observer locations.
For Earth-based geodetic locations, the module provides a set of functions to convert between various (local) coordinate systems. For example, calculating parallactic angles (a.k.a. vertical position angle or VPA), for converting between horizontal and equatorial offsets, and various functions for converting equatorial vectors to a line-of-sight, or u, v, w coordinate bases, and vice-versa.
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'.
jd_tt | [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2]) | |
coord_sys | The astrometric reference system of the coordinates. If 'coord_sys' is NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to J2000 ecliptic coordinates. | |
accuracy | NOVAS_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. |
References 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'.
jd_tt | [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2]) | |
coord_sys | The astrometric reference system type of the coordinates | |
accuracy | NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) | |
in | Position vector, referred to specified ecliptic and equinox of date. | |
[out] | out | Position vector, referred to specified equator and equinox of date. It can be the same vector as the input. |
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.
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'.
jd_tt | [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2]) | |
coord_sys | The astrometric reference system of the coordinates. If 'coord_sys' is NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to J2000 ecliptic coordinates. | |
accuracy | NOVAS_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. |
References 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'.
jd_tt | [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' is NOVAS_GCRS_EQUATOR[2]) | |
coord_sys | The astrometric reference system type of the coordinates. | |
accuracy | NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) | |
in | Position vector, referred to specified equator and equinox of date. | |
[out] | out | Position 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. |
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.
int equ2gal | ( | double | ra, |
double | dec, | ||
double *restrict | glon, | ||
double *restrict | glat ) |
Converts ICRS right ascension and declination to galactic longitude and latitude.
REFERENCES:
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. |
int equ2hor | ( | double | jd_ut1, |
double | ut1_to_tt, | ||
enum novas_accuracy | accuracy, | ||
double | xp, | ||
double | yp, | ||
const on_surface *restrict | location, | ||
double | ra, | ||
double | dec, | ||
enum novas_refraction_model | ref_option, | ||
double *restrict | zd, | ||
double *restrict | az, | ||
double *restrict | rar, | ||
double *restrict | decr ) |
novas_app_to_hor()
instead, or else the more explicit (less ambiguous) tod_to_itrs()
followed by itrs_to_hor()
, and possibly following it with an atmospheric refraction correction if appropriate.Transforms topocentric (TOD) apparent right ascension and declination to zenith distance and azimuth. This method should not be used to convert CIRS apparent coordinates (IAU 2000 standard) – for those you should use cirs_to_itrs() followed by itrs_to_hor() instead.
It uses a method that properly accounts for polar motion, which is significant at the sub-arcsecond level. This function can also adjust coordinates for atmospheric refraction.
NOTES:
REFERENCES:
jd_ut1 | [day] UT1 based Julian date | |
ut1_to_tt | [s] TT - UT1 Time difference in seconds | |
accuracy | NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) | |
xp | [arcsec] Conventionally-defined x coordinate of celestial intermediate pole with respect to ITRS reference pole, e.g. from IERS Bulletin A. If you have defined pole offsets to be incorporated into the TOD input coordinates (pre-IAU2000 method) via cel_pole() , then you should set this to 0. | |
yp | [arcsec] Conventionally-defined y coordinate of celestial intermediate pole with respect to ITRS reference pole, e.g. from IERS Bulletin A. If you have defined pole offsets to be incorporated into the TOD input coordinates (pre-IAU2000 method) via cel_pole() , then you should set this to 0. | |
location | Geodetic (ITRF / GRS80) observer location on Earth | |
ra | [h] Topocentric apparent (TOD) right ascension of object of interest, referred to true equator and equinox of date. | |
dec | [deg] Topocentric apparent (TOD) declination of object of interest, referred to true equator and equinox of date. | |
ref_option | Refraction model to use. E.g., NOVAS_STANDARD_ATMOSPHERE (1), or NOVAS_WEATHER_AT_LOCATION (2) if to use the weather. | |
[out] | zd | [deg] Topocentric zenith distance in degrees (unrefracted). |
[out] | az | [deg] Topocentric azimuth (measured east from north) in degrees. |
[out] | rar | [h] Topocentric right ascension of object of interest, in hours, referred to true equator and equinox of date, affected by refraction if 'ref_option' is non-zero. (It may be NULL if not required) |
[out] | decr | [deg] Topocentric declination of object of interest, in degrees, referred to true equator and equinox of date. (It may be NULL if not required) |
References EROT_GST, NOVAS_DYNAMICAL_CLASS, refract_astro(), and ter2cel().
int gal2equ | ( | double | glon, |
double | glat, | ||
double *restrict | ra, | ||
double *restrict | dec ) |
Converts galactic longitude and latitude to ICRS right ascension and declination.
REFERENCES:
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. |
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.
location | Geodetic (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] | itrs | Unit 3-vector direction in 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.
location | Geodetic (ITRF / GRS80) observer location on Earth. | |
itrs | 3-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. |
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:
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. |
References novas_h2e_offset().
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.
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 |
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:
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. |
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.
az | [deg] Azimuth angle |
el | [deg] Elevation angle |
lat | [deg] Geodetic latitude of observer |
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.
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. |
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.
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. |
References novas_los_to_xyz().
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.
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. |
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.
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. |
References novas_xyz_to_los().