This document describes version 0.1d6 of ATSOFA, the International Astronomical Union's SOFA Collection of libraries repackaged using the GNU Autotools.
ATSOFA is neither distributed, nor supported, nor endorsed by the International Astronomical Union. Any use of this pacakge should comply with SOFA's license and terms of use. Especially, but not exclusively, any published work or commercial products which includes results achieved by using ATSOFA shall acknowledge that the SOFA software was used in obtaining those results.
The package is distributed under the terms of the International Astronomical Union's SOFA license and its development takes place at:
and as backup at:
proper release tarballs are available from:
the International Astronomical Union's SOFA project home page is at:
Copyright © International Astronomical Union Standards of Fundamental Astronomy (http://www.iausofa.org).
The repackaging into ATSOFA is by Marco Maggi marco.maggi-ipsu@poste.it.
Appendices
Indexes
The International Astronomical Union's SOFA Collection consists of two libraries of routines, one coded in Fortran 77 the other in ANSI C. There is a suite of vector/matrix routines and various utilities that underpin the astronomy algorithms, which include routines for the following:
This package makes use of the GNU Autotools to allow easy distribution and installation of SOFA.
ATSOFA is neither distributed, nor supported, nor endorsed by the International Astronomical Union. Any use of this pacakge should comply with SOFA's license and terms of use. Especially, but not exclusively, any published work or commercial products which includes results achieved by using ATSOFA shall acknowledge that the SOFA software was used in obtaining those results.
SOFA stands for “Standards Of Fundamental Astronomy”. The SOFA software libraries are a collection of subprograms, in source code form, which implement official IAU algorithms for fundamental astronomy computations. The subprograms at present comprise 131 ”astronomy” routines supported by 55 “vector/matrix” routines, available in both Fortran77 and C implementations.
SOFA is an IAU Service which operates under Division 1 (Fundamental Astronomy) and reports through Commission 19 (Rotation of the Earth).
The IAU set up the SOFA initiative at the 1994 General Assembly, to promulgate an authoritative set of fundamental astronomy constants and algorithms. At the subsequent General Assembly, in 1997, the appointment of a SOFA Review Board and the selection of a site for the SOFA Center (the outlet for SOFA products) were announced.
The SOFA initiative was originally proposed by the IAU Working Group on Astronomical Standards (WGAS), under the chairmanship of Toshio Fukushima. The proposal was for “...new arrangements to establish and maintain an accessible and authoritative set of constants, algorithms and procedures that implement standard models used in fundamental astronomy”. The SOFA Software Libraries implement the “algorithms” part of the SOFA initiative. They were developed under the supervision of an international panel called the SOFA Review Board. The current membership of this panel is listed in an appendix.
A feature of the original SOFA software proposals was that the products would be self–contained and not depend on other software. This includes basic documentation, which, like the present file, will mostly be plain ASCII text. It should also be noted that there is no assumption that the software will be used on a particular computer and Operating System. Although OS–related facilities may be present (Unix make files for instance, use by the SOFA Center of automatic code management systems, HTML versions of some documentation), the routines themselves will be visible as individual text files and will run on a variety of platforms.
The SOFA Review Board's initial goal has been to create a set of callable subprograms. Whether “subroutines” or “functions”, they are all referred to simply as “routines”. They are designed for use by software developers wishing to write complete applications; no runnable, free–standing applications are included in SOFA's present plans.
The algorithms are drawn from a variety of sources. Because most of the routines so far developed have either been standard “text–book” operations or implement well–documented standard algorithms, it has not been necessary to invite the whole community to submit algorithms, though consultation with authorities has occurred where necessary. It should also be noted that consistency with the conventions published by the International Earth Rotation Service was a stipulation in the original SOFA proposals, further constraining the software designs. This state of affairs will continue to exist for some time, as there is a large backlog of agreed extensions to work on. However, in the future the Board may decide to call for proposals, and is in the meantime willing to look into any suggestions that are received by the SOFA Center.
The routines currently available are listed in the next two chapters of this document.
The “astronomy” library comprises 131 routines (plus one obsolete Fortran routine that now appears under a revised name). The areas addressed include calendars, time scales, ephemerides, precession–nutation, star space–motion, star catalog transformations and geodetic/geocentric transformations.
The “vector–matrix” library, comprising 55 routines, contains a collection of simple tools for manipulating the vectors, matrices and angles used by the astronomy routines.
There is no explicit commitment by SOFA to support historical models, though as time goes on a legacy of superseded models will naturally accumulate. There is, for example, no support of B1950/FK4 star coordinates, or pre-1976 precession models, though these capabilities could be added were there significant demand.
Though the SOFA software libraries are rather limited in scope, and are likely to remain so for a considerable time, they do offer distinct advantages to prospective users. In particular, the routines are:
Once it has been published, an issue is never revised or updated, and remains accessible indefinitely. Subsequent issues may, however, include corrected versions under the original routine name and filenames. However, where a different model is introduced, it will have a different name.
The issues will be referred to by the date when they were announced. The frequency of re–issue will be decided by the Board, taking into account the importance of the changes and the impact on the user community.
At present there is little free–standing documentation about individual routines. However, each routine has preamble comments which specify in detail what the routine does and how it is used.
The file sofa_pn.pdf describes the SOFA tools for precession–nutation and other aspects of Earth attitude and includes example code and (see the appendix) diagrams showing the interrelationships between the routines supporting the latest (IAU 2006/2000A) models.
The SOFA routines are available in two programming languages at present: Fortran77 and ANSI C. Related software in other languages is under consideration.
The Fortran code conforms to ANSI X3.9-1978 in all but two
minor respects: each has an IMPLICIT NONE
declaration, and its
name has a prefix of iau_
and may be longer than 6 characters. A
global edit to erase both of these will produce ANSI
–compliant
code with no change in its function.
Coding style, and restrictions on the range of language features, have been much debated by the Board, and the results comply with the majority view. There is (at present) no document that defines the standards, but the code itself offers a wide range of examples of what is acceptable.
The Fortran routines contain explicit numerical constants (the
INCLUDE
statement is not part of ANSI Fortran77).
These are drawn from the file consts.lis, which is listed in an
appendix. Constants for the SOFA/C functions are defined in a header
file sofam.h.
The naming convention is such that a SOFA routine referred to
generically as EXAMPL
exists as a Fortran subprogram
iau_EXAMPL
and a C function iauExampl
. The calls for the
two versions are very similar, with the same arguments in the same
order. In a few cases, the C equivalent of a Fortran SUBROUTINE
subprogram uses a return value rather than an argument.
Each language version includes a “testbed” main–program that can be used to verify that the SOFA routines have been correctly compiled on the end user's system. The Fortran and C versions are called t_sofa_f.for and t_sofa_c.c respectively. The testbeds execute every SOFA routine and check that the results are within expected accuracy margins. It is not possible to guarantee that all platforms will meet the rather stringent criteria that have been used, and an occasional warning message may be encountered on some systems.
Copyright for all of the SOFA software and documentation is owned by the IAU SOFA Review Board. The Software is made available free of charge for all classes of user, including commercial. However, there are strict rules designed to avoid unauthorized variants coming into circulation. It is permissible to distribute derived works and other modifications, but they must be clearly marked to avoid confusion with the SOFA originals.
Further details are included in the block of comments which concludes every routine. The text is also set out in an appendix to the present document.
The SOFA policy is to organize the calculations so that the machine accuracy is fully exploited. The gap between the precision of the underlying model or theory and the computational resolution has to be kept as large as possible, hopefully leaving several orders of magnitude of headroom.
The SOFA routines in some cases involve design compromises between rigor and ease of use (and also speed, though nowadays this is seldom a major concern).
The Board is indebted to a number of contributors, who are acknowledged in the preamble comments of the routines concerned.
The Board's effort is provided by the members' individual institutes.
Resources for operating the SOFA Center are provided by Her Majesty's Nautical Almanac Office, operated by the United Kingdom Hydrographic Office.
The routines described here comprise the SOFA astronomy library. Their general appearance and coding style conforms to conventions agreed by the SOFA Review Board, and their functions, names and algorithms have been ratified by the Board. Procedures for soliciting and agreeing additions to the library are still evolving.
The SOFA routines are available in two programming languages at present: Fortran 77 and ANSI C.
Except for a single obsolete Fortran routine, which has no C equivalent,
there is a one–to–one relationship between the two language versions.
The naming convention is such that a SOFA routine referred to
generically as EXAMPL
exists as a Fortran subprogram
iau_EXAMPL
and a C function iauExampl
. The calls for the
two versions are very similar, with the same arguments in the same
order. In a few cases, the C equivalent of a Fortran SUBROUTINE
subprogram uses a return value rather than an argument.
The principal function of the SOFA Astronomy Library is to provide definitive algorithms. A secondary function is to provide software suitable for convenient direct use by writers of astronomical applications.
The astronomy routines call on the SOFA vector/matrix library routines, which are separately listed.
The routines are designed to exploit the full floating–point accuracy of the machines on which they run, and not to rely on compiler optimizations. Within these constraints, the intention is that the code corresponds to the published formulation (if any).
Dates are always Julian Dates (except in calendar conversion routines) and are expressed as two double precision numbers which sum to the required value.
A distinction is made between routines that implement IAU–approved models and those that use those models to create other results. The former are referred to as “canonical models” in the preamble comments; the latter are described as “support routines”.
Using the library requires knowledge of positional astronomy and time–scales. These topics are covered in “Explanatory Supplement to the Astronomical Almanac”, P. Kenneth Seidelmann (ed.), University Science Books, 1992. Recent developments are documented in the journals, and references to the relevant papers are given in the SOFA code as required. The IERS Conventions are also an essential reference. The routines concerned with Earth attitude (precession–nutation etc.) are described in the SOFA document sofa_pn.pdf.
Calendars:
CAL2JD
EPB
EPB2JD
EPJ
EPJ2JD
JD2CAL
JDCALF
Time scales:
D2DTF
DAT
DTDB
DTF2D
TAITT
TAIUT1
TAIUTC
TCBTDB
TCGTT
TDBTCB
TDBTT
TTTAI
TTTCG
TTTDB
TTUT1
UT1TAI
UT1TT
UT1UTC
UTCTAI
UTCUT1
Earth rotation angle and sidereal time:
EE00
EE00A
EE00B
EE06A
EECT00
EQEQ94
ERA00
GMST00
GMST06
GMST82
GST00A
GST00B
GST06
GST06A
GST94
Ephemerides (limited precision):
EPV00
PLAN94
Precession, nutation, polar motion:
BI00
BP00
BP06
BPN2XY
C2I00A
C2I00B
C2I06A
C2IBPN
C2IXY
C2IXYS
C2T00A
C2T00B
C2T06A
C2TCIO
C2TEQX
C2TPE
C2TXY
EO06A
EORS
FW2M
FW2XY
NUM00A
NUM00B
NUM06A
NUMAT
NUT00A
NUT00B
NUT06A
NUT80
NUTM80
OBL06
OBL80
PB06
PFW06
PMAT00
PMAT06
PMAT76
PN00
PN00A
PN00B
PN06
PN06A
PNM00A
PNM00B
PNM06A
PNM80
P06E
POM00
PR00
PREC76
S00
S00A
S00B
S06
S06A
SP00
XY06
XYS00A
XYS00B
XYS06A
Fundamental arguments for nutation etc.
FAD03
FAE03
FAF03
FAJU03
FAL03
FALP03
FAMA03
FAME03
FANE03
FAOM03
FAPA03
FASA03
FAUR03
FAVE03
Star space motion:
PVSTAR
STARPV
Star catalog conversions:
FK52H
FK5HIP
FK5HZ
H2FK5
HFK5Z
STARPM
Geodetic/geocentric:
EFORM
GC2GD
GC2GDE
GD2GC
GD2GCE
Obsolete:
C2TCEO
Calls: Fortran version:
CALL iau_BI00 ( DPSIBI, DEPSBI, DRA ) CALL iau_BP00 ( DATE1, DATE2, RB, RP, RBP ) CALL iau_BP06 ( DATE1, DATE2, RB, RP, RBP ) CALL iau_BPN2XY ( RBPN, X, Y ) CALL iau_C2I00A ( DATE1, DATE2, RC2I ) CALL iau_C2I00B ( DATE1, DATE2, RC2I ) CALL iau_C2I06A ( DATE1, DATE2, RC2I ) CALL iau_C2IBPN ( DATE1, DATE2, RBPN, RC2I ) CALL iau_C2IXY ( DATE1, DATE2, X, Y, RC2I ) CALL iau_C2IXYS ( X, Y, S, RC2I ) CALL iau_C2T00A ( TTA, TTB, UTA, UTB, XP, YP, RC2T ) CALL iau_C2T00B ( TTA, TTB, UTA, UTB, XP, YP, RC2T ) CALL iau_C2T06A ( TTA, TTB, UTA, UTB, XP, YP, RC2T ) CALL iau_C2TCEO ( RC2I, ERA, RPOM, RC2T ) CALL iau_C2TCIO ( RC2I, ERA, RPOM, RC2T ) CALL iau_C2TEQX ( RBPN, GST, RPOM, RC2T ) CALL iau_C2TPE ( TTA, TTB, UTA, UTB, DPSI, DEPS, XP, YP, RC2T ) CALL iau_C2TXY ( TTA, TTB, UTA, UTB, X, Y, XP, YP, RC2T ) CALL iau_CAL2JD ( IY, IM, ID, DJM0, DJM, J ) CALL iau_D2DTF ( SCALE, NDP, D1, D2, IY, IM, ID, IHMSF, J ) CALL iau_DAT ( IY, IM, ID, FD, DELTAT, J ) D = iau_DTDB ( DATE1, DATE2, UT, ELONG, U, V ) CALL iau_DTF2D ( SCALE, IY, IM, ID, IHR, IMN, SEC, D1, D2, J ) D = iau_EE00 ( DATE1, DATE2, EPSA, DPSI ) D = iau_EE00A ( DATE1, DATE2 ) D = iau_EE00B ( DATE1, DATE2 ) D = iau_EE06A ( DATE1, DATE2 ) D = iau_EECT00 ( DATE1, DATE2 ) CALL iau_EFORM ( N, A, F, J ) D = iau_EO06A ( DATE1, DATE2 ) D = iau_EORS ( RNPB, S ) D = iau_EPB ( DJ1, DJ2 ) CALL iau_EPB2JD ( EPB, DJM0, DJM ) D = iau_EPJ ( DJ1, DJ2 ) CALL iau_EPJ2JD ( EPJ, DJM0, DJM ) CALL iau_EPV00 ( DJ1, DJ2, PVH, PVB, J ) D = iau_EQEQ94 ( DATE1, DATE2 ) D = iau_ERA00 ( DJ1, DJ2 ) D = iau_FAD03 ( T ) D = iau_FAE03 ( T ) D = iau_FAF03 ( T ) D = iau_FAJU03 ( T ) D = iau_FAL03 ( T ) D = iau_FALP03 ( T ) D = iau_FAMA03 ( T ) D = iau_FAME03 ( T ) D = iau_FANE03 ( T ) D = iau_FAOM03 ( T ) D = iau_FAPA03 ( T ) D = iau_FASA03 ( T ) D = iau_FAUR03 ( T ) D = iau_FAVE03 ( T ) CALL iau_FK52H ( R5, D5, DR5, DD5, PX5, RV5, RH, DH, DRH, DDH, PXH, RVH ) CALL iau_FK5HIP ( R5H, S5H ) CALL iau_FK5HZ ( R5, D5, DATE1, DATE2, RH, DH ) CALL iau_FW2M ( GAMB, PHIB, PSI, EPS, R ) CALL iau_FW2XY ( GAMB, PHIB, PSI, EPS, X, Y ) CALL iau_GC2GD ( N, XYZ, ELONG, PHI, HEIGHT, J ) CALL iau_GC2GDE ( A, F, XYZ, ELONG, PHI, HEIGHT, J ) CALL iau_GD2GC ( N, ELONG, PHI, HEIGHT, XYZ, J ) CALL iau_GD2GCE ( A, F, ELONG, PHI, HEIGHT, XYZ, J ) D = iau_GMST00 ( UTA, UTB, TTA, TTB ) D = iau_GMST06 ( UTA, UTB, TTA, TTB ) D = iau_GMST82 ( UTA, UTB ) D = iau_GST00A ( UTA, UTB, TTA, TTB ) D = iau_GST00B ( UTA, UTB ) D = iau_GST06 ( UTA, UTB, TTA, TTB, RNPB ) D = iau_GST06A ( UTA, UTB, TTA, TTB ) D = iau_GST94 ( UTA, UTB ) CALL iau_H2FK5 ( RH, DH, DRH, DDH, PXH, RVH, R5, D5, DR5, DD5, PX5, RV5 ) CALL iau_HFK5Z ( RH, DH, DATE1, DATE2, R5, D5, DR5, DD5 ) CALL iau_JD2CAL ( DJ1, DJ2, IY, IM, ID, FD, J ) CALL iau_JDCALF ( NDP, DJ1, DJ2, IYMDF, J ) CALL iau_NUM00A ( DATE1, DATE2, RMATN ) CALL iau_NUM00B ( DATE1, DATE2, RMATN ) CALL iau_NUM06A ( DATE1, DATE2, RMATN ) CALL iau_NUMAT ( EPSA, DPSI, DEPS, RMATN ) CALL iau_NUT00A ( DATE1, DATE2, DPSI, DEPS ) CALL iau_NUT00B ( DATE1, DATE2, DPSI, DEPS ) CALL iau_NUT06A ( DATE1, DATE2, DPSI, DEPS ) CALL iau_NUT80 ( DATE1, DATE2, DPSI, DEPS ) CALL iau_NUTM80 ( DATE1, DATE2, RMATN ) D = iau_OBL06 ( DATE1, DATE2 ) D = iau_OBL80 ( DATE1, DATE2 ) CALL iau_PB06 ( DATE1, DATE2, BZETA, BZ, BTHETA ) CALL iau_PFW06 ( DATE1, DATE2, GAMB, PHIB, PSIB, EPSA ) CALL iau_PLAN94 ( DATE1, DATE2, NP, PV, J ) CALL iau_PMAT00 ( DATE1, DATE2, RBP ) CALL iau_PMAT06 ( DATE1, DATE2, RBP ) CALL iau_PMAT76 ( DATE1, DATE2, RMATP ) CALL iau_PN00 ( DATE1, DATE2, DPSI, DEPS, EPSA, RB, RP, RBP, RN, RBPN ) CALL iau_PN00A ( DATE1, DATE2, DPSI, DEPS, EPSA, RB, RP, RBP, RN, RBPN ) CALL iau_PN00B ( DATE1, DATE2, DPSI, DEPS, EPSA, RB, RP, RBP, RN, RBPN ) CALL iau_PN06 ( DATE1, DATE2, DPSI, DEPS, EPSA, RB, RP, RBP, RN, RBPN ) CALL iau_PN06A ( DATE1, DATE2, DPSI, DEPS, RB, RP, RBP, RN, RBPN ) CALL iau_PNM00A ( DATE1, DATE2, RBPN ) CALL iau_PNM00B ( DATE1, DATE2, RBPN ) CALL iau_PNM06A ( DATE1, DATE2, RNPB ) CALL iau_PNM80 ( DATE1, DATE2, RMATPN ) CALL iau_P06E ( DATE1, DATE2, EPS0, PSIA, OMA, BPA, BQA, PIA, BPIA, EPSA, CHIA, ZA, ZETAA, THETAA, PA, GAM, PHI, PSI ) CALL iau_POM00 ( XP, YP, SP, RPOM ) CALL iau_PR00 ( DATE1, DATE2, DPSIPR, DEPSPR ) CALL iau_PREC76 ( EP01, EP02, EP11, EP12, ZETA, Z, THETA ) CALL iau_PVSTAR ( PV, RA, DEC, PMR, PMD, PX, RV, J ) D = iau_S00 ( DATE1, DATE2, X, Y ) D = iau_S00A ( DATE1, DATE2 ) D = iau_S00B ( DATE1, DATE2 ) D = iau_S06 ( DATE1, DATE2, X, Y ) D = iau_S06A ( DATE1, DATE2 ) D = iau_SP00 ( DATE1, DATE2 ) CALL iau_STARPM ( RA1, DEC1, PMR1, PMD1, PX1, RV1, EP1A, EP1B, EP2A, EP2B, RA2, DEC2, PMR2, PMD2, PX2, RV2, J ) CALL iau_STARPV ( RA, DEC, PMR, PMD, PX, RV, PV, J ) CALL iau_TAITT ( TAI1, TAI2, TT1, TT2, J ) CALL iau_TAIUT1 ( TAI1, TAI2, DTA, UT11, UT12, J ) CALL iau_TAIUTC ( TAI1, TAI2, UTC1, UTC2, J ) CALL iau_TCBTDB ( TCB1, TCB2, TDB1, TDB2, J ) CALL iau_TCGTT ( TCG1, TCG2, TT1, TT2, J ) CALL iau_TDBTCB ( TDB1, TDB2, TCB1, TCB2, J ) CALL iau_TDBTT ( TDB1, TDB2, DTR, TT1, TT2, J ) CALL iau_TTTAI ( TT1, TT2, TAI1, TAI2, J ) CALL iau_TTTCG ( TT1, TT2, TCG1, TCG2, J ) CALL iau_TTTDB ( TT1, TT2, DTR, TDB1, TDB2, J ) CALL iau_TTUT1 ( TT1, TT2, DT, UT11, UT12, J ) CALL iau_UT1TAI ( UT11, UT12, TAI1, TAI2, J ) CALL iau_UT1TT ( UT11, UT12, DT, TT1, TT2, J ) CALL iau_UT1UTC ( UT11, UT12, DUT, UTC1, UTC2, J ) CALL iau_UTCTAI ( UTC1, UTC2, DTA, TAI1, TAI2, J ) CALL iau_UTCUT1 ( UTC1, UTC2, DUT, UT11, UT12, J ) CALL iau_XY06 ( DATE1, DATE2, X, Y ) CALL iau_XYS00A ( DATE1, DATE2, X, Y, S ) CALL iau_XYS00B ( DATE1, DATE2, X, Y, S ) CALL iau_XYS06A ( DATE1, DATE2, X, Y, S )
Calls: C version:
iauBi00 ( &dpsibi, &depsbi, &dra ); iauBp00 ( date1, date2, rb, rp, rbp ); iauBp06 ( date1, date2, rb, rp, rbp ); iauBpn2xy ( rbpn, &x, &y ); iauC2i00a ( date1, date2, rc2i ); iauC2i00b ( date1, date2, rc2i ); iauC2i06a ( date1, date2, rc2i ); iauC2ibpn ( date1, date2, rbpn, rc2i ); iauC2ixy ( date1, date2, x, y, rc2i ); iauC2ixys ( x, y, s, rc2i ); iauC2t00a ( tta, ttb, uta, utb, xp, yp, rc2t ); iauC2t00b ( tta, ttb, uta, utb, xp, yp, rc2t ); iauC2t06a ( tta, ttb, uta, utb, xp, yp, rc2t ); iauC2tcio ( rc2i, era, rpom, rc2t ); iauC2teqx ( rbpn, gst, rpom, rc2t ); iauC2tpe ( tta, ttb, uta, utb, dpsi, deps, xp, yp, rc2t ); iauC2txy ( tta, ttb, uta, utb, x, y, xp, yp, rc2t ); i = iauCal2jd ( iy, im, id, &djm0, &djm ); i = iauD2dtf ( scale, ndp, d1, d2, &iy, &im, &id, ihmsf ); i = iauDat ( iy, im, id, fd, &deltat ); d = iauDtdb ( date1, date2, ut, elong, u, v ); i = iauDtf2d ( scale, iy, im, id, ihr, imn, sec, &d1, &d2 ); d = iauEe00 ( date1, date2, epsa, dpsi ); d = iauEe00a ( date1, date2 ); d = iauEe00b ( date1, date2 ); d = iauEe06 ( date1, date2 ); d = iauEect00 ( date1, date2 ); i = iauEform ( n, &a, &f ); d = iauEo06 ( date1, date2 ); d = iauEors ( rnpb, s ); d = iauEpb ( dj1, dj2 ); iauEpb2jd ( epb, &djm0, &djm ); d = iauEpj ( dj1, dj2 ); iauEpj2jd ( epj, &djm0, &djm ); i = iauEpv00 ( dj1, dj2, pvh, pvb ); d = iauEqeq94 ( date1, date2 ); d = iauEra00 ( dj1, dj2 ); d = iauFad03 ( t ); d = iauFae03 ( t ); d = iauFaf03 ( t ); d = iauFaju03 ( t ); d = iauFal03 ( t ); d = iauFalp03 ( t ); d = iauFama03 ( t ); d = iauFame03 ( t ); d = iauFane03 ( t ); d = iauFaom03 ( t ); d = iauFapa03 ( t ); d = iauFasa03 ( t ); d = iauFaur03 ( t ); d = iauFave03 ( t ); iauFk52h ( r5, d5, dr5, dd5, px5, rv5, &rh, &dh, &drh, &ddh, &pxh, &rvh ); iauFk5hip ( r5h, s5h ); iauFk5hz ( r5, d5, date1, date2, &rh, &dh ); iauFw2m ( gamb, phib, psi, eps, r ); iauFw2xy ( gamb, phib, psi, eps, &x, &y ); i = iauGc2gd ( n, xyz, &elong, &phi, &height ); i = iauGc2gde ( a, f, xyz, &elong, &phi, &height ); i = iauGd2gc ( n, elong, phi, height, xyz ); i = iauGd2gce ( a, f, elong, phi, height, xyz ); d = iauGmst00 ( uta, utb, tta, ttb ); d = iauGmst06 ( uta, utb, tta, ttb ); d = iauGmst82 ( uta, utb ); d = iauGst00a ( uta, utb, tta, ttb ); d = iauGst00b ( uta, utb ); d = iauGst06 ( uta, utb, tta, ttb, rnpb ); d = iauGst06a ( uta, utb, tta, ttb ); d = iauGst94 ( uta, utb ); iauH2fk5 ( rh, dh, drh, ddh, pxh, rvh, &r5, &d5, &dr5, &dd5, &px5, &rv5 ); iauHfk5z ( rh, dh, date1, date2, &r5, &d5, &dr5, &dd5 ); i = iauJd2cal ( dj1, dj2, &iy, &im, &id, &fd ); i = iauJdcalf ( ndp, dj1, dj2, iymdf ); iauNum00a ( date1, date2, rmatn ); iauNum00b ( date1, date2, rmatn ); iauNum06a ( date1, date2, rmatn ); iauNumat ( epsa, dpsi, deps, rmatn ); iauNut00a ( date1, date2, &dpsi, &deps ); iauNut00b ( date1, date2, &dpsi, &deps ); iauNut06a ( date1, date2, &dpsi, &deps ); iauNut80 ( date1, date2, &dpsi, &deps ); iauNutm80 ( date1, date2, rmatn ); d = iauObl06 ( date1, date2 ); d = iauObl80 ( date1, date2 ); iauPb06 ( date1, date2, &bzeta, &bz, &btheta ); iauPfw06 ( date1, date2, &gamb, &phib, &psib, &epsa ); i = iauPlan94 ( date1, date2, np, pv ); iauPmat00 ( date1, date2, rbp ); iauPmat06 ( date1, date2, rbp ); iauPmat76 ( date1, date2, rmatp ); iauPn00 ( date1, date2, dpsi, deps, &epsa, rb, rp, rbp, rn, rbpn ); iauPn00a ( date1, date2, &dpsi, &deps, &epsa, rb, rp, rbp, rn, rbpn ); iauPn00b ( date1, date2, &dpsi, &deps, &epsa, rb, rp, rbp, rn, rbpn ); iauPn06 ( date1, date2, dpsi, deps, &epsa, rb, rp, rbp, rn, rbpn ); iauPn06a ( date1, date2, &dpsi, &deps, &epsa, rb, rp, rbp, rn, rbpn ); iauPnm00a ( date1, date2, rbpn ); iauPnm00b ( date1, date2, rbpn ); iauPnm06a ( date1, date2, rnpb ); iauPnm80 ( date1, date2, rmatpn ); iauP06e ( date1, date2, &eps0, &psia, &oma, &bpa, &bqa, &pia, &bpia, &epsa, &chia, &za, &zetaa, &thetaa, &pa, &gam, &phi, &psi ); iauPom00 ( xp, yp, sp, rpom ); iauPr00 ( date1, date2, &dpsipr, &depspr ); iauPrec76 ( ep01, ep02, ep11, ep12, &zeta, &z, &theta ); i = iauPvstar ( pv, &ra, &dec, &pmr, &pmd, &px, &rv ); d = iauS00 ( date1, date2, x, y ); d = iauS00a ( date1, date2 ); d = iauS00b ( date1, date2 ); d = iauS06 ( date1, date2, x, y ); d = iauS06a ( date1, date2 ); d = iauSp00 ( date1, date2 ); i = iauStarpm ( ra1, dec1, pmr1, pmd1, px1, rv1, ep1a, ep1b, ep2a, ep2b, &ra2, &dec2, &pmr2, &pmd2, &px2, &rv2 ); i = iauStarpv ( ra, dec, pmr, pmd, px, rv, pv ); i = iauTaitt ( tai1, tai2, &tt1, &tt2 ); i = iauTaiut1 ( tai1, tai2, dta, &ut11, &ut12 ); i = iauTaiutc ( tai1, tai2, &utc1, &utc2 ); i = iauTcbtdb ( tcb1, tcb2, &tdb1, &tdb2 ); i = iauTcgtt ( tcg1, tcg2, &tt1, &tt2 ); i = iauTdbtcb ( tdb1, tdb2, &tcb1, &tcb2 ); i = iauTdbtt ( tdb1, tdb2, dtr, &tt1, &tt2 ); i = iauTttai ( tt1, tt2, &tai1, &tai2 ); i = iauTttcg ( tt1, tt2, &tcg1, &tcg2 ); i = iauTttdb ( tt1, tt2, dtr, &tdb1, &tdb2 ); i = iauTtut1 ( tt1, tt2, dt, &ut11, &ut12 ); i = iauUt1tai ( ut11, ut12, &tai1, &tai2 ); i = iauUt1tt ( ut11, ut12, dt, &tt1, &tt2 ); i = iauUt1utc ( ut11, ut12, dut, &utc1, &utc2 ); i = iauUtctai ( utc1, utc2, dta, &tai1, &tai2 ); i = iauUtcut1 ( utc1, utc2, dut, &ut11, &ut12 ); iauXy06 ( date1, date2, &x, &y ); iauXys00a ( date1, date2, &x, &y, &s ); iauXys00b ( date1, date2, &x, &y, &s ); iauXys06a ( date1, date2, &x, &y, &s );
The routines described here comprise the SOFA vector/matrix library. Their general appearance and coding style conforms to conventions agreed by the SOFA Review Board, and their functions, names and algorithms have been ratified by the Board. Procedures for soliciting and agreeing additions to the library are still evolving.
The SOFA routines are available in two programming languages at present: Fortran 77 and ANSI C.
There is a one–to–one relationship between the two language versions.
The naming convention is such that a SOFA routine referred to
generically as EXAMPL
exists as a Fortran subprogram
iau_EXAMPL
and a C function iauExampl
. The calls for the
two versions are very similar, with the same arguments in the same
order. In a few cases, the C equivalent of a Fortran SUBROUTINE
subprogram uses a return value rather than an argument.
The library consists mostly of routines which operate on ordinary
Cartesian vectors (x, y, z) and 3x3 rotation matrices. However,
there is also support for vectors which represent velocity as well as
position and vectors which represent rotation instead of position. The
vectors which represent both position and velocity may be considered
still to have dimensions (3)
, but to comprise elements each of
which is two numbers, representing the value itself and the time
derivative. Thus:
p
vectors (or just plain 3-vectors) have
dimension (3)
in Fortran and [3]
in C.
pv
vectors have dimensions (3,2)
in Fortran and [2][3]
in C.
r
matrices have dimensions (3,3)
in
Fortran and [3][3]
in C. When used for rotation, they are
“orthogonal”; the inverse of such a matrix is equal to the transpose.
Most of the routines in this library do not assume that
r
-matrices are necessarily orthogonal and in fact work on any 3x3
matrix.
r
vectors have dimensions (3)
in Fortran
and [3]
in C. Such vectors are a combination of the Euler axis
and angle and are convertible to and from r
-matrices. The
direction is the axis of rotation and the magnitude is the angle of
rotation, in radians. Because the amount of rotation can be scaled up
and down simply by multiplying the vector by a scalar, r
-vectors
are useful for representing spins about an axis which is fixed.
pv
-vector follow the three position components.
Application code is permitted to exploit this and all other knowledge of
the internal layouts: that x, y and z appear in
that order and are in a right–handed Cartesian coordinate system etc.
For example, the cp
function (copy a p
-vector) can be used
to copy the velocity component of a pv
-vector (indeed, this is
how the CPV
routine is coded).
In addition to the vector/matrix routines, the library contains some routines related to spherical angles, including conversions to and from sexagesimal format.
Using the library requires knowledge of vector/matrix methods, spherical trigonometry, and methods of attitude representation. These topics are covered in many textbooks, including “Spacecraft Attitude Determination and Control”, James R. Wertz (ed.), Astrophysics and Space Science Library, Vol. 73, D. Reidel Publishing Company, 1986.
p
-vectors and r
-matricesInitialize:
ZP
p
-vector.
ZR
r
-matrix to null.
IR
r
-matrix to identity.
Copy/extend/extract:
CP
p
-vector.
CR
r
-matrix.
Build rotations:
RX
r
-matrix about x.
RY
r
-matrix about y.
RZ
r
-matrix about z.
Spherical/Cartesian conversions:
S2C
C2S
S2P
p
-vector.
P2S
p
-vector to spherical.
Operations on vectors:
PPP
p
-vector plus p
-vector.
PMP
p
-vector minus p
-vector.
PPSP
p
-vector plus scaled p
-vector.
PDP
p
-vectors.
PXP
p
-vectors.
PM
p
-vector.
PN
p
-vector returning modulus.
SXP
p
-vector by scalar.
Operations on matrices:
RXR
r
-matrix multiply.
TR
r
-matrix.
Matrix–vector products:
RXP
r
-matrix and p
-vector.
TRXP
r
-matrix and p
-vector.
Separation and position–angle:
SEPP
p
-vectors.
SEPS
PAP
p
-vectors.
PAS
Rotation vectors:
RV2M
r
-vector to r
-matrix.
RM2V
r
-matrix to r
-vector.
pv
-vectorsInitialize:
ZPV
pv
-vector.
Copy/extend/extract:
CPV
pv
-vector.
P2PV
p
-vector.
PV2P
pv
-vector.
Spherical/Cartesian conversions:
S2PV
pv
-vector.
PV2S
pv
-vector to spherical.
Operations on vectors:
PVPPV
pv
-vector plus pv
-vector.
PVMPV
pv
-vector minus pv
-vector.
PVDPV
pv
-vectors.
PVXPV
pv
-vectors.
PVM
pv
-vector.
SXPV
pv
-vector by scalar.
S2XPV
pv
-vector by two scalars.
PVU
pv
-vector.
PVUP
pv
-vector discarding velocity.
Matrix–vector products:
RXPV
r
-matrix and pv
-vector.
TRXPV
r
-matrix and pv
-vector.
Operations on angles:
ANP
0
to 2pi
.
ANPM
-pi
to +pi
.
A2TF
A2AF
AF2A
D2TF
TF2A
TF2D
Calls: Fortran version
CALL iau_A2AF ( NDP, ANGLE, SIGN, IDMSF ) CALL iau_A2TF ( NDP, ANGLE, SIGN, IHMSF ) CALL iau_AF2A ( S, IDEG, IAMIN, ASEC, RAD, J ) D = iau_ANP ( A ) D = iau_ANPM ( A ) CALL iau_C2S ( P, THETA, PHI ) CALL iau_CP ( P, C ) CALL iau_CPV ( PV, C ) CALL iau_CR ( R, C ) CALL iau_D2TF ( NDP, DAYS, SIGN, IHMSF ) CALL iau_IR ( R ) CALL iau_P2PV ( P, PV ) CALL iau_P2S ( P, THETA, PHI, R ) CALL iau_PAP ( A, B, THETA ) CALL iau_PAS ( AL, AP, BL, BP, THETA ) CALL iau_PDP ( A, B, ADB ) CALL iau_PM ( P, R ) CALL iau_PMP ( A, B, AMB ) CALL iau_PN ( P, R, U ) CALL iau_PPP ( A, B, APB ) CALL iau_PPSP ( A, S, B, APSB ) CALL iau_PV2P ( PV, P ) CALL iau_PV2S ( PV, THETA, PHI, R, TD, PD, RD ) CALL iau_PVDPV ( A, B, ADB ) CALL iau_PVM ( PV, R, S ) CALL iau_PVMPV ( A, B, AMB ) CALL iau_PVPPV ( A, B, APB ) CALL iau_PVU ( DT, PV, UPV ) CALL iau_PVUP ( DT, PV, P ) CALL iau_PVXPV ( A, B, AXB ) CALL iau_PXP ( A, B, AXB ) CALL iau_RM2V ( R, P ) CALL iau_RV2M ( P, R ) CALL iau_RX ( PHI, R ) CALL iau_RXP ( R, P, RP ) CALL iau_RXPV ( R, PV, RPV ) CALL iau_RXR ( A, B, ATB ) CALL iau_RY ( THETA, R ) CALL iau_RZ ( PSI, R ) CALL iau_S2C ( THETA, PHI, C ) CALL iau_S2P ( THETA, PHI, R, P ) CALL iau_S2PV ( THETA, PHI, R, TD, PD, RD, PV ) CALL iau_S2XPV ( S1, S2, PV ) CALL iau_SEPP ( A, B, S ) CALL iau_SEPS ( AL, AP, BL, BP, S ) CALL iau_SXP ( S, P, SP ) CALL iau_SXPV ( S, PV, SPV ) CALL iau_TF2A ( S, IHOUR, IMIN, SEC, RAD, J ) CALL iau_TF2D ( S, IHOUR, IMIN, SEC, DAYS, J ) CALL iau_TR ( R, RT ) CALL iau_TRXP ( R, P, TRP ) CALL iau_TRXPV ( R, PV, TRPV ) CALL iau_ZP ( P ) CALL iau_ZPV ( PV ) CALL iau_ZR ( R )
Calls: C version:
iauA2af ( ndp, angle, &sign, idmsf ); iauA2tf ( ndp, angle, &sign, ihmsf ); i = iauAf2a ( s, ideg, iamin, asec, &rad ); d = iauAnp ( a ); d = iauAnpm ( a ); iauC2s ( p, &theta, &phi ); iauCp ( p, c ); iauCpv ( pv, c ); iauCr ( r, c ); iauD2tf ( ndp, days, &sign, ihmsf ); iauIr ( r ); iauP2pv ( p, pv ); iauP2s ( p, &theta, &phi, &r ); d = iauPap ( a, b ); d = iauPas ( al, ap, bl, bp ); d = iauPdp ( a, b ); d = iauPm ( p ); iauPmp ( a, b, amb ); iauPn ( p, &r, u ); iauPpp ( a, b, apb ); iauPpsp ( a, s, b, apsb ); iauPv2p ( pv, p ); iauPv2s ( pv, &theta, &phi, &r, &td, &pd, &rd ); iauPvdpv ( a, b, adb ); iauPvm ( pv, &r, &s ); iauPvmpv ( a, b, amb ); iauPvppv ( a, b, apb ); iauPvu ( dt, pv, upv ); iauPvup ( dt, pv, p ); iauPvxpv ( a, b, axb ); iauPxp ( a, b, axb ); iauRm2v ( r, p ); iauRv2m ( p, r ); iauRx ( phi, r ); iauRxp ( r, p, rp ); iauRxpv ( r, pv, rpv ); iauRxr ( a, b, atb ); iauRy ( theta, r ); iauRz ( psi, r ); iauS2c ( theta, phi, c ); iauS2p ( theta, phi, r, p ); iauS2pv ( theta, phi, r, td, pd, rd, pV ); iauS2xpv ( s1, s2, pv ); d = iauSepp ( a, b ); d = iauSeps ( al, ap, bl, bp ); iauSxp ( s, p, sp ); iauSxpv ( s, pv, spv ); i = iauTf2a ( s, ihour, imin, sec, &rad ); i = iauTf2d ( s, ihour, imin, sec, &days ); iauTr ( r, rt ); iauTrxp ( r, p, trp ); iauTrxpv ( r, pv, trpv ); iauZp ( p ); iauZpv ( pv ); iauZr ( r );
Gregorian Calendar to Julian Date.
Status: support function.
Given:
iy,im,id int year, month, day in Gregorian calendar (Note 1)Returned:
djm0 double MJD zero-point: always 2400000.5 djm double Modified Julian Date for 0 hrsReturned (function value):
int status: 0 = OK -1 = bad year (Note 3: JD not computed) -2 = bad month (JD not computed) -3 = bad day (JD computed)Notes:
- The algorithm used is valid from -4800 March 1, but this implementation rejects dates before -4799 January 1.
- The Julian Date is returned in two pieces, in the usual SOFA manner, which is designed to preserve time resolution. The Julian Date is available as a single number by adding djm0 and djm.
- In early eras the conversion is from the “Proleptic Gregorian Calendar”; no account is taken of the date(s) of adoption of the Gregorian Calendar, nor is the AD/BC numbering convention observed.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 12.92 (p604).
Julian Date to Besselian Epoch.
Status: support function.
Given:
dj1,dj2 double Julian Date (see note)Returned (function value):
double Besselian Epoch.Note:
- The Julian Date is supplied in two pieces, in the usual SOFA manner, which is designed to preserve time resolution. The Julian Date is available as a single number by adding dj1 and dj2. The maximum resolution is achieved if dj1 is 2451545D0 (J2000.0).
Reference:
- Lieske,J.H., 1979. Astron.Astrophys.,73,282.
Besselian Epoch to Julian Date.
Status: support function.
Given:
epb double Besselian Epoch (e.g. 1957.3D0)Returned:
djm0 double MJD zero-point: always 2400000.5 djm double Modified Julian DateNote:
- The Julian Date is returned in two pieces, in the usual SOFA manner, which is designed to preserve time resolution. The Julian Date is available as a single number by adding djm0 and djm.
Reference:
- Lieske, J.H., 1979, Astron.Astrophys. 73, 282.
Julian Date to Julian Epoch.
Status: support function.
Given:
dj1,dj2 double Julian Date (see note)Returned (function value):
double Julian EpochNote:
- The Julian Date is supplied in two pieces, in the usual SOFA manner, which is designed to preserve time resolution. The Julian Date is available as a single number by adding dj1 and dj2. The maximum resolution is achieved if dj1 is 2451545D0 (J2000.0).
Reference:
- Lieske, J.H., 1979, Astron.Astrophys. 73, 282.
Julian Epoch to Julian Date.
Status: support function.
Given:
epj double Julian Epoch (e.g. 1996.8D0)Returned:
djm0 double MJD zero-point: always 2400000.5 djm double Modified Julian DateNote:
- The Julian Date is returned in two pieces, in the usual SOFA manner, which is designed to preserve time resolution. The Julian Date is available as a single number by adding djm0 and djm.
Reference:
- Lieske, J.H., 1979, Astron.Astrophys. 73, 282.
Julian Date to Gregorian year, month, day, and fraction of a day.
Status: support function.
Given:
dj1,dj2 double Julian Date (Notes 1, 2)Returned (arguments):
iy int year im int month id int day fd double fraction of dayReturned (function value):
int status: 0 = OK -1 = unacceptable date (Note 3)Notes:
- The earliest valid date is -68569.5 (-4900 March 1). The largest value accepted is 10^9.
- The Julian Date is apportioned in any convenient way between the arguments dj1 and dj2. For example,
JD = 2450123.7
could be expressed in any of these ways, among others:dj1 dj2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)- In early eras the conversion is from the “proleptic Gregorian calendar”; no account is taken of the date(s) of adoption of the Gregorian calendar, nor is the AD/BC numbering convention observed.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 12.92 (p604).
Julian Date to Gregorian Calendar, expressed in a form convenient for formatting messages: rounded to a specified precision.
Status: support function.
Given:
ndp int number of decimal places of days in fraction dj1,dj2 double dj1+dj2 = Julian Date (Note 1)Returned:
iymdf int[4] year, month, day, fraction in Gregorian calendarReturned (function value):
int status: -1 = date out of range 0 = OK +1 = NDP not 0-9 (interpreted as 0)Notes:
- The Julian Date is apportioned in any convenient way between the arguments dj1 and dj2. For example,
JD = 2450123.7
could be expressed in any of these ways, among others:dj1 dj2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)- In early eras the conversion is from the “Proleptic Gregorian Calendar”; no account is taken of the date(s) of adoption of the Gregorian Calendar, nor is the AD/BC numbering convention observed.
- Refer to the function iauJd2cal.
- NDP should be 4 or less if internal overflows are to be avoided on machines which use 16-bit integers.
Called:
iauJd2cal
- JD to Gregorian calendar.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 12.92 (p604).
Earth position and velocity, heliocentric and barycentric, with respect to the Barycentric Celestial Reference System.
Status: support function.
Given:
date1,date2 double TDB date (Note 1)Returned:
pvh double[2][3] heliocentric Earth position/velocity pvb double[2][3] barycentric Earth position/velocityReturned (function value):
int status: 0 = OK +1 = warning: date outside the range 1900-2100 ADNotes:
- The TDB date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TDB) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience. However, the accuracy of the result is more likely to be limited by the algorithm itself than the way the date has been expressed.
n.b. TT can be used instead of TDB in most applications.- On return, the arrays pvh and pvb contain the following:
pvh[0][0] x } pvh[0][1] y } heliocentric position, AU pvh[0][2] z } pvh[1][0] xdot } pvh[1][1] ydot } heliocentric velocity, AU/d pvh[1][2] zdot } pvb[0][0] x } pvb[0][1] y } barycentric position, AU pvb[0][2] z } pvb[1][0] xdot } pvb[1][1] ydot } barycentric velocity, AU/d pvb[1][2] zdot }The vectors are with respect to the Barycentric Celestial Reference System. The time unit is one day in TDB.
- The function is a simplified solution from the planetary theory VSOP2000 (X. Moisson, P. Bretagnon, 2001, Celes. Mechanics & Dyn. Astron., 80, 3/4, 205-213) and is an adaptation of original Fortran code supplied by P. Bretagnon (private comm., 2000).
- Comparisons over the time span 1900-2100 with this simplified solution and the JPL DE405 ephemeris give the following results:
RMS max Heliocentric: position error 3.7 11.2 km velocity error 1.4 5.0 mm/s Barycentric: position error 4.6 13.4 km velocity error 1.4 4.9 mm/sComparisons with the JPL DE406 ephemeris show that by 1800 and 2200 the position errors are approximately double their 1900-2100 size. By 1500 and 2500 the deterioration is a factor of 10 and by 1000 and 3000 a factor of 60. The velocity accuracy falls off at about half that rate.
- It is permissible to use the same array for pvh and pvb, which will receive the barycentric values.
Approximate heliocentric position and velocity of a nominated major planet: Mercury, Venus, EMB, Mars, Jupiter, Saturn, Uranus or Neptune (but not the Earth itself).
Status: support function.
Given:
date1 double TDB date part A (Note 1) date2 double TDB date part B (Note 1) np int planet (1=Mercury, 2=Venus, 3=EMB, 4=Mars, 5=Jupiter, 6=Saturn, 7=Uranus, 8=Neptune)Returned (argument):
pv double[2][3] planet p,v (heliocentric, J2000.0, AU,AU/d)Returned (function value):
int status: -1 = illegal NP (outside 1-8) 0 = OK +1 = warning: year outside 1000-3000 +2 = warning: failed to convergeNotes:
- The date date1
+
date2 is in the TDB time scale (in practice TT can be used) and is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TDB) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience. The limited accuracy of the present algorithm is such that any of the methods is satisfactory.
- If an np value outside the range 1-8 is supplied, an error status (function value -1) is returned and the pv vector set to zeroes.
- For np
= 3
the result is for the Earth–Moon Barycenter. To obtain the heliocentric position and velocity of the Earth, use instead the SOFA functioniauEpv00
.- On successful return, the array pv contains the following:
pv[0][0] x } pv[0][1] y } heliocentric position, AU pv[0][2] z } pv[1][0] xdot } pv[1][1] ydot } heliocentric velocity, AU/d pv[1][2] zdot }The reference frame is equatorial and is with respect to the mean equator and equinox of epoch J2000.0.
- The algorithm is due to J.L. Simon, P. Bretagnon, J. Chapront, M. Chapront-Touze, G. Francou and J. Laskar (Bureau des Longitudes, Paris, France). From comparisons with JPL ephemeris DE102, they quote the following maximum errors over the interval 1800-2050:
L (arcsec) B (arcsec) R (km) Mercury 4 1 300 Venus 5 1 800 EMB 6 1 1000 Mars 17 1 7700 Jupiter 71 5 76000 Saturn 81 13 267000 Uranus 86 7 712000 Neptune 11 1 253000Over the interval 1000-3000, they report that the accuracy is no worse than 1.5 times that over 1800-2050. Outside 1000-3000 the accuracy declines.
Comparisons of the present function with the JPL DE200 ephemeris give the following RMS errors over the interval 1960-2025:
position (km) velocity (m/s) Mercury 334 0.437 Venus 1060 0.855 EMB 2010 0.815 Mars 7690 1.98 Jupiter 71700 7.70 Saturn 199000 19.4 Uranus 564000 16.4 Neptune 158000 14.4Comparisons against DE200 over the interval 1800-2100 gave the following maximum absolute differences. (The results using DE406 were essentially the same.)
L (arcsec) B (arcsec) R (km) Rdot (m/s) Mercury 7 1 500 0.7 Venus 7 1 1100 0.9 EMB 9 1 1300 1.0 Mars 26 1 9000 2.5 Jupiter 78 6 82000 8.2 Saturn 87 14 263000 24.6 Uranus 86 7 661000 27.4 Neptune 11 2 248000 21.4- The present SOFA re-implementation of the original Simon et al. Fortran code differs from the original in the following respects:
- C instead of Fortran.
- The date is supplied in two parts.
- The result is returned only in equatorial Cartesian form; the ecliptic longitude, latitude and radius vector are not returned.
- The result is in the J2000.0 equatorial frame, not ecliptic.
- More is done in–line: there are fewer calls to subroutines.
- Different error/warning status values are used.
- A different Kepler's–equation–solver is used (avoiding use of double precision complex).
- Polynomials in t are nested to minimize rounding errors.
- Explicit double constants are used to avoid mixed–mode expressions.
None of the above changes affects the result significantly.
- The returned status indicates the most serious condition encountered during execution of the function. Illegal np is considered the most serious, overriding failure to converge, which in turn takes precedence over the remote date warning.
Called:
iauAnp
- Normalize angle into range 0 to 2pi.
Reference:
Simon, J.L, Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., and Laskar, J., Astron. Astrophys. 282, 663 (1994).
Fundamental argument, IERS Conventions (2003): mean elongation of the Moon from the Sun.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double D, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and is from Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
Fundamental argument, IERS Conventions (2003): mean longitude of Earth.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Earth, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and comes from Souchay et al. (1999) after Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
Fundamental argument, IERS Conventions (2003): mean longitude of the Moon minus mean longitude of the ascending node.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double F, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and is from Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
Fundamental argument, IERS Conventions (2003): mean longitude of Jupiter.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Jupiter, radians (Note 2)Notes:
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
- The expression used is as adopted in IERS Conventions (2003) and comes from Souchay et al. (1999) after Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
Fundamental argument, IERS Conventions (2003): mean anomaly of the Moon.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double l, radians (Note 2)Notes:
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
Fundamental argument, IERS Conventions (2003): mean anomaly of the Sun.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double l', radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and is from Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
Fundamental argument, IERS Conventions (2003): mean longitude of Mars.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Mars, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and comes from Souchay et al. (1999) after Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
Fundamental argument, IERS Conventions (2003): mean longitude of Mercury.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Mercury, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and comes from Souchay et al. (1999) after Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
Fundamental argument, IERS Conventions (2003): mean longitude of Neptune.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Neptune, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and is adapted from Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
Fundamental argument, IERS Conventions (2003): mean longitude of the Moon's ascending node.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double Omega, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and is from Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
Fundamental argument, IERS Conventions (2003): general accumulated precession in longitude.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double general precession in longitude, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003). It is taken from Kinoshita & Souchay (1990) and comes originally from Lieske et al. (1977).
References:
- Kinoshita, H. and Souchay J. 1990, Celest.Mech. and Dyn.Astron. 48, 187.
- Lieske, J.H., Lederle, T., Fricke, W. & Morando, B. 1977, Astron.Astrophys. 58, 1-16.
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Fundamental argument, IERS Conventions (2003): mean longitude of Saturn.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Saturn, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and comes from Souchay et al. (1999) after Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
Fundamental argument, IERS Conventions (2003): mean longitude of Uranus.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Uranus, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and is adapted from Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
Fundamental argument, IERS Conventions (2003): mean longitude of Venus.
Status: canonical model.
Given:
t double TDB, Julian centuries since J2000.0 (Note 1)Returned (function value):
double mean longitude of Venus, radians (Note 2)Notes:
- Though t is strictly TDB, it is usually more convenient to use TT, which makes no significant difference.
- The expression used is as adopted in IERS Conventions (2003) and comes from Souchay et al. (1999) after Simon et al. (1994).
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
Frame bias components of IAU 2000 precession–nutation models (part of MHB2000 with additions).
Status: canonical model.
Returned:
dpsibi,depsbi double longitude and obliquity corrections dra double the ICRS RA of the J2000.0 mean equinoxNotes:
- The frame bias corrections in longitude and obliquity (radians) are required in order to correct for the offset between the GCRS pole and the mean J2000.0 pole. They define, with respect to the GCRS frame, a J2000.0 mean pole that is consistent with the rest of the IAU 2000A precession–nutation model.
- In addition to the displacement of the pole, the complete description of the frame bias requires also an offset in right ascension. This is not part of the IAU 2000A model, and is from Chapront et al. (2002). It is returned in radians.
- This is a supplemented implementation of one aspect of the IAU 2000A nutation model, formally adopted by the IAU General Assembly in 2000, namely MHB2000 (Mathews et al. 2002).
References:
- Chapront, J., Chapront-Touze, M. & Francou, G., Astron. Astrophys., 387, 700, 2002.
- Mathews, P.M., Herring, T.A., Buffet, B.A., “Modeling of nutation and precession New nutation series for nonrigid Earth and insights into the Earth's interior”, J.Geophys.Res., 107, B4, 2002. The MHB2000 code itself was obtained on 9th September 2002 from:
Frame bias and precession, IAU 2000.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rb double[3][3] frame bias matrix (Note 2) rp double[3][3] precession matrix (Note 3) rbp double[3][3] bias-precession matrix (Note 4)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix rb transforms vectors from GCRS to mean J2000.0 by applying frame bias.
- The matrix rp transforms vectors from J2000.0 mean equator and equinox to mean equator and equinox of date by applying precession.
- The matrix rbp transforms vectors from GCRS to mean equator and equinox of date by applying frame bias then precession. It is the product rp
*
rb.- It is permissible to reuse the same array in the returned arguments. The arrays are filled in the order given.
Called:
iauBi00
- Frame bias components, IAU 2000.
iauPr00
- IAU 2000 precession adjustments.
iauIr
- Initialize
r
-matrix to identity.iauRx
- Rotate around x-axis.
iauRy
- Rotate around y-axis.
iauRz
- Rotate around z-axis.
iauCr
- Copy
r
-matrix.iauRxr
- Product of two
r
-matrices.Reference:
- “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.
Frame bias and precession, IAU 2006.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rb double[3][3] frame bias matrix (Note 2) rp double[3][3] precession matrix (Note 3) rbp double[3][3] bias-precession matrix (Note 4)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix rb transforms vectors from GCRS to mean J2000.0 by applying frame bias.
- The matrix rp transforms vectors from mean J2000.0 to mean of date by applying precession.
- The matrix rbp transforms vectors from GCRS to mean of date by applying frame bias then precession. It is the product rp
*
rb.Called:
iauPfw06
- Bias–precession F-W angles, IAU 2006.
iauFw2m
- F-W angles to
r
-matrix.iauPmat06
- PB matrix, IAU 2006.
iauTr
- Transpose
r
-matrix.iauRxr
- Product of two
r
-matrices.References:
Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855.
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
Extract from the bias–precession–nutation matrix the x, y coordinates of the Celestial Intermediate Pole.
Status: support function.
Given:
rbpn double[3][3] celestial-to-true matrix (Note 1)Returned:
x,y double Celestial Intermediate Pole (Note 2)Notes:
- The matrix rbpn transforms vectors from GCRS to true equator (and CIO or equinox) of date, and therefore the Celestial Intermediate Pole unit vector is the bottom row of the matrix.
- The arguments x, y are components of the Celestial Intermediate Pole unit vector in the Geocentric Celestial Reference System.
Reference:
- “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003)
n.b. The celestial ephemeris origin (CEO) was renamed "celestial intermediate origin" (CIO) by IAU 2006 Resolution 2.
Form the celestial–to–intermediate matrix for a given date using the IAU 2000A precession-nutation model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rc2i double[3][3] celestial-to-intermediate matrix (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix rc2i is the first stage in the transformation from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * rc2i * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003), ERA is the Earth Rotation Angle and RPOM is the polar motion matrix.- A faster, but slightly less accurate result (about 1 mas), can be obtained by using instead the
iauC2i00b
function.Called:
iauPnm00a
- Classical NPB matrix, IAU 2000A.
iauC2ibpn
- Celestial–to–intermediate matrix, given NPB matrix.
References:
- “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003)
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate originW” (CIO) by IAU 2006 Resolution 2.- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Form the celestial–to–intermediate matrix for a given date using the IAU 2000B precession-nutation model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rc2i double[3][3] celestial-to-intermediate matrix (Note 2)Notes:
- The TT date date1+date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,
JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix rc2i is the first stage in the transformation from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * rc2i * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003), ERA is the Earth Rotation Angle and RPOM is the polar motion matrix.- The present function is faster, but slightly less accurate (about 1 mas), than the
iauC2i00a
function.Called:
iauPnm00b
- Classical NPB matrix, IAU 2000B.
iauC2ibpn
- Celestial–to–intermediate matrix, given NPB matrix.
References:
- “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession-nutation model”, Astron.Astrophys. 400, 1145-1154 (2003)
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004)
Form the celestial–to–intermediate matrix for a given date using the IAU 2006 precession and IAU 2000A nutation models.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rc2i double[3][3] celestial-to-intermediate matrix (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix rc2i is the first stage in the transformation from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * rc2i * [CRS] = RC2T * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003), ERA is the Earth Rotation Angle and RPOM is the polar motion matrix.Called:
iauPnm06a
- Classical NPB matrix, IAU 2006/2000A.
iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
iauS06
- The CIO locator s, Given X,Y, IAU 2006.
iauC2ixys
- Celestial–to–intermediate matrix, Given X,Y and s.
References:
- McCarthy, D. D., Petit, G. (eds.), 2004, IERS Conventions (2003), IERS Technical Note No. 32, BKG
Form the celestial–to–intermediate matrix for a given date given the bias–precession–nutation matrix. IAU 2000.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1) rbpn double[3][3] celestial-to-true matrix (Note 2)Returned:
rc2i double[3][3] celestial-to-intermediate matrix (Note 3)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix rbpn transforms vectors from GCRS to true equator (and CIO or equinox) of date. Only the CIP (bottom row) is used.
- The matrix rc2i is the first stage in the transformation from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * rc2i * [CRS] = RC2T * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),ERA
is the Earth Rotation Angle andRPOM
is the polar motion matrix.- Although its name does not include "00", This function is in fact specific to the IAU 2000 models.
Called:
iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
iauC2ixy
- Celestial–to–intermediate matrix, given X,Y.
References:
- “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession-nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004)
Form the celestial to intermediate–frame–of–date matrix for a given date when the CIP X,Y coordinates are known. IAU 2000.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1) x,y double Celestial Intermediate Pole (Note 2)Returned:
rc2i double[3][3] celestial-to-intermediate matrix (Note 3)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The Celestial Intermediate Pole coordinates are the x, y components of the unit vector in the Geocentric Celestial Reference System.
- The matrix rc2i is the first stage in the transformation from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * rc2i * [CRS] = RC2T * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),ERA
is the Earth Rotation Angle andRPOM
is the polar motion matrix.- Although its name does not include
00
, this function is in fact specific to the IAU 2000 models.Called:
iauC2ixys
- Celestial–to–intermediate matrix, given X,Y and s.
iauS00
- The CIO locator s, given X,Y, IAU 2000A.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Form the celestial to intermediate–frame–of–date matrix given the CIP x, y and the CIO locator s.
Status: support function.
Given:
x,y double Celestial Intermediate Pole (Note 1) s double the CIO locator s (Note 2)Returned:
rc2i double[3][3] celestial-to-intermediate matrix (Note 3)Notes:
- The Celestial Intermediate Pole coordinates are the x, y components of the unit vector in the Geocentric Celestial Reference System.
- The CIO locator s (in radians) positions the Celestial Intermediate Origin on the equator of the CIP.
- The matrix rc2i is the first stage in the transformation from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * rc2i * [CRS] = RC2T * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),ERA
is the Earth Rotation Angle andRPOM
is the polar motion matrix.Called:
iauIr
- Initialize
r
-matrix to identity.iauRz
- Rotate around Z-axis.
iauRy
- Rotate around Y-axis.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004)
Form the celestial to terrestrial matrix given the date, the UT1 and the polar motion, using the IAU 2000A nutation model.
Status: support function.
Given:
tta,ttb double TT as a 2-part Julian Date (Note 1) uta,utb double UT1 as a 2-part Julian Date (Note 1) xp,yp double coordinates of the pole (radians, Note 2)Returned:
rc2t double[3][3] celestial-to-terrestrial matrix (Note 3)Notes:
- The TT and UT1 dates tta
+
ttb and uta+
utb are Julian Dates, apportioned in any convenient way between the arguments uta and utb. For example,JD(UT1) = 2450123.7
could be expressed in any of these ways, among others:uta utb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. In the case of uta, utb, the date & time method is best matched to the Earth rotation angle algorithm used: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- The arguments xp and yp are the coordinates (in radians) of the Celestial Intermediate Pole with respect to the International Terrestrial Reference System (see IERS Conventions 2003), measured along the meridians to 0 and 90 deg west respectively.
- The matrix rc2t transforms from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * RC2I * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),RC2I
is the celestial–to–intermediate matrix,ERA
is the Earth rotation angle andRPOM
is the polar motion matrix.- A faster, but slightly less accurate result (about 1 mas), can be obtained by using instead the
iauC2t00b
function.Called:
iauC2i00a
- Celestial–to–intermediate matrix, IAU 2000A.
iauEra00
- Earth rotation angle, IAU 2000.
iauSp00
- The TIO locator s', IERS 2000.
iauPom00
- Polar motion matrix.
iauC2tcio
- Form CIO-based celestial–to–terrestrial matrix.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004)
Form the celestial to terrestrial matrix given the date, the UT1 and the polar motion, using the IAU 2000B nutation model.
Status: support function.
Given:
tta,ttb double TT as a 2-part Julian Date (Note 1) uta,utb double UT1 as a 2-part Julian Date (Note 1) xp,yp double coordinates of the pole (radians, Note 2)Returned:
rc2t double[3][3] celestial-to-terrestrial matrix (Note 3)Notes:
- The TT and UT1 dates tta
+
ttb and uta+
utb are Julian Dates, apportioned in any convenient way between the arguments uta and utb. For example,JD(UT1) = 2450123.7
could be expressed in any of these ways, among others:uta utb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. In the case of uta, utb, the date & time method is best matched to the Earth rotation angle algorithm used: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- The arguments xp and yp are the coordinates (in radians) of the Celestial Intermediate Pole with respect to the International Terrestrial Reference System (see IERS Conventions 2003), measured along the meridians to 0 and 90 deg west respectively.
- The matrix rc2t transforms from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * RC2I * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),RC2I
is the celestial–to–intermediate matrix,ERA
is the Earth rotation angle andRPOM
is the polar motion matrix.- The present function is faster, but slightly less accurate (about 1 mas), than the
iauC2t00a
function.Called:
iauC2i00b
- Celestial–to–intermediate matrix, IAU 2000B.
iauEra00
- Earth rotation angle, IAU 2000.
iauPom00
- Polar motion matrix.
iauC2tcio
- Form CIO-based celestial–to–terrestrial matrix.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004)
Form the celestial to terrestrial matrix given the date, the UT1 and the polar motion, using the IAU 2006 precession and IAU 2000A nutation models.
Status: support function.
Given:
tta,ttb double TT as a 2-part Julian Date (Note 1) uta,utb double UT1 as a 2-part Julian Date (Note 1) xp,yp double coordinates of the pole (radians, Note 2)Returned:
rc2t double[3][3] celestial-to-terrestrial matrix (Note 3)Notes:
- The TT and UT1 dates tta
+
ttb and uta+
utb are Julian Dates, apportioned in any convenient way between the arguments uta and utb. For example,JD(UT1) = 2450123.7
could be expressed in any of these ways, among others:uta utb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. In the case of uta, utb, the date & time method is best matched to the Earth rotation angle algorithm used: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- The arguments xp and yp are the coordinates (in radians) of the Celestial Intermediate Pole with respect to the International Terrestrial Reference System (see IERS Conventions 2003), measured along the meridians to 0 and 90 deg west respectively.
- The matrix rc2t transforms from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * RC2I * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),RC2I
is the celestial–to–intermediate matrix,ERA
is the Earth rotation angle andRPOM
is the polar motion matrix.Called:
iauC2i06a
- Celestial–to–intermediate matrix, IAU 2006/2000A.
iauEra00
- Earth rotation angle, IAU 2000.
iauSp00
- The TIO locator s', IERS 2000.
iauPom00
- Polar motion matrix.
iauC2tcio
- Form CIO-based celestial–to–terrestrial matrix.
Reference:
- McCarthy, D. D., Petit, G. (eds.), 2004, IERS Conventions (2003), IERS Technical Note No. 32, BKG
Assemble the celestial to terrestrial matrix from CIO-based components (the celestial–to–intermediate matrix, the Earth Rotation Angle and the polar motion matrix).
Status: support function.
Given:
rc2i double[3][3] celestial-to-intermediate matrix era double Earth rotation angle rpom double[3][3] polar-motion matrixReturned:
rc2t double[3][3] celestial-to-terrestrial matrixNotes:
- This function constructs the rotation matrix that transforms vectors in the celestial system into vectors in the terrestrial system. It does so starting from precomputed components, namely the matrix which rotates from celestial coordinates to the intermediate frame, the Earth rotation angle and the polar motion matrix. One use of the present function is when generating a series of celestial–to–terrestrial matrices where only the Earth Rotation Angle changes, avoiding the considerable overhead of recomputing the precession–nutation more often than necessary to achieve given accuracy objectives.
- The relationship between the arguments is as follows:
[TRS] = RPOM * R_3(ERA) * rc2i * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003).Called:
iauCr
- Copy
r
-matrix.iauRz
- Rotate around Z-axis.
iauRxr
- Product of two
r
-matrices.Reference:
- McCarthy, D. D., Petit, G. (eds.), 2004, IERS Conventions (2003), IERS Technical Note No. 32, BKG
Assemble the celestial to terrestrial matrix from equinox-based components (the celestial–to–true matrix, the Greenwich Apparent Sidereal Time and the polar motion matrix).
Status: support function.
Given:
rbpn double[3][3] celestial-to-true matrix gst double Greenwich (apparent) Sidereal Time rpom double[3][3] polar-motion matrixReturned:
rc2t double[3][3] celestial-to-terrestrial matrix (Note 2)Notes:
- This function constructs the rotation matrix that transforms vectors in the celestial system into vectors in the terrestrial system. It does so starting from precomputed components, namely the matrix which rotates from celestial coordinates to the true equator and equinox of date, the Greenwich Apparent Sidereal Time and the polar motion matrix. One use of the present function is when generating a series of celestial–to–terrestrial matrices where only the Sidereal Time changes, avoiding the considerable overhead of recomputing the precession-nutation more often than necessary to achieve given accuracy objectives.
- The relationship between the arguments is as follows:
[TRS] = rpom * R_3(gst) * rbpn * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003).Called:
iauCr
- Copy
r
-matrix.iauRz
- Rotate around Z-axis.
iauRxr
- Product of two
r
-matrices.Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004)
Form the celestial to terrestrial matrix given the date, the UT1, the nutation and the polar motion. IAU 2000.
Status: support function.
Given:
tta,ttb double TT as a 2-part Julian Date (Note 1) uta,utb double UT1 as a 2-part Julian Date (Note 1) dpsi,deps double nutation (Note 2) xp,yp double coordinates of the pole (radians, Note 3)Returned:
rc2t double[3][3] celestial-to-terrestrial matrix (Note 4)Notes:
- The TT and UT1 dates tta
+
ttb and uta+
utb are Julian Dates, apportioned in any convenient way between the arguments uta and utb. For example,JD(UT1) = 2450123.7
could be expressed in any of these ways, among others:uta utb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. In the case of uta, utb, the date & time method is best matched to the Earth rotation angle algorithm used: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- The caller is responsible for providing the nutation components; they are in longitude and obliquity, in radians and are with respect to the equinox and ecliptic of date. For high–accuracy applications, free core nutation should be included as well as any other relevant corrections to the position of the CIP.
- The arguments xp and yp are the coordinates (in radians) of the Celestial Intermediate Pole with respect to the International Terrestrial Reference System (see IERS Conventions 2003), measured along the meridians to 0 and 90 deg west respectively.
- The matrix rc2t transforms from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(GST) * RBPN * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),RBPN
is the bias–precession–nutation matrix,GST
is the Greenwich (apparent) Sidereal Time andRPOM
is the polar motion matrix.- Although its name does not include
00
, this function is in fact specific to the IAU 2000 models.Called:
iauPn00
- Bias/precession/nutation results, IAU 2000.
iauGmst00
- Greenwich mean sidereal time, IAU 2000.
iauSp00
- The TIO locator s', IERS 2000.
iauEe00
- Equation of the equinoxes, IAU 2000.
iauPom00
- Polar motion matrix.
iauC2teqx
- Form equinox–based celestial–to–terrestrial matrix.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Form the celestial to terrestrial matrix given the date, the UT1, the CIP coordinates and the polar motion. IAU 2000.
Status: support function.
Given:
tta,ttb double TT as a 2-part Julian Date (Note 1) uta,utb double UT1 as a 2-part Julian Date (Note 1) x,y double Celestial Intermediate Pole (Note 2) xp,yp double coordinates of the pole (radians, Note 3)Returned:
rc2t double[3][3] celestial-to-terrestrial matrix (Note 4)Notes:
- The TT and UT1 dates tta
+
ttb and uta+
utb are Julian Dates, apportioned in any convenient way between the arguments uta and utb. For example,JD(UT1) = 2450123.7
could be expressed in any o these ways, among others:uta utb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. In the case of uta, utb, the date & time method is best matched to the Earth rotation angle algorithm used: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- The Celestial Intermediate Pole coordinates are the x, y components of the unit vector in the Geocentric Celestial Reference System.
- The arguments xp and yp are the coordinates (in radians) of the Celestial Intermediate Pole with respect to the International Terrestrial Reference System (see IERS Conventions 2003), measured along the meridians to 0 and 90 deg west respectively.
- The matrix rc2t transforms from celestial to terrestrial coordinates:
[TRS] = RPOM * R_3(ERA) * RC2I * [CRS] = rc2t * [CRS]where
[CRS]
is a vector in the Geocentric Celestial Reference System and[TRS]
is a vector in the International Terrestrial Reference System (see IERS Conventions 2003),ERA
is the Earth Rotation Angle andRPOM
is the polar motion matrix.- Although its name does not include
00
, this function is in fact specific to the IAU 2000 models.Called:
iauC2ixy
- Celestial–to–intermediate matrix, given X, Y.
iauEra00
- Earth rotation angle, IAU 2000.
iauSp00
- The TIO locator s', IERS 2000.
iauPom00
- Polar motion matrix.
iauC2tcio
- Form CIO-based celestial–to–terrestrial matrix.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Equation of the origins, IAU 2006 precession and IAU 2000A nutation.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double equation of the origins in radiansNotes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The equation of the origins is the distance between the true equinox and the celestial intermediate origin and, equivalently, the difference between Earth rotation angle and Greenwich apparent sidereal time (ERA-GST). It comprises the precession (since J2000.0) in right ascension plus the equation of the equinoxes (including the small correction terms).
Called:
iauPnm06a
- Classical NPB matrix, IAU 2006/2000A.
iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
iauS06
- The CIO locator s, given X,Y, IAU 2006.
iauEors
- Equation of the origins, Given NPB matrix and s.
References:
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
Equation of the origins, given the classical NPB matrix and the quantity s.
Status: support function.
Given:
rnpb double[3][3] classical nutation x precession x bias matrix s double the quantity s (the CIO locator)Returned (function value):
double the equation of the origins in radians.Notes:
- The equation of the origins is the distance between the true equinox and the celestial intermediate origin and, equivalently, the difference between Earth rotation angle and Greenwich apparent sidereal time (ERA-GST). It comprises the precession (since J2000.0) in right ascension plus the equation of the equinoxes (including the small correction terms).
- The algorithm is from Wallace & Capitaine (2006).
References:
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855
- Wallace, P. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
Form rotation matrix given the Fukushima–Williams angles.
Status: support function.
Given:
gamb double F-W angle gamma_bar (radians) phib double F-W angle phi_bar (radians) psi double F-W angle psi (radians) eps double F-W angle epsilon (radians)Returned:
r double[3][3] rotation matrixNotes:
- Naming the following points:
e = J2000.0 ecliptic pole, p = GCRS pole, E = ecliptic pole of date, P = CIP,the four Fukushima–Williams angles are as follows:
gamb = gamma = epE phib = phi = pE psi = psi = pEP eps = epsilon = EP- The matrix representing the combined effects of frame bias, precession and nutation is:
NxPxB = R_1(-eps).R_3(-psi).R_1(phib).R_3(gamb)- Three different matrices can be constructed, depending on the supplied angles:
- To obtain the nutation x precession x frame bias matrix, generate the four precession angles, generate the nutation components and add them to the psi_bar and epsilon_A angles, and call the present function.
- To obtain the precession x frame bias matrix, generate the four precession angles and call the present function.
- To obtain the frame bias matrix, generate the four precession angles for date J2000.0 and call the present function.
The nutation–only and precession–only matrices can if necessary be obtained by combining these three appropriately.
Called:
iauIr
- Initialize
r
-matrix to identity.iauRz
- Rotate around Z-axis.
iauRx
- Rotate around X-axis.
Reference:
- Hilton, J. et al., 2006, Celest.Mech.Dyn.Astron. 94, 351.
CIP X,Y given Fukushima–Williams bias–precession–nutation angles.
Status: support function.
Given:
gamb double F-W angle gamma_bar (radians) phib double F-W angle phi_bar (radians) psi double F-W angle psi (radians) eps double F-W angle epsilon (radians)Returned:
x,y double CIP X,Y ("radians")Notes:
- Naming the following points:
e = J2000.0 ecliptic pole, p = GCRS pole E = ecliptic pole of date, P = CIP,- the four Fukushima–Williams angles are as follows:
gamb = gamma = epE phib = phi = pE psi = psi = pEP eps = epsilon = EP- The matrix representing the combined effects of frame bias, precession and nutation is:
NxPxB = R_1(-epsA).R_3(-psi).R_1(phib).R_3(gamb)X,Y are elements (3,1) and (3,2) of the matrix.
Called:
iauFw2m
- F-W angles to
r
-matrix.iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
Reference:
- Hilton, J. et al., 2006, Celest.Mech.Dyn.Astron. 94, 351.
Form the matrix of nutation for a given date, IAU 2000A model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rmatn double[3][3] nutation matrixNotes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(true) = rmatn * V(mean)
, where thep
-vectorV(true)
is with respect to the true equatorial triad of date and thep
-vectorV(mean)
is with respect to the mean equatorial triad of date.- A faster, but slightly less accurate result (about 1 mas), can be obtained by using instead the iauNum00b function.
Called:
iauPn00a
- Bias/precession/nutation, IAU 2000A.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 3.222-3 (p114).
Form the matrix of nutation for a given date, IAU 2000B model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rmatn double[3][3] nutation matrixNotes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(true) = rmatn * V(mean)
, where thep
-vectorV(true)
is with respect to the true equatorial triad of date and thep
-vectorV(mean)
is with respect to the mean equatorial triad of date.- The present function is faster, but slightly less accurate (about 1 mas), than the iauNum00a function.
Called:
iauPn00b
- Bias/precession/nutation, IAU 2000B.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 3.222-3 (p114).
Form the matrix of nutation for a given date, IAU 2006/2000A model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rmatn double[3][3] nutation matrixNotes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(true) = rmatn * V(mean)
, where thep
-vectorV(true)
is with respect to the true equatorial triad of date and thep
-vectorV(mean)
is with respect to the mean equatorial triad of date.Called:
iauObl06
- Mean obliquity, IAU 2006.
iauNut06a
- Nutation, IAU 2006/2000A.
iauNumat
- Form nutation matrix.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 3.222-3 (p114).
Form the matrix of nutation.
Status: support function.
Given:
epsa double mean obliquity of date (Note 1) dpsi,deps double nutation (Note 2)Returned:
rmatn double[3][3] nutation matrix (Note 3)Notes:
- The supplied mean obliquity epsa, must be consistent with the precession–nutation models from which dpsi and deps were obtained.
- The caller is responsible for providing the nutation components; they are in longitude and obliquity, in radians and are with respect to the equinox and ecliptic of date.
- The matrix operates in the sense
V(true) = rmatn * V(mean)
, where thep
-vectorV(true)
is with respect to the true equatorial triad of date and thep
-vectorV(mean)
is with respect to the mean equatorial triad of date.Called:
iauIr
- Initialize
r
-matrix to identity.iauRx
- Rotate around X-axis.
iauRz
- Rotate around Z-axis.
Reference:
Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 3.222-3 (p114).
Nutation, IAU 2000A model (MHB2000 luni-solar and planetary nutation with free core nutation omitted).
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsi,deps double nutation, luni-solar + planetary (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The nutation components in longitude and obliquity are in radians and with respect to the equinox and ecliptic of date. The obliquity at J2000.0 is assumed to be the Lieske et al. (1977) value of 84381.448 arcsec.
Both the luni–solar and planetary nutations are included. The latter are due to direct planetary nutations and the perturbations of the lunar and terrestrial orbits.
- The function computes the MHB2000 nutation series with the associated corrections for planetary nutations. It is an implementation of the nutation part of the IAU 2000A precession–nutation model, formally adopted by the IAU General Assembly in 2000, namely MHB2000 (Mathews et al. 2002), but with the free core nutation (FCN — see Note 4) omitted.
- The full MHB2000 model also contains contributions to the nutations in longitude and obliquity due to the free-excitation of the free–core–nutation during the period 1979-2000. These FCN terms, which are time–dependent and unpredictable, are NOT included in the present function and, if required, must be independently computed. With the FCN corrections included, the present function delivers a pole which is at current epochs accurate to a few hundred microarcseconds. The omission of FCN introduces further errors of about that size.
- The present function provides classical nutation. The MHB2000 algorithm, from which it is adapted, deals also with (i) the offsets between the GCRS and mean poles and (ii) the adjustments in longitude and obliquity due to the changed precession rates. These additional functions, namely frame bias and precession adjustments, are supported by the SOFA functions iauBi00 and iauPr00.
- The MHB2000 algorithm also provides “total” nutations, comprising the arithmetic sum of the frame bias, precession adjustments, luni-solar nutation and planetary nutation. These total nutations can be used in combination with an existing IAU 1976 precession implementation, such as iauPmat76, to deliver GCRS- to-true predictions of sub-mas accuracy at current dates. However, there are three shortcomings in the MHB2000 model that must be taken into account if more accurate or definitive results are required (see Wallace 2002):
- The MHB2000 total nutations are simply arithmetic sums, yet in reality the various components are successive Euler rotations. This slight lack of rigor leads to cross terms that exceed 1 mas after a century. The rigorous procedure is to form the GCRS-to-true rotation matrix by applying the bias, precession and nutation in that order.
- Although the precession adjustments are stated to be with respect to Lieske et al. (1977), the MHB2000 model does not specify which set of Euler angles are to be used and how the adjustments are to be applied. The most literal and straightforward procedure is to adopt the 4-rotation epsilon_0, psi_A, omega_A, xi_A option, and to add DPSIPR to psi_A and DEPSPR to both omega_A and eps_A.
- The MHB2000 model predates the determination by Chapront et al. (2002) of a 14.6 mas displacement between the J2000.0 mean equinox and the origin of the ICRS frame. It should, however, be noted that neglecting this displacement when calculating star coordinates does not lead to a 14.6 mas change in right ascension, only a small second- order distortion in the pattern of the precession-nutation effect.
For these reasons, the SOFA functions do not generate the “total nutations” directly, though they can of course easily be generated by calling iauBi00, iauPr00 and the present function and adding the results.
- The MHB2000 model contains 41 instances where the same frequency appears multiple times, of which 38 are duplicates and three are triplicates. To keep the present code close to the original MHB algorithm, this small inefficiency has not been corrected.
Called:
iauFal03
- Mean anomaly of the Moon.
iauFaf03
- Mean argument of the latitude of the Moon.
iauFaom03
- Mean longitude of the Moon's ascending node.
iauFame03
- Mean longitude of Mercury.
iauFave03
- Mean longitude of Venus.
iauFae03
- Mean longitude of Earth.
iauFama03
- Mean longitude of Mars.
iauFaju03
- Mean longitude of Jupiter.
iauFasa03
- Mean longitude of Saturn.
iauFaur03
- Mean longitude of Uranus.
iauFapa03
- General accumulated precession in longitude.
References:
- Chapront, J., Chapront-Touze, M. & Francou, G. 2002, Astron.Astrophys. 387, 700
- Lieske, J.H., Lederle, T., Fricke, W. & Morando, B. 1977, Astron.Astrophys. 58, 1-16
- Mathews, P.M., Herring, T.A., Buffet, B.A. 2002, J.Geophys.Res. 107, B4. The MHB_2000 code itself was obtained on 9th September 2002 from:
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
- Wallace, P.T., “Software for Implementing the IAU 2000 Resolutions”, in IERS Workshop 5.1 (2002).
Nutation, IAU 2000B model.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsi,deps double nutation, luni-solar + planetary (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The nutation components in longitude and obliquity are in radians and with respect to the equinox and ecliptic of date. The obliquity at J2000.0 is assumed to be the Lieske et al. (1977) value of 84381.448 arcsec. (The errors that result from using this function with the IAU 2006 value of 84381.406 arcsec can be neglected.)
The nutation model consists only of luni–solar terms, but includes also a fixed offset which compensates for certain long–period planetary terms (Note 7).
- This function is an implementation of the IAU 2000B abridged nutation model formally adopted by the IAU General Assembly in 2000. The function computes the MHB_2000_SHORT luni–solar nutation series (Luzum 2001), but without the associated corrections for the precession rate adjustments and the offset between the GCRS and J2000.0 mean poles.
- The full IAU 2000A (MHB2000) nutation model contains nearly 1400 terms. The IAU 2000B model (McCarthy & Luzum 2003) contains only 77 terms, plus additional simplifications, yet still delivers results of 1 mas accuracy at present epochs. This combination of accuracy and size makes the IAU 2000B abridged nutation model suitable for most practical applications.
The function delivers a pole accurate to 1 mas from 1900 to 2100 (usually better than 1 mas, very occasionally just outside 1 mas). The full IAU 2000A model, which is implemented in the function iauNut00a (q.v.), delivers considerably greater accuracy at current dates; however, to realize this improved accuracy, corrections for the essentially unpredictable free-core-nutation (FCN) must also be included.
- The present function provides classical nutation. The MHB_2000_SHORT algorithm, from which it is adapted, deals also with (i) the offsets between the GCRS and mean poles and (ii) the adjustments in longitude and obliquity due to the changed precession rates. These additional functions, namely frame bias and precession adjustments, are supported by the SOFA functions iauBi00 and iauPr00.
- The MHB_2000_SHORT algorithm also provides “total” nutations, comprising the arithmetic sum of the frame bias, precession adjustments, and nutation (luni-solar + planetary). These total nutations can be used in combination with an existing IAU 1976 precession implementation, such as
iauPmat76
, to deliver GCRS–to–true predictions of mas accuracy at current epochs. However, for symmetry with the iauNut00a function (q.v. for the reasons), the SOFA functions do not generate the “total nutations” directly. Should they be required, they could of course easily be generated by callingiauBi00
,iauPr00
and the present function and adding the results.- The IAU 2000B model includes “planetary bias” terms that are fixed in size but compensate for long-period nutations. The amplitudes quoted in McCarthy & Luzum (2003), namely
Dpsi = -1.5835 mas
andDepsilon = +1.6339 mas
, are optimized for the “total nutations” method described in Note 6. The Luzum (2001) values used in this SOFA implementation, namely -0.135 mas and +0.388 mas, are optimized for the “rigorous” method, where frame bias, precession and nutation are applied separately and in that order. During the interval 1995-2050, the SOFA implementation delivers a maximum error of 1.001 mas (not including FCN).References:
- Lieske, J.H., Lederle, T., Fricke, W., Morando, B., “Expressions for the precession quantities based upon the IAU /1976/ system of astronomical constants”, Astron.Astrophys. 58, 1-2, 1-16. (1977).
- Luzum, B., private communication, 2001 (Fortran code MHB_2000_SHORT).
- McCarthy, D.D. & Luzum, B.J., “An abridged model of the precession–nutation of the celestial pole”, Cel.Mech.Dyn.Astron. 85, 37-49 (2003).
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J., Astron.Astrophys. 282, 663-683 (1994).
IAU 2000A nutation with adjustments to match the IAU 2006 precession.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsi,deps double nutation, luni-solar + planetary (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The nutation components in longitude and obliquity are in radians and with respect to the mean equinox and ecliptic of date, IAU 2006 precession model (Hilton et al. 2006, Capitaine et al. 2005).
- The function first computes the IAU 2000A nutation, then applies adjustments for (i) the consequences of the change in obliquity from the IAU 1980 ecliptic to the IAU 2006 ecliptic and (ii) the secular variation in the Earth's dynamical form factor J2.
- The present function provides classical nutation, complementing the IAU 2000 frame bias and IAU 2006 precession. It delivers a pole which is at current epochs accurate to a few tens of microarcseconds, apart from the free core nutation.
Called:
iauNut00a
- Nutation, IAU 2000A.
References:
- Chapront, J., Chapront-Touze, M. & Francou, G. 2002, Astron.Astrophys. 387, 700
- Lieske, J.H., Lederle, T., Fricke, W. & Morando, B. 1977, Astron.Astrophys. 58, 1-16
- Mathews, P.M., Herring, T.A., Buffet, B.A. 2002, J.Geophys.Res. 107, B4. The MHB_2000 code itself was obtained on 9th September 2002 from:
- Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G., Laskar, J. 1994, Astron.Astrophys. 282, 663-683.
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M. 1999, Astron.Astrophys.Supp.Ser. 135, 111.
- Wallace, P.T., “Software for Implementing the IAU 2000 Resolutions”, in IERS Workshop 5.1 (2002).
Nutation, IAU 1980 model.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsi double nutation in longitude (radians) deps double nutation in obliquity (radians)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The nutation components are with respect to the ecliptic of date.
Called:
iauAnpm
- Normalize angle into range +/- pi.
Reference:
Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 3.222 (p111).
Form the matrix of nutation for a given date, IAU 1980 model.
Status: support function.
Given:
date1,date2 double TDB date (Note 1)Returned:
rmatn double[3][3] nutation matrixNotes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(true) = rmatn * V(mean)
, where thep
-vectorV(true)
is with respect to the true equatorial triad of date and thep
-vectorV(mean)
is with respect to the mean equatorial triad of date.Called:
iauNut80
- Nutation, IAU 1980.
iauObl80
- Mean obliquity, IAU 1980.
iauNumat
- Form nutation matrix.
Mean obliquity of the ecliptic, IAU 2006 precession model.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double obliquity of the ecliptic (radians, Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The result is the angle between the ecliptic and mean equator of date date1
+
date2.Reference:
- Hilton, J. et al., 2006, Celest.Mech.Dyn.Astron. 94, 351.
Mean obliquity of the ecliptic, IAU 1980 model.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double obliquity of the ecliptic (radians, Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The result is the angle between the ecliptic and mean equator of date date1
+
date2.Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Expression 3.222-1 (p114).
Precession angles, IAU 2006, equinox based.
Status: canonical models.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (see Note 2):
eps0 double epsilon_0 psia double psi_A oma double omega_A bpa double P_A bqa double Q_A pia double pi_A bpia double Pi_A epsa double obliquity epsilon_A chia double chi_A za double z_A zetaa double zeta_A thetaa double theta_A pa double p_A gam double F-W angle gamma_J2000 phi double F-W angle phi_J2000 psi double F-W angle psi_J2000Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- This function returns the set of equinox based angles for the Capitaine et al. “P03” precession theory, adopted by the IAU in 2006. The angles are set out in Table 1 of Hilton et al. (2006):
eps0 epsilon_0 obliquity at J2000.0 psia psi_A luni-solar precession oma omega_A inclination of equator wrt J2000.0 ecliptic bpa P_A ecliptic pole x, J2000.0 ecliptic triad bqa Q_A ecliptic pole -y, J2000.0 ecliptic triad pia pi_A angle between moving and J2000.0 ecliptics bpia Pi_A longitude of ascending node of the ecliptic epsa epsilon_A obliquity of the ecliptic chia chi_A planetary precession za z_A equatorial precession: -3rd 323 Euler angle zetaa zeta_A equatorial precession: -1st 323 Euler angle thetaa theta_A equatorial precession: 2nd 323 Euler angle pa p_A general precession gam gamma_J2000 J2000.0 RA difference of ecliptic poles phi phi_J2000 J2000.0 codeclination of ecliptic pole psi psi_J2000 longitude difference of equator poles, J2000.0The returned values are all radians.
- Hilton et al. (2006) Table 1 also contains angles that depend on models distinct from the P03 precession theory itself, namely the IAU 2000A frame bias and nutation. The quoted polynomials are used in other SOFA functions:
iauXy06
- Contains the polynomial parts of the X and Y series.
iauS06
- Contains the polynomial part of the s+XY/2 series.
iauPfw06
- Implements the series for the Fukushima–Williams angles that are with respect to the GCRS pole (i.e. the variants that include frame bias).
- The IAU resolution stipulated that the choice of parameterization was left to the user, and so an IAU compliant precession implementation can be constructed using various combinations of the angles returned by the present function.
- The parameterization used by SOFA is the version of the Fukushima- Williams angles that refers directly to the GCRS pole. These angles may be calculated by calling the function iauPfw06. SOFA also supports the direct computation of the CIP GCRS X,Y by series, available by calling iauXy06.
- The agreement between the different parameterizations is at the 1 microarcsecond level in the present era.
- When constructing a precession formulation that refers to the GCRS pole rather than the dynamical pole, it may (depending on the choice of angles) be necessary to introduce the frame bias explicitly.
- It is permissible to re-use the same variable in the returned arguments. The quantities are stored in the stated order.
Called:
iauObl06
- Mean obliquity, IAU 2006.
Reference:
- Hilton, J. et al., 2006, Celest.Mech.Dyn.Astron. 94, 351.
This function forms three Euler angles which implement general precession from epoch J2000.0, using the IAU 2006 model. Frame bias (the offset between ICRS and mean J2000.0) is included.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
bzeta double 1st rotation: radians cw around z bz double 3rd rotation: radians cw around z btheta double 2nd rotation: radians ccw around yNotes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The traditional accumulated precession angles zeta_A, z_A, theta_A cannot be obtained in the usual way, namely through polynomial expressions, because of the frame bias. The latter means that two of the angles undergo rapid changes near this date. They are instead the results of decomposing the precession–bias matrix obtained by using the Fukushima–Williams method, which does not suffer from the problem. The decomposition returns values which can be used in the conventional formulation and which include frame bias.
- The three angles are returned in the conventional order, which is not the same as the order of the corresponding Euler rotations. The precession–bias matrix is:
R_3(-z) x R_2(+theta) x R_3(-zeta).- Should zeta_A, z_A, theta_A angles be required that do not contain frame bias, they are available by calling the SOFA function iauP06e.
Called:
iauPmat06
- PB matrix, IAU 2006.
iauRz
- Rotate around Z-axis.
Precession angles, IAU 2006 (Fukushima–Williams 4-angle formulation).
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
gamb double F-W angle gamma_bar (radians) phib double F-W angle phi_bar (radians) psib double F-W angle psi_bar (radians) epsa double F-W angle epsilon_A (radians)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- Naming the following points:
e = J2000.0 ecliptic pole, p = GCRS pole, E = mean ecliptic pole of date, P = mean pole of date,the four Fukushima–Williams angles are as follows:
gamb = gamma_bar = epE phib = phi_bar = pE psib = psi_bar = pEP epsa = epsilon_A = EP- The matrix representing the combined effects of frame bias and precession is:
PxB = R_1(-epsa).R_3(-psib).R_1(phib).R_3(gamb)- The matrix representing the combined effects of frame bias, precession and nutation is simply:
NxPxB = R_1(-epsa-dE).R_3(-psib-dP).R_1(phib).R_3(gamb)where dP and dE are the nutation components with respect to the ecliptic of date.
Called:
iauObl06
- Mean obliquity, IAU 2006.
Reference:
- Hilton, J. et al., 2006, Celest.Mech.Dyn.Astron. 94, 351.
Precession matrix (including frame bias) from GCRS to a specified date, IAU 2000 model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rbp double[3][3] bias-precession matrix (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(date) = rbp * V(GCRS)
, where thep
-vectorV(GCRS)
is with respect to the Geocentric Celestial Reference System (IAU, 2000) and thep
-vectorV(date)
is with respect to the mean equatorial triad of the given date.Called:
iauBp00
- Frame bias and precession matrices, IAU 2000.
Reference:
- IAU: Trans. International Astronomical Union, Vol. XXIVB; Proc. 24th General Assembly, Manchester, UK. Resolutions B1.3, B1.6. (2000).
Precession matrix (including frame bias) from GCRS to a specified date, IAU 2006 model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rbp double[3][3] bias-precession matrix (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(date) = rbp * V(GCRS)
, where thep
-vectorV(GCRS)
is with respect to the Geocentric Celestial Reference System (IAU, 2000) and thep
-vector V(date) is with respect to the mean equatorial triad of the given date.Called:
iauPfw06
- Bias-precession F-W angles, IAU 2006.
iauFw2m
- F-W angles to
r
-matrix.References:
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
Precession matrix from J2000.0 to a specified date, IAU 1976 model.
Status: support function.
Given:
date1,date2 double ending date, TT (Note 1)Returned:
rmatp double[3][3] precession matrix, J2000.0 -> date1+
date2Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(date) = RMATP * V(J2000)
, where thep
-vectorV(J2000)
is with respect to the mean equatorial triad of epoch J2000.0 and thep
-vectorV(date)
is with respect to the mean equatorial triad of the given date.- Though the matrix method itself is rigorous, the precession angles are expressed through canonical polynomials which are valid only for a limited time span. In addition, the IAU 1976 precession rate is known to be imperfect. The absolute accuracy of the present formulation is better than 0.1 arcsec from 1960AD to 2040AD, better than 1 arcsec from 1640AD to 2360AD, and remains below 3 arcsec for the whole of the period 500BC to 3000AD. The errors exceed 10 arcsec outside the range 1200BC to 3900AD, exceed 100 arcsec outside 4200BC to 5600AD and exceed 1000 arcsec outside 6800BC to 8200AD.
Called:
iauPrec76
- Accumulated precession angles, IAU 1976.
iauIr
- Initialize
r
-matrix to identity.iauRz
- Rotate around Z-axis.
iauRy
- Rotate around Y-axis.
iauCr
- Copy
r
-matrix.References:
- Lieske, J.H., 1979, Astron.Astrophys. 73, 282.
- Equations (6) & (7), p283.
- Kaplan,G.H., 1981. USNO circular no. 163, pA2.
Precession–nutation, IAU 2000 model: a multi–purpose function, supporting classical (equinox–based) use directly and CIO–based use indirectly.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1) dpsi,deps double nutation (Note 2)Returned:
epsa double mean obliquity (Note 3) rb double[3][3] frame bias matrix (Note 4) rp double[3][3] precession matrix (Note 5) rbp double[3][3] bias-precession matrix (Note 6) rn double[3][3] nutation matrix (Note 7) rbpn double[3][3] GCRS-to-true matrix (Note 8)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The caller is responsible for providing the nutation components; they are in longitude and obliquity, in radians and are with respect to the equinox and ecliptic of date. For high–accuracy applications, free core nutation should be included as well as any other relevant corrections to the position of the CIP.
- The returned mean obliquity is consistent with the IAU 2000 precession–nutation models.
- The matrix rb transforms vectors from GCRS to J2000.0 mean equator and equinox by applying frame bias.
- The matrix rp transforms vectors from J2000.0 mean equator and equinox to mean equator and equinox of date by applying precession.
- The matrix rbp transforms vectors from GCRS to mean equator and equinox of date by applying frame bias then precession. It is the product rp
*
rb.- The matrix rn transforms vectors from mean equator and equinox of date to true equator and equinox of date by applying the nutation (luni–solar + planetary).
- The matrix rbpn transforms vectors from GCRS to true equator and equinox of date. It is the product rn
*
rbp, applying frame bias, precession and nutation in that order.- It is permissible to re-use the same array in the returned arguments. The arrays are filled in the order given.
Called:
iauPr00
- IAU 2000 precession adjustments.
iauObl80
- Mean obliquity, IAU 1980.
iauBp00
- Frame bias and precession matrices, IAU 2000.
iauCr
- Copy
r
-matrix.iauNumat
- Form nutation matrix.
iauRxr
- Product of two
r
-matrices.Reference:
- Capitaine, N., Chapront, J., Lambert, S. and Wallace, P., “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.
Precession–nutation, IAU 2000A model: a multi–purpose function, supporting classical (equinox–based) use directly and CIO–based use indirectly.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsi,deps double nutation (Note 2) epsa double mean obliquity (Note 3) rb double[3][3] frame bias matrix (Note 4) rp double[3][3] precession matrix (Note 5) rbp double[3][3] bias-precession matrix (Note 6) rn double[3][3] nutation matrix (Note 7) rbpn double[3][3] GCRS-to-true matrix (Notes 8,9)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The nutation components (luni–solar + planetary, IAU 2000A) in longitude and obliquity are in radians and with respect to the equinox and ecliptic of date. Free core nutation is omitted; for the utmost accuracy, use the
iauPn00
function, where the nutation components are caller–specified. For faster but slightly less accurate results, use theiauPn00b
function.- The mean obliquity is consistent with the IAU 2000 precession.
- The matrix rb transforms vectors from GCRS to J2000.0 mean equator and equinox by applying frame bias.
- The matrix rp transforms vectors from J2000.0 mean equator and equinox to mean equator and equinox of date by applying precession.
- The matrix rbp transforms vectors from GCRS to mean equator and equinox of date by applying frame bias then precession. It is the product rp
*
rb.- The matrix rn transforms vectors from mean equator and equinox of date to true equator and equinox of date by applying the nutation (luni–solar + planetary).
- The matrix rbpn transforms vectors from GCRS to true equator and equinox of date. It is the product rn
*
rbp, applying frame bias, precession and nutation in that order.- The X, Y, Z coordinates of the IAU 2000B Celestial Intermediate Pole are elements (3,1-3) of the matrix rbpn.
- It is permissible to re-use the same array in the returned arguments. The arrays are filled in the order given.
Called:
iauNut00a
- Nutation, IAU 2000A.
iauPn00
- Bias/precession/nutation results, IAU 2000.
Reference:
- Capitaine, N., Chapront, J., Lambert, S. and Wallace, P., “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.
Precession–nutation, IAU 2000B model: a multi–purpose function, supporting classical (equinox–based) use directly and CIO–based use indirectly.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsi,deps double nutation (Note 2) epsa double mean obliquity (Note 3) rb double[3][3] frame bias matrix (Note 4) rp double[3][3] precession matrix (Note 5) rbp double[3][3] bias-precession matrix (Note 6) rn double[3][3] nutation matrix (Note 7) rbpn double[3][3] GCRS-to-true matrix (Notes 8,9)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The nutation components (luni–solar + planetary, IAU 2000B) in longitude and obliquity are in radians and with respect to the equinox and ecliptic of date. For more accurate results, but at the cost of increased computation, use the
iauPn00a
function. For the utmost accuracy, use theiauPn00
function, where the nutation components are caller–specified.- The mean obliquity is consistent with the IAU 2000 precession.
- The matrix rb transforms vectors from GCRS to J2000.0 mean equator and equinox by applying frame bias.
- The matrix rp transforms vectors from J2000.0 mean equator and equinox to mean equator and equinox of date by applying precession.
- The matrix rbp transforms vectors from GCRS to mean equator and equinox of date by applying frame bias then precession. It is the product rp
*
rb.- The matrix rn transforms vectors from mean equator and equinox of date to true equator and equinox of date by applying the nutation (luni–solar + planetary).
- The matrix rbpn transforms vectors from GCRS to true equator and equinox of date. It is the product rn
*
rbp, applying frame bias, precession and nutation in that order.- The X, Y, Z coordinates of the IAU 2000B Celestial Intermediate Pole are elements (3,1-3) of the matrix rbpn.
- It is permissible to re-use the same array in the returned arguments. The arrays are filled in the stated order.
Called:
iauNut00b
- Nutation, IAU 2000B.
iauPn00
- Bias/precession/nutation results, IAU 2000.
Reference:
- Capitaine, N., Chapront, J., Lambert, S. and Wallace, P., “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.
Precession–nutation, IAU 2006 model: a multi–purpose function, supporting classical (equinox–based) use directly and CIO–based use indirectly.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1) dpsi,deps double nutation (Note 2)Returned:
epsa double mean obliquity (Note 3) rb double[3][3] frame bias matrix (Note 4) rp double[3][3] precession matrix (Note 5) rbp double[3][3] bias-precession matrix (Note 6) rn double[3][3] nutation matrix (Note 7) rbpn double[3][3] GCRS-to-true matrix (Note 8)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The caller is responsible for providing the nutation components; they are in longitude and obliquity, in radians and are with respect to the equinox and ecliptic of date. For high–accuracy applications, free core nutation should be included as well as any other relevant corrections to the position of the CIP.
- The returned mean obliquity is consistent with the IAU 2006 precession.
- The matrix rb transforms vectors from GCRS to J2000.0 mean equator and equinox by applying frame bias.
- The matrix rp transforms vectors from J2000.0 mean equator and equinox to mean equator and equinox of date by applying precession.
- The matrix rbp transforms vectors from GCRS to mean equator and equinox of date by applying frame bias then precession. It is the product rp
*
rb.- The matrix rn transforms vectors from mean equator and equinox of date to true equator and equinox of date by applying the nutation (luni–solar + planetary).
- The matrix rbpn transforms vectors from GCRS to true equator and equinox of date. It is the product rn
*
rbp, applying frame bias, precession and nutation in that order.- The X, Y, Z coordinates of the IAU 2000B Celestial Intermediate Pole are elements (3,1-3) of the matrix rbpn.
- It is permissible to re-use the same array in the returned arguments. The arrays are filled in the stated order.
Called:
iauPfw06
- Bias-precession F-W angles, IAU 2006.
iauFw2m
- F-W angles to
r
-matrix.iauCr
- Copy
r
-matrix.iauTr
- Transpose
r
-matrix.iauRxr
- Product of two
r
-matrices.References:
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
Precession–nutation, IAU 2006/2000A models: a multi–purpose function, supporting classical (equinox–based) use directly and CIO–based use indirectly.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsi,deps double nutation (Note 2) epsa double mean obliquity (Note 3) rb double[3][3] frame bias matrix (Note 4) rp double[3][3] precession matrix (Note 5) rbp double[3][3] bias-precession matrix (Note 6) rn double[3][3] nutation matrix (Note 7) rbpn double[3][3] GCRS-to-true matrix (Notes 8,9)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The nutation components (luni–solar + planetary, IAU 2000A) in longitude and obliquity are in radians and with respect to the equinox and ecliptic of date. Free core nutation is omitted; for the utmost accuracy, use the
iauPn06
function, where the nutation components are caller–specified.- The mean obliquity is consistent with the IAU 2006 precession.
- The matrix rb transforms vectors from GCRS to mean J2000.0 by applying frame bias.
- The matrix rp transforms vectors from mean J2000.0 to mean of date by applying precession.
- The matrix rbp transforms vectors from GCRS to mean of date by applying frame bias then precession. It is the product rp
*
rb.- The matrix rn transforms vectors from mean of date to true of date by applying the nutation (luni–solar + planetary).
- The matrix rbpn transforms vectors from GCRS to true of date (CIP/equinox). It is the product rn
*
rbp, applying frame bias, precession and nutation in that order.- The X, Y, Z coordinates of the IAU 2006/2000A Celestial Intermediate Pole are elements (1,1-3) of the matrix rbpn.
- It is permissible to re-use the same array in the returned arguments. The arrays are filled in the stated order.
Called:
iauNut06a
- Nutation, IAU 2006/2000A.
iauPn06
- Bias/precession/nutation results, IAU 2006.
Reference:
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855.
Form the matrix of precession–nutation for a given date (including frame bias), equinox–based, IAU 2000A model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rbpn double[3][3] classical NPB matrix (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(date) = rbpn * V(GCRS)
, where thep
-vectorV(date)
is with respect to the true equatorial triad of date date1+
date2 and thep
-vectorV(GCRS)
is with respect to the Geocentric Celestial Reference System (IAU, 2000).- A faster, but slightly less accurate result (about 1 mas), can be obtained by using instead the
iauPnm00b
function.Called:
iauPn00a
- Bias/precession/nutation, IAU 2000A.
Reference:
- IAU: Trans. International Astronomical Union, Vol. XXIVB; Proc. 24th General Assembly, Manchester, UK. Resolutions B1.3, B1.6. (2000).
Form the matrix of precession–nutation for a given date (including frame bias), equinox–based, IAU 2000B model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rbpn double[3][3] bias-precession-nutation matrix (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(date) = rbpn * V(GCRS)
, where thep
-vectorV(date)
is with respect to the true equatorial triad of date date1+
date2 and thep
-vector V(GCRS) is with respect to the Geocentric Celestial Reference System (IAU, 2000).- The present function is faster, but slightly less accurate (about 1 mas), than the
iauPnm00a
function.Called:
iauPn00b
- Bias/precession/nutation, IAU 2000B.
Reference:
- IAU: Trans. International Astronomical Union, Vol. XXIVB; Proc. 24th General Assembly, Manchester, UK. Resolutions B1.3, B1.6. (2000).
Form the matrix of precession–nutation for a given date (including frame bias), IAU 2006 precession and IAU 2000A nutation models.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
rnpb double[3][3] bias-precession-nutation matrix (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(date) = rnpb * V(GCRS)
, where thep
-vectorV(date)
is with respect to the true equatorial triad of date date1+
date2 and thep
-vectorV(GCRS)
is with respect to the Geocentric Celestial Reference System (IAU, 2000).Called:
iauPfw06
- Bias-precession F-W angles, IAU 2006.
iauNut06a
- Nutation, IAU 2006/2000A.
iauFw2m
- F-W angles to
r
-matrix.Reference:
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855.
Form the matrix of precession/nutation for a given date, IAU 1976 precession model, IAU 1980 nutation model.
Status: support function.
Given:
date1,date2 double TDB date (Note 1)Returned:
rmatpn double[3][3] combined precession/nutation matrixNotes:
- The TDB date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TDB) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 450123.7 0.0 (JD method) 451545.0 -1421.3 (J2000 method) 400000.5 50123.2 (MJD method) 450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The matrix operates in the sense
V(date) = rmatpn * V(J2000)
, where thep
-vector V(date) is with respect to the true equatorial triad of date date1+
date2 and thep
-vectorV(J2000)
is with respect to the mean equatorial triad of epoch J2000.0.Called:
iauPmat76
- Precession matrix, IAU 1976.
iauNutm80
- Nutation matrix, IAU 1980.
iauRxr
- Product of two
r
-matrices.Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 3.3 (p145).
Form the matrix of polar motion for a given date, IAU 2000.
Status: support function.
Given:
xp,yp double coordinates of the pole (radians, Note 1) sp double the TIO locator s' (radians, Note 2)Returned:
rpom double[3][3] polar-motion matrix (Note 3)Notes:
- The arguments xp and yp are the coordinates (in radians) of the Celestial Intermediate Pole with respect to the International Terrestrial Reference System (see IERS Conventions 2003), measured along the meridians to 0 and 90 deg west respectively.
- The argument sp is the TIO locator s', in radians, which positions the Terrestrial Intermediate Origin on the equator. It is obtained from polar motion observations by numerical integration, and so is in essence unpredictable. However, it is dominated by a secular drift of about 47 microarcseconds per century, and so can be taken into account by using s' = -47*t, where t is centuries since J2000.0. The function
iauSp00
implements this approximation.- The matrix operates in the sense
V(TRS) = rpom * V(CIP)
, meaning that it is the final rotation when computing the pointing direction to a celestial source.Called:
iauIr
- Initialize
r
-matrix to identity.iauRz
- Rotate around Z-axis.
iauRy
- Rotate around Y-axis.
iauRx
- Rotate around X-axis.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Precession–rate part of the IAU 2000 precession–nutation models (part of MHB2000).
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
dpsipr,depspr double precession corrections (Notes 2,3)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The precession adjustments are expressed as “nutation components”, corrections in longitude and obliquity with respect to the J2000.0 equinox and ecliptic.
- Although the precession adjustments are stated to be with respect to Lieske et al. (1977), the MHB2000 model does not specify which set of Euler angles are to be used and how the adjustments are to be applied. The most literal and straightforward procedure is to adopt the 4-rotation epsilon_0, psi_A, omega_A, xi_A option, and to add dpsipr to psi_A and depspr to both omega_A and eps_A.
- This is an implementation of one aspect of the IAU 2000A nutation model, formally adopted by the IAU General Assembly in 2000, namely MHB2000 (Mathews et al. 2002).
References:
- Lieske, J.H., Lederle, T., Fricke, W. & Morando, B., “Expressions for the precession quantities based upon the IAU (1976) System of Astronomical Constants”, Astron.Astrophys., 58, 1-16 (1977).
- Mathews, P.M., Herring, T.A., Buffet, B.A., “Modeling of nutation and precession New nutation series for nonrigid Earth and insights into the Earth's interior”, J.Geophys.Res., 107, B4, 2002. The MHB2000 code itself was obtained on 9th September 2002 from:
- Wallace, P.T., “Software for Implementing the IAU 2000 Resolutions”, in IERS Workshop 5.1 (2002).
IAU 1976 precession model. This function forms the three Euler angles which implement general precession between two epochs, using the IAU 1976 model (as for the FK5 catalog).
Status: canonical model.
Given:
ep01,ep02 double TDB starting epoch (Note 1) ep11,ep12 double TDB ending epoch (Note 1)Returned:
zeta double 1st rotation: radians cw around z z double 3rd rotation: radians cw around z theta double 2nd rotation: radians ccw around yNotes:
- The epochs ep01
+
ep02 and ep11+
ep12 are Julian Dates, apportioned in any convenient way between the arguments epn1 and epn2. For example,JD(TDB) = 2450123.7
could be expressed in any of these ways, among others:epn1 epn2 450123.7 0.0 (JD method) 451545.0 -1421.3 (J2000 method) 400000.5 50123.2 (MJD method) 450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience. The two epochs may be expressed using different methods, but at the risk of losing some resolution.
- The accumulated precession angles zeta, z, theta are expressed through canonical polynomials which are valid only for a limited time span. In addition, the IAU 1976 precession rate is known to be imperfect. The absolute accuracy of the present formulation is better than 0.1 arcsec from 1960AD to 2040AD, better than 1 arcsec from 1640AD to 2360AD, and remains below 3 arcsec for the whole of the period 500BC to 3000AD. The errors exceed 10 arcsec outside the range 1200BC to 3900AD, exceed 100 arcsec outside 4200BC to 5600AD and exceed 1000 arcsec outside 6800BC to 8200AD.
- The three angles are returned in the conventional order, which is not the same as the order of the corresponding Euler rotations. The precession matrix is:
R_3(-z) * R_2(+theta) * R_3(-zeta)Reference:
- Lieske, J.H., 1979, Astron.Astrophys. 73, 282, equations (6) & (7), p283.
The CIO locator s, positioning the Celestial Intermediate Origin on the equator of the Celestial Intermediate Pole, given the CIP's X, Y coordinates. Compatible with IAU 2000A precession–nutation.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1) x,y double CIP coordinates (Note 3)Returned (function value):
double the CIO locator s in radians (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The CIO locator s is the difference between the right ascensions of the same point in two systems: the two systems are the GCRS and the CIP, CIO, and the point is the ascending node of the CIP equator. The quantity s remains below 0.1 arcsecond throughout 1900-2100.
- The series used to compute s is in fact for s + X Y/2, where X and Y are the x and y components of the CIP unit vector; this series is more compact than a direct series for s would be. This function requires X, Y to be supplied by the caller, who is responsible for providing values that are consistent with the supplied date.
- The model is consistent with the IAU 2000A precession–nutation.
Called:
iauFal03
- Mean anomaly of the Moon.
iauFalp03
- Mean anomaly of the Sun.
iauFaf03
- Mean argument of the latitude of the Moon.
iauFad03
- Mean elongation of the Moon from the Sun.
iauFaom03
- Mean longitude of the Moon's ascending node.
iauFave03
- Mean longitude of Venus.
iauFae03
- Mean longitude of Earth.
iauFapa03
- General accumulated precession in longitude.
References:
- Capitaine, N., Chapront, J., Lambert, S. and Wallace, P., “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
The CIO locator s, positioning the Celestial Intermediate Origin on the equator of the Celestial Intermediate Pole, using the IAU 2000A precession–nutation model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double the CIO locator s in radians (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The CIO locator s is the difference between the right ascensions of the same point in two systems. The two systems are the GCRS and the CIP, CIO, and the point is the ascending node of the CIP equator. The CIO locator s remains a small fraction of 1 arcsecond throughout 1900-2100.
- The series used to compute s is in fact for s + X Y/2, where X and Y are the x and y components of the CIP unit vector; this series is more compact than a direct series for s would be. The present function uses the full IAU 2000A nutation model when predicting the CIP position. Faster results, with no significant loss of accuracy, can be obtained via the function
iauS00b
, which uses instead the IAU 2000B truncated model.Called:
iauPnm00a
- Classical NPB matrix, IAU 2000A.
iauBnp2xy
- Extract CIP X, Y from the BPN matrix.
iauS00
- The CIO locator s, given X, Y, IAU 2000A.
References:
- Capitaine, N., Chapront, J., Lambert, S. and Wallace, P., “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
The CIO locator s, positioning the Celestial Intermediate Origin on the equator of the Celestial Intermediate Pole, using the IAU 2000B precession–nutation model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double the CIO locator s in radians (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The CIO locator s is the difference between the right ascensions of the same point in two systems. The two systems are the GCRS and the CIP, CIO, and the point is the ascending node of the CIP equator. The CIO locator s remains a small fraction of 1 arcsecond throughout 1900-2100.
- The series used to compute s is in fact for s + X Y/2, where X and Y are the x and y components of the CIP unit vector; this series is more compact than a direct series for s would be. The present function uses the IAU 2000B truncated nutation model when predicting the CIP position. The function iauS00a uses instead the full IAU 2000A model, but with no significant increase in accuracy and at some cost in speed.
Called:
iauPnm00b
- Classical NPB matrix, IAU 2000B.
iauBnp2xy
- Extract CIP X, Y from the BPN matrix.
iauS00
- The CIO locator s, given X, Y, IAU 2000A.
References:
- Capitaine, N., Chapront, J., Lambert, S. and Wallace, P., “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
The CIO locator s, positioning the Celestial Intermediate Origin on the equator of the Celestial Intermediate Pole, given the CIP's X, Y coordinates. Compatible with IAU 2006/2000A precession–nutation.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1) x,y double CIP coordinates (Note 3)Returned (function value):
double the CIO locator s in radians (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The CIO locator s is the difference between the right ascensions of the same point in two systems: the two systems are the GCRS and the CIP, CIO, and the point is the ascending node of the CIP equator. The quantity s remains below 0.1 arcsecond throughout 1900-2100.
- The series used to compute s is in fact for s + X Y/2, where X and Y are the x and y components of the CIP unit vector; this series is more compact than a direct series for s would be. This function requires X, Y to be supplied by the caller, who is responsible for providing values that are consistent with the supplied date.
- The model is consistent with the “P03” precession (Capitaine et al. 2003), adopted by IAU 2006 Resolution 1, 2006, and the IAU 2000A nutation (with P03 adjustments).
Called:
iauFal03
- Mean anomaly of the Moon.
iauFalp03
- Mean anomaly of the Sun.
iauFaf03
- Mean argument of the latitude of the Moon.
iauFad03
- Mean elongation of the Moon from the Sun.
iauFaom03
- Mean longitude of the Moon's ascending node.
iauFave03
- Mean longitude of Venus.
iauFae03
- Mean longitude of Earth.
iauFapa03
- General accumulated precession in longitude.
References:
- Capitaine, N., Wallace, P.T. & Chapront, J., 2003, Astron. Astrophys. 432, 355.
- McCarthy, D.D., Petit, G. (eds.) 2004, IERS Conventions (2003), IERS Technical Note No. 32, BKG.
The CIO locator s, positioning the Celestial Intermediate Origin on the equator of the Celestial Intermediate Pole, using the IAU 2006 precession and IAU 2000A nutation models.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double the CIO locator s in radians (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The CIO locator s is the difference between the right ascensions of the same point in two systems. The two systems are the GCRS and the CIP, CIO, and the point is the ascending node of the CIP equator. The CIO locator s remains a small fraction of 1 arcsecond throughout 1900-2100.
- The series used to compute s is in fact for s + X Y/2, where X and Y are the x and y components of the CIP unit vector; this series is more compact than a direct series for s would be. The present function uses the full IAU 2000A nutation model when predicting the CIP position.
Called:
iauPnm06a
- Classical NPB matrix, IAU 2006/2000A.
iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
iauS06
- The CIO locator s, given X,Y, IAU 2006.
References:
- Capitaine, N., Chapront, J., Lambert, S. and Wallace, P., “Expressions for the Celestial Intermediate Pole and Celestial Ephemeris Origin consistent with the IAU 2000A precession–nutation model”, Astron.Astrophys. 400, 1145-1154 (2003).
n.b. The celestial ephemeris origin (CEO) was renamed “celestial intermediate origin” (CIO) by IAU 2006 Resolution 2.- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855
- McCarthy, D. D., Petit, G. (eds.), 2004, IERS Conventions (2003), IERS Technical Note No. 32, BKG.
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
The TIO locator s', positioning the Terrestrial Intermediate Origin on the equator of the Celestial Intermediate Pole.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double the TIO locator s' in radians (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The TIO locator s' is obtained from polar motion observations by numerical integration, and so is in essence unpredictable. However, it is dominated by a secular drift of about 47 microarcseconds per century, which is the approximation evaluated by the present function.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
X, Y coordinates of celestial intermediate pole from series based on IAU 2006 precession and IAU 2000A nutation.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
x,y double CIP X,Y coordinates (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The X, Y coordinates are those of the unit vector towards the celestial intermediate pole. They represent the combined effects of frame bias, precession and nutation.
- The fundamental arguments used are as adopted in IERS Conventions (2003) and are from Simon et al. (1994) and Souchay et al. (1999).
- This is an alternative to the angles–based method, via the SOFA function
iauFw2xy
and as used iniauXys06a
for example. The two methods agree at the 1 microarcsecond level (at present), a negligible amount compared with the intrinsic accuracy of the models. However, it would be unwise to mix the two methods (angles–based and series–based) in a single application.Called:
iauFal03
- Mean anomaly of the Moon.
iauFalp03
- Mean anomaly of the Sun.
iauFaf03
- Mean argument of the latitude of the Moon.
iauFad03
- Mean elongation of the Moon from the Sun.
iauFaom03
- Mean longitude of the Moon's ascending node.
iauFame03
- Mean longitude of Mercury.
iauFave03
- Mean longitude of Venus.
iauFae03
- Mean longitude of Earth.
iauFama03
- Mean longitude of Mars.
iauFaju03
- Mean longitude of Jupiter.
iauFasa03
- Mean longitude of Saturn.
iauFaur03
- Mean longitude of Uranus.
iauFane03
- Mean longitude of Neptune.
iauFapa03
- General accumulated precession in longitude.
References:
- Capitaine, N., Wallace, P.T. & Chapront, J., 2003, Astron.Astrophys., 412, 567.
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855
- McCarthy, D. D., Petit, G. (eds.), 2004, IERS Conventions (2003), IERS Technical Note No. 32, BKG.
- Simon, J.L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G. & Laskar, J., Astron.Astrophys., 1994, 282, 663
- Souchay, J., Loysel, B., Kinoshita, H., Folgueira, M., 1999, Astron.Astrophys.Supp.Ser. 135, 111
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
For a given TT date, compute the X, Y coordinates of the Celestial Intermediate Pole and the CIO locator s, using the IAU 2000A precession–nutation model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
x,y double Celestial Intermediate Pole (Note 2) s double the CIO locator s (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The Celestial Intermediate Pole coordinates are the x, y components of the unit vector in the Geocentric Celestial Reference System.
- The CIO locator s (in radians) positions the Celestial Intermediate Origin on the equator of the CIP.
- A faster, but slightly less accurate result (about 1 mas for X, Y), can be obtained by using instead the
iauXys00b
function.Called:
iauPnm00a
- Classical NPB matrix, IAU 2000A.
iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
iauS00
- The CIO locator s, given X,Y, IAU 2000A.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
For a given TT date, compute the X, Y coordinates of the Celestial Intermediate Pole and the CIO locator s, using the IAU 2000B precession–nutation model.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
x,y double Celestial Intermediate Pole (Note 2) s double the CIO locator s (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The Celestial Intermediate Pole coordinates are the X, Y components of the unit vector in the Geocentric Celestial Reference System.
- The CIO locator s (in radians) positions the Celestial Intermediate Origin on the equator of the CIP.
- The present function is faster, but slightly less accurate (about 1 mas in X, Y), than the
iauXys00a
function.Called:
iauPnm00b
- Classical NPB matrix, IAU 2000B.
iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
iauS00
- The CIO locator s, given X,Y, IAU 2000A.
Reference:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
For a given TT date, compute the X,Y coordinates of the Celestial Intermediate Pole and the CIO locator s, using the IAU 2006 precession and IAU 2000A nutation models.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned:
x,y double Celestial Intermediate Pole (Note 2) s double the CIO locator s (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and
- The Celestial Intermediate Pole coordinates are the x, y components of the unit vector in the Geocentric Celestial Reference System.
- The CIO locator s (in radians) positions the Celestial Intermediate Origin on the equator of the CIP.
- Series–based solutions for generating X and Y are also available: see Capitaine & Wallace (2006) and
iauXy06
.Called:
iauPnm06a
- Classical NPB matrix, IAU 2006/2000A.
iauBpn2xy
- Extract CIP X, Y coordinates from NPB matrix.
iauS06
- The CIO locator s, given X, Y, IAU 2006.
References:
- Capitaine, N. & Wallace, P.T., 2006, Astron.Astrophys. 450, 855
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
The equation of the equinoxes, compatible with IAU 2000 resolutions, given the nutation in longitude and the mean obliquity.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1) epsa double mean obliquity (Note 2) dpsi double nutation in longitude (Note 3)Returned (function value):
double equation of the equinoxes (Note 4)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The obliquity, in radians, is mean of date.
- The result, which is in radians, operates in the following sense:
Greenwich apparent ST = GMST + equation of the equinoxes- The result is compatible with the IAU 2000 resolutions. For further details, see IERS Conventions 2003 and Capitaine et al. (2002).
Called:
iauEect00
- Equation of the equinoxes complementary terms.
References:
- Capitaine, N., Wallace, P.T. and McCarthy, D.D., “Expressions to implement the IAU 2000 definition of UT1”, Astronomy & Astrophysics, 406, 1135-1149 (2003)
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Equation of the equinoxes, compatible with IAU 2000 resolutions.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double equation of the equinoxes (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The result, which is in radians, operates in the following sense:
Greenwich apparent ST = GMST + equation of the equinoxes- The result is compatible with the IAU 2000 resolutions. For further details, see IERS Conventions 2003 and Capitaine et al. (2002).
Called:
iauPr00
- IAU 2000 precession adjustments.
iauObl80
- Mean obliquity, IAU 1980.
iauNut00a
- Nutation, IAU 2000A.
iauEe00
- Equation of the equinoxes, IAU 2000.
References:
- Capitaine, N., Wallace, P.T. and McCarthy, D.D., “Expressions to implement the IAU 2000 definition of UT1”, Astronomy & Astrophysics, 406, 1135-1149 (2003).
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Equation of the equinoxes, compatible with IAU 2000 resolutions but using the truncated nutation model IAU 2000B.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double equation of the equinoxes (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The result, which is in radians, operates in the following sense:
Greenwich apparent ST = GMST + equation of the equinoxes- The result is compatible with the IAU 2000 resolutions except that accuracy has been compromised for the sake of speed. For further details, see McCarthy & Luzum (2001), IERS Conventions 2003 and Capitaine et al. (2003).
Called:
iauPr00
- IAU 2000 precession adjustments.
iauObl80
- Mean obliquity, IAU 1980.
iauNut00b
- Nutation, IAU 2000B.
iauEe00
- Equation of the equinoxes, IAU 2000.
References:
- Capitaine, N., Wallace, P.T. and McCarthy, D.D., “Expressions to implement the IAU 2000 definition of UT1”, Astronomy & Astrophysics, 406, 1135-1149 (2003).
- McCarthy, D.D. & Luzum, B.J., “An abridged model of the precession–nutation of the celestial pole”, Celestial Mechanics & Dynamical Astronomy, 85, 37-49 (2003).
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Equation of the equinoxes, compatible with IAU 2000 resolutions and IAU 2006/2000A precession–nutation.
Status: support function.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double equation of the equinoxes (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The result, which is in radians, operates in the following sense:
Greenwich apparent ST = GMST + equation of the equinoxesCalled:
iauAnpm
- Normalize angle into range +/- pi.
iauGst06a
- Greenwich apparent sidereal time, IAU 2006/2000A.
iauGmst06
- Greenwich mean sidereal time, IAU 2006.
Reference:
- McCarthy, D. D., Petit, G. (eds.), 2004, IERS Conventions (2003), IERS Technical Note No. 32, BKG
Equation of the equinoxes complementary terms, consistent with IAU 2000 resolutions.
Status: canonical model.
Given:
date1,date2 double TT as a 2-part Julian Date (Note 1)Returned (function value):
double complementary terms (Note 2)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The “complementary terms” are part of the equation of the equinoxes (EE), classically the difference between apparent and mean Sidereal Time:
GAST = GMST + EEwith:
EE = dpsi * cos(eps)where
dpsi
is the nutation in longitude and eps is the obliquity of date. However, if the rotation of the Earth were constant in an inertial frame the classical formulation would lead to apparent irregularities in the UT1 timescale traceable to side–effects of precession–nutation. In order to eliminate these effects from UT1, “complementary terms” were introduced in 1994 (IAU, 1994) and took effect from 1997 (Capitaine and Gontier, 1993):GAST = GMST + CT + EEBy convention, the complementary terms are included as part of the equation of the equinoxes rather than as part of the mean Sidereal Time. This slightly compromises the “geometrical” interpretation of mean sidereal time but is otherwise inconsequential.
The present function computes
CT
in the above expression, compatible with IAU 2000 resolutions (Capitaine et al., 2002, and IERS Conventions 2003).Called:
iauFal03
- Mean anomaly of the Moon.
iauFalp03
- Mean anomaly of the Sun.
iauFaf03
- Mean argument of the latitude of the Moon.
iauFad03
- Mean elongation of the Moon from the Sun.
iauFaom03
- Mean longitude of the Moon's ascending node.
iauFave03
- Mean longitude of Venus.
iauFae03
- Mean longitude of Earth.
iauFapa03
- General accumulated precession in longitude.
References:
- Capitaine, N. & Gontier, A.-M., Astron. Astrophys., 275, 645-650 (1993).
- Capitaine, N., Wallace, P.T. and McCarthy, D.D., “Expressions to implement the IAU 2000 definition of UT1”, Astronomy & Astrophysics, 406, 1135-1149 (2003).
- IAU Resolution C7, Recommendation 3 (1994).
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Equation of the equinoxes, IAU 1994 model.
Status: canonical model.
Given:
date1,date2 double TDB date (Note 1)Returned (function value):
double equation of the equinoxes (Note 2)Notes:
- The date date1 date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,
JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The result, which is in radians, operates in the following sense:
Greenwich apparent ST = GMST + equation of the equinoxesCalled:
iauNut80
- Nutation, IAU 1980.
iauObl80
- Mean obliquity, IAU 1980.
References:
- IAU Resolution C7, Recommendation 3 (1994).
- Capitaine, N. & Gontier, A.-M., 1993, Astron. Astrophys., 275, 645-650.
Earth rotation angle (IAU 2000 model).
Status: canonical model.
Given:
dj1,dj2 double UT1 as a 2-part Julian Date (see note)Returned (function value):
double Earth rotation angle (radians), range 0-2piNotes:
- The UT1 date dj1
+
dj2 is a Julian Date, apportioned in any convenient way between the arguments dj1 and dj2. For example,JD(UT1) = 2450123.7
could be expressed in any of these ways, among others:dj1 dj2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. The date & time method is best matched to the algorithm used: maximum precision is delivered when the dj1 argument is for 0hrs UT1 on the day in question and the dj2 argument lies in the range 0 to 1, or vice versa.
- The algorithm is adapted from Expression 22 of Capitaine et al. 2000. The time argument has been expressed in days directly, and, to retain precision, integer contributions have been eliminated. The same formulation is given in IERS Conventions (2003), Chap. 5, Eq. 14.
Called:
iauAnp
- Normalize angle into range 0 to 2pi.
References:
- Capitaine N., Guinot B. and McCarthy D.D, 2000, Astron. Astrophys., 355, 398-405.
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Greenwich mean sidereal time (model consistent with IAU 2000 resolutions).
Status: canonical model.
Given:
uta,utb double UT1 as a 2-part Julian Date (Notes 1,2) tta,ttb double TT as a 2-part Julian Date (Notes 1,2)Returned (function value):
double Greenwich mean sidereal time (radians)Notes:
- The UT1 and TT dates uta
+
utb and tta+
ttb respectively, are both Julian Dates, apportioned in any convenient way between the argument pairs. For example,JD = 2450123.7
could be expressed in any of these ways, among others:Part A Part B 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable (in the case of UT; the TT is not at all critical in this respect). The J2000 and MJD methods are good compromises between resolution and convenience. For UT, the date & time method is best matched to the algorithm that is used by the Earth Rotation Angle function, called internally: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- Both UT1 and TT are required, UT1 to predict the Earth rotation and TT to predict the effects of precession. If UT1 is used for both purposes, errors of order 100 microarcseconds result.
- This GMST is compatible with the IAU 2000 resolutions and must be used only in conjunction with other IAU 2000 compatible components such as precession-nutation and equation of the equinoxes.
- The result is returned in the range 0 to 2pi.
- The algorithm is from Capitaine et al. (2003) and IERS Conventions 2003.
Called:
iauEra00
- Earth rotation angle, IAU 2000.
iauAnp
- Normalize angle into range 0 to 2pi.
References:
- Capitaine, N., Wallace, P.T. and McCarthy, D.D., “Expressions to implement the IAU 2000 definition of UT1”, Astronomy & Astrophysics, 406, 1135-1149 (2003)
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Greenwich mean sidereal time (consistent with IAU 2006 precession).
Status: canonical model.
Given:
uta,utb double UT1 as a 2-part Julian Date (Notes 1,2) tta,ttb double TT as a 2-part Julian Date (Notes 1,2)Returned (function value):
double Greenwich mean sidereal time (radians)Notes:
- The UT1 and TT dates uta
+
utb and tta+
ttb respectively, are both Julian Dates, apportioned in any convenient way between the argument pairs. For example,JD = 2450123.7
could be expressed in any of these ways, among others:Part A Part B 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable (in the case of UT; the TT is not at all critical in this respect). The J2000 and MJD methods are good compromises between resolution and convenience. For UT, the date & time method is best matched to the algorithm that is used by the Earth rotation angle function, called internally: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- Both UT1 and TT are required, UT1 to predict the Earth rotation and TT to predict the effects of precession. If UT1 is used for both purposes, errors of order 100 microarcseconds result.
- This GMST is compatible with the IAU 2006 precession and must not be used with other precession models.
- The result is returned in the range 0 to 2pi.
Called:
iauEra00
- Earth rotation angle, IAU 2000.
iauAnp
- Normalize angle into range 0 to 2pi.
Reference:
- Capitaine, N., Wallace, P.T. & Chapront, J., 2005, Astron.Astrophys. 432, 355.
Universal Time to Greenwich mean sidereal time (IAU 1982 model).
Status: canonical model.
Given:
dj1,dj2 double UT1 Julian Date (see note)Returned (function value):
double Greenwich mean sidereal time (radians)Notes:
- The UT1 date dj1
+
dj2 is a Julian Date, apportioned in any convenient way between the arguments dj1 and dj2. For example,JD(UT1) = 2450123.7
could be expressed in any of these ways, among others:dj1 dj2 2450123.7D0 0D0 (JD method) 2451545D0 -1421.3D0 (J2000 method) 2400000.5D0 50123.2D0 (MJD method) 2450123.5D0 0.2D0 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. The date & time method is best matched to the algorithm used: maximum accuracy (or, at least, minimum noise) is delivered when the dj1 argument is for 0hrs UT1 on the day in question and the dj2 argument lies in the range 0 to 1, or vice versa.
- The algorithm is based on the IAU 1982 expression. This is always described as giving the GMST at 0 hours UT1. In fact, it gives the difference between the GMST and the UT, the steady 4–minutes–per–day drawing–ahead of ST with respect to UT. When whole days are ignored, the expression happens to equal the GMST at 0 hours UT1 each day.
- In this function, the entire UT1 (the sum of the two arguments dj1 and dj2) is used directly as the argument for the standard formula, the constant term of which is adjusted by 12 hours to take account of the noon phasing of Julian Date. The UT1 is then added, but omitting whole days to conserve accuracy.
Called:
iauAnp
- Normalize angle into range 0 to 2pi.
References:
- Transactions of the International Astronomical Union, XVIII B, 67 (1983).
- Aoki et al., Astron. Astrophys. 105, 359-361 (1982).
Greenwich apparent sidereal time (consistent with IAU 2000 resolutions).
Status: canonical model.
Given:
uta,utb double UT1 as a 2-part Julian Date (Notes 1,2) tta,ttb double TT as a 2-part Julian Date (Notes 1,2)Returned (function value):
double Greenwich apparent sidereal time (radians)Notes:
- The UT1 and TT dates uta
+
utb and tta+
ttb respectively, are both Julian Dates, apportioned in any convenient way between the argument pairs. For example,JD = 2450123.7
could be expressed in any of these ways, among others:Part A Part B 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable (in the case of UT; the TT is not at all critical in this respect). The J2000 and MJD methods are good compromises between resolution and convenience. For UT, the date & time method is best matched to the algorithm that is used by the Earth Rotation Angle function, called internally: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- Both UT1 and TT are required, UT1 to predict the Earth rotation and TT to predict the effects of precession–nutation. If UT1 is used for both purposes, errors of order 100 microarcseconds result.
- This GAST is compatible with the IAU 2000 resolutions and must be used only in conjunction with other IAU 2000 compatible components such as precession–nutation.
- The result is returned in the range 0 to 2pi.
- The algorithm is from Capitaine et al. (2003) and IERS Conventions 2003.
Called:
iauGmst00
- Greenwich mean sidereal time, IAU 2000.
iauEe00a
- Equation of the equinoxes, IAU 2000A.
iauAnp
- Normalize angle into range 0 to 2pi.
References:
- Capitaine, N., Wallace, P.T. and McCarthy, D.D., “Expressions to implement the IAU 2000 definition of UT1”, Astronomy & Astrophysics, 406, 1135-1149 (2003)
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Greenwich apparent sidereal time (consistent with IAU 2000 resolutions but using the truncated nutation model IAU 2000B).
Status: support function.
Given:
uta,utb double UT1 as a 2-part Julian Date (Notes 1,2)Returned (function value):
double Greenwich apparent sidereal time (radians)Notes:
- The UT1 date uta
+
utb is a Julian Date, apportioned in any convenient way between the argument pair. For example,JD = 2450123f.7
could be expressed in any of these ways, among others:uta utb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. For UT, the date & time method is best matched to the algorithm that is used by the Earth Rotation Angle function, called internally: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- The result is compatible with the IAU 2000 resolutions, except that accuracy has been compromised for the sake of speed and convenience in two respects:
- UT is used instead of TDB (or TT) to compute the precession component of GMST and the equation of the equinoxes. This results in errors of order 0.1 mas at present.
- The IAU 2000B abridged nutation model (McCarthy & Luzum, 2001) is used, introducing errors of up to 1 mas.
- This GAST is compatible with the IAU 2000 resolutions and must be used only in conjunction with other IAU 2000 compatible components such as precession-nutation.
- The result is returned in the range 0 to 2pi.
- The algorithm is from Capitaine et al. (2003) and IERS Conventions 2003.
Called:
iauGmst00
- Greenwich mean sidereal time, IAU 2000.
iauEe00b
- Equation of the equinoxes, IAU 2000B.
iauAnp
- Normalize angle into range 0 to 2pi.
References:
- Capitaine, N., Wallace, P.T. and McCarthy, D.D., “Expressions to implement the IAU 2000 definition of UT1”, Astronomy & Astrophysics, 406, 1135-1149 (2003)
- McCarthy, D.D. & Luzum, B.J., “An abridged model of the precession–nutation of the celestial pole”, Celestial Mechanics & Dynamical Astronomy, 85, 37-49 (2003).
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
Greenwich apparent sidereal time, IAU 2006, given the NPB matrix.
Status: support function.
Given:
uta,utb double UT1 as a 2-part Julian Date (Notes 1,2) tta,ttb double TT as a 2-part Julian Date (Notes 1,2) rnpb double[3][3] nutation x precession x bias matrixReturned (function value):
double Greenwich apparent sidereal time (radians)Notes:
- The UT1 and TT dates uta
+
utb and tta+
ttb respectively, are both Julian Dates, apportioned in any convenient way between the argument pairs. For example,JD = 2450123.7
could be expressed in any of these ways, among others:Part A Part B 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable (in the case of UT; the TT is not at all critical in this respect). The J2000 and MJD methods are good compromises between resolution and convenience. For UT, the date & time method is best matched to the algorithm that is used by the Earth rotation angle function, called internally: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- Both UT1 and TT are required, UT1 to predict the Earth rotation and TT to predict the effects of precession–nutation. If UT1 is used for both purposes, errors of order 100 microarcseconds result.
- Although the function uses the IAU 2006 series for
s + XY/2
, it is otherwise independent of the precession–nutation model and can in practice be used with any equinox–based NPB matrix.- The result is returned in the range 0 to 2pi.
Called:
iauBpn2xy
- Extract CIP X,Y coordinates from NPB matrix.
iauS06
- The CIO locator s, given X,Y, IAU 2006.
iauAnp
- Normalize angle into range 0 to 2pi.
iauEra00
- Earth rotation angle, IAU 2000.
iauEors
- Equation of the origins, given NPB matrix and s.
Reference:
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
Greenwich apparent sidereal time (consistent with IAU 2000 and 2006 resolutions).
Status: canonical model.
Given:
uta,utb double UT1 as a 2-part Julian Date (Notes 1,2) tta,ttb double TT as a 2-part Julian Date (Notes 1,2)Returned (function value):
double Greenwich apparent sidereal time (radians)Notes:
- The UT1 and TT dates uta
+
utb and tta+
ttb respectively, are both Julian Dates, apportioned in any convenient way between the argument pairs. For example,JD = 2450123.7
could be expressed in any of these ways, among others:Part A Part B 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable (in the case of UT; the TT is not at all critical in this respect). The J2000 and MJD methods are good compromises between resolution and convenience. For UT, the date & time method is best matched to the algorithm that is used by the Earth rotation angle function, called internally: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- Both UT1 and TT are required, UT1 to predict the Earth rotation and TT to predict the effects of precession-nutation. If UT1 is used for both purposes, errors of order 100 microarcseconds result.
- This GAST is compatible with the IAU 2000/2006 resolutions and must be used only in conjunction with IAU 2006 precession and IAU 2000A nutation.
- The result is returned in the range 0 to 2pi.
Called:
iauPnm06a
- Classical NPB matrix, IAU 2006/2000A.
iauGst06
- Greenwich apparent ST, IAU 2006, given NPB matrix.
Reference:
- Wallace, P.T. & Capitaine, N., 2006, Astron.Astrophys. 459, 981.
Greenwich apparent sidereal time (consistent with IAU 1982/94 resolutions).
Status: support function.
Given:
uta,utb double UT1 as a 2-part Julian Date (Notes 1,2)Returned (function value):
double Greenwich apparent sidereal time (radians)Notes:
- The UT1 date uta
+
utb is a Julian Date, apportioned in any convenient way between the argument pair. For example,JD = 2450123.7
could be expressed in any of these ways, among others:uta utb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 and MJD methods are good compromises between resolution and convenience. For UT, the date & time method is best matched to the algorithm that is used by the Earth Rotation Angle function, called internally: maximum precision is delivered when the uta argument is for 0hrs UT1 on the day in question and the utb argument lies in the range 0 to 1, or vice versa.
- The result is compatible with the IAU 1982 and 1994 resolutions, except that accuracy has been compromised for the sake of convenience in that UT is used instead of TDB (or TT) to compute the equation of the equinoxes.
- This GAST must be used only in conjunction with contemporaneous IAU standards such as 1976 precession, 1980 obliquity and 1982 nutation. It is not compatible with the IAU 2000 resolutions.
- The result is returned in the range 0 to 2pi.
Called:
iauGmst82
- Greenwich mean sidereal time, IAU 1982.
iauEqeq94
- Equation of the equinoxes, IAU 1994.
iauAnp
- Normalize angle into range 0 to 2pi.
References:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
- IAU Resolution C7, Recommendation 3 (1994).
Convert star position+velocity vector to catalog coordinates.
Status: support function.
Given (Note 1):
pv double[2][3] pv-vector (AU, AU/day)Returned (Note 2):
ra double right ascension (radians) dec double declination (radians) pmr double RA proper motion (radians/year) pmd double Dec proper motion (radians/year) px double parallax (arcsec) rv double radial velocity (km/s, positive = receding)Returned (function value):
int status: 0 = OK -1 = superluminal speed (Note 5) -2 = null position vectorNotes:
- The specified pv-vector is the coordinate direction (and its rate of change) for the date at which the light leaving the star reached the solar–system barycenter.
- The star data returned by this function are “observables” for an imaginary observer at the solar–system barycenter. Proper motion and radial velocity are, strictly, in terms of barycentric coordinate time, TCB. For most practical applications, it is permissible to neglect the distinction between TCB and ordinary “proper” time on Earth (TT/TAI). The result will, as a rule, be limited by the intrinsic accuracy of the proper–motion and radial–velocity data; moreover, the supplied pv-vector is likely to be merely an intermediate result (for example generated by the function iauStarpv), so that a change of time unit will cancel out overall.
In accordance with normal star–catalog conventions, the object's right ascension and declination are freed from the effects of secular aberration. The frame, which is aligned to the catalog equator and equinox, is Lorentzian and centered on the SSB.
Summarizing, the specified pv-vector is for most stars almost identical to the result of applying the standard geometrical “space motion” transformation to the catalog data. The differences, which are the subject of the Stumpff paper cited below, are:
- In stars with significant radial velocity and proper motion, the constantly changing light–time distorts the apparent proper motion. Note that this is a classical, not a relativistic, effect.
- The transformation complies with special relativity.
- Care is needed with units. The star coordinates are in radians and the proper motions in radians per Julian year, but the parallax is in arcseconds; the radial velocity is in km/s, but the pv-vector result is in AU and AU/day.
- The proper motions are the rate of change of the right ascension and declination at the catalog epoch and are in radians per Julian year. The RA proper motion is in terms of coordinate angle, not true angle, and will thus be numerically larger at high declinations.
- Straight–line motion at constant speed in the inertial frame is assumed. If the speed is greater than or equal to the speed of light, the function aborts with an error status.
- The inverse transformation is performed by the function iauStarpv.
Called:
iauPn
- Decompose
p
-vector into modulus and direction.iauPdp
- Scalar product of two
p
-vectors.iauSxp
- Multiply
p
-vector by scalar.iauPmp
p
-vector minusp
-vector.iauPm
- Modulus of
p
-vector.iauPpp
p
-vector plusp
-vector.iauPv2s
pv
-vector to spherical.iauAnp
- Normalize angle into range 0 to 2pi.
Reference:
- Stumpff, P., 1985, Astron.Astrophys. 144, 232-240.
Convert star catalog coordinates to position+velocity vector.
Status: support function.
Given (Note 1):
ra double right ascension (radians) dec double declination (radians) pmr double RA proper motion (radians/year) pmd double Dec proper motion (radians/year) px double parallax (arcseconds) rv double radial velocity (km/s, positive = receding)Returned (Note 2):
pv double[2][3] pv-vector (AU, AU/day)Returned (function value):
int status: 0 = no warnings 1 = distance overridden (Note 6) 2 = excessive speed (Note 7) 4 = solution didn't converge (Note 8) else = binary logical OR of the aboveNotes:
- The star data accepted by this function are “observables” for an imaginary observer at the solar–system barycenter. Proper motion and radial velocity are, strictly, in terms of barycentric coordinate time, TCB. For most practical applications, it is permissible to neglect the distinction between TCB and ordinary “proper” time on Earth (TT/TAI). The result will, as a rule, be limited by the intrinsic accuracy of the proper–motion and radial–velocity data; moreover, the pv-vector is likely to be merely an intermediate result, so that a change of time unit would cancel out overall.
In accordance with normal star–catalog conventions, the object's right ascension and declination are freed from the effects of secular aberration. The frame, which is aligned to the catalog equator and equinox, is Lorentzian and centered on the SSB.
- The resulting position and velocity pv-vector is with respect to the same frame and, like the catalog coordinates, is freed from the effects of secular aberration. Should the “coordinate direction”, where the object was located at the catalog epoch, be required, it may be obtained by calculating the magnitude of the position vector pv
[0][0-2]
dividing by the speed of light in AU/day to give the light–time, and then multiplying the space velocity pv[1][0-2]
by this light–time and adding the result to pv[0][0-2]
.Summarizing, the pv-vector returned is for most stars almost identical to the result of applying the standard geometrical “space motion” transformation. The differences, which are the subject of the Stumpff paper referenced below, are:
- In stars with significant radial velocity and proper motion, the constantly changing light–time distorts the apparent proper motion. Note that this is a classical, not a relativistic, effect.
- The transformation complies with special relativity.
- Care is needed with units. The star coordinates are in radians and the proper motions in radians per Julian year, but the parallax is in arcseconds; the radial velocity is in km/s, but the pv-vector result is in AU and AU/day.
- The RA proper motion is in terms of coordinate angle, not true angle. If the catalog uses arcseconds for both RA and Dec proper motions, the RA proper motion will need to be divided by
cos(Dec)
before use.- Straight–line motion at constant speed, in the inertial frame, is assumed.
- An extremely small (or zero or negative) parallax is interpreted to mean that the object is on the “celestial sphere”, the radius of which is an arbitrary (large) value (see the constant
PXMIN
). When the distance is overridden in this way, the status, initially zero, has 1 added to it.- If the space velocity is a significant fraction of c (see the constant
VMAX
), it is arbitrarily set to zero. When this action occurs, 2 is added to the status.- The relativistic adjustment involves an iterative calculation. If the process fails to converge within a set number (
IMAX
) of iterations, 4 is added to the status.- The inverse transformation is performed by the function
iauPvstar
.Called:
iauS2pv
- Spherical coordinates to pv-vector.
iauPm
- Modulus of
p
-vector.iauZp
- Zero
p
-vector.iauPn
- Decompose
p
-vector into modulus and direction.iauPdp
- Scalar product of two
p
-vectors.iauSxp
- Multiply
p
-vector by scalar.iauPmp
p
-vector minusp
-vector.iauPpp
p
-vector plusp
-vector.Reference:
- Stumpff, P., 1985, Astron.Astrophys. 144, 232-240.
Transform FK5 (J2000.0) star data into the Hipparcos system.
Status: support function.
Given (all FK5, equinox J2000.0, epoch J2000.0):
r5 double RA (radians) d5 double Dec (radians) dr5 double proper motion in RA (dRA/dt, rad/Jyear) dd5 double proper motion in Dec (dDec/dt, rad/Jyear) px5 double parallax (arcsec) rv5 double radial velocity (km/s, positive = receding)Returned (all Hipparcos, epoch J2000.0):
rh double RA (radians) dh double Dec (radians) drh double proper motion in RA (dRA/dt, rad/Jyear) ddh double proper motion in Dec (dDec/dt, rad/Jyear) pxh double parallax (arcsec) rvh double radial velocity (km/s, positive = receding)Notes:
- This function transforms FK5 star positions and proper motions into the system of the Hipparcos catalog.
- The proper motions in RA are
dRA/dt
rather thancos(Dec) * dRA/dt
, and are per year rather than per century.- The FK5 to Hipparcos transformation is modeled as a pure rotation and spin; zonal errors in the FK5 catalog are not taken into account.
- See also
iauH2fk5
,iauFk5hz
,iauHfk5z
.Called:
iauStarpv
- Star catalog data to space motion
pv
-vector.iauFk5hip
- FK5 to Hipparcos rotation and spin.
iauRxp
- Product of
r
-matrix andp
-vector.iauPxp
- Vector product of two
p
-vectors.iauPpp
p
-vector plusp
-vector.iauPvstar
- Space motion
pv
-vector to star catalog data.Reference:
- F.Mignard & M.Froeschle, Astron. Astrophys. 354, 732-739 (2000).
FK5 to Hipparcos rotation and spin.
Status: support function.
Returned:
r5h double[3][3] r-matrix: FK5 rotation wrt Hipparcos (Note 2) s5h double[3] r-vector: FK5 spin wrt Hipparcos (Note 3)Notes:
- This function models the FK5 to Hipparcos transformation as a pure rotation and spin; zonal errors in the FK5 catalogue are not taken into account.
- The
r
-matrix r5h operates in the sense:P_Hipparcos = r5h x P_FK5where
P_FK5
is ap
-vector in the FK5 frame, andP_Hipparcos
is the equivalent Hipparcosp
-vector.- The
r
-vector s5h represents the time derivative of the FK5 to Hipparcos rotation. The units are radians per year (Julian, TDB).Called:
iauRv2m
r
-vector tor
-matrix.Reference:
- F.Mignard & M.Froeschle, Astron. Astrophys. 354, 732-739 (2000).
Transform an FK5 (J2000.0) star position into the system of the Hipparcos catalogue, assuming zero Hipparcos proper motion.
Status: support function.
Given:
r5 double FK5 RA (radians), equinox J2000.0, at date d5 double FK5 Dec (radians), equinox J2000.0, at date date1,date2 double TDB date (Notes 1,2)Returned:
rh double Hipparcos RA (radians) dh double Hipparcos Dec (radians)Notes:
- This function converts a star position from the FK5 system to the Hipparcos system, in such a way that the Hipparcos proper motion is zero. Because such a star has, in general, a non-zero proper motion in the FK5 system, the function requires the date at which the position in the FK5 system was determined.
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The FK5 to Hipparcos transformation is modeled as a pure rotation and spin; zonal errors in the FK5 catalogue are not taken into account.
- The position returned by this function is in the Hipparcos reference system but at date date1+date2.
- See also
iauFk52h
,iauH2fk5
,iauHfk5z
.Called:
iauS2c
- Spherical coordinates to unit vector.
iauFk5hip
- FK5 to Hipparcos rotation and spin.
iauSxp
- Multiply
p
-vector by scalar.iauRv2m
r
-vector tor
-matrix.iauTrxp
- Product of transpose of
r
-matrix andp
-vector.iauPxp
- Vector product of two
p
-vectors.iauC2s
p
-vector to spherical.iauAnp
- Normalize angle into range 0 to 2pi.
Reference:
- F.Mignard & M.Froeschle, 2000, Astron.Astrophys. 354, 732-739.
Transform Hipparcos star data into the FK5 (J2000.0) system.
Status: support function.
Given (all Hipparcos, epoch J2000.0):
rh double RA (radians) dh double Dec (radians) drh double proper motion in RA (dRA/dt, rad/Jyear) ddh double proper motion in Dec (dDec/dt, rad/Jyear) pxh double parallax (arcsec) rvh double radial velocity (km/s, positive = receding)Returned (all FK5, equinox J2000.0, epoch J2000.0):
r5 double RA (radians) d5 double Dec (radians) dr5 double proper motion in RA (dRA/dt, rad/Jyear) dd5 double proper motion in Dec (dDec/dt, rad/Jyear) px5 double parallax (arcsec) rv5 double radial velocity (km/s, positive = receding)Notes:
- This function transforms Hipparcos star positions and proper motions into FK5 J2000.0.
- The proper motions in RA are
dRA/dt
rather thancos(Dec) * dRA/dt
, and are per year rather than per century.- The FK5 to Hipparcos transformation is modeled as a pure rotation and spin; zonal errors in the FK5 catalog are not taken into account.
- See also
iauFk52h
,iauFk5hz
,iauHfk5z
.Called:
iauStarpv
- Star catalog data to space motion
pv
-vector.iauFk5hip
- FK5 to Hipparcos rotation and spin.
iauRv2m
r
-vector tor
-matrix.iauRxp
- Product of
r
-matrix andp
-vector.iauTrxp
- Product of transpose of
r
-matrix andp
-vector.iauPxp
- Vector product of two
p
-vectors.iauPmp
p
-vector minusp
-vector.iauPvstar
- Space motion
pv
-vector to star catalog data.Reference:
- F.Mignard & M.Froeschle, Astron. Astrophys. 354, 732-739 (2000).
Transform a Hipparcos star position into FK5 J2000.0, assuming zero Hipparcos proper motion.
Status: support function.
Given:
rh double Hipparcos RA (radians) dh double Hipparcos Dec (radians) date1,date2 double TDB date (Note 1)Returned (all FK5, equinox J2000.0, date date1
+
date2):r5 double RA (radians) d5 double Dec (radians) dr5 double FK5 RA proper motion (rad/year, Note 4) dd5 double Dec proper motion (rad/year, Note 4)Notes:
- The TT date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- The proper motion in RA is
dRA/dt
rather thancos(Dec) * dRA/dt
.- The FK5 to Hipparcos transformation is modeled as a pure rotation and spin; zonal errors in the FK5 catalogue are not taken into account.
- It was the intention that Hipparcos should be a close approximation to an inertial frame, so that distant objects have zero proper motion; such objects have (in general) non-zero proper motion in FK5, and this function returns those fictitious proper motions.
- The position returned by this function is in the FK5 J2000.0 reference system but at date date1
+
date2.- See also
iauFk52h
,iauH2fk5
,iauFk5zhz
.Called:
iauS2c
- Spherical coordinates to unit vector.
iauFk5hip
- FK5 to Hipparcos rotation and spin.
iauRxp
- Product of
r
-matrix andp
-vector.iauSxp
- Multiply
p
-vector by scalar.iauRxr
- Product of two
r
-matrices.iauTrxp
- Product of transpose of
r
-matrix andp
-vector.iauPxp
- Vector product of two
p
-vectors.iauPv2s
pv
-vector to spherical.iauAnp
- Normalize angle into range 0 to 2pi.
Reference:
- F.Mignard & M.Froeschle, 2000, Astron.Astrophys. 354, 732-739.
Star proper motion: update star catalog data for space motion.
Status: support function.
Given:
ra1 double right ascension (radians), before dec1 double declination (radians), before pmr1 double RA proper motion (radians/year), before pmd1 double Dec proper motion (radians/year), before px1 double parallax (arcseconds), before rv1 double radial velocity (km/s, +ve = receding), before ep1a double "before" epoch, part A (Note 1) ep1b double "before" epoch, part B (Note 1) ep2a double "after" epoch, part A (Note 1) ep2b double "after" epoch, part B (Note 1)Returned:
ra2 double right ascension (radians), after dec2 double declination (radians), after pmr2 double RA proper motion (radians/year), after pmd2 double Dec proper motion (radians/year), after px2 double parallax (arcseconds), after rv2 double radial velocity (km/s, +ve = receding), afterReturned (function value):
int status: -1 = system error (should not occur) 0 = no warnings or errors 1 = distance overridden (Note 6) 2 = excessive velocity (Note 7) 4 = solution didn't converge (Note 8) else = binary logical OR of the above warningsNotes:
- The starting and ending TDB dates ep1a
+
ep1b and ep2a+
ep2b are Julian Dates, apportioned in any convenient way between the two parts (A and B). For example,JD(TDB) = 2450123.7
could be expressed in any of these ways, among others:epna epnb 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
- In accordance with normal star–catalog conventions, the object's right ascension and declination are freed from the effects of secular aberration. The frame, which is aligned to the catalog equator and equinox, is Lorentzian and centered on the SSB.
The proper motions are the rate of change of the right ascension and declination at the catalog epoch and are in radians per TDB Julian year.
The parallax and radial velocity are in the same frame.
- Care is needed with units. The star coordinates are in radians and the proper motions in radians per Julian year, but the parallax is in arcseconds.
- The RA proper motion is in terms of coordinate angle, not true angle. If the catalog uses arcseconds for both RA and Dec proper motions, the RA proper motion will need to be divided by
cos(Dec)
before use.- Straight–line motion at constant speed, in the inertial frame, is assumed.
- An extremely small (or zero or negative) parallax is interpreted to mean that the object is on the “celestial sphere”, the radius of which is an arbitrary (large) value (see the
iauStarpv
function for the value used). When the distance is overridden in this way, the status, initially zero, has 1 added to it.- If the space velocity is a significant fraction of c (see the constant
VMAX
in the functioniauStarpv
), it is arbitrarily set to zero. When this action occurs, 2 is added to the status.- The relativistic adjustment carried out in the
iauStarpv
function involves an iterative calculation. If the process fails to converge within a set number of iterations, 4 is added to the status.Called:
iauStarpv
- Star catalog data to space motion
pv
-vector.iauPvu
- Update a pv-vector.
iauPdp
- Scalar product of two
p
-vectors.iauPvstar
- Space motion
pv
-vector to star catalog data.
Earth reference ellipsoids.
Status: canonical.
Given:
n int ellipsoid identifier (Note 1)Returned:
a double equatorial radius (meters, Note 2) f double flattening (Note 2)Returned (function value):
int status: 0 = OK -1 = illegal identifier (Note 3)Notes:
- The identifier n is a number that specifies the choice of reference ellipsoid. The following are supported:
n ellipsoid 1 WGS84 2 GRS80 3 WGS72The n value has no significance outside the SOFA software. For convenience, symbols WGS84 etc. are defined in sofam.h.
- The ellipsoid parameters are returned in the form of equatorial radius in meters (a) and flattening (f). The latter is a number around 0.00335, i.e. around 1/298.
- For the case where an unsupported n value is supplied, zero a and f are returned, as well as error status.
References:
- Department of Defense World Geodetic System 1984, National Imagery and Mapping Agency Technical Report 8350.2, Third Edition, p3-2.
- Moritz, H., Bull. Geodesique 66-2, 187 (1992).
- The Department of Defense World Geodetic System 1972, World Geodetic System Committee, May 1974.
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), p220.
Transform geocentric coordinates to geodetic using the specified reference ellipsoid.
Status: canonical transformation.
Given:
n int ellipsoid identifier (Note 1) xyz double[3] geocentric vector (Note 2)Returned:
elong double longitude (radians, east +ve) phi double latitude (geodetic, radians, Note 3) height double height above ellipsoid (geodetic, Notes 2,3)Returned (function value):
int status: 0 = OK -1 = illegal identifier (Note 3) -2 = internal error (Note 3)Notes:
- The identifier n is a number that specifies the choice of reference ellipsoid. The following are supported:
n ellipsoid 1 WGS84 2 GRS80 3 WGS72The n value has no significance outside the SOFA software. For convenience, symbols WGS84 etc. are defined in sofam.h.
- The geocentric vector (xyz, given) and height (height, returned) are in meters.
- An error status -1 means that the identifier n is illegal. An error status -2 is theoretically impossible. In all error cases, phi and height are both set to -1e9.
- The inverse transformation is performed in the function
iauGd2gc
.Called:
iauEform
- Earth reference ellipsoids.
iauGc2gde
- Geocentric to geodetic transformation, general.
Transform geocentric coordinates to geodetic for a reference ellipsoid of specified form.
Status: support function.
Given:
a double equatorial radius (Notes 2,4) f double flattening (Note 3) xyz double[3] geocentric vector (Note 4)Returned:
elong double longitude (radians, east +ve) phi double latitude (geodetic, radians) height double height above ellipsoid (geodetic, Note 4)Returned (function value):
int status: 0 = OK -1 = illegal f -2 = illegal aNotes:
- This function is based on the GCONV2H Fortran subroutine by Toshio Fukushima (see reference).
- The equatorial radius, a, can be in any units, but meters is the conventional choice.
- The flattening, f, is (for the Earth) a value around 0.00335, i.e. around 1/298.
- The equatorial radius, a, and the geocentric vector, xyz, must be given in the same units, and determine the units of the returned height, height.
- If an error occurs (status < 0), elong, phi and height are unchanged.
- The inverse transformation is performed in the function
iauGd2gce
.- The transformation for a standard ellipsoid (such as WGS84) can more conveniently be performed by calling iauGc2gd, which uses a numerical code to identify the required A and F values.
Reference:
- Fukushima, T., “Transformation from Cartesian to geodetic coordinates accelerated by Halley's method”, J.Geodesy (2006) 79: 689-693
Transform geodetic coordinates to geocentric using the specified reference ellipsoid.
Status: canonical transformation.
Given:
n int ellipsoid identifier (Note 1) elong double longitude (radians, east +ve) phi double latitude (geodetic, radians, Note 3) height double height above ellipsoid (geodetic, Notes 2,3)Returned:
xyz double[3] geocentric vector (Note 2)Returned (function value):
int status: 0 = OK -1 = illegal identifier (Note 3) -2 = illegal case (Note 3)Notes:
- The identifier n is a number that specifies the choice of reference ellipsoid. The following are supported:
n ellipsoid 1 WGS84 2 GRS80 3 WGS72The n value has no significance outside the SOFA software. For convenience, symbols WGS84 etc. are defined in sofam.h.
- The height (height, given) and the geocentric vector (xyz, returned) are in meters.
- No validation is performed on the arguments elong, phi and height. An error status -1 means that the identifier n is illegal. An error status -2 protects against cases that would lead to arithmetic exceptions. In all error cases, xyz is set to zeros.
- The inverse transformation is performed in the function
iauGc2gd
.Called:
iauEform
- Earth reference ellipsoids.
iauGd2gce
- Geodetic to geocentric transformation, general.
iauZp
- Zero
p
-vector.
Transform geodetic coordinates to geocentric for a reference ellipsoid of specified form.
Status: support function.
Given:
a double equatorial radius (Notes 1,4) f double flattening (Notes 2,4) elong double longitude (radians, east +ve) phi double latitude (geodetic, radians, Note 4) height double height above ellipsoid (geodetic, Notes 3,4)Returned:
xyz double[3] geocentric vector (Note 3)Returned (function value):
int status: 0 = OK -1 = illegal case (Note 4)Notes:
- The equatorial radius, a, can be in any units, but meters is the conventional choice.
- The flattening, f, is (for the Earth) a value around 0.00335, i.e. around 1/298.
- The equatorial radius, a, and the height, height, must be given in the same units, and determine the units of the returned geocentric vector, xyz.
- No validation is performed on individual arguments. The error status -1 protects against (unrealistic) cases that would lead to arithmetic exceptions. If an error occurs, xyz is unchanged.
- The inverse transformation is performed in the function
iauGc2gde
.- The transformation for a standard ellipsoid (such as WGS84) can more conveniently be performed by calling
iauGd2gc
, which uses a numerical code to identify the required a and f values.References:
- Green, R.M., Spherical Astronomy, Cambridge University Press, (1985) Section 4.5, p96.
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992), Section 4.22, p202.
Format for output a 2-part Julian Date (or in the case of UTC a quasi-JD form that includes special provision for leap seconds).
Status: support function.
Given:
scale char[] time scale ID (Note 1) ndp int resolution (Note 2) d1,d2 double time as a 2-part Julian Date (Notes 3,4)Returned:
iy,im,id int year, month, day in Gregorian calendar (Note 5) ihmsf int[4] hours, minutes, seconds, fraction (Note 1)Returned (function value):
int status: +1 = dubious year (Note 5) 0 = OK -1 = unacceptable date (Note 6)Notes:
- scale identifies the time scale. Only the value
UTC
(in upper case) is significant, and enables handling of leap seconds (see Note 4).- ndp is the number of decimal places in the seconds field, and can have negative as well as positive values, such as:
ndp resolution -4 1 00 00 -3 0 10 00 -2 0 01 00 -1 0 00 10 0 0 00 01 1 0 00 00.1 2 0 00 00.01 3 0 00 00.001The limits are platform dependent, but a safe range is -5 to +9.
- d1
+
d2 is Julian Date, apportioned in any convenient way between the two arguments, for example where d1 is the Julian Day Number and d2 is the fraction of a day. In the case of UTC, where the use of JD is problematical, special conventions apply: see the next note.- JD cannot unambiguously represent UTC during a leap second unless special measures are taken. The SOFA internal convention is that the quasi-JD day represents UTC days whether the length is 86399, 86400 or 86401 SI seconds.
- The warning status “dubious year” flags UTCs that predate the introduction of the time scale and that are too far in the future to be trusted. See
iauDat
for further details.- For calendar conventions and limitations, see
iauCal2jd
.Called:
iauJd2cal
- JD to Gregorian calendar.
iauD2tf
- Decompose days to hms.
iauDat
delta(AT) = TAI-UTC
.
For a given UTC date, calculate
delta(AT) = TAI-UTC
.IMPORTANT A new version of this function must be produced whenever a new leap second is announced. There are four items to change on each such occasion:
- A new line must be added to the set of statements that initialize the array “changes”.
- The parameter IYV must be set to the current year.
- The “Latest leap second” comment below must be set to the new leap second date.
- The “This revision” comment, later, must be set to the current date.
Change (2) must also be carried out whenever the function is re–issued, even if no leap seconds have been added.
Latest leap second: 2012 June 30
Status: support function.
Given:
iy int UTC: year (Notes 1 and 2) im int month (Note 2) id int day (Notes 2 and 3) fd double fraction of day (Note 4)Returned:
deltat double TAI minus UTC, secondsReturned (function value):
int status (Note 5): 1 = dubious year (Note 1) 0 = OK -1 = bad year -2 = bad month -3 = bad day (Note 3) -4 = bad fraction (Note 4)Notes:
- UTC began at 1960 January 1.0 (JD 2436934.5) and it is improper to call the function with an earlier date. If this is attempted, zero is returned together with a warning status.
Because leap seconds cannot, in principle, be predicted in advance, a reliable check for dates beyond the valid range is impossible. To guard against gross errors, a year five or more after the release year of the present function (see parameter IYV) is considered dubious. In this case a warning status is returned but the result is computed in the normal way.
For both too–early and too–late years, the warning status is j = +1. This is distinct from the error status j = -1, which signifies a year so early that JD could not be computed.
- If the specified date is for a day which ends with a leap second, the UTC-TAI value returned is for the period leading up to the leap second. If the date is for a day which begins as a leap second ends, the UTC-TAI returned is for the period following the leap second.
- The day number must be in the normal calendar range, for example 1 through 30 for April. The “almanac” convention of allowing such dates as January 0 and December 32 is not supported in this function, in order to avoid confusion near leap seconds.
- The fraction of day is used only for dates before the introduction of leap seconds, the first of which occurred at the end of 1971. It is tested for validity (0 to 1 is the valid range) even if not used; if invalid, zero is used and status j = -4 is returned. For many applications, setting fd to zero is acceptable; the resulting error is always less than 3 ms (and occurs only pre-1972).
- The status value returned in the case where there are multiple errors refers to the first error detected. For example, if the month and day are 13 and 32 respectively, j = -2 (bad month) will be returned.
- In cases where a valid result is not available, zero is returned.
References:
- For dates from 1961 January 1 onwards, the expressions from the file ftp://maia.usno.navy.mil/ser7/tai-utc.dat are used.
- The 5ms timestep at 1961 January 1 is taken from 2.58.1 (p87) of the 1992 Explanatory Supplement.
Called:
iauCal2jd
- Gregorian calendar to Julian Day number.
An approximation to TDB-TT, the difference between barycentric dynamical time and terrestrial time, for an observer on the Earth.
The different time scales — proper, coordinate and realized — are related to each other:
TAI <- physically realized : offset <- observed (nominally +32.184s) : TT <- terrestrial time : rate adjustment (L_G) <- definition of TT : TCG <- time scale for GCRS : "periodic" terms <- iauDtdb is an implementation : rate adjustment (L_C) <- function of solar-system ephemeris : TCB <- time scale for BCRS : rate adjustment (-L_B) <- definition of TDB : TDB <- TCB scaled to track TT : "periodic" terms <- -iauDtdb is an approximation : TT <- terrestrial timeAdopted values for the various constants can be found in the IERS Conventions (McCarthy & Petit 2003).
Status: support routine.
Given:
date1,date2 double date, TDB (Notes 1-3) ut double universal time (UT1, fraction of one day) elong double longitude (east positive, radians) u double distance from Earth spin axis (km) v double distance north of equatorial plane (km)Returned (function value):
double TDB-TT (seconds)Notes:
- The date date1
+
date2 is a Julian Date, apportioned in any convenient way between the two arguments. For example,JD(TT) = 2450123.7
could be expressed in any of these ways, among others:date1 date2 2450123.7 0.0 (JD method) 2451545.0 -1421.3 (J2000 method) 2400000.5 50123.2 (MJD method) 2450123.5 0.2 (date & time method)The JD method is the most natural and convenient to use in cases where the loss of several decimal digits of resolution is acceptable. The J2000 method is best matched to the way the argument is handled internally and will deliver the optimum resolution. The MJD method and the date & time methods are both good compromises between resolution and convenience.
Although the date is, formally, barycentric dynamical time (TDB), the terrestrial dynamical time (TT) can be used with no practical effect on the accuracy of the prediction.
- TT can be regarded as a coordinate time that is realized as an offset of 32.184s from International Atomic Time, TAI. TT is a specific linear transformation of geocentric coordinate time TCG, which is the time scale for the Geocentric Celestial Reference System, GCRS.
- TDB is a coordinate time, and is a specific linear transformation of barycentric coordinate time TCB, which is the time scale for the Barycentric Celestial Reference System, BCRS.
- The difference TCG-TCB depends on the masses and positions of the bodies of the solar system and the velocity of the Earth. It is dominated by a rate difference, the residual being of a periodic character. The latter, which is modeled by the present function, comprises a main (annual) sinusoidal term of amplitude approximately 0.00166 seconds, plus planetary terms up to about 20 microseconds, and lunar and diurnal terms up to 2 microseconds. These effects come from the changing transverse Doppler effect and gravitational red-shift as the observer (on the Earth's surface) experiences variations in speed (with respect to the BCRS) and gravitational potential.
- TDB can be regarded as the same as TCB but with a rate adjustment to keep it close to TT, which is convenient for many applications. The history of successive attempts to define TDB is set out in Resolution 3 adopted by the IAU General Assembly in 2006, which defines a fixed TDB(TCB) transformation that is consistent with contemporary solar-system ephemerides. Future ephemerides will imply slightly changed transformations between TCG and TCB, which could introduce a linear drift between TDB and TT; however, any such drift is unlikely to exceed 1 nanosecond per century.
- The geocentric TDB-TT model used in the present function is that of Fairhead & Bretagnon (1990), in its full form. It was originally supplied by Fairhead (private communications with P.T.Wallace, 1990) as a Fortran subroutine. The present C function contains an adaptation of the Fairhead code. The numerical results are essentially unaffected by the changes, the differences with respect to the Fairhead & Bretagnon original being at the 1e-20 s level.
The topocentric part of the model is from Moyer (1981) and Murray (1983), with fundamental arguments adapted from Simon et al. 1994. It is an approximation to the expression
( v / c ) . ( r / c )
, wherev
is the barycentric velocity of the Earth,r
is the geocentric position of the observer andc
is the speed of light.By supplying zeroes for u and v, the topocentric part of the model can be nullified, and the function will return the Fairhead & Bretagnon result alone.
- During the interval 1950-2050, the absolute accuracy is better than +/- 3 nanoseconds relative to time ephemerides obtained by direct numerical integrations based on the JPL DE405 solar system ephemeris.
- It must be stressed that the present function is merely a model, and that numerical integration of solar-system ephemerides is the definitive method for predicting the relationship between TCG and TCB and hence between TT and TDB.
References:
- Fairhead, L., & Bretagnon, P., Astron.Astrophys., 229, 240-247 (1990).
- IAU 2006 Resolution 3.
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Moyer, T.D., Cel.Mech., 23, 33 (1981).
- Murray, C.A., Vectorial Astrometry, Adam Hilger (1983).
- Seidelmann, P.K. et al., Explanatory Supplement to the Astronomical Almanac, Chapter 2, University Science Books (1992).
- Simon, J.L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G. & Laskar, J., Astron.Astrophys., 282, 663-683 (1994).
Encode date and time fields into 2-part Julian Date (or in the case of UTC a quasi-JD form that includes special provision for leap seconds).
Status: support function.
Given:
scale char[] time scale ID (Note 1) iy,im,id int year, month, day in Gregorian calendar (Note 2) ihr,imn int hour, minute sec double secondsReturned:
d1,d2 double 2-part Julian Date (Notes 3,4)Returned (function value):
int status: +3 = both of next two +2 = time is after end of day (Note 5) +1 = dubious year (Note 6) 0 = OK -1 = bad year -2 = bad month -3 = bad day -4 = bad hour -5 = bad minute -6 = bad second (<0)Notes:
- scale identifies the time scale. Only the value
UTC
(in upper case) is significant, and enables handling of leap seconds (see Note 4).- For calendar conventions and limitations, see iauCal2jd.
- The sum of the results, d1
+
d2, is Julian Date, where normally d1 is the Julian Day Number and d2 is the fraction of a day. In the case of UTC, where the use of JD is problematical, special conventions apply: see the next note.- JD cannot unambiguously represent UTC during a leap second unless special measures are taken. The SOFA internal convention is that the quasi-JD day represents UTC days whether the length is 86399, 86400 or 86401 SI seconds.
- The warning status “time is after end of day” usually means that the sec argument is greater than 60.0. However, in a day ending in a leap second the limit changes to 61.0 (or 59.0 in the case of a negative leap second).
- The warning status “dubious year” flags UTCs that predate the introduction of the time scale and that are too far in the future to be trusted. See
iauDat
for further details.- Only in the case of continuous and regular time scales (TAI, TT, TCG, TCB and TDB) is the result d1+d2 a Julian Date, strictly speaking. In the other cases (UT1 and UTC) the result must be used with circumspection; in particular the difference between two such results cannot be interpreted as a precise time interval.
Called:
iauCal2jd
- Gregorian calendar to JD.
iauDat
delta(AT) = TAI-UTC
.iauJd2cal
- JD to Gregorian calendar.
Time scale transformation: International Atomic Time, TAI, to Terrestrial Time, TT.
Status: canonical.
Given:
tai1,tai2 double TAI as a 2-part Julian DateReturned:
tt1,tt2 double TT as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tai1
+
tai2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tai1 is the Julian Day Number and tai2 is the fraction of a day. The returned tt1, tt2 follow suit.References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: International Atomic Time, TAI, to Universal Time, UT1.
Status: canonical.
Given:
tai1,tai2 double TAI as a 2-part Julian Date dta double UT1-TAI in secondsReturned:
ut11,ut12 double UT1 as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tai1
+
tai2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tai1 is the Julian Day Number and tai2 is the fraction of a day. The returned ut11, ut12 follow suit.- The argument dta, i.e. UT1-TAI, is an observed quantity, and is available from IERS tabulations.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: International Atomic Time, TAI, to Coordinated Universal Time, UTC.
Status: canonical.
Given:
tai1,tai2 double TAI as a 2-part Julian Date (Note 1)Returned:
utc1,utc2 double UTC as a 2-part quasi Julian Date (Notes 1-3)Returned (function value):
int status: +1 = dubious year (Note 4) 0 = OK -1 = unacceptable dateNotes:
- tai1
+
tai2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tai1 is the Julian Day Number and tai2 is the fraction of a day. The returned utc1 and utc2 form an analogous pair, except that a special convention is used, to deal with the problem of leap seconds—see the next note.- JD cannot unambiguously represent UTC during a leap second unless special measures are taken. The convention in the present function is that the JD day represents UTC days whether the length is 86399, 86400 or 86401 SI seconds.
- The function
iauD2dtf
can be used to transform the UTC quasi-JD into calendar date and clock time, including UTC leap second handling.- The warning status “dubious year” flags UTCs that predate the introduction of the time scale and that are too far in the future to be trusted. See
iauDat
for further details.Called:
iauJd2cal
- JD to Gregorian calendar.
iauDat
delta(AT) = TAI-UTC
.iauCal2jd
- Gregorian calendar to JD.
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: Barycentric Coordinate Time, TCB, to Barycentric Dynamical Time, TDB.
Status: canonical.
Given:
tcb1,tcb2 double TCB as a 2-part Julian DateReturned:
tdb1,tdb2 double TDB as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tcb1
+
tcb2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tcb1 is the Julian Day Number and tcb2 is the fraction of a day. The returned tdb1, tdb2 follow suit.- The 2006 IAU General Assembly introduced a conventional linear transformation between TDB and TCB. This transformation compensates for the drift between TCB and terrestrial time TT, and keeps TDB approximately centered on TT. Because the relationship between TT and TCB depends on the adopted solar system ephemeris, the degree of alignment between TDB and TT over long intervals will vary according to which ephemeris is used. Former definitions of TDB attempted to avoid this problem by stipulating that TDB and TT should differ only by periodic effects. This is a good description of the nature of the relationship but eluded precise mathematical formulation. The conventional linear relationship adopted in 2006 sidestepped these difficulties whilst delivering a TDB that in practice was consistent with values before that date.
- TDB is essentially the same as Teph, the time argument for the JPL solar system ephemerides.
Reference:
- IAU 2006 Resolution B3.
Time scale transformation: Geocentric Coordinate Time, TCG, to Terrestrial Time, TT.
Status: canonical.
Given:
tcg1,tcg2 double TCG as a 2-part Julian DateReturned:
tt1,tt2 double TT as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tcg1
+
tcg2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tcg1 is the Julian Day Number and tcg2 is the fraction of a day. The returned tt1, tt2 follow suit.References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003),. IERS Technical Note No. 32, BKG (2004).
- IAU 2000 Resolution B1.9.
Time scale transformation: Barycentric Dynamical Time, TDB, to Barycentric Coordinate Time, TCB.
Status: canonical.
Given:
tdb1,tdb2 double TDB as a 2-part Julian DateReturned:
tcb1,tcb2 double TCB as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tdb1
+
tdb2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tdb1 is the Julian Day Number and tdb2 is the fraction of a day. The returned tcb1, tcb2 follow suit.- The 2006 IAU General Assembly introduced a conventional linear transformation between TDB and TCB. This transformation compensates for the drift between TCB and terrestrial time TT, and keeps TDB approximately centered on TT. Because the relationship between TT and TCB depends on the adopted solar system ephemeris, the degree of alignment between TDB and TT over long intervals will vary according to which ephemeris is used. Former definitions of TDB attempted to avoid this problem by stipulating that TDB and TT should differ only by periodic effects. This is a good description of the nature of the relationship but eluded precise mathematical formulation. The conventional linear relationship adopted in 2006 sidestepped these difficulties whilst delivering a TDB that in practice was consistent with values before that date.
- TDB is essentially the same as Teph, the time argument for the JPL solar system ephemerides.
Reference:
- IAU 2006 Resolution B3.
Time scale transformation: Barycentric Dynamical Time, TDB, to Terrestrial Time, TT.
Status: canonical.
Given:
tdb1,tdb2 double TDB as a 2-part Julian Date dtr double TDB-TT in secondsReturned:
tt1,tt2 double TT as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tdb1
+
tdb2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tdb1 is the Julian Day Number and tdb2 is the fraction of a day. The returned tt1, tt2 follow suit.- The argument dtr represents the quasi–periodic component of the GR transformation between TT and TCB. It is dependent upon the adopted solar–system ephemeris, and can be obtained by numerical integration, by interrogating a precomputed time ephemeris or by evaluating a model such as that implemented in the SOFA function
iauDtdb
. The quantity is dominated by an annual term of 1.7 ms amplitude.- TDB is essentially the same as Teph, the time argument for the JPL solar system ephemerides.
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- IAU 2006 Resolution 3.
Time scale transformation: Terrestrial Time, TT, to International Atomic Time, TAI.
Status: canonical.
Given:
tt1,tt2 double TT as a 2-part Julian DateReturned:
tai1,tai2 double TAI as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tt1
+
tt2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tt1 is the Julian Day Number and tt2 is the fraction of a day. The returned tai1, tai2 follow suit.References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: Terrestrial Time, TT, to Geocentric Coordinate Time, TCG.
Status: canonical.
Given:
tt1,tt2 double TT as a 2-part Julian DateReturned:
tcg1,tcg2 double TCG as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tt1
+
tt2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tt1 is the Julian Day Number and tt2 is the fraction of a day. The returned tcg1, tcg2 follow suit.References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- IAU 2000 Resolution B1.9.
Time scale transformation: Terrestrial Time, TT, to Barycentric Dynamical Time, TDB.
Status: canonical.
Given:
tt1,tt2 double TT as a 2-part Julian Date dtr double TDB-TT in secondsReturned:
tdb1,tdb2 double TDB as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tt1
+
tt2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tt1 is the Julian Day Number and tt2 is the fraction of a day. The returned tdb1, tdb2 follow suit.- The argument dtr represents the quasi–periodic component of the GR transformation between TT and TCB. It is dependent upon the adopted solar–system ephemeris, and can be obtained by numerical integration, by interrogating a precomputed time ephemeris or by evaluating a model such as that implemented in the SOFA function
iauDtdb
. The quantity is dominated by an annual term of 1.7 ms amplitude.- TDB is essentially the same as Teph, the time argument for the JPL solar system ephemerides.
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- IAU 2006 Resolution 3.
Time scale transformation: Terrestrial Time, TT, to Universal Time, UT1.
Status: canonical.
Given:
tt1,tt2 double TT as a 2-part Julian Date dt double TT-UT1 in secondsReturned:
ut11,ut12 double UT1 as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- tt1
+
tt2 is Julian Date, apportioned in any convenient way between the two arguments, for example where tt1 is the Julian Day Number and tt2 is the fraction of a day. The returned ut11, ut12 follow suit.- The argument dt is classical Delta T.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: Universal Time, UT1, to International Atomic Time, TAI.
Status: canonical.
Given:
ut11,ut12 double UT1 as a 2-part Julian Date dta double UT1-TAI in secondsReturned:
tai1,tai2 double TAI as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- ut11
+
ut12 is Julian Date, apportioned in any convenient way between the two arguments, for example where ut11 is the Julian Day Number and ut12 is the fraction of a day. The returned tai1, tai2 follow suit.- The argument dta, i.e. UT1-TAI, is an observed quantity, and is available from IERS tabulations.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: Universal Time, UT1, to Terrestrial Time, TT.
Status: canonical.
Given:
ut11,ut12 double UT1 as a 2-part Julian Date dt double TT-UT1 in secondsReturned:
tt1,tt2 double TT as a 2-part Julian DateReturned (function value):
int status: 0 = OKNotes:
- ut11
+
ut12 is Julian Date, apportioned in any convenient way between the two arguments, for example where ut11 is the Julian Day Number and ut12 is the fraction of a day. The returned tt1, tt2 follow suit.- The argument dt is classical Delta T.
Reference:
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: Universal Time, UT1, to Coordinated Universal Time, UTC.
Status: canonical.
Given:
ut11,ut12 double UT1 as a 2-part Julian Date (Note 1) dut1 double Delta UT1: UT1-UTC in seconds (Note 2)Returned:
utc1,utc2 double UTC as a 2-part quasi Julian Date (Notes 3,4)Returned (function value):
int status: +1 = dubious year (Note 5) 0 = OK -1 = unacceptable dateNotes:
- ut11
+
ut12 is Julian Date, apportioned in any convenient way between the two arguments, for example where ut11 is the Julian Day Number and ut12 is the fraction of a day. The returned utc1 and utc2 form an analogous pair, except that a special convention is used, to deal with the problem of leap seconds—see Note 3.- Delta UT1 can be obtained from tabulations provided by the International Earth Rotation and Reference Systems Service. The value changes abruptly by 1 second at a leap second; however, close to a leap second the algorithm used here is tolerant of the “wrong” choice of value being made.
- JD cannot unambiguously represent UTC during a leap second unless special measures are taken. The convention in the present function is that the returned quasi JD day UTC1
+
utc2 represents UTC days whether the length is 86399, 86400 or 86401 SI seconds.- The function
iauD2dtf
can be used to transform the UTC quasi-JD into calendar date and clock time, including UTC leap second handling.- The warning status “dubious year” flags UTCs that predate the introduction of the time scale and that are too far in the future to be trusted. See
iauDat
for further details.Called:
iauJd2cal
- JD to Gregorian calendar.
iauDat
delta(AT) = TAI-UTC
.iauCal2jd
- Gregorian calendar to JD.
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: Coordinated Universal Time, UTC, to International Atomic Time, TAI.
Status: canonical.
Given:
utc1,utc2 double UTC as a 2-part quasi Julian Date (Notes 1-4)Returned:
tai1,tai2 double TAI as a 2-part Julian Date (Note 5)Returned (function value):
int status: +1 = dubious year (Note 3) 0 = OK -1 = unacceptable dateNotes:
- utc1
+
utc2 is quasi Julian Date (see Note 2), apportioned in any convenient way between the two arguments, for example where utc1 is the Julian Day Number and utc2 is the fraction of a day.- JD cannot unambiguously represent UTC during a leap second unless special measures are taken. The convention in the present function is that the JD day represents UTC days whether the length is 86399, 86400 or 86401 SI seconds.
- The warning status “dubious year” flags UTCs that predate the introduction of the time scale and that are too far in the future to be trusted. See
iauDat
for further details.- The function
iauDtf2d
converts from calendar date and time of day into 2-part Julian Date, and in the case of UTC implements the leap–second–ambiguity convention described above.- The returned TAI1, TAI2 are such that their sum is the TAI Julian Date.
Called:
iauJd2cal
- JD to Gregorian calendar.
iauDat
delta(AT) = TAI-UTC
.iauCal2jd
- Gregorian calendar to JD.
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Time scale transformation: Coordinated Universal Time, UTC, to Universal Time, UT1.
Status: canonical.
Given:
utc1,utc2 double UTC as a 2-part quasi Julian Date (Notes 1-4) dut1 double Delta UT1 = UT1-UTC in seconds (Note 5)Returned:
ut11,ut12 double UT1 as a 2-part Julian Date (Note 6)Returned (function value):
int status: +1 = dubious year (Note 7) 0 = OK -1 = unacceptable dateNotes:
- utc1
+
utc2 is quasi Julian Date (see Note 2), apportioned in any convenient way between the two arguments, for example where utc1 is the Julian Day Number and utc2 is the fraction of a day.- JD cannot unambiguously represent UTC during a leap second unless special measures are taken. The convention in the present function is that the JD day represents UTC days whether the length is 86399, 86400 or 86401 SI seconds.
- The warning status “dubious year” flags UTCs that predate the introduction of the time scale and that are too far in the future to be trusted. See
iauDat
for further details.- The function
iauDtf2d
converts from calendar date and time of day into 2-part Julian Date, and in the case of UTC implements the leap–second–ambiguity convention described above.- Delta UT1 can be obtained from tabulations provided by the International Earth Rotation and Reference Systems Service. It It is the caller's responsibility to supply a DUT argument containing the UT1-UTC value that matches the given UTC.
- The returned ut11, ut12 are such that their sum is the UT1 Julian Date.
- The warning status “dubious year” flags UTCs that predate the introduction of the time scale and that are too far in the future to be trusted. See
iauDat
for further details.Called:
iauJd2cal
- JD to Gregorian calendar.
iauDat
delta(AT) = TAI-UTC
.iauUtctai
- UTC to TAI.
iauTaiut1
- TAI to UT1.
References:
- McCarthy, D. D., Petit, G. (eds.), IERS Conventions (2003), IERS Technical Note No. 32, BKG (2004).
- Explanatory Supplement to the Astronomical Almanac, P. Kenneth Seidelmann (ed), University Science Books (1992).
Decompose radians into degrees, arcminutes, arcseconds, fraction.
Status: vector/matrix support function.
Given:
ndp int resolution (Note 1) angle double angle in radiansReturned:
sign char '+' or '-' idmsf int[4] degrees, arcminutes, arcseconds, fractionCalled:
iauD2tf
- Decompose days to hms.
Notes:
- The argument ndp is interpreted as follows:
ndp resolution : ...0000 00 00 -7 1000 00 00 -6 100 00 00 -5 10 00 00 -4 1 00 00 -3 0 10 00 -2 0 01 00 -1 0 00 10 0 0 00 01 1 0 00 00.1 2 0 00 00.01 3 0 00 00.001 : 0 00 00.000...- The largest positive useful value for ndp is determined by the size of angle, the format of doubles on the target platform, and the risk of overflowing
idmsf[3]
. On a typical platform, for angle up to2pi
, the available floating–point precision might correspond tondp=12
. However, the practical limit is typicallyndp=9
, set by the capacity of a 32-bit int, orndp=4
ifint
is only 16 bits.- The absolute value of angle may exceed
2pi
. In cases where it does not, it is up to the caller to test for and handle the case where angle is very nearly2pi
and rounds up to 360 degrees, by testing foridmsf[0]=360
and settingidmsf[0-3]
to zero.
Decompose radians into hours, minutes, seconds, fraction.
Status: vector/matrix support function.
Given:
ndp int resolution (Note 1) angle double angle in radiansReturned:
sign char '+' or '-' ihmsf int[4] hours, minutes, seconds, fractionCalled:
iauD2tf
- Decompose days to hms.
Notes:
- The argument ndp is interpreted as follows:
ndp resolution : ...0000 00 00 -7 1000 00 00 -6 100 00 00 -5 10 00 00 -4 1 00 00 -3 0 10 00 -2 0 01 00 -1 0 00 10 0 0 00 01 1 0 00 00.1 2 0 00 00.01 3 0 00 00.001 : 0 00 00.000...- The largest positive useful value for ndp is determined by the size of angle, the format of doubles on the target platform, and the risk of overflowing
ihmsf[3]
. On a typical platform, for angle up to2pi
, the available floating–point precision might correspond tondp=12
. However, the practical limit is typicallyndp=9
, set by the capacity of a 32-bit int, orndp=4
ifint
is only 16 bits.- The absolute value of angle may exceed
2pi
. In cases where it does not, it is up to the caller to test for and handle the case where angle is very nearly2pi
and rounds up to 24 hours, by testing forihmsf[0]=24
and settingihmsf(0-3)
to zero.
Convert degrees, arcminutes, arcseconds to radians.
Status: support function.
Given:
s char sign: '-' = negative, otherwise positive ideg int degrees iamin int arcminutes asec double arcsecondsReturned:
rad double angle in radiansReturned (function value):
int status: 0 = OK 1 = ideg outside range 0-359 2 = iamin outside range 0-59 3 = asec outside range 0-59.999...Notes:
- The result is computed even if any of the range checks fail.
- Negative ideg, iamin and/or asec produce a warning status, but the absolute value is used in the conversion.
- If there are multiple errors, the status value reflects only the first, the smallest taking precedence.
Normalize angle into the range
0 <=
a< 2pi
.Status: vector/matrix support function.
Given:
a double angle (radians)Returned (function value):
double angle in range 0-2pi
Normalize angle into the range
-pi <=
a< +pi
.Status: vector/matrix support function.
Given:
a double angle (radians)Returned (function value):
double angle in range +/-pi
Decompose days to hours, minutes, seconds, fraction.
Status: vector/matrix support function.
Given:
ndp int resolution (Note 1) days double interval in daysReturned:
sign char '+' or '-' ihmsf int[4] hours, minutes, seconds, fractionNotes:
- The argument ndp is interpreted as follows:
ndp resolution : ...0000 00 00 -7 1000 00 00 -6 100 00 00 -5 10 00 00 -4 1 00 00 -3 0 10 00 -2 0 01 00 -1 0 00 10 0 0 00 01 1 0 00 00.1 2 0 00 00.01 3 0 00 00.001 : 0 00 00.000...- The largest positive useful value for ndp is determined by the size of days, the format of double on the target platform, and the risk of overflowing
ihmsf[3]
. On a typical platform, for days up to 1.0, the available floating–point precision might correspond to ndp= 12
. However, the practical limit is typically ndp= 9
, set by the capacity of a 32-bit int, or ndp= 4
if int is only 16 bits.- The absolute value of days may exceed 1.0. In cases where it does not, it is up to the caller to test for and handle the case where days is very nearly 1.0 and rounds up to 24 hours, by testing for ihmsf
[0] = 24
and setting ihmsf[0-3]
to zero.
Convert hours, minutes, seconds to radians.
Status: support function.
Given:
s char sign: '-' = negative, otherwise positive ihour int hours imin int minutes sec double secondsReturned:
rad double angle in radiansReturned (function value):
int status: 0 = OK 1 = ihour outside range 0-23 2 = imin outside range 0-59 3 = sec outside range 0-59.999...Notes:
- The result is computed even if any of the range checks fail.
- Negative ihour, imin and/or sec produce a warning status, but the absolute value is used in the conversion.
- If there are multiple errors, the status value reflects only the first, the smallest taking precedence.
Convert hours, minutes, seconds to days.
Status: support function.
Given:
s char sign: '-' = negative, otherwise positive ihour int hours imin int minutes sec double secondsReturned:
days double interval in daysReturned (function value):
int status: 0 = OK 1 = ihour outside range 0-23 2 = imin outside range 0-59 3 = sec outside range 0-59.999...Notes:
- The result is computed even if any of the range checks fail.
- Negative ihour, imin and/or sec produce a warning status, but the absolute value is used in the conversion.
- If there are multiple errors, the status value reflects only the first, the smallest taking precedence.
Rotate an
r
-matrix about the X-axis.Status: vector/matrix support function.
Given:
phi double angle (radians)Given and returned:
r double[3][3] r-matrix, rotatedNotes:
- Calling this function with positive phi incorporates in the supplied
r
-matrix r an additional rotation, about the X-axis, anticlockwise as seen looking towards the origin from positive x.- The additional rotation can be represented by this matrix:
( 1 0 0 ) ( 0 +cos(phi) +sin(phi) ) ( 0 -sin(phi) +cos(phi) )
Rotate an
r
-matrix about the Y-axis.Status: vector/matrix support function.
Given:
theta double angle (radians)Given and returned:
r double[3][3] r-matrix, rotatedNotes:
- Calling this function with positive theta incorporates in the supplied
r
-matrix r an additional rotation, about the Y-axis, anticlockwise as seen looking towards the origin from positive Y.- The additional rotation can be represented by this matrix:
( +cos(theta) 0 -sin(theta) ) ( 0 1 0 ) ( +sin(theta) 0 +cos(theta) )
Rotate an
r
-matrix about the Z-axis.Status: vector/matrix support function.
Given:
psi double angle (radians)Given and returned:
r double[3][3] r-matrix, rotatedNotes:
- Calling this function with positive psi incorporates in the supplied
r
-matrix r an additional rotation, about the Z-axis, anticlockwise as seen looking towards the origin from positive Z.- The additional rotation can be represented by this matrix:
( +cos(psi) +sin(psi) 0 ) ( -sin(psi) +cos(psi) 0 ) ( 0 0 1 )
Copy a
p
-vector.Status: vector/matrix support function.
Given:
p double[3] p-vector to be copiedReturned:
c double[3] copy
Copy a position/velocity vector.
Status: vector/matrix support function.
Given:
pv double[2][3] position/velocity vector to be copiedReturned:
c double[2][3] copyCalled:
iauCp
- Copy
p
-vector.
Copy an
r
-matrix.Status: vector/matrix support function.
Given:
r double[3][3] r-matrix to be copiedReturned:
char[] double[3][3] copyCalled:
iauCp
- Copy
p
-vector.
Extend a
p
-vector to apv
-vector by appending a zero velocity.Status: vector/matrix support function.
Given:
p double[3] p-vectorReturned:
pv double[2][3] pv-vectorCalled:
iauCp
- Copy
p
-vector.iauZp
- Zero
p
-vector.
Discard velocity component of a pv-vector.
Status: vector/matrix support function.
Given:
pv double[2][3] pv-vectorReturned:
p double[3] p-vectorCalled:
iauCp
- Copy
p
-vector.
Initialize an
r
-matrix to the identity matrix.Status: vector/matrix support function.
Returned:
r double[3][3] r-matrix
Zero a
p
-vector.Status: vector/matrix support function.
Returned:
p double[3] p-vector
Zero a
pv
-vector.Status: vector/matrix support function.
Returned:
pv double[2][3] pv-vectorCalled:
iauZp
- Zero
p
-vector.
Initialize an
r
-matrix to the null matrix.Status: vector/matrix support function.
Returned:
r double[3][3] r-matrix
Multiply two
r
-matrices.Status: vector/matrix support function.
Given:
a double[3][3] first r-matrix b double[3][3] second r-matrixReturned:
atb double[3][3] a * bNotes:
- It is permissible to re-use the same array for any of the arguments.
Called:
iauCr
- Copy
r
-matrix.
Transpose an
r
-matrix.Status: vector/matrix support function.
Given:
r double[3][3] r-matrixReturned:
rt double[3][3] transposeNotes:
- It is permissible for r and rt to be the same array.
Called:
iauCr
- Copy
r
-matrix.
Multiply a
p
-vector by anr
-matrix.Status: vector/matrix support function.
Given:
r double[3][3] r-matrix p double[3] p-vectorReturned:
rp double[3] r * pNotes:
It is permissible for p and rp to be the same array.
Called:
iauCp
- Copy
p
-vector.
Multiply a pv-vector by an
r
-matrix.Status: vector/matrix support function.
Given:
r double[3][3] r-matrix pv double[2][3] pv-vectorReturned:
rpv double[2][3] r * pvNotes:
- It is permissible for pv and rpv to be the same array.
Called:
iauRxp
- Product of
r
-matrix andp
-vector.
Multiply a
p
-vector by the transpose of anr
-matrix.Status: vector/matrix support function.
Given:
r double[3][3] r-matrix p double[3] p-vectorReturned:
trp double[3] r * pNotes:
- It is permissible for p and trp to be the same array.
Called:
iauTr
- Transpose
r
-matrix.iauRxp
- Product of
r
-matrix andp
-vector.
Multiply a pv-vector by the transpose of an
r
-matrix.Status: vector/matrix support function.
Given:
r double[3][3] r-matrix pv double[2][3] pv-vectorReturned:
trpv double[2][3] r * pvNotes:
- It is permissible for pv and trpv to be the same array.
Called:
iauTr
- Transpose
r
-matrix.iauRxpv
- Product of
r
-matrix and pv-vector.
Express an
r
-matrix as anr
-vector.Status: vector/matrix support function.
Given:
r double[3][3] rotation matrixReturned:
w double[3] rotation vector (Note 1)Notes:
- A rotation matrix describes a rotation through some angle about some arbitrary axis called the Euler axis. The “rotation vector” returned by this function has the same direction as the Euler axis, and its magnitude is the angle in radians. (The magnitude and direction can be separated by means of the function
iauPn
.)- If r is null, so is the result. If r is not a rotation matrix the result is undefined; r must be proper (i.e. have a positive determinant) and real orthogonal (inverse = transpose).
- The reference frame rotates clockwise as seen looking along the rotation vector from the origin.
Form the
r
-matrix corresponding to a givenr
-vector.Status: vector/matrix support function.
Given:
w double[3] rotation vector (Note 1)Returned:
r double[3][3] rotation matrixNotes:
- A rotation matrix describes a rotation through some angle about some arbitrary axis called the Euler axis. The “rotation vector” supplied to This function has the same direction as the Euler axis, and its magnitude is the angle in radians.
- If w is null, the unit matrix is returned.
- The reference frame rotates clockwise as seen looking along the rotation vector from the origin.
Position–angle from two
p
-vectors.Status: vector/matrix support function.
Given:
a double[3] direction of reference point b double[3] direction of point whose PA is requiredReturned (function value):
double position angle of b with respect to a (radians)Notes:
- The result is the position angle, in radians, of direction b with respect to direction a. It is in the range -pi to +pi. The sense is such that if b is a small distance “north” of a the position angle is approximately zero, and if b is a small distance “east” of a the position angle is approximately +pi/2.
- The vectors a and b need not be of unit length.
- Zero is returned if the two directions are the same or if either vector is null.
- If vector a is at a pole, the result is ill–defined.
Called:
iauPn
- Decompose
p
-vector into modulus and direction.iauPm
- Modulus of
p
-vector.iauPxp
- Vector product of two
p
-vectors.iauPmp
p
-vector minusp
-vector.iauPdp
- Scalar product of two
p
-vectors.
Position–angle from spherical coordinates.
Status: vector/matrix support function.
Given:
al double longitude of point A (e.g. RA) in radians ap double latitude of point A (e.g. Dec) in radians bl double longitude of point B bp double latitude of point BReturned (function value):
double position angle of B with respect to ANotes:
- The result is the bearing (position angle), in radians, of point B with respect to point A. It is in the range -pi to +pi. The sense is such that if B is a small distance “east” of point A, the bearing is approximately +pi/2.
- Zero is returned if the two points are coincident.
Angular separation between two
p
-vectors.Status: vector/matrix support function.
Given:
a double[3] first p-vector (not necessarily unit length) b double[3] second p-vector (not necessarily unit length)Returned (function value):
double angular separation (radians, always positive)Notes:
- If either vector is null, a zero result is returned.
- The angular separation is most simply formulated in terms of scalar product. However, this gives poor accuracy for angles near zero and \pi. The present algorithm uses both cross product and dot product, to deliver full accuracy whatever the size of the angle.
Called:
iauPxp
- Vector product of two
p
-vectors.iauPm
- Modulus of
p
-vector.iauPdp
- Scalar product of two
p
-vectors.
Angular separation between two sets of spherical coordinates.
Status: vector/matrix support function.
Given:
al double first longitude (radians) ap double first latitude (radians) bl double second longitude (radians) bp double second latitude (radians)Returned (function value):
double angular separation (radians)Called:
iauS2c
- Spherical coordinates to unit vector.
iauSepp
- Angular separation between two
p
-vectors.
p
-vector to spherical coordinates.Status: vector/matrix support function.
Given:
p double[3] p-vectorReturned:
theta double longitude angle (radians) phi double latitude angle (radians)
Notes:
- The vector p can have any magnitude; only its direction is used.
- If p is null, zero theta and phi are returned.
- At either pole, zero theta is returned.
p
-vector to spherical polar coordinates.Status: vector/matrix support function.
Given:
p double[3] p-vectorReturned:
theta double longitude angle (radians) phi double latitude angle (radians) r double radial distanceNotes:
- If p is null, zero theta, phi and r are returned.
- At either pole, zero theta is returned.
Called:
iauC2s
p
-vector to spherical.iauPm
- Modulus of
p
-vector.
Convert position/velocity from Cartesian to spherical coordinates.
Status: vector/matrix support function.
Given:
pv double[2][3] pv-vectorReturned:
theta double longitude angle (radians) phi double latitude angle (radians) r double radial distance td double rate of change of theta pd double rate of change of phi rd double rate of change of rNotes:
- If the position part of pv is null, theta, phi, td and pd are indeterminate. This is handled by extrapolating the position through unit time by using the velocity part of pv. This moves the origin without changing the direction of the velocity component. If the position and velocity components of pv are both null, zeroes are returned for all six results.
- If the position is a pole, theta, td and pd are indeterminate. In such cases zeroes are returned for all three.
Convert spherical coordinates to Cartesian.
Status: vector/matrix support function.
Given:
theta double longitude angle (radians) phi double latitude angle (radians)Returned:
c double[3] direction cosines
Convert spherical polar coordinates to
p
-vector.Status: vector/matrix support function.
Given:
theta double longitude angle (radians) phi double latitude angle (radians) r double radial distanceReturned:
p double[3] Cartesian coordinatesCalled:
iauS2c
- Spherical coordinates to unit vector.
iauSxp
- Multiply
p
-vector by scalar.
Convert position/velocity from spherical to Cartesian coordinates.
Status: vector/matrix support function.
Given:
theta double longitude angle (radians) phi double latitude angle (radians) r double radial distance td double rate of change of theta pd double rate of change of phi rd double rate of change of rReturned:
pv double[2][3] pv-vector
p
-vector inner (=scalar=dot) product.Status: vector/matrix support function.
Given:
a double[3] first p-vector b double[3] second p-vectorReturned (function value):
double a . b
Modulus of
p
-vector.Status: vector/matrix support function.
Given:
p double[3] p-vectorReturned (function value):
double modulus
p
-vector subtraction.Status: vector/matrix support function.
Given:
a double[3] first p-vector b double[3] second p-vectorReturned:
amb double[3] a - bNotes:
- It is permissible to re-use the same array for any of the arguments.
Convert a
p
-vector into modulus and unit vector.Status: vector/matrix support function.
Given:
p double[3] p-vectorReturned:
r double modulus u double[3] unit vectorNotes:
- If p is null, the result is null. Otherwise the result is a unit vector.
- It is permissible to re-use the same array for any of the arguments.
Called:
iauPm
- Modulus of
p
-vector.iauZp
- Zero
p
-vector.iauSxp
- Multiply
p
-vector by scalar.
p
-vector addition.Status: vector/matrix support function.
Given:
a double[3] first p-vector b double[3] second p-vectorReturned:
apb double[3] a + bNotes:
It is permissible to re-use the same array for any of the arguments.
p
-vector plus scaledp
-vector.Status: vector/matrix support function.
Given:
a double[3] first p-vector s double scalar (multiplier for b) b double[3] second p-vectorReturned:
apsb double[3] a + s*bNotes:
- It is permissible for any of a, b and apsb to be the same array.
Called:
iauSxp
- Multiply
p
-vector by scalar.iauPpp
p
-vector plusp
-vector.
Inner (=scalar=dot) product of two pv-vectors.
Status: vector/matrix support function.
Given:
a double[2][3] first pv-vector b double[2][3] second pv-vectorReturned:
adb double[2] a . b (see note)Notes:
- If the position and velocity components of the two pv-vectors are
(ap, av)
and(bp, bv)
, the result, a*
b, is the pair of numbers(ap * bp, ap * bv + av bp)
. The two numbers are the dot–product of the twop
-vectors and its derivative.Called:
iauPdp
- Scalar product of two
p
-vectors.
Modulus of pv-vector.
Status: vector/matrix support function.
Given:
pv double[2][3] pv-vectorReturned:
r double modulus of position component s double modulus of velocity componentCalled:
iauPm
- Modulus of
p
-vector.
Subtract one pv-vector from another.
Status: vector/matrix support function.
Given:
a double[2][3] first pv-vector b double[2][3] second pv-vectorReturned:
amb double[2][3] a - bNotes:
- It is permissible to re-use the same array for any of the arguments.
Called:
iauPmp
p
-vector minusp
-vector.
Add one pv-vector to another.
Status: vector/matrix support function.
Given:
a double[2][3] first pv-vector b double[2][3] second pv-vectorReturned:
apb double[2][3] a + bNotes:
- It is permissible to re-use the same array for any of the arguments.
Called:
iauPpp
p
-vector plusp
-vector.
Update a pv-vector.
Status: vector/matrix support function.
Given:
dt double time interval pv double[2][3] pv-vectorReturned:
upv double[2][3] p updated, v unchangedNotes:
- “Update” means “refer the position component of the vector to a new date dt time units from the existing date”.
- The time units of dt must match those of the velocity.
- It is permissible for pv and upv to be the same array.
Called:
iauPpsp
p
-vector plus scaledp
-vector.iauCp
- Copy
p
-vector.
Update a pv-vector, discarding the velocity component.
Status: vector/matrix support function.
Given:
dt double time interval pv double[2][3] pv-vectorReturned:
p double[3] p-vectorNotes:
- “Update” means “refer the position component of the vector to a new date dt time units from the existing date”.
- The time units of dt must match those of the velocity.
Outer (=vector=cross) product of two pv-vectors.
Status: vector/matrix support function.
Given:
a double[2][3] first pv-vector b double[2][3] second pv-vectorReturned:
axb double[2][3] a x bNotes:
- If the position and velocity components of the two pv-vectors are
(ap, av)
and(bp, bv)
, the result,a x b
, is the pair of vectors(ap x bp, ap x bv + av x bp)
. The two vectors are the cross–product of the twop
-vectors and its derivative.- It is permissible to re-use the same array for any of the arguments.
Called:
iauCpv
- Copy
pv
-vector.iauPxp
- Vector product of two
p
-vectors.iauPpp
p
-vector plusp
-vector.
p
-vector outer (=vector=cross) product.Status: vector/matrix support function.
Given:
a double[3] first p-vector b double[3] second p-vectorReturned:
axb double[3] a x bNotes:
- It is permissible to re-use the same array for any of the arguments.
Multiply a pv-vector by two scalars.
Status: vector/matrix support function.
Given:
s1 double scalar to multiply position component by s2 double scalar to multiply velocity component by pv double[2][3] pv-vectorReturned:
spv double[2][3] pv-vector: p scaled by s1, v scaled by s2Notes:
- It is permissible for pv and spv to be the same array.
Called:
iauSxp
- Multiply
p
-vector by scalar.
Multiply a
p
-vector by a scalar.Status: vector/matrix support function.
Given:
s double scalar p double[3] p-vectorReturned:
sp double[3] s * pNotes:
- It is permissible for p and sp to be the same array.
Multiply a pv-vector by a scalar.
Status: vector/matrix support function.
Given:
s double scalar pv double[2][3] pv-vectorReturned:
spv double[2][3] s * pvNotes:
- It is permissible for pv and psv to be the same array.
Called:
iauS2xpv
- Multiply pv-vector by two scalars.
These must be used exactly as presented below.
DOUBLE PRECISION DPI PARAMETER ( DPI = 3.141592653589793238462643D0 )
DOUBLE PRECISION D2PI PARAMETER ( D2PI = 6.283185307179586476925287D0 )
DOUBLE PRECISION DR2H PARAMETER ( DR2H = 3.819718634205488058453210D0 )
DOUBLE PRECISION DR2S PARAMETER ( DR2S = 13750.98708313975701043156D0 )
DOUBLE PRECISION DR2D PARAMETER ( DR2D = 57.29577951308232087679815D0 )
DOUBLE PRECISION DR2AS PARAMETER ( DR2AS = 206264.8062470963551564734D0 )
DOUBLE PRECISION DH2R PARAMETER ( DH2R = 0.2617993877991494365385536D0 )
DOUBLE PRECISION DS2R PARAMETER ( DS2R = 7.272205216643039903848712D-5 )
DOUBLE PRECISION DD2R PARAMETER ( DD2R = 1.745329251994329576923691D-2 )
DOUBLE PRECISION DAS2R PARAMETER ( DAS2R = 4.848136811095359935899141D-6 )
The constants used by the C version of SOFA are defined in the header file sofam.h.
Truncate to nearest whole number towards zero (double).
Current membership:
Past members:
The e-mail for the Board chair is Catherine.Hohenkerk@ukho.gov.uk.
The copyright of the SOFA Software belongs to the Standards Of Fundamental Astronomy Board of the International Astronomical Union.
By using this software you accept the following six terms and conditions which apply to its use.
iau
or sofa
or trivial modifications thereof such
as changes of case.
Note that, as originally distributed, the SOFA software is intended to be a definitive implementation of the IAU standards, and consequently third–party modifications are discouraged. All variations, no matter how minor, must be explicitly marked as such, as explained above.
Correspondence concerning SOFA software should be addressed as follows:
IAU SOFA Center
HM Nautical Almanac Office
UK Hydrographic Office
Admiralty Way
Taunton, TA1 2DN
United Kingdom
If your use of SOFA results in a publication, presentation or product, please include a citation. Evidence of your use of the SOFA libraries is necessary to ensure continuing support of the initiative.
The following is a suitable form of words: “Software Routines from the IAU SOFA Collection were used. Copyright © International Astronomical Union Standards of Fundamental Astronomy (http://www.iausofa.org)”.
dint
: api consts cdnint
: api consts cdsign
: api consts ciauA2af
: api AngleOpsiauA2tf
: api AngleOpsiauAf2a
: api AngleOpsiauAnp
: api AngleOpsiauAnpm
: api AngleOpsiauBi00
: api PrecNutPolariauBp00
: api PrecNutPolariauBp06
: api PrecNutPolariauBpn2xy
: api PrecNutPolariauC2i00a
: api PrecNutPolariauC2i00b
: api PrecNutPolariauC2i06a
: api PrecNutPolariauC2ibpn
: api PrecNutPolariauC2ixy
: api PrecNutPolariauC2ixys
: api PrecNutPolariauC2s
: api SphericalCartesianiauC2t00a
: api PrecNutPolariauC2t00b
: api PrecNutPolariauC2t06a
: api PrecNutPolariauC2tcio
: api PrecNutPolariauC2teqx
: api PrecNutPolariauC2tpe
: api PrecNutPolariauC2txy
: api PrecNutPolariauCal2jd
: api CalendarsiauCp
: api CopyExtendExtractiauCpv
: api CopyExtendExtractiauCr
: api CopyExtendExtractiauD2dtf
: api TimescalesiauD2tf
: api AngleOpsiauDat
: api TimescalesiauDtdb
: api TimescalesiauDtf2d
: api TimescalesiauEe00
: api RotationAndTimeiauEe00a
: api RotationAndTimeiauEe00b
: api RotationAndTimeiauEe06a
: api RotationAndTimeiauEect00
: api RotationAndTimeiauEform
: api Geodetic/GeocentriciauEo06a
: api PrecNutPolariauEors
: api PrecNutPolariauEpb
: api CalendarsiauEpb2jd
: api CalendarsiauEpj
: api CalendarsiauEpj2jd
: api CalendarsiauEpv00
: api EphemeridesiauEqeq94
: api RotationAndTimeiauEra00
: api RotationAndTimeiauFad03
: api FundamentalArgsiauFae03
: api FundamentalArgsiauFaf03
: api FundamentalArgsiauFaju03
: api FundamentalArgsiauFal03
: api FundamentalArgsiauFalp03
: api FundamentalArgsiauFama03
: api FundamentalArgsiauFame03
: api FundamentalArgsiauFane03
: api FundamentalArgsiauFaom03
: api FundamentalArgsiauFapa03
: api FundamentalArgsiauFasa03
: api FundamentalArgsiauFaur03
: api FundamentalArgsiauFave03
: api FundamentalArgsiauFk52h
: api StarCatalogsiauFk5hip
: api StarCatalogsiauFk5hz
: api StarCatalogsiauFw2m
: api PrecNutPolariauFw2xy
: api PrecNutPolariauGc2gd
: api Geodetic/GeocentriciauGc2gde
: api Geodetic/GeocentriciauGd2gc
: api Geodetic/GeocentriciauGd2gce
: api Geodetic/GeocentriciauGmst00
: api RotationAndTimeiauGmst06
: api RotationAndTimeiauGmst82
: api RotationAndTimeiauGst00a
: api RotationAndTimeiauGst00b
: api RotationAndTimeiauGst06
: api RotationAndTimeiauGst06a
: api RotationAndTimeiauGst94
: api RotationAndTimeiauH2fk5
: api StarCatalogsiauHfk5z
: api StarCatalogsiauIr
: api InitializationiauJd2cal
: api CalendarsiauJdcalf
: api CalendarsiauNum00a
: api PrecNutPolariauNum00b
: api PrecNutPolariauNum06a
: api PrecNutPolariauNumat
: api PrecNutPolariauNut00a
: api PrecNutPolariauNut00b
: api PrecNutPolariauNut06a
: api PrecNutPolariauNut80
: api PrecNutPolariauNutm80
: api PrecNutPolariauObl06
: api PrecNutPolariauObl80
: api PrecNutPolariauP06e
: api PrecNutPolariauP2pv
: api CopyExtendExtractiauP2s
: api SphericalCartesianiauPap
: api SeparationAndAngleiauPas
: api SeparationAndAngleiauPb06
: api PrecNutPolariauPdp
: api VectorOpsiauPfw06
: api PrecNutPolariauPlan94
: api EphemeridesiauPm
: api VectorOpsiauPmat00
: api PrecNutPolariauPmat06
: api PrecNutPolariauPmat76
: api PrecNutPolariauPmp
: api VectorOpsiauPn
: api VectorOpsiauPn00
: api PrecNutPolariauPn00a
: api PrecNutPolariauPn00b
: api PrecNutPolariauPn06
: api PrecNutPolariauPn06a
: api PrecNutPolariauPnm00a
: api PrecNutPolariauPnm00b
: api PrecNutPolariauPnm06a
: api PrecNutPolariauPnm80
: api PrecNutPolariauPom00
: api PrecNutPolariauPpp
: api VectorOpsiauPpsp
: api VectorOpsiauPr00
: api PrecNutPolariauPrec76
: api PrecNutPolariauPv2p
: api CopyExtendExtractiauPv2s
: api SphericalCartesianiauPvdpv
: api VectorOpsiauPvm
: api VectorOpsiauPvmpv
: api VectorOpsiauPvppv
: api VectorOpsiauPvstar
: api SpaceMotioniauPvu
: api VectorOpsiauPvup
: api VectorOpsiauPvxpv
: api VectorOpsiauPxp
: api VectorOpsiauRm2v
: api RotationVectorsiauRv2m
: api RotationVectorsiauRx
: api BuildRotationsiauRxp
: api MatrixVectorProductsiauRxpv
: api MatrixVectorProductsiauRxr
: api MatrixOpsiauRy
: api BuildRotationsiauRz
: api BuildRotationsiauS00
: api PrecNutPolariauS00a
: api PrecNutPolariauS00b
: api PrecNutPolariauS06
: api PrecNutPolariauS06a
: api PrecNutPolariauS2c
: api SphericalCartesianiauS2p
: api SphericalCartesianiauS2pv
: api SphericalCartesianiauS2xpv
: api VectorOpsiauSepp
: api SeparationAndAngleiauSeps
: api SeparationAndAngleiauSp00
: api PrecNutPolariauStarpm
: api StarCatalogsiauStarpv
: api SpaceMotioniauSxp
: api VectorOpsiauSxpv
: api VectorOpsiauTaitt
: api TimescalesiauTaiut1
: api TimescalesiauTaiutc
: api TimescalesiauTcbtdb
: api TimescalesiauTcgtt
: api TimescalesiauTdbtcb
: api TimescalesiauTdbtt
: api TimescalesiauTf2a
: api AngleOpsiauTf2d
: api AngleOpsiauTr
: api MatrixOpsiauTrxp
: api MatrixVectorProductsiauTrxpv
: api MatrixVectorProductsiauTttai
: api TimescalesiauTttcg
: api TimescalesiauTttdb
: api TimescalesiauTtut1
: api TimescalesiauUt1tai
: api TimescalesiauUt1tt
: api TimescalesiauUt1utc
: api TimescalesiauUtctai
: api TimescalesiauUtcut1
: api TimescalesiauXy06
: api PrecNutPolariauXys00a
: api PrecNutPolariauXys00b
: api PrecNutPolariauXys06a
: api PrecNutPolariauZp
: api InitializationiauZpv
: api InitializationiauZr
: api InitializationD2PI
: api consts cDAS2R
: api consts cDAU
: api consts cDAYSEC
: api consts cDC
: api consts cDD2R
: api consts cDJ00
: api consts cDJC
: api consts cDJM
: api consts cDJM0
: api consts cDJM00
: api consts cDJM77
: api consts cDJY
: api consts cDMAS2R
: api consts cDPI
: api consts cDR2AS
: api consts cDS2R
: api consts cDTY
: api consts cELB
: api consts cELG
: api consts cGRS80
: api consts cTDB0
: api consts cTTMTAI
: api consts cTURNAS
: api consts cWGS72
: api consts cWGS84
: api consts c