SuperNOVAS v1.5
The NOVAS C library, made better
Loading...
Searching...
No Matches
Earth orientation

Data Structures

struct  novas_delaunay_args
 Fundamental Delaunay arguments of the Sun and Moon, from Simon section 3.4(b.3). More...
 

Macros

#define NOVAS_EARTH_ANGVEL   7.2921150e-5
 [rad/s] Rotational angular velocity of Earth from IERS Conventions (2003).
 

Typedefs

typedef int(* novas_nutation_provider) (double jd_tt_high, double jd_tt_low, double *restrict dpsi, double *restrict deps)
 Function type definition for the IAU 2000 nutation series calculation.
 

Enumerations

enum  novas_equator_type { NOVAS_MEAN_EQUATOR = 0 , NOVAS_TRUE_EQUATOR , NOVAS_GCRS_EQUATOR }
 Constants that determine the type of equator to be used for the coordinate system. More...
 
enum  novas_equinox_type { NOVAS_MEAN_EQUINOX = 0 , NOVAS_TRUE_EQUINOX }
 The type of equinox used in the pre IAU 2006 (old) methodology. More...
 

Functions

double accum_prec (double t)
 Returns the general precession in longitude (Simon et al.
 
int e_tilt (double jd_tdb, enum novas_accuracy accuracy, double *restrict mobl, double *restrict tobl, double *restrict ee, double *restrict dpsi, double *restrict deps)
 (primarily for internal use) Computes quantities related to the orientation of the Earth's rotation axis at the specified Julian date.
 
double era (double jd_ut1_high, double jd_ut1_low)
 Returns the value of the Earth Rotation Angle (θ) for a given UT1 Julian date.
 
int fund_args (double t, novas_delaunay_args *restrict a)
 Compute the fundamental (a.k.a.
 
novas_nutation_provider get_nutation_lp_provider ()
 Returns the function configured for low-precision IAU 2000 nutation calculations instead of the default nu2000k().
 
int iau2000a (double jd_tt_high, double jd_tt_low, double *restrict dpsi, double *restrict deps)
 Computes the IAU 2000A high-precision nutation series for the specified date, to 0.1 μas accuracy.
 
int iau2000b (double jd_tt_high, double jd_tt_low, double *restrict dpsi, double *restrict deps)
 Computes the forced nutation of the non-rigid Earth based at reduced precision.
 
double mean_obliq (double jd_tdb)
 Computes the mean obliquity of the ecliptic.
 
int novas_diurnal_eop (double gmst, const novas_delaunay_args *restrict delaunay, double *restrict xp, double *restrict yp, double *restrict dut1)
 Calculate corrections to the Earth orientation parameters (EOP) due to short term (diurnal and semidiurnal) libration and the ocean tides.
 
int novas_diurnal_eop_at_time (const novas_timespec *restrict time, double *restrict dxp, double *restrict dyp, double *restrict dut1)
 Calculate corrections to the Earth orientation parameters (EOP) due to short term (diurnal and semidiurnal) libration and the ocean tides at a given astromtric time.
 
int novas_diurnal_libration (double gmst, const novas_delaunay_args *restrict delaunay, double *restrict xp, double *restrict yp, double *restrict dut1)
 Calculate diurnal and semi-diurnal libration corrections to the Earth orientation parameters (EOP) for the non-rigid Earth.
 
int novas_diurnal_ocean_tides (double gmst, const novas_delaunay_args *restrict delaunay, double *restrict xp, double *restrict yp, double *restrict dut1)
 Calculate corrections to the Earth orientation parameters (EOP) due to the ocean tides.
 
double novas_gast (double jd_ut1, double ut1_to_tt, enum novas_accuracy accuracy)
 Returns the Greenwich Apparent Sidereal Time (GAST) for a given UT1 date.
 
double novas_gmst (double jd_ut1, double ut1_to_tt)
 Returns the Greenwich Mean Sidereal Time (GMST) for a given UT1 date, using eq.
 
int novas_itrf_transform_eop (int from_year, double from_xp, double from_yp, double from_dut1, int to_year, double *restrict to_xp, double *restrict to_yp, double *restrict to_dut1)
 Transforms Earth orientation parameters (xp, yp, dUT1) from one ITRF realization to another.
 
double novas_time_gst (const novas_timespec *restrict time, enum novas_accuracy accuracy)
 Returns the Greenwich (apparent) Sidereal Time (GST/GaST) for a given astronomical time specification.
 
int nu2000k (double jd_tt_high, double jd_tt_low, double *restrict dpsi, double *restrict deps)
 Computes the forced nutation of the non-rigid Earth: Model NU2000K.
 
int nutation (double jd_tdb, enum novas_nutation_direction direction, enum novas_accuracy accuracy, const double *in, double *out)
 Nutates equatorial rectangular coordinates from mean equator and equinox of epoch to true equator and equinox of epoch.
 
int nutation_angles (double t, enum novas_accuracy accuracy, double *restrict dpsi, double *restrict deps)
 Returns the IAU2000 / 2006 values for nutation in longitude and nutation in obliquity for a given TDB Julian date and the desired level of accuracy.
 
short precession (double jd_tdb_in, const double *in, double jd_tdb_out, double *out)
 Precesses equatorial rectangular coordinates from one epoch to another using the IAU2006 (P03) precession model of Capitaine et al.
 
int set_nutation_lp_provider (novas_nutation_provider func)
 Set the function to use for low-precision IAU 2000 nutation calculations instead of the default nu2000k().
 
int terra (const on_surface *restrict location, double gast, double *restrict pos, double *restrict vel)
 Computes the position and velocity vectors of a terrestrial observer with respect to the center of the Earth, based on the GRS80 reference ellipsoid, used for the International Terrestrial Reference Frame (ITRF) and its realizations.
 

Detailed Description

Typedef Documentation

◆ novas_nutation_provider

typedef int(* novas_nutation_provider) (double jd_tt_high, double jd_tt_low, double *restrict dpsi, double *restrict deps)

Function type definition for the IAU 2000 nutation series calculation.

Parameters
jd_tt_high[day] High-order part of the Terrestrial Time (TT) based Julian date. Typically it may be the integer part of a split date for the highest precision, or the full date for normal (reduced) precision.
jd_tt_low[day] Low-order part of the Terrestrial Time (TT) based Julian date. Typically it may be the fractional part of a split date for the highest precision, or 0.0 for normal (reduced) precision.
[out]dpsi[rad] δψ Nutation (luni-solar + planetary) in longitude, in radians.
[out]deps[rad] δε Nutation (luni-solar + planetary) in obliquity, in radians.
Returns
0 if successful, or else -1 (errno should be set to indicate the type of error).
See also
nutation(), nutation_angles(), iau2000a(), iau2000b(), iau2000k()
Author
Attila Kovacs
Since
1.0

Enumeration Type Documentation

◆ novas_equator_type

Constants that determine the type of equator to be used for the coordinate system.

See also
NOVAS_EQUATOR_TYPES
equ2ecl(), ecl2equ()
Enumerator
NOVAS_MEAN_EQUATOR 

Mean celestial equator of date without nutation (pre IAU 2006 system).

NOVAS_TRUE_EQUATOR 

True celestial equator of date (pre IAU 2006 system).

NOVAS_GCRS_EQUATOR 

Geocentric Celestial Reference System (GCRS) equator.

◆ novas_equinox_type

The type of equinox used in the pre IAU 2006 (old) methodology.

See also
ira_equinox()
Enumerator
NOVAS_TRUE_EQUINOX 

Mean equinox: includes precession but not nutation.

True apparent equinox: includes both precession and nutation

Function Documentation

◆ accum_prec()

double accum_prec ( double t)

Returns the general precession in longitude (Simon et al.

1994), equivalent to 5028.8200 arcsec/cy at J2000.

Parameters
t[cy] Julian centuries since J2000
Returns
[rad] the approximate precession angle [-π:π].
See also
planet_lon(), nutation_angles(), e_tilt(), NOVAS_JD_J2000
Since
1.0
Author
Attila Kovacs

References TWOPI.

◆ e_tilt()

int e_tilt ( double jd_tdb,
enum novas_accuracy accuracy,
double *restrict mobl,
double *restrict tobl,
double *restrict ee,
double *restrict dpsi,
double *restrict deps )

(primarily for internal use) Computes quantities related to the orientation of the Earth's rotation axis at the specified Julian date.

In the pre-IAU2000 method, unmodelled corrections to earth orientation can be defined via cel_pole() prior to this call. However, we strongly recommend against that approach, and suggest you apply Earth orientation corrections only in novas_make_frame() or wobble().

NOTES:

  1. This function caches the results of the last calculation in case it may be re-used at no extra computational cost for the next call.
Parameters
jd_tdb[day] Barycentric Dynamical Time (TDB) based Julian date.
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
[out]mobl[deg] Mean obliquity of the ecliptic. It may be NULL if not required.
[out]tobl[deg] True obliquity of the ecliptic. It may be NULL if not required.
[out]ee[s] Equation of the equinoxes in seconds of time. It may be NULL if not required.
[out]dpsi[arcsec] Nutation in longitude. It may be NULL if not required.
[out]deps[arcsec] Nutation in obliquity. It may be NULL if not required.
Returns
0 if successful, or -1 if the accuracy argument is invalid
See also
novas_gast(), nutation(), ira_equinox(), equ2ecl(), ecl2equ()

References mean_obliq(), NOVAS_FULL_ACCURACY, NOVAS_REDUCED_ACCURACY, and nutation_angles().

◆ era()

double era ( double jd_ut1_high,
double jd_ut1_low )

Returns the value of the Earth Rotation Angle (θ) for a given UT1 Julian date.

The expression used is taken from the note to IAU Resolution B1.8 of 2000. The input Julian date can be split into an into high and low order parts (e.g. integer and fractional parts) for improved accuracy, or else one of the components (e.g. the low part) can be set to zero if no split is desired.

The algorithm used here is equivalent to the canonical θ = 0.7790572732640 + 1.00273781191135448 * t, where t is the time in days from J2000 (t = jd_high + jd_low - JD_J2000), but it avoids many two-PI 'wraps' that decrease precision (adopted from SOFA Fortran routine iau_era00).

REFERENCES:

  1. IERS Conventions 2010, Chapter 5, Eq. 5.15
  2. IAU Resolution B1.8, adopted at the 2000 IAU General Assembly, Manchester, UK.
  3. Kaplan, G. (2005), US Naval Observatory Circular 179.
Parameters
jd_ut1_high[day] High-order part of UT1 Julian date.
jd_ut1_low[day] Low-order part of UT1 Julian date.
Returns
[deg] The Earth Rotation Angle (theta) in degrees [0:360].
See also
novas_gast(), novas_time_gst()
cirs_to_itrs(), itrs_to_cirs()

◆ fund_args()

int fund_args ( double t,
novas_delaunay_args *restrict a )

Compute the fundamental (a.k.a.

Delaunay) arguments (mean elements) of the Sun and Moon.

REFERENCES:

  1. IERS Conventions 2010, Chapter 5, Eq. 5.43.
  2. Simon et al. (1994) Astronomy and Astrophysics 282, 663-683, esp. Sections 3.4-3.5.
Parameters
t[cy] TDB time in Julian centuries since J2000.0
[out]a[rad] Fundamental arguments data to populate (5 doubles) [0:2π]
Returns
0 if successful, or -1 if the output pointer argument is NULL.
See also
nutation_angles(), e_tilt(), NOVAS_JD_J2000

References novas_norm_ang().

◆ get_nutation_lp_provider()

novas_nutation_provider get_nutation_lp_provider ( )

Returns the function configured for low-precision IAU 2000 nutation calculations instead of the default nu2000k().

Returns
the function to use for low-precision IAU 2000 nutation calculations
Since
1.3
Author
Attila Kovacs
See also
set_nutation_lp_provider(), nutation_angles(), nutation()

◆ iau2000a()

int iau2000a ( double jd_tt_high,
double jd_tt_low,
double *restrict dpsi,
double *restrict deps )

Computes the IAU 2000A high-precision nutation series for the specified date, to 0.1 μas accuracy.

It is rather expensive computationally.

NOTES:

  1. As of SuperNOVAS v1.4.2, this function has been modified to replace the original IAU2000A model coefficients with the IAU2006 (a.k.a. IAU2000A R06) model coefficients, to provide an updated nutation model, which is dynamically consistent with the IAU2006 (P03) precesion model of Capitaine et al. 2003. This is now the same model with respect to which the IERS Earth orinetation parameters are computed and published.

REFERENCES:

  1. IERS Conventions (2003), Chapter 5.
  2. Simon et al. (1994) Astronomy and Astrophysics 282, 663-683, esp. Sections 3.4-3.5.
  3. Capitaine, N. et al. (2003), Astronomy And Astrophysics 412, pp. 567-586.
  4. https://hpiers.obspm.fr/eop-pc/models/nutations/nut.html
Parameters
jd_tt_high[day] High-order part of the Terrestrial Time (TT) based Julian date. Typically it may be the integer part of a split date for the highest precision, or the full date for normal (reduced) precision.
jd_tt_low[day] Low-order part of the Terrestrial Time (TT) based Julian date. Typically, it may be the fractional part of a split date for the highest precision, or 0.0 for normal (reduced) precision.
[out]dpsi[rad] δψ Nutation (luni-solar + planetary) in longitude. It may be NULL if not required.
[out]deps[rad] δε Nutation (luni-solar + planetary) in obliquity. It may be NULL if not required.
Returns
0
See also
iau2000b(), nu2000k(), nutation_angles(), nutation()

◆ iau2000b()

int iau2000b ( double jd_tt_high,
double jd_tt_low,
double *restrict dpsi,
double *restrict deps )

Computes the forced nutation of the non-rigid Earth based at reduced precision.

It reproduces the IAU 2000A (R06) model to a precision of 1 milliarcsecond in the interval 1995-2020, while being about 15x faster than iau2000a().

NOTES

  1. Originally this was the IAU2000B series of McCarthy & Luzum (2003), consistent with the original IAU2000 precession model

  2. As of SuperNOVAS v1.4.2, this function has been modified to use a truncated series for the IAU2006 (a.k.a. IAU 2000A R06) nutation model, whereby terms with amplitudes larger than 100 μas are omitted, resulting in 102 terms for the longitude and 57 terms the obliquity. This results in similar, or slightly better, precision than the original IAU2000B series of McCarthy & Luzum (2003), and it is now dynamically consistent with the IAU2006 (P03) precession model (Capitaine et al. 2005).

REFERENCES:

  1. McCarthy, D. and Luzum, B. (2003). "An Abridged Model of the Precession & Nutation of the Celestial Pole," Celestial Mechanics and Dynamical Astronomy, Volume 85, Issue 1, Jan. 2003, p. 37. (IAU 2000B) IERS Conventions (2003), Chapter 5.
  2. Capitaine, N. et al. (2003), Astronomy And Astrophysics 412, pp. 567-586.
Parameters
jd_tt_high[day] High-order part of the Terrestrial Time (TT) based Julian date. Typically it may be the integer part of a split date for the highest precision, or the full date for normal (reduced) precision.
jd_tt_low[day] Low-order part of the Terrestrial Time (TT) based Julian date. Typically it may be the fractional part of a split date for the highest precision, or 0.0 for normal (reduced) precision.
[out]dpsi[rad] δψ Nutation (luni-solar + planetary) in longitude, It may be NULL if not required.
[out]deps[rad] δε Nutation (luni-solar + planetary) in obliquity. It may be NULL if not required.
Returns
0
See also
iau2000a(), nu2000k(), nutation_angles(), set_nutation_lp_provider(), nutation()

◆ mean_obliq()

double mean_obliq ( double jd_tdb)

Computes the mean obliquity of the ecliptic.

REFERENCES:

  1. Capitaine et al. (2003), Astronomy and Astrophysics 412, 567-586.
Parameters
jd_tdb[day] Barycentric Dynamic Time (TDB) based Julian date
Returns
[arcsec] Mean obliquity of the ecliptic in arcseconds.
See also
e_tilt(), equ2ecl(), ecl2equ(), tt2tdb(), novas_get_time()

◆ novas_diurnal_eop()

int novas_diurnal_eop ( double gmst,
const novas_delaunay_args *restrict delaunay,
double *restrict dxp,
double *restrict dyp,
double *restrict dut1 )

Calculate corrections to the Earth orientation parameters (EOP) due to short term (diurnal and semidiurnal) libration and the ocean tides.

See Chapters 5 and 8 of the IERS Conventions. This one is faster than novas_diurnal_eop_at_time() if the GMST and/or fundamental arguments are readily avaialbe, e.g. from a previous calculation.

These corrections are typically at the sub-milliarcsecond level in polar offsets x and y and sub-millisecond level in UT1. Thus, these diurnal and semi-diurnal variations are important for the highest precision astrometry only, when converting between ITRS and TIRS frames.

REFERENCES:

  1. IERS Conventions 2010, Chapter 5, https://iers-conventions.obspm.fr/content/chapter5/icc5.pdf
  2. IERS Conventions 2010, Chapter 8, https://iers-conventions.obspm.fr/content/chapter8/icc8.pdf
Parameters
gmst[h] Greenwich Mean Sidereal Time of observation, e.g. as obtained by novas_gmst().
delaunay[rad] Delaunay arguments, e.g. as returned by fund_args().
[out]dxp[arcsec] x-pole correction for libration and ocean tides, or NULL if not required.
[out]dyp[arcsec] y-pole corrections for libration and ocean tides, or NULL if not required.
[out]dut1[s] UT1 correction for libration and ocean tides, or NULL if not required.
Returns
0 if successful, or else -1 if the delaunay arguments are NULL (errno set to EINVAL).
Since
1.5
Author
Attila Kovacs
See also
novas_diurnal_eop_at_time()
novas_diurnal_librations(), novas_diurnal_ocean_tides()

◆ novas_diurnal_eop_at_time()

int novas_diurnal_eop_at_time ( const novas_timespec *restrict time,
double *restrict dxp,
double *restrict dyp,
double *restrict dut1 )

Calculate corrections to the Earth orientation parameters (EOP) due to short term (diurnal and semidiurnal) libration and the ocean tides at a given astromtric time.

See Chapters 5 and 8 of the IERS Conventions.

REFERENCES:

  1. IERS Conventions 2010, Chapter 5, https://iers-conventions.obspm.fr/content/chapter5/icc5.pdf
  2. IERS Conventions 2010, Chapter 8, https://iers-conventions.obspm.fr/content/chapter8/icc8.pdf
Parameters
timeAstrometric time specification
dxp[arcsec] x-pole correction for libration and ocean tides, or NULL if not required.
dyp[arcsec] y-pole corrections for libration and ocean tides, or NULL if not required.
dut1[s] UT1 correction for libration and ocean tides, or NULL if not required.
Returns
0 if successful, or else -1 if the delaunay arguments are NULL (errno set to EINVAL).
Since
1.5
Author
Attila Kovacs
See also
novas_diurnal_eop()

References fund_args(), novas_diurnal_eop(), NOVAS_JD_J2000, NOVAS_REDUCED_ACCURACY, and novas_time_gst().

◆ novas_diurnal_libration()

int novas_diurnal_libration ( double gmst,
const novas_delaunay_args *restrict delaunay,
double *restrict dxp,
double *restrict dyp,
double *restrict dut1 )

Calculate diurnal and semi-diurnal libration corrections to the Earth orientation parameters (EOP) for the non-rigid Earth.

Only terms with periods less than 2 days are considered, since the longer period terms are included in the IAU2000A nutation model. See Chapter 5 of the IERS Conventions 2010. The typical amplitude of librations is tens of micro-arcseconds (μas) in polar offsets x and y and tens of micro-seconds in UT1.

Normally, you would need the combined effect from librations and ocean tides together to apply diurnal corrections to the Earth orientation parameters (EOP) published by IERS. For that, novas_diurnal_eop() or novas_diurnal_eop_at_time() is more directly suited.

REFERENCES:

  1. IERS Conventions 2010, Chapter 5, Tablea 5.1a and 5.1b, see https://iers-conventions.obspm.fr/content/chapter5/icc5.pdf
Parameters
gmst[h] Greenwich Mean Sidereal Time of observation, e.g. as obtained by novas_gmst().
delaunay[rad] Delaunay arguments, e.g. as returned by fund_args().
[out]dxp[arcsec] x-pole correction due to librations, or NULL if not required.
[out]dyp[arcsec] y-pole correction due to librations, or NULL if not required.
[out]dut1[s] UT1 correction due to librations, or NULL if not required.
Returns
0 if successful, or else -1 if the delaunay arguments are NULL (errno set to EINVAL).
Since
1.5
Author
Attila Kovacs
See also
novas_diurnal_eop(), novas_diurnal_ocean_tides()

◆ novas_diurnal_ocean_tides()

int novas_diurnal_ocean_tides ( double gmst,
const novas_delaunay_args *restrict delaunay,
double *restrict dxp,
double *restrict dyp,
double *restrict dut1 )

Calculate corrections to the Earth orientation parameters (EOP) due to the ocean tides.

See Chapter 8 of the IERS Conventions 2010. Ocean tides manifest at the sub-milliarcsecond level in polar offsets x and y and sub-millisecond level in UT1.

Normally, you would need the combined effect from librations and ocean tides together to apply diurnal corrections to the Earth orientation parameters (EOP) published by IERS. For that, novas_diurnal_eop() or novas_diurnal_eop_at_time() is more directly suited.

REFERENCES:

  1. IERS Conventions 2010, Chapter 8, Tables 8.2a, 8.2b, 8.3a, and 8.3b, see https://iers-conventions.obspm.fr/content/chapter8/icc8.pdf
Parameters
gmst[h] Greenwich Mean Sidereal Time of observation, e.g. as obtained by novas_gmst().
delaunay[rad] Delaunay arguments, e.g. as returned by fund_args().
[out]dxp[arcsec] x-pole correction due to ocean tides, or NULL if not required.
[out]dyp[arcsec] y-pole correction due to ocean tides, or NULL if not required.
[out]dut1[s] UT1 correction due to ocean tides, or NULL if not required.
Returns
0 if successful, or else -1 if the delaunay arguments are NULL (errno set to EINVAL).
Since
1.5
Author
Attila Kovacs
See also
novas_diurnal_eop(), novas_diurnal_librations()

◆ novas_gast()

double novas_gast ( double jd_ut1,
double ut1_to_tt,
enum novas_accuracy accuracy )

Returns the Greenwich Apparent Sidereal Time (GAST) for a given UT1 date.

Parameters
jd_ut1[day] UT1-based Julian Date.
ut1_to_tt[s] UT1 - UTC time difference.
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
Returns
[h] The Greenwich Apparent Sidereal Time (GAST) in the 0-24 range.

REFERENCES:

  1. Capitaine, N. et al. (2003), Astronomy and Astrophysics 412, 567-586, eq. (42).
  2. Kaplan, G. (2005), US Naval Observatory Circular 179.
  3. IERS Conventions 2010, Chapter 5, Eq. 5.30 and Table 5.2e, see at https://iers-conventions.obspm.fr/content/chapter5/icc5.pdf and https://iers-conventions.obspm.fr/content/chapter5/additional_info/tab5.2e.txt
Since
1.5
Author
Attila Kovacs
See also
novas_gmst()

References e_tilt(), and novas_gmst().

◆ novas_gmst()

double novas_gmst ( double jd_ut1,
double ut1_to_tt )

Returns the Greenwich Mean Sidereal Time (GMST) for a given UT1 date, using eq.

42 from Capitaine et al. (2003).

REFERENCES:

  1. Capitaine, N. et al. (2003), Astronomy and Astrophysics 412, 567-586, eq. (42).
  2. IERS Conventions 2010, Chapter 5, Eq. 5.32, see at https://iers-conventions.obspm.fr/content/chapter5/icc5.pdf
Parameters
jd_ut1[day] UT1-based Julian Date.
ut1_to_tt[s] UT1 - UTC time difference.
Returns
[h] The Greenwich Mean Sidereal Time (GMST) in the 0-24 range.
Since
1.5
Author
Attila Kovacs
See also
novas_gast(), novas_time_gst()

References era().

◆ novas_itrf_transform_eop()

int novas_itrf_transform_eop ( int from_year,
double from_xp,
double from_yp,
double from_dut1,
int to_year,
double *restrict to_xp,
double *restrict to_yp,
double *restrict to_dut1 )

Transforms Earth orientation parameters (xp, yp, dUT1) from one ITRF realization to another.

For the highest precision applications, observing sites should be defined in the same ITRF realization as the IERS Earth orientation parameters (EOP). To reconcile you may transform either the site location or the EOP between different realizations to match.

REFERENCES:

  1. IERS Conventions, Chapter 4, Eq. 4.14 and Table 4.1. See https://iers-conventions.obspm.fr/content/chapter4/icc4.pdf
Parameters
from_year[yr] ITRF realization year of input coordinates / rates. E.g. 1992 for ITRF92.
[in]from_xp[arcsec] x-pole Earth orientation parameter (angle) in the input ITRF realization.
[in]from_yp[arcsec] y-pole Earth orientation parameter (angle) in the input ITRF realization.
[in]from_dut1[s] UT1-UTC time difference in the input ITRF realization.
to_year[yr] ITRF realization year of input coordinates / rates. E.g. 2000 for ITRF2000.
[out]to_xp[arcsec] x-pole Earth orientation parameter (angle) in the output ITRF realization, or NULL if not required.
[out]to_yp[arcsec] y-pole Earth orientation parameter (angle) in the output ITRF realization, or NULL if not required.
[out]to_dut1[s] UT1-UTC time difference in the output ITRF realization, or NULL if not required.
Returns
0
Since
1.5
Author
Attila Kovacs
See also
novas_itrf_transform(), novas_itrf_transform_site()
novas_make_frame(), novas_timespec, wobble()

References NOVAS_DAY, NOVAS_EARTH_FLATTENING, and TWOPI.

◆ novas_time_gst()

double novas_time_gst ( const novas_timespec *restrict time,
enum novas_accuracy accuracy )

Returns the Greenwich (apparent) Sidereal Time (GST/GaST) for a given astronomical time specification.

If you need mean sidereal time (GMST), you should use sidereal_time() instead.

Parameters
timeThe astronomoical time specification.
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1).
Returns
[h] The Greenwich apparent Sidereal Time (GST / GaST) in the [0:24) hour range, or else NAN if there was an error (errno will indicate the type of error).
Since
1.4
Author
Attila Kovacs
See also
novas_time_lst(), novas_set_time()

References novas_gast(), novas_get_time(), and NOVAS_UT1.

◆ nu2000k()

int nu2000k ( double jd_tt_high,
double jd_tt_low,
double *restrict dpsi,
double *restrict deps )

Computes the forced nutation of the non-rigid Earth: Model NU2000K.

This model is a modified version of the original IAU 2000A, which has been truncated for speed of execution. NU2000K agrees with IAU 2000A at the 0.1 milliarcsecond level from 1700 to 2300, while being is about 5x faster than the more precise iau2000a().

NU2000K was compared to IAU 2000A over six centuries (1700-2300). The average error in dψ is 20 microarcseconds, with 98% of the errors < 60 microarcseconds; the average error in dεis 8 microarcseconds, with 100% of the errors < 60 microarcseconds.

NU2000K was developed by G. Kaplan (USNO) in March 2004

NOTES:

  1. As of SuperNOVAS v1.4.2, this function has been modified to include a rescaling of the original IAU2000 nutation compatible values to 'IAU2006' (i.e. IAU2000A R06) compatible values, according to Coppola, Seago, and Vallado (2009). The rescaling makes this model more dynamically consistent with the IAU2006 (P03) precession model of Capitaine et al. (2003).

REFERENCES:

  1. IERS Conventions (2003), Chapter 5.
  2. Simon et al. (1994) Astronomy and Astrophysics 282, 663-683, esp. Sections 3.4-3.5.
  3. Capitaine, N. et al. (2003), Astronomy And Astrophysics 412, pp. 567-586.
  4. Coppola, V., Seago, G.H., & Vallado, D.A. (2009), AAS 09-159
Parameters
jd_tt_high[day] High-order part of the Terrestrial Time (TT) based Julian date. Typically it may be the integer part of a split date for the highest precision, or the full date for normal (reduced) precision.
jd_tt_low[day] Low-order part of the Terrestrial Time (TT) based Julian date. Typically it may be the fractional part of a split date for the highest precision, or 0.0 for normal (reduced) precision.
[out]dpsi[rad] δψ Nutation (luni-solar + planetary) in longitude, It may be NULL if not required.
[out]deps[rad] δε Nutation (luni-solar + planetary) in obliquity. It may be NULL if not required.
Returns
0
See also
iau2000a(), iau2000b(), nutation_angles(), set_nutation_lp_provider(), nutation()

References accum_prec(), novas_delaunay_args::D, novas_delaunay_args::F, fund_args(), novas_delaunay_args::l, novas_delaunay_args::l1, NOVAS_EARTH, NOVAS_JUPITER, NOVAS_MARS, NOVAS_SATURN, NOVAS_VENUS, novas_delaunay_args::Omega, and planet_lon().

◆ nutation()

int nutation ( double jd_tdb,
enum novas_nutation_direction direction,
enum novas_accuracy accuracy,
const double * in,
double * out )

Nutates equatorial rectangular coordinates from mean equator and equinox of epoch to true equator and equinox of epoch.

Inverse transformation may be applied by setting flag 'direction'.

This is the old (pre IAU 2006) method of nutation calculation. If you follow the now standard IAU 2000 / 2006 methodology you will want to use nutation_angles() instead.

REFERENCES:

  1. Explanatory Supplement To The Astronomical Almanac, pp. 114-115.
Parameters
jd_tdb[day] Barycentric Dynamic Time (TDB) based Julian date
directionNUTATE_MEAN_TO_TRUE (0) or NUTATE_TRUE_TO_MEAN (-1; or non-zero)
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
inPosition 3-vector, geocentric equatorial rectangular coordinates, referred to mean equator and equinox of epoch.
[out]outPosition vector, geocentric equatorial rectangular coordinates, referred to true equator and equinox of epoch. It can be the same as the input position.
Returns
0 if successful, or -1 if one of the vector arguments is NULL.
See also
nutation_angles()
tt2tdb(), novas_get_time(), NOVAS_MOD, NOVAS_TOD

References e_tilt(), and NUTATE_MEAN_TO_TRUE.

◆ nutation_angles()

int nutation_angles ( double t,
enum novas_accuracy accuracy,
double *restrict dpsi,
double *restrict deps )

Returns the IAU2000 / 2006 values for nutation in longitude and nutation in obliquity for a given TDB Julian date and the desired level of accuracy.

For NOVAS_FULL_ACCURACY (0), the IAU 2000A R06 nutation model is used. Otherwise, the model set by set_nutation_lp_provider() is used, or else the iau2000b() by default.

REFERENCES:

  1. Kaplan, G. (2005), US Naval Observatory Circular 179.
  2. Capitaine, N., P.T. Wallace and J. Chapront (2005), “Improvement of the IAU 2000 precession model.” Astronomy & Astrophysics, Vol. 432, pp. 355–67.
  3. Coppola, V., Seago, G.H., & Vallado, D.A. (2009), AAS 09-159
Parameters
t[cy] TDB time in Julian centuries since J2000.0
accuracyNOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
[out]dpsi[arcsec] Nutation in longitude in arcseconds.
[out]deps[arcsec] Nutation in obliquity in arcseconds.
Returns
0 if successful, or -1 if the output pointer arguments are NULL
See also
nutation(), set_nutation_lp_provider()
iau2000a(), iau2000b(), nu2000k()
NOVAS_JD_J2000

References get_nutation_lp_provider(), iau2000a(), and NOVAS_FULL_ACCURACY.

◆ precession()

short precession ( double jd_tdb_in,
const double * in,
double jd_tdb_out,
double * out )

Precesses equatorial rectangular coordinates from one epoch to another using the IAU2006 (P03) precession model of Capitaine et al.

2003.

NOTE:

  1. Unlike the original NOVAS C 3.1 version, this one does not require that one of the time arguments must be J2000. You can precess from any date to any other date, and the intermediate epoch of J2000 will be handled internally as needed.

  2. This function caches the results of the last calculation in case it may be re-used at no extra computational cost for the next call.

REFERENCES:

  1. Explanatory Supplement To The Astronomical Almanac, pp. 103-104.
  2. Capitaine, N. et al. (2003), Astronomy And Astrophysics 412, pp. 567-586.
  3. Hilton, J. L. et al. (2006), IAU WG report, Celest. Mech., 94, pp. 351-367.
  4. Capitaine, N., P.T. Wallace and J. Chapront (2005), “Improvement of the IAU 2000 precession model.” Astronomy & Astrophysics, Vol. 432, pp. 355–67.
  5. Liu, J.-C., & Capitaine, N. (2017), A&A 597, A83
Parameters
jd_tdb_in[day] Barycentric Dynamic Time (TDB) based Julian date of the input epoch
inPosition 3-vector, geocentric equatorial rectangular coordinates, referred to mean dynamical equator and equinox of the initial epoch.
jd_tdb_out[day] Barycentric Dynamic Time (TDB) based Julian date of the output epoch
[out]outPosition 3-vector, geocentric equatorial rectangular coordinates, referred to mean dynamical equator and equinox of the final epoch. It can be the same vector as the input.
Returns
0 if successful, or -1 if either of the position vectors is NULL.
See also
nutation()
tt2tdb(), novas_get_time(), frame_tie(), novas_epoch(), NOVAS_MOD, NOVAS_JD_J2000, NOVAS_JD_B1950, NOVAS_JD_B1900

References precession().

◆ set_nutation_lp_provider()

int set_nutation_lp_provider ( novas_nutation_provider func)

Set the function to use for low-precision IAU 2000 nutation calculations instead of the default nu2000k().

Parameters
functhe new function to use for low-precision IAU 2000 nutation calculations
Returns
0 if successful, or -1 if the input argument is NULL
Since
1.0
Author
Attila Kovacs
See also
get_nutation_lp_provider(), nutation_angles(), nutation()

◆ terra()

int terra ( const on_surface *restrict location,
double gast,
double *restrict pos,
double *restrict vel )

Computes the position and velocity vectors of a terrestrial observer with respect to the center of the Earth, based on the GRS80 reference ellipsoid, used for the International Terrestrial Reference Frame (ITRF) and its realizations.

This function ignores polar motion, unless the observer's longitude and latitude have been corrected for it, and variation in the length of day (angular velocity of earth).

The true equator and equinox of date do not form an inertial system. Therefore, with respect to an inertial system, the very small velocity component (several meters/day) due to the precession and nutation of the Earth's axis is not accounted for here.

REFERENCES:

  1. Kaplan, G. H. et. al. (1989). Astron. Journ. 97, 1197-1210.
Parameters
locationITRF / GRS80 geodetic location of observer in Earth's rotating frame.
gast[h] Greenwich apparent sidereal time (GAST).
[out]pos[AU] Position vector of observer with respect to center of Earth, equatorial rectangular coordinates, referred to true equator and equinox of date, components in AU. If reference meridian is Greenwich and 'lst' = 0, 'pos' is effectively referred to equator and Greenwich. (It may be NULL if no position data is required).
[out]vel[AU/day] Velocity vector of observer with respect to center of Earth, equatorial rectangular coordinates, referred to true equator and equinox of date, components in AU/day. (It must be distinct from the pos output vector, and may be NULL if no velocity data is required).
Returns
0 if successful, or -1 if location is NULL or if the pos and vel output arguments are identical pointers.
See also
geo_posvel()
make_gps_site(), make_itrf_site(), make_xyz_site(), novas_gast(), novas_time_gst()

References NOVAS_KM.