Library

Documentation for SatelliteToolboxTransformations.jl.

SatelliteToolboxTransformations.EopIau1980Type
EopIau1980{T}

Earth orientation parameters for the model IAU 1980.

Note

Each field will be an AbstractInterpolation indexed by the Julian Day. Hence, if one wants to obtain, for example, the X component of the polar motion with respect to the crust at 19 June 2018, the following can be used:

x[DateTime(2018, 6, 19, 0, 0, 0) |> datetime2julian]

Fields

  • x, y: Polar motion with respect to the crust [arcsec].
  • Δut1_utc: Irregularities of the rotation angle [s].
  • lod: Length of day offset [ms].
  • δΔψ, δΔϵ: Celestial pole offsets referred to the model IAU1980 [milliarcsec].
  • *_error: Errors in the components [same unit as the component].
source
SatelliteToolboxTransformations.EopIau2000AType
EopIau2000A{T}

Earth orientation parameters for the model IAU 2000A.

Note

Each field will be an AbstractInterpolation indexed by the Julian Day. Hence, if one want to obtain, for example, the X component of the polar motion with respect to the crust at 19 June 2018, the following can be used:

x[DateTime(2018, 6, 19, 0, 0, 0) |> datetime2julian]

Fields

  • x, y: Polar motion with respect to the crust [arcsec].
  • Δut1_utc: Irregularities of the rotation angle [s].
  • lod: Length of day offset [ms].
  • δx, δy: Celestial pole offsets referred to the model IAU2000A [milliarcsec].
  • *_error: Errors in the components [same unit as the component].
source
SatelliteToolboxTransformations.cio_iau2006Method
cio_iau2006(jd_tt::Number) -> Float64, Float64, Float64

Compute the coordinates X and Y of the Celestial Intermediate Pole (CIP) with respect to the Geocentric Celestial Reference Frame (GCRF), and the CIO locator s. The algorithm is based on the IAU-2006 theory.

The CIO locator s provides the position of the CIO on the Equator of the CIP corresponding to the kinematical definition of the non-rotation origin in the GCRS when the CIP is moving with respect to the GCRS between the reference epoch and the epoch due to precession and nutation [1](p. 214).

Returns

  • Float64: The coordinate X of the CIP w.r.t. the GCRF.
  • Float64: The coordinate Y of the CIP w.r.t. the GCRF.
  • Float64: The CIO locator s.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.compute_δΔϵ_δΔψMethod
compute_δϵ_δψ(eop_iau2000a::EopIau2000A, JD::Number) -> (Float64, Float64)

Compute the celestial pole offsets in obliquity (δΔϵ_2000) and longitude (δΔΨ_2000) [arcsec] given the IERS EOP IAU 2000A eop_iau2000a.

This function obtains those values by converting the celestial pole offsets with respect to the GCRS (δx and δy). These values are necessary in the equinox-based IAU-2006 theory.

The algorithm was obtained from [1](eq. 5.25) and [2](DPSIDEPS2000_DXDY2000).

References

  • [1]: IERS (2010). Transformation between the International Terrestrial Reference System and the Geocentric Celestial Reference System. IERS Technical Note No. 36, Chapter 5.
  • [2]: ftp://hpiers.obspm.fr/eop-pc/models/uai2000.package
source
SatelliteToolboxTransformations.ecef_to_geocentricMethod
ecef_to_geocentric(r_e::AbstractVector{T}) -> NTuple{3, float(T)}

Convert the vector r_e represented in the Earth-Centered, Earth-Fixed (ECEF) reference frame into geocentric coordinates (geocentric latitude, longitude, and distance from Earth's center).

Returns

  • float(T): Geocentric latitude [rad] ∈ [-π / 2, π / 2].
  • float(T): Longitude [rad] ∈ [-π , π].
  • float(T): Distance from Earth's center [m].
source
SatelliteToolboxTransformations.ecef_to_geodeticMethod
ecef_to_geodetic(r_e::AbstractVector; ellipsoid::Ellipsoid{T} = WGS84_ELLIPSOID) where T<:Number -> NTuple{3, T}

Convert the vector r_e [m] represented in the Earth-Centered, Earth-Fixed (ECEF) reference frame into Geodetic coordinates for a custom target ellipsoid (defaults to WGS-84).

Info

The algorithm is based in [1].

Returns

  • T: Latitude [rad].
  • T: Longitude [rad].
  • T: Altitude [m].

Reference

  • [1]: mu-blox ag (1999). Datum Transformations of GPS Positions. Application Note.
source
SatelliteToolboxTransformations.ecef_to_nedMethod
ecef_to_ned(r_ecef::AbstractVector{T1}, lat::T2, lon::T3, h::T4; translate::Bool = false) -> SVector{3, T}

Convert a vector r_ecef represented in the Earth-Centered, Earth-Fixed (ECEF) frame to the local reference frame NED (North, East, Down) at the geodetic position lat [rad], lon [rad], and h [m].

If translate is false, this function computes only the rotation between ECEF and NED. Otherwise, it will also translate the vector considering the distance between the Earth's center and NED origin.

The element type T of the returned vector is obtained by promoting T1, T2, T3, and T4 to a float.

Remarks

This algorithm was based on the information in [1].

References

source
SatelliteToolboxTransformations.fetch_iers_eopMethod
fetch_iers_eop([data_type]; kwargs...) -> EopIau1980 | EopIau2000A

Download and parse the IERS EOP C04 data. The data type is specified by data_type. Supported values are:

  • Val(:IAU1980): Get IERS EOP C04 IAU1980 data.
  • Val(:IAU2000A): Get IERS EOP C04 IAU2000A data.

If data_type is omitted, then it defaults to Val(:IAU1980).

This function returns a structure (EopIau1980 or EopIau2000A, depending on data_type) with the interpolations of the EOP parameters. Notice that the interpolation indexing is set to the Julian Day.

Keywords

  • force_download::Bool: If the EOP file exists and is less than 7 days old, it will not be downloaded again. A new download can be forced by setting this keyword to true. (Default = false)
  • url::String: URL of the EOP file.
julia> eop = fetch_iers_eop()
[ Info: Downloading file 'EOP_IAU1980.TXT' from 'https://datacenter.iers.org/data/csv/finals.all.csv'...
EopIau1980:
     Data │ Timespan
 ─────────┼──────────────────────────────────────────────
        x │ 1973-01-02T00:00:00 -- 2022-05-28T00:00:00
        y │ 1973-01-02T00:00:00 -- 2022-05-28T00:00:00
  UT1-UTC │ 1973-01-02T00:00:00 -- 2022-05-28T00:00:00
      LOD │ 1973-01-02T00:00:00 -- 2021-05-19T00:00:00
       dψ │ 1973-01-02T00:00:00 -- 2021-08-05T00:00:00
       dϵ │ 1973-01-02T00:00:00 -- 2021-08-05T00:00:00

julia> eop = fetch_iers_eop(Val(:IAU2000A))
[ Info: Downloading file 'EOP_IAU2000A.TXT' from 'https://datacenter.iers.org/data/csv/finals2000A.all.csv'...
EopIau2000A:
     Data │ Timespan
 ─────────┼──────────────────────────────────────────────
        x │ 1973-01-02T00:00:00 -- 2022-05-28T00:00:00
        y │ 1973-01-02T00:00:00 -- 2022-05-28T00:00:00
  UT1-UTC │ 1973-01-02T00:00:00 -- 2022-05-28T00:00:00
      LOD │ 1973-01-02T00:00:00 -- 2021-05-19T00:00:00
       dX │ 1973-01-02T00:00:00 -- 2021-08-05T00:00:00
       dY │ 1973-01-02T00:00:00 -- 2021-08-05T00:00:00
source
SatelliteToolboxTransformations.geocentric_to_ecefMethod
geocentric_to_ecef(lat::Number, lon::Number, r::Number) -> SVector{3, T}

Convert the geocentric coordinates (latitude lat [rad], longitude lon [rad], and distance from Earth's center r [m]) into a Earth-Centered, Earth-Fixed vector [m].

Note

The output type T is obtained by promoting the input types T1, T2, and T3 to float.

source
SatelliteToolboxTransformations.geocentric_to_geodeticMethod
geocentric_to_geodetic(ϕ_gc::Number, r::Number; ellipsoid::Ellipsoid{T} = WGS84_ELLIPSOID) where T<:Number -> T, T

Compute the geodetic latitude and altitude above the reference ellipsoid (defaults to WGS-84) from the geocentric latitude ϕ_gc (-π/2, π/2) [rad] and radius r [m]. Notice that the longitude is the same in both geocentric and geodetic coordinates.

Info

The algorithm is based in [1].

Returns

  • T: Geodetic latitude [rad].
  • T: Altitude above the reference ellipsoid (defaults to WGS-84) [m].

References

  • [1] Borkowski, K. M (1987). Transformation of geocentric to geodetic coordinates without approximations. Astrophysics and Space Science, vol. 139, pp. 1-4.
source
SatelliteToolboxTransformations.geodetic_to_ecefMethod
geodetic_to_ecef(lat::Number, lon::Number, h::Number; ellipsoid::Ellipsoid{T} = wgs84_ellipsoid) where T<:Number -> SVector{3, T}

Convert the latitude lat [rad], longitude lon [rad], and altitude h [m] above the reference ellipsoid (defaults to WGS-84) into a vector represented on the Earth-Centered, Earth-Fixed (ECEF) reference frame.

Info

The algorithm is based in [1].

Reference

  • [1]: mu-blox ag (1999). Datum Transformations of GPS Positions. Application Note.
source
SatelliteToolboxTransformations.geodetic_to_geocentricMethod
geodetic_to_geocentric(ϕ_gd::Number, h::Number; ellipsoid::Ellipsoid{T} = WGS84_ELLIPSOID) where T<:Number -> T, T

Compute the geocentric latitude and radius from the geodetic latitude ϕ_gd (-π/2, π/2) [rad] and height above the reference ellipsoid h [m] (defaults to WGS-84). Notice that the longitude is the same in both geocentric and geodetic coordinates.

Info

The algorithm is based in [1](p. 3).

Returns

  • T: Geocentric latitude [rad].
  • T: Radius from the center of the Earth [m].

References

  • [1] ISO TC 20/SC 14 N (2011). Geomagnetic Reference Models.
source
SatelliteToolboxTransformations.get_ΔatMethod
get_Δat(JD::Number) -> Float64

Get the accumulated leap seconds (ΔAT) [s] between UTC and International Atomic Time (TAI) in the given JD. This function search for ΔAT in the array ΔAT_Data.

Remarks

If JD is before _ΔAT[1, 1], then 10 will be returned. Notice that this can lead to errors.

If JD is after _ΔAT[end, 1], then _ΔAT[end, 2] will be returned, because it is not possible yet to predict when leap seconds will be added.

source
SatelliteToolboxTransformations.jd_tt_to_utcMethod
jd_tt_to_utc(JD_TT::Number, ΔAT::Number) -> Float64

Convert the Julian Day in TT JD_TT (Terrestrial Time) to the Julian Day in UTC (Terrestrial Time) using the accumulated difference ΔAT between UTC and the International Atomic Time (TAI). If no value is provided, then the leap seconds will be obtained from the table ΔAT_Data. Notice that, in this case, if a date previous to 1973 is provided, then a fixed value of 10 will be used, leading to wrong computations.

source
SatelliteToolboxTransformations.jd_ut1_to_utcMethod
jd_utc_to_ut1(JD_UTC::Number, eop::Union{EopIau1980, EopIau2000A}) -> Float64

Convert the Julian Day in UT1 JD_UT1 to the Julian Day in UTC using the accumulated difference given by the EOP Data eop (see fetch_iers_eop). Notice that the accumulated difference will be interpolated.

source
SatelliteToolboxTransformations.jd_utc_to_ttMethod
jd_utc_to_tt(JD_UTC::Number[, ΔAT::Number]) -> Float64

Convert the Julian Day in UTC JD_UTC to the Julian Day in TT (Terrestrial Time) using the accumulated difference ΔAT between UTC and the International Atomic Time (TAI). If no value is provided, then the leap seconds will be obtained from the table ΔAT_Data. Notice that, in this case, if a date previous to 1973 is provided, then a fixed value of 10 will be used, leading to wrong computations.

source
SatelliteToolboxTransformations.jd_utc_to_ut1Method
jd_utc_to_ut1(JD_UTC::Number, eop::Union{EopIau1980, EopIau2000A}) -> Float64

Convert the Julian Day in UTC JD_UTC to the Julian Day in UT1 using the accumulated difference given by the EOP Data eop (see fetch_iers_eop). Notice that the accumulated difference will be interpolated.

source
SatelliteToolboxTransformations.luni_solar_args_iau2006Method
luni_solar_args_iau2006(jd_tt::Number) -> NTuple{5, Float64}

Compute the fundamental arguments related to the luni-solar effect for the IAU-2006 theory [1](p. 211).

The returned values are in [rad].

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.mean_obliquity_iau2006Method
mean_obliquity_iau2006(jd_tt::Number) -> Float64

Compute the mean obliquity of the ecliptic [rad] using the equinox-based IAU-2006 theory in the Julian day jd_tt [Terrestiral Time].

The algorithm was obtained in [1].

Reference

  • [1]: Wallace, P. T., Capitaine, N (2006). Precession-nutation procedures consistent with IAU 2006 resolutions. Astronomy & Astrophysics.
source
SatelliteToolboxTransformations.ned_to_ecefMethod
ned_to_ecef(r_ned::AbstractVector{T1}, lat::T2, lon::T3, h::T4; translate::Bool = false) -> SVector{3, T}

Convert a vector r_ned represented in the local reference frame NED (North, East, Down) at the geodetic position lat [rad], lon [rad], and h [m] to the Earth-Centered, Earth-Fixed (ECEF) frame.

If translate is false, then this function computes only the rotation between NED and ECEF. Otherwise, it will also translate the vector considering the distance between the Earth's center and NED origin.

The element type T of the returned vector is obtained by promoting T1, T2, T3, and T4 to a float.

Remarks

This algorithm was based on the information in [1].

References

source
SatelliteToolboxTransformations.nutation_eo_iau2006Function
nutation_eo_iau2006(jd_tt::Number) -> NTuple{4, Float64}

Compute the nutation parameters and the Equation of Origins (EO) at the Julian Day jd_tt [TT] using the equinox-based 2006 IAU Theory of Nutation. Notice that one can provide corrections for the nutation in obliquity (δΔϵ_2000) [rad] and in longitude (δΔψ_2000) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

Returns

  • Float64: The mean obliquity of the ecliptic [rad].
  • Float64: The nutation in obliquity of the ecliptic [rad].
  • Float64: The nutation in longitude [rad].
  • Float64: The Equation of Origins (EO) [rad].
source
SatelliteToolboxTransformations.nutation_fk5Function
nutation_fk5(jd_tt::Number, n_max::Number = 106, nut_coefs_1980::Matrix = _IAU_1980_NUTATION_COEFFICIENTS)

Compute the nutation parameters at the Julian Day jd_tt [Terrestrial Time] using the 1980 IAU Theory of Nutation. The coefficients are nut_coefs_1980 that must be a matrix in which each line has the following syntax [1](p. 1043):

an1  an2  an3  an4  an5  Ai  Bi  Ci  Di

where the units of Ai and Ci are [0.0001"] and the units of Bi and Di are [0.0001"/JC]. The user can also specify the number of coefficients n_max that will be used when computing the nutation. If n_max is omitted, the it defaults to 106.

Returns

  • Float64: The mean obliquity of the ecliptic [rad].
  • Float64: The nutation in obliquity of the ecliptic [rad].
  • Float64: The nutation in longitude [rad].

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.planetary_args_iau2006Method
planetary_args_iau2006(jd_tt::Number) -> NTuple{9, Float64}

Compute the fundamental arguments related to the planetary effects for the IAU-2006 theory [1](p. 211).

The returned values are in [rad].

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.precession_fk5Method
precession_fk5(jd_tt::Number) -> (Float64, Float64, Float64)

Compute the angles related to the precession movement in the Julian Day jd_tt [Terrestrial Time] using the theory IAU-76/FK5.

Returns

  • NTuple{3, Float64}: The angles (ζ, Θ, z) as described in [1](p. 226-228).

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.precession_iau2006Method
precession_iau2006(jd_tt::Number) -> NTuple{3, Float64}

Compute the precession angles [rad] according to equinox-based IAU-2006 theory in the Julia day jd_tt [Terrestrial Time].

This algorithm was obtained from [1](p. 49).

References

  • [1]: IERS (2010). Transformation between the International Terrestrial Reference System and the Geocentric Celestial Reference System. IERS Technical Note No. 36, Chapter 5.
source
SatelliteToolboxTransformations.r_cirs_to_gcrf_iau2006Function
r_cirs_to_gcrf_iau2006([T, ]jd_tt::Number, δx::Number = 0, δy::Number = 0) -> T

Compute the rotation that aligns the Celestial Intermediate Reference System (CIRS) with the Geocentric Celestial Reference Frame (GCRF) at the Julian Day jd_tt [TT] and considering the IERS EOP Data δx [rad] and δy [rad] (see fetch_iers_eop). This algorithm uses the IAU-2006 theory.

The IERS EOP Data δx and δy accounts for the free-core nutation and time dependent effects of the Celestial Intermediate Pole (CIP) position with respect to the GCRF.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the CIRS frame with the GCRF frame.
source
SatelliteToolboxTransformations.r_cirs_to_tirs_iau2006Method
r_cirs_to_tirs_iau2006([T, ]jd_ut1::Number) -> T

Compute the rotation that aligns the Celestial Intermediate Reference System (CIRS) with the Terrestrial Intermediate Reference System (TIRS) at the Julian Day jd_ut1 [UT1]. This algorithm uses the IAU-2006 theory.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

The rotation that aligns the CIRS frame with the TIRS frame. The rotation representation is selected by the optional parameter T.

Remarks

The reference frames TIRS and CIRS are separated by a rotation about the Z-axis of the Earth Rotation Angle, which is the angle between the Conventional International Origin (CIO) and the Terrestrial Intermediate Origin (TIO) [1]. The latter is a reference meridian on Earth that is located about 100m away from Greenwich meridian along the equator of the Celestial Intermediate Pole (CIP) [1].

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_ecef_to_ecefMethod
r_ecef_to_ecef([T, ]ECEFo, ECEFf, jd_utc::Number, eop) -> T

Compute the rotation from an Earth-Centered, Earth-Fixed (ECEF) reference frame to another ECEF reference frame at the Julian Day [UTC] jd_utc. The rotation description that will be used is given by T, which can be DCM or Quaternion. The origin ECEF frame is selected by the input ECEFo and the destination ECEF frame is selected by the input ECEFf. The model used to compute the rotation is specified by the selection of the origin and destination frames. Currently, there are two supported models: IAU-76/FK5 and IAU-2006 with 2010 conventions.

Rotation description

The rotations that aligns the origin ECEF frame with the destination ECEF frame can be described by Direction Cosine Matrices or Quaternions. This is selected by the parameter T.

The possible values are:

  • DCM: The rotation will be described by a Direction Cosine Matrix.
  • Quaternion: The rotation will be described by a Quaternion.

If no value is specified, then it falls back to DCM.

Conversion model

The model that will be used to compute the rotation is automatically inferred given the selection of the origin and destination frames. Notice that mixing IAU-76/FK5 and IAU-2006/2010 frames is not supported.

ECEF Frame

The supported ECEF frames for both origin ECEFo and destination ECEFf are:

  • ITRF(): ECEF will be selected as the International Terrestrial Reference Frame (ITRF).
  • PEF(): ECEF will be selected as the Pseudo-Earth Fixed (PEF) reference frame.
  • TIRS(): ECEF will be selected as the Terrestrial Intermediate Reference System (TIRS).

Earth orientation parameters (EOP)

The conversion between the supported ECEF frames always depends on EOP (see fetch_iers_eop and read_iers_eop). If IAU-76/FK5 model is used, then the type of eop must be EopIau1980. Otherwise, if IAU-2006/2010 model is used, then the type of eop must be EopIau2000A.

Returns

  • T: The rotation that aligns the ECEF reference frame with the ECI reference frame.

Examples

julia> eop_iau1980 = fetch_iers_eop();

julia> r_ecef_to_ecef(PEF(), ITRF(), date_to_jd(1986, 6, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
  1.0          0.0         -4.34677e-7
 -6.29476e-13  1.0         -1.44815e-6
  4.34677e-7   1.44815e-6   1.0

julia> r_ecef_to_ecef(Quaternion, PEF(), ITRF(), date_to_jd(1986, 6, 19, 21, 35, 0), eop_IAU1980)
Quaternion{Float64}:
  + 1.0 - 7.24073e-7⋅i + 2.17339e-7⋅j + 2.17339e-7⋅k

julia> eop_IAU2000A = fetch_iers_eop(Val(:IAU2000A));

julia> r_ecef_to_ecef(TIRS(), ITRF(), date_to_jd(1986, 6, 19, 21, 35, 0), eop_IAU2000A)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
  1.0          3.08408e-11  -4.34677e-7
 -3.14703e-11  1.0          -1.44815e-6
  4.34677e-7   1.44815e-6    1.0

julia> r_ecef_to_ecef(Quaternion, TIRS(), ITRF(), date_to_jd(1986, 6, 19, 21, 35, 0), eop_IAU2000A)
Quaternion{Float64}:
  + 1.0 - 7.24073e-7⋅i + 2.17339e-7⋅j + 2.17339e-7⋅k
source
SatelliteToolboxTransformations.r_ecef_to_eciMethod
r_ecef_to_eci([T, ]ECEF, ECI, jd_utc::Number[, eop]) -> T

Compute the rotation from an Earth-Centered, Earth-Fixed (ECEF) reference frame to an Earth-Centered Inertial (ECI) reference frame at the Julian Day [UTC] jd_utc. The rotation description that will be used is given by T, which can be DCM or Quaternion. The ECEF frame is selected by the input ECEF and the ECI frame is selected by the input ECI. The possible values are listed below. The model used to compute the rotation is specified by the selection of the origin and destination frames. Currently, there are two supported models: IAU-76/FK5 and IAU-2006 with 2010 conventions (CIO and equinox approaches).

Rotation description

The rotations that aligns the ECEF with ECI can be described by Direction Cosine Matrices or Quaternions. This is selected by the parameter T. The possible values are:

  • DCM: The rotation will be described by a Direction Cosine Matrix.
  • Quaternion: The rotation will be described by a Quaternion.

If no value is specified, then it falls back to DCM.

Conversion model

The model that will be used to compute the rotation is automatically inferred given the selection of the origin and destination frames. Notice that mixing IAU-76/FK5 and IAU-2006/2010 frames is not supported.

ECEF Frame

The ECEF frame is selected by the parameter ECEF. The possible values are:

  • ITRF(): ECEF will be selected as the International Terrestrial Reference Frame (ITRF).
  • PEF(): ECEF will be selected as the Pseudo-Earth Fixed (PEF) reference frame.
  • TIRS(): ECEF will be selected as the Terrestrial Intermediate Reference System (TIRS).

ECI Frame

The ECI frame is selected by the parameter ECI. The possible values are:

  • TEME(): ECI will be selected as the True Equator Mean Equinox (TEME) reference frame.
  • TOD(): ECI will be selected as the True of Date (TOD).
  • MOD(): ECI will be selected as the Mean of Date (MOD).
  • J2000(): ECI will be selected as the J2000 reference frame.
  • GCRF(): ECI will be selected as the Geocentric Celestial Reference Frame (GCRF).
  • CIRS(): ECI will be selected as the Celestial Intermediate Reference System (CIRS).
  • ERS(): ECI will be selected as the Earth Reference System (ERS).
  • MOD06(): ECI will be selected as the Mean of Date (MOD) according to the definition in IAU-2006/2010 theory.
  • MJ2000(): ECI will be selected as the J2000 mean equatorial frame (MJ2000).
Note

The frames MOD() and MOD06() are virtually the same. However, we selected different names to make clear which theory are being used since mixing transformation between frames from IAU-76/FK5 and IAU-2006/2010 must be performed with caution.

Earth orientation parameters (EOP)

The conversion between the frames depends on EOP (see fetch_iers_eop and read_iers_eop). If IAU-76/FK5 model is used, then the type of eop must be EopIau1980. Otherwise, if IAU-2006/2010 model is used, then the type of eop must be EopIau2000A. The following table shows the requirements for EOP data given the selected frames.

ModelECEFECIEOP Data
IAU-76/FK5ITRFGCRFEOP IAU1980
IAU-76/FK5ITRFJ2000EOP IAU1980
IAU-76/FK5ITRFMODEOP IAU1980
IAU-76/FK5ITRFTODEOP IAU1980
IAU-76/FK5ITRFTEMEEOP IAU1980
IAU-76/FK5PEFGCRFEOP IAU1980
IAU-76/FK5PEFJ2000Not required¹
IAU-76/FK5PEFMODNot required¹
IAU-76/FK5PEFTODNot required¹
IAU-76/FK5PEFTEMENot required¹
IAU-2006/2010 CIO-basedITRFCIRSEOP IAU2000A
IAU-2006/2010 CIO-basedITRFGCRFEOP IAU2000A
IAU-2006/2010 CIO-basedTIRSCIRSNot required¹
IAU-2006/2010 CIO-basedTIRSGCRFNot required¹ ²
IAU-2006/2010 Equinox-basedITRFERSEOP IAU2000A
IAU-2006/2010 Equinox-basedITRFMOD06EOP IAU2000A
IAU-2006/2010 Equinox-basedITRFMJ2000EOP IAU2000A
IAU-2006/2010 Equinox-basedTIRSERSNot required¹ ³
IAU-2006/2010 Equinox-basedTIRSMOD06Not required¹ ³
IAU-2006/2010 Equinox-basedTIRSMJ2000Not required¹ ³

¹: In this case, UTC will be assumed equal to UT1 to compute the Greenwich Mean Sidereal Time. This is an approximation, but should be sufficiently accurate for some applications. Notice that, if EOP Data is provided, UT1 will be accurately computed.

²: In this case, the terms that account for the free core nutation and time dependent effects of the Celestial Intermediate Pole (CIP) position with respect to the GCRF will not be available, reducing the precision.

³: In this case, the terms that corrects the nutation in obliquity and in longitude due to the free core nutation will not be available, reducing the precision.

MOD and TOD

In this function, if EOP corrections are not provided, then MOD and TOD frames will be computed considering the original IAU-76/FK5 theory. Otherwise, the corrected frame will be used.

Returns

  • T: The rotation that aligns the ECEF reference frame with the ECI reference frame.

Examples

julia> eop_iau1980 = fetch_iers_eop(Val(:IAU1980));

julia> r_ecef_to_eci(DCM, ITRF(), GCRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267      0.78518     -0.00132979
 -0.78518      -0.619267     3.33509e-5
 -0.000797312   0.00106478   0.999999

julia> r_ecef_to_eci(ITRF(), GCRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267      0.78518     -0.00132979
 -0.78518      -0.619267     3.33509e-5
 -0.000797312   0.00106478   0.999999

julia> r_ecef_to_eci(PEF(), J2000(), date_to_jd(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619271      0.785176    -0.00133066
 -0.785177     -0.619272     3.45854e-5
 -0.000796885   0.00106622   0.999999

julia> r_ecef_to_eci(PEF(), J2000(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267      0.78518     -0.00133066
 -0.78518      -0.619267     3.45854e-5
 -0.000796879   0.00106623   0.999999

julia> r_ecef_to_eci(Quaternion, ITRF(), GCRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
Quaternion{Float64}:
  + 0.43631 - 0.000590997⋅i + 0.000305106⋅j + 0.000305106⋅k

julia> eop_iau2000a = fetch_iers_eop(Val(:IAU2000A));

julia> r_ecef_to_eci(ITRF(), GCRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau2000a)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267      0.78518     -0.00132979
 -0.78518      -0.619267     3.33516e-5
 -0.000797311   0.00106478   0.999999

julia> r_ecef_to_eci(TIRS(), GCRF(), date_to_jd(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619271      0.785176    -0.00133066
 -0.785177     -0.619272     3.45884e-5
 -0.000796885   0.00106623   0.999999

julia> r_ecef_to_eci(Quaternion, ITRF(), GCRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau2000a)
Quaternion{Float64}:
  + 0.43631 - 0.000590997⋅i + 0.000305106⋅j + 0.000305106⋅k
source
SatelliteToolboxTransformations.r_eci_to_ecefMethod
r_eci_to_ecef([T, ]ECI, ECEF, jd_utc::Number[, eop]) -> T

Compute the rotation from an Earth-Centered Inertial (ECI) reference frame to an Earth-Centered, Earth-Fixed (ECEF) reference frame at the Julian Day [UTC] jd_utc. The rotation description that will be used is given by T, which can be DCM or Quaternion. The ECI frame is selected by the input ECI and the ECEF frame is selected by the input ECEF. The possible values are listed below. The model used to compute the rotation is specified by the selection of the origin and destination frames. Currently, there are two supported models: IAU-76/FK5 and IAU-2006 with 2010 conventions (CIO and equinox approaches).

Rotation description

The rotations that aligns the ECI with ECEF can be described by Direction Cosine Matrices or Quaternions. This is selected by the parameter T. The possible values are:

  • DCM: The rotation will be described by a Direction Cosine Matrix.
  • Quaternion: The rotation will be described by a Quaternion.

If no value is specified, then it falls back to DCM.

Conversion model

The model that will be used to compute the rotation is automatically inferred given the selection of the origin and destination frames. Notice that mixing IAU-76/FK5 and IAU-2006/2010 frames is not supported.

ECI Frame

The ECI frame is selected by the parameter ECI. The possible values are:

  • TEME(): ECI will be selected as the True Equator Mean Equinox (TEME) reference frame.
  • TOD(): ECI will be selected as the True of Date (TOD).
  • MOD(): ECI will be selected as the Mean of Date (MOD).
  • J2000(): ECI will be selected as the J2000 reference frame.
  • GCRF(): ECI will be selected as the Geocentric Celestial Reference Frame (GCRF).
  • CIRS(): ECEF will be selected as the Celestial Intermediate Reference System (CIRS).
  • ERS(): ECI will be selected as the Earth Reference System (ERS).
  • MOD06(): ECI will be selected as the Mean of Date (MOD) according to the definition in IAU-2006/2010 theory.
  • MJ2000(): ECI will be selected as the J2000 mean equatorial frame (MJ2000).
Note

The frames MOD() and MOD06() are virtually the same. However, we selected different names to make clear which theory are being used since mixing transformation between frames from IAU-76/FK5 and IAU-2006/2010 must be performed with caution.

ECEF Frame

The ECEF frame is selected by the parameter ECEF. The possible values are:

  • ITRF(): ECEF will be selected as the International Terrestrial Reference Frame (ITRF).
  • PEF(): ECEF will be selected as the Pseudo-Earth Fixed (PEF) reference frame.
  • TIRS(): ECEF will be selected as the Terrestrial Intermediate Reference System (TIRS).

Earth orientation parameters (EOP)

The conversion between the frames depends on EOP (see fetch_iers_eop and read_iers_eop). If IAU-76/FK5 model is used, then the type of eop must be EopIau1980. Otherwise, if IAU-2006/2010 model is used, then the type of eop must be EopIau2000A. The following table shows the requirements for EOP data given the selected frames.

ModelECIECEFEOP Data
IAU-76/FK5GCRFITRFEOP IAU1980
IAU-76/FK5J2000ITRFEOP IAU1980
IAU-76/FK5MODITRFEOP IAU1980
IAU-76/FK5TODITRFEOP IAU1980
IAU-76/FK5TEMEITRFEOP IAU1980
IAU-76/FK5GCRFPEFEOP IAU1980
IAU-76/FK5J2000PEFNot required¹
IAU-76/FK5MODPEFNot required¹
IAU-76/FK5TODPEFNot required¹
IAU-76/FK5TEMEPEFNot required¹
IAU-2006/2010 CIO-basedCIRSITRFEOP IAU2000A
IAU-2006/2010 CIO-basedGCRFITRFEOP IAU2000A
IAU-2006/2010 CIO-basedCIRSTIRSNot required¹
IAU-2006/2010 CIO-basedGCRFTIRSNot required¹ ²
IAU-2006/2010 Equinox-basedERSTIRSEOP IAU2000A
IAU-2006/2010 Equinox-basedMOD06ITRFEOP IAU2000A
IAU-2006/2010 Equinox-basedMJ2000ITRFEOP IAU2000A
IAU-2006/2010 Equinox-basedERSTIRSNot required¹ ³
IAU-2006/2010 Equinox-basedMOD06TIRSNot required¹ ³
IAU-2006/2010 Equinox-basedMJ2000TIRSNot required¹ ³

¹: In this case, UTC will be assumed equal to UT1 to compute the Greenwich Mean Sidereal Time. This is an approximation, but should be sufficiently accurate for some applications. Notice that, if EOP Data is provided, UT1 will be accurately computed.

²: In this case, the terms that account for the free-core nutation and time dependent effects of the Celestial Intermediate Pole (CIP) position with respect to the GCRF will not be available, reducing the precision.

MOD and TOD

In this function, if EOP corrections are not provided, then MOD and TOD frames will be computed considering the original IAU-76/FK5 theory. Otherwise, the corrected frame will be used.

Returns

  • T: The rotation that aligns the ECI reference with the ECEF reference frame.

Examples

julia> eop_iau1980 = fetch_iers_eop(Val(:IAU1980));

julia> r_eci_to_ecef(DCM, GCRF(), ITRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267    -0.78518     -0.000797312
  0.78518     -0.619267     0.00106478
 -0.00132979   3.33509e-5   0.999999

julia> r_eci_to_ecef(GCRF(), ITRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267    -0.78518     -0.000797312
  0.78518     -0.619267     0.00106478
 -0.00132979   3.33509e-5   0.999999

julia> r_eci_to_ecef(J2000(), PEF(), date_to_jd(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619271    -0.785177    -0.000796885
  0.785176    -0.619272     0.00106622
 -0.00133066   3.45854e-5   0.999999

julia> r_eci_to_ecef(J2000(), PEF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267    -0.78518     -0.000796879
  0.78518     -0.619267     0.00106623
 -0.00133066   3.45854e-5   0.999999

julia> r_eci_to_ecef(Quaternion, GCRF(), ITRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau1980)
Quaternion{Float64}:
  + 0.43631 + 0.000590997⋅i - 0.000305106⋅j - 0.000305106⋅k

julia> eop_iau2000a = fetch_iers_eop(Val(:IAU2000A));

julia> r_eci_to_ecef(GCRF(), ITRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau2000a)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619267    -0.78518     -0.000797311
  0.78518     -0.619267     0.00106478
 -0.00132979   3.33516e-5   0.999999

julia> r_eci_to_ecef(GCRF(), TIRS(), date_to_jd(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 -0.619271    -0.785177    -0.000796885
  0.785176    -0.619272     0.00106623
 -0.00133066   3.45884e-5   0.999999

julia> r_eci_to_ecef(Quaternion, GCRF(), ITRF(), date_to_jd(1986, 06, 19, 21, 35, 0), eop_iau2000a)
Quaternion{Float64}:
  + 0.43631 + 0.000590997⋅i - 0.000305106⋅j - 0.000305106⋅k
source
SatelliteToolboxTransformations.r_eci_to_eciMethod
r_eci_to_eci([T, ]ECIo, ECIf, jd_utc::Number[, eop]) -> T
r_eci_to_eci([T, ]ECIo, jd_utco::Number, ECIf, jd_utcf::Number[, eop])

Compute the rotation from an Earth-Centered Inertial (ECI) reference frame to another ECI reference frame. If the origin and destination frame contain only one of date frame, then the first signature is used and jd_utc is the epoch of this frame. On the other hand, if the origin and destination frame contain two of date frame¹, e.g. TOD => MOD, then the second signature must be used in which jd_utco is the epoch of the origin frame and jd_utcf is the epoch of the destination frame.

The rotation description that will be used is given by T, which can be DCM or Quaternion. The origin ECI frame is selected by the input ECIo and the destination ECI frame is selected by the input ECIf. The model used to compute the rotation is specified by the selection of the origin and destination frames. Currently, there are two supported models: IAU-76/FK5 and IAU-2006 with 2010 conventions (CIO and equinox approaches).

¹: TEME is an of date frame.

Rotation description

The rotations that aligns the origin ECI frame with the destination ECI frame can be described by Direction Cosine Matrices or Quaternions. This is selected by the parameter T.

The possible values are:

  • DCM: The rotation will be described by a Direction Cosine Matrix.
  • Quaternion: The rotation will be described by a Quaternion.

If no value is specified, then it falls back to DCM.

Conversion model

The model that will be used to compute the rotation is automatically inferred given the selection of the origin and destination frames. Notice that mixing IAU-76/FK5 and IAU-2006/2010 frames is not supported.

ECI Frame

The supported ECI frames for both origin ECIo and destination ECIf are:

  • TEME(): ECI will be selected as the True Equator Mean Equinox (TEME) reference frame.
  • TOD(): ECI will be selected as the True of Date (TOD).
  • MOD(): ECI will be selected as the Mean of Date (MOD).
  • J2000(): ECI will be selected as the J2000 reference frame.
  • GCRF(): ECI will be selected as the Geocentric Celestial Reference Frame (GCRF).
  • CIRS(): ECEF will be selected as the Celestial Intermediate Reference System (CIRS).
  • ERS(): ECI will be selected as the Earth Reference System (ERS).
  • MOD06(): ECI will be selected as the Mean of Date (MOD) according to the definition in IAU-2006/2010 theory.
  • MJ2000(): ECI will be selected as the J2000 mean equatorial frame (MJ2000).
Note

The frames MOD() and MOD06() are virtually the same. However, we selected different names to make clear which theory are being used since mixing transformation between frames from IAU-76/FK5 and IAU-2006/2010 must be performed with caution.

Earth orientation parameters (EOP)

The conversion between the frames depends on EOP Data (see fetch_iers_eop and read_iers_eop). If IAU-76/FK5 model is used, then the type of eop must be EopIau1980. Otherwise, if IAU-2006/2010 model is used, then the type of eop must be EopIau2000A. The following table shows the requirements for EOP data given the selected frames.

ModelECIoECIfEOP DataFunction Signature
IAU-76/FK5GCRFJ2000EOP IAU1980First
IAU-76/FK5GCRFMODEOP IAU1980First
IAU-76/FK5GCRFTODEOP IAU1980First
IAU-76/FK5GCRFTEMEEOP IAU1980First
IAU-76/FK5J2000GCRFEOP IAU1980First
IAU-76/FK5J2000MODNot requiredFirst
IAU-76/FK5J2000TODNot requiredFirst
IAU-76/FK5J2000TEMENot requiredFirst
IAU-76/FK5MODGCRFEOP IAU1980First
IAU-76/FK5MODJ2000Not requiredFirst
IAU-76/FK5MODTODNot requiredSecond
IAU-76/FK5MODTEMENot requiredSecond
IAU-76/FK5TODGCRFEOP IAU1980First
IAU-76/FK5TODJ2000Not requiredFirst
IAU-76/FK5TODMODNot requiredSecond
IAU-76/FK5TODTEMENot requiredSecond
IAU-76/FK5TEMEGCRFEOP IAU1980First
IAU-76/FK5TEMEJ2000Not requiredFirst
IAU-76/FK5TEMEMODNot requiredSecond
IAU-76/FK5TEMETODNot requiredSecond
IAU-2006/2010 CIO-basedGCRFCIRSNot required¹First
IAU-2006/2010 CIO-basedCIRSCIRSNot required¹Second
IAU-2006/2010 Equinox-basedGCRFMJ2000Not requiredFirst²
IAU-2006/2010 Equinox-basedGCRFMOD06Not requiredFirst
IAU-2006/2010 Equinox-basedGCRFERSNot required³First
IAU-2006/2010 Equinox-basedMJ2000GCRFNot requiredFirst²
IAU-2006/2010 Equinox-basedMJ2000MOD06Not requiredFirst
IAU-2006/2010 Equinox-basedMJ2000ERSNot required³First
IAU-2006/2010 Equinox-basedMOD06GCRFNot requiredFirst
IAU-2006/2010 Equinox-basedMOD06MJ2000Not requiredFirst
IAU-2006/2010 Equinox-basedMOD06ERSNot required³First
IAU-2006/2010 Equinox-basedERSGCRFNot required³First
IAU-2006/2010 Equinox-basedERSMJ2000Not required³First
IAU-2006/2010 Equinox-basedERSMOD06Not required³First

¹: In this case, the terms that account for the free-core nutation and time dependent effects of the Celestial Intermediate Pole (CIP) position with respect to the GCRF will not be available, reducing the precision.

²: The transformation between GCRF and MJ2000 is a constant rotation matrix called bias. Hence, the date does not modify it. However, this signature was kept to avoid complications in the API.

³: In this case, the terms that corrects the nutation in obliquity and in longitude due to the free core nutation will not be available, reducing the precision.

MOD and TOD

In this function, if EOP corrections are not provided, then MOD and TOD frames will be computed considering the original IAU-76/FK5 theory. Otherwise, the corrected frame will be used.

Returns

  • T: The rotation that aligns the origin ECI reference frame with the destination ECI reference frame.

Examples

julia> eop_iau1980 = fetch_iers_eop(Val(:IAU1980));

julia> r_eci_to_eci(DCM, GCRF(), J2000(), date_to_jd(1986, 6, 19, 21, 35, 0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
  1.0          -4.71326e-12   1.53474e-9
  4.71332e-12   1.0          -3.53979e-9
 -1.53474e-9    3.53979e-9    1.0

julia> r_eci_to_eci(Quaternion, TEME(), GCRF(), date_to_jd(1986, 6, 19, 21, 35, 0), eop_iau1980)
Quaternion{Float64}:
  + 0.999999 + 1.83013e-5⋅i + 0.000665304⋅j + 0.000665304⋅k

julia> r_eci_to_eci(TOD(), date_to_jd(1986,6,19,21,35,0), TOD(), date_to_jd(1987,5,19,3,0,0), eop_iau1980)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 1.0          -0.000224088  -9.73787e-5
 0.000224087   1.0          -5.80065e-6
 9.738e-5      5.77883e-6    1.0

julia> r_eci_to_eci(Quaternion, TOD(), JD_J2000, MOD(), JD_J2000, eop_iau1980)
Quaternion{Float64}:
  + 1.0 - 1.40025e-5⋅i + 1.34736e-5⋅j + 1.34736e-5⋅k

julia> r_eci_to_eci(J2000(), TEME(), date_to_jd(1986,6,19,21,35,0))
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
  0.999995    0.0030265    0.00133055
 -0.00302645  0.999995    -3.86125e-5
 -0.00133066  3.45854e-5   0.999999

julia> eop_iau2000a = fetch_iers_eop(Val(:IAU2000A));

julia> r_eci_to_eci(CIRS(), GCRF(), date_to_jd(1986,6,19,21,35,0), eop_iau2000a)
3×3 StaticArrays.SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 0.999999     3.88389e-8  -0.00133066
 7.18837e-9   1.0          3.45897e-5
 0.00133066  -3.45897e-5   0.999999

julia> r_eci_to_eci(Quaternion, CIRS(), GCRF(), date_to_jd(1986,6,19,21,35,0), eop_iau2000a)
Quaternion{Float64}:
  + 1.0 + 1.72949e-5⋅i + 0.000665332⋅j + 0.000665332⋅k
source
SatelliteToolboxTransformations.r_ers_to_mod_iau2006Function
r_ers_to_mod_iau2006([T, ]jd_tt::Number, δΔϵ_2000::Number = 0, δΔΨ_2000::Number = 0) -> T

Compute the rotation that aligns the Earth Reference System (ERS) with the Mean of Date (MOD) reference frame at Julian day jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

Notice that one can provide corrections for the nutation in obliquity (δΔϵ_2000) and in longitude (δΔψ_2000) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop and compute_δΔϵ_δΔψ). This corrections are related to Free Core Nutation (FCN) that models the effect of a liquid Earth core.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Info

The reference systems ERS and MOD are separated by the nutation of the pole.

Returns

  • T: The rotation that aligns the ERS frame with the MOD frame.
source
SatelliteToolboxTransformations.r_ers_to_tirs_iau2006Function
r_ers_to_tirs_iau2006(jd_ut1::Number, jd_tt::Number, δΔΨ_2000::Number = 0) -> T

Compute the rotation that aligns the Earth Reference System (ERS) with the Terrestrial Intermediate Reference System (TIRS) at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

Notice that one can provide corrections for the nutation in longitude (δΔψ_2000) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop and compute_δΔϵ_δΔψ). This corrections are related to Free Core Nutation (FCN) that models the effect of a liquid Earth core.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the ERS frame with the TIRS frame.

Remarks

The reference frames TIRS and ERS are separated by a rotation about the Z-axis of the Greenwhich apparent sidereal angle (GAST). This angle is computed using the IAU-2006 theory, which consist of obtaining the Earth Rotation Angle (ERA) and subtracting the result of the Equation of Origins (EO).

source
SatelliteToolboxTransformations.r_gcrf_to_cirs_iau2006Function
r_gcrf_to_cirs_iau2006([T, ]jd_tt::Number, δx::Number = 0, δy::Number = 0) -> T

Compute the rotation that aligns the Geocentric Celestial Reference Frame (GCRF) with the Celestial Intermediate Reference System (CIRS) at the Julian Day jd_tt [TT] and considering the IERS EOP Data δx [rad] and δy [rad] (see fetch_iers_eop). This algorithm uses the IAU-2006 theory.

The IERS EOP Data δx and δy accounts for the free-core nutation and time dependent effects of the Celestial Intermediate Pole (CIP) position with respect to the GCRF.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the GCRF frame with the CIRS frame.
source
SatelliteToolboxTransformations.r_gcrf_to_itrf_fk5Function
r_gcrf_to_itrf_fk5([T, ]jd_ut1::Number, jd_tt::Number, x_p::Number, y_p::Number[, δΔϵ_1980::Number, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the Geocentric Celestial Reference Frame (GCRF) with the International Terrestrial Reference Frame (ITRF) at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time], and considering the IERS EOP Data x_p [rad], y_p [rad], δΔϵ_1980 [rad], and δΔψ_1980 [rad] (see fetch_iers_eop). This algorithm uses the IAU-76/FK5 theory.

x_p is the polar motion displacement about X-axis, which is the IERS Reference Meridian direction (positive south along the 0˚ longitude meridian). y_p is the polar motion displacement about Y-axis (90˚W or 270˚E meridian). δΔϵ_1980 is the nutation in obliquity. δΔψ_1980 is the nutation in longitude.

The Julian Day in UT1 is used to compute the Greenwich Mean Sidereal Time (GMST) (see jd_to_gmst), whereas the Julian Day in Terrestrial Time is used to compute the nutation in the longitude. Notice that the Julian Day in UT1 and in Terrestrial Time must be equivalent, i.e. must be related to the same instant. This function does not check this.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the GCRF frame with the ITRF frame.

Remarks

The EOP data related to the polar motion (x_p and y_p) is required, since this is the only way available to compute the conversion ITRF <=> PEF (the models are highly imprecise since the motion is still not very well understood [1]). However, the EOP data related to the nutation of the obliquity (δΔϵ_1980) and the nutation of the longitude (δΔψ_1980) can be omitted. In this case, the GCRF frame is what is usually called J2000 reference frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_gcrf_to_mj2000_iau2006Function
r_gcrf_to_mj2000_iau2006([T, ]jd_tt::Number = 0) -> T

Compute the rotation that aligns the Geocentric Celestial Reference Frame (GCRF) with the J2000 mean equatorial frame. This algorithm uses the IAU-2006 theory. Notice that this rotation is just a bias matrix that does not depend on the date. However, this function receives the argument jd_tt just to keep the API compatibility.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Info

According to [1], the frame bias that converts MJ2000 <=> GCRF is not a precise transformation for all the times.

Returns

  • T: The rotation that aligns the MJ2000 frame with the MOD frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_gcrf_to_mod_fk5Method
r_gcrf_to_mod_fk5([T, ]jd_tt::Number) -> T

Compute the rotation that aligns the Geocentric Celestial Reference Frame (GCRF) with the Mean of Date (MOD) frame at the Julian Day [Terrestrial Time] jd_tt. This algorithm uses the IAU-76/FK5 theory.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the GCRF frame with the MOD frame.

Remarks

The Geocentric Celestial Reference Frame (GCRF) is rotated into the Mean of Date (MOD) frame considering the IAU 1976 Precession model.

Notice that if the conversion MOD => TOD is performed without considering the EOP corrections, then the GCRF in this rotation is what is usually called the J2000 reference frame.

source
SatelliteToolboxTransformations.r_gcrf_to_temeFunction
r_gcrf_to_teme([T, ]jd_tt::Number[, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0]) -> T

Compute the rotation that aligns the GCRF frame with the True Equator Mean Equinox (TEME) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233). Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Info

The EOP data related to the nutation of the obliquity (δΔϵ_1980) and the nutation of the longitude (δΔψ_1980) can be omitted. In this case, the GCRF frame is what is usually called J2000 reference frame.

Returns

  • T: The rotation that aligns the GCRF frame with the TEME frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_itrf_to_gcrf_fk5Function
r_itrf_to_gcrf_fk5([T, ]jd_ut1::Number, jd_tt::Number, x_p::Number, y_p::Number[, δΔϵ_1980::Number, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the International Terrestrial Reference Frame (ITRF) with the Geocentric Celestial Reference Frame (GCRF) at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time], and considering the IERS EOP Data x_p [rad], y_p [rad], δΔϵ_1980 [rad], and δΔψ_1980 [rad] (see fetch_iers_eop). This algorithm uses the IAU-76/FK5 theory.

x_p is the polar motion displacement about X-axis, which is the IERS Reference Meridian direction (positive south along the 0˚ longitude meridian). y_p is the polar motion displacement about Y-axis (90˚W or 270˚E meridian). δΔϵ_1980 is the nutation in obliquity. δΔψ_1980 is the nutation in longitude.

The Julian Day in UT1 is used to compute the Greenwich Mean Sidereal Time (GMST) (see jd_to_gmst), whereas the Julian Day in Terrestrial Time is used to compute the nutation in the longitude. Notice that the Julian Day in UT1 and in Terrestrial Time must be equivalent, i.e. must be related to the same instant. This function does not check this.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the ITRF frame with the GCRF frame.

Remarks

The EOP data related to the polar motion (x_p and y_p) is required, since this is the only way available to compute the conversion ITRF <=> PEF (the models are highly imprecise since the motion is still not very well understood [1]). However, the EOP data related to the nutation of the obliquity (δΔϵ_1980) and the nutation of the longitude (δΔψ_1980) can be omitted. In this case, the GCRF frame is what is usually called J2000 reference frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_itrf_to_pef_fk5Method
r_itrf_to_pef_fk5([T, ]x_p::Number, y_p::Number) -> T

Compute the rotation that aligns the International Terrestrial Reference Frame (ITRF) with the Pseudo-Earth Fixed (PEF) frame considering the polar motion represented by the angles x_p [rad] and y_p [rad] that are obtained from IERS EOP Data (see fetch_iers_eop).

x_p is the polar motion displacement about X-axis, which is the IERS Reference Meridian direction (positive south along the 0˚ longitude meridian). y_p is the polar motion displacement about Y-axis (90˚W or 270˚E meridian).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the ITRF frame with the PEF frame.

Remarks

The ITRF is defined based on the International Reference Pole (IRP), which is the location of the terrestrial pole agreed by international committees [1]. The Pseudo-Earth Fixed, on the other hand, is defined based on the Earth axis of rotation, or the Celestial Intermediate Pole (CIP). Hence, PEF XY-plane contains the True Equator. Furthermore, since the recovered latitude and longitude are sensitive to the CIP, then it should be computed considering the PEF frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_itrf_to_tirs_iau2006Method
r_itrf_to_tirs_iau2006([T, ]jd_tt::Number, x_p::Number, y_p::Number) -> T

Compute the rotation that aligns the International Terrestrial Reference Frame (ITRF) with the Terrestrial Intermediate Reference System (TIRS) considering the polar motion represented by the angles x_p [rad] and y_p [rad] that are obtained from IERS EOP Data (see fetch_iers_eop).

x_p is the polar motion displacement about X-axis, which is the IERS Reference Meridian direction (positive south along the 0˚ longitude meridian). y_p is the polar motion displacement about Y-axis (90˚W or 270˚E meridian).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the ITRF frame with the TIRS frame.

Remarks

The ITRF is defined based on the International Reference Pole (IRP), which is the location of the terrestrial pole agreed by international committees [1]. The Terrestrial Intermediate Reference Frame (TIRS), on the other hand, is defined based on the Earth axis of rotation, or the Celestial Intermediate Pole (CIP). Hence, TIRS XY-plane contains the True Equator. Furthermore, since the recovered latitude and longitude are sensitive to the CIP, then it should be computed considering the TIRS frame.

The TIRS and PEF (IAU-76/FK5) are virtually the same reference frame, but according to [1] it is convenient to separate the names as the exact formulae differ.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_mj2000_to_gcrf_iau2006Function
r_mj2000_to_gcrf_iau2006([T, ]jd_tt::Number = 0) -> T

Compute the rotation that aligns the J2000 mean equatorial frame with the Geocentric Celestial Reference Frame (GCRF). This algorithm uses the IAU-2006 theory. Notice that this rotation is just a bias matrix that does not depend on the date. However, this function receives the argument jd_tt just to keep the API compatibility.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Info

According to [1], the frame bias that converts MJ2000 <=> GCRF is not a precise transformation for all the times.

Returns

  • T: The rotation that aligns the MJ2000 frame with the MOD frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_mj2000_to_mod_iau2006Method
r_mj2000_to_mod_iau2006([T, ]jd_tt::Number) -> T

Compute the rotation that aligns the J2000 mean equatorial frame with the Mean of Date (MOD) reference frame with the at Julian day jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the MJ2000 frame with the MOD frame.

Remarks

The J2000 reference frame here is not equal to the previous definition in FK5 theory. It is the reason why it is internally called MJ2000. According to [1]:

The mean equinox of J2000.0 to be considered is not the “rotational dynamical mean equinox of J2000.0” as used in the past, but the “inertial dynamical mean equinox of J2000.0” to which the recent numerical or analytical solutions refer. The latter is associated with the ecliptic in the inertial sense, which is the plane perpendicular to the angular momentum vector of the orbital motion of the Earth-Moon barycenter as computed from the velocity of the barycenter relative to an inertial system. The rotational equinox is associated with the ecliptic in the rotational sense, which is perpendicular to the angular momentum vector computed from the velocity referred to the rotating orbital plane of the Earth-Moon barycenter. (The difference between the two angular momenta is the angular momentum associated with the rotation of the orbital plane.)

References

  • [1]: IERS (2010). Transformation between the International Terrestrial Reference System and the Geocentric Celestial Reference System. IERS Technical Note No. 36, Chapter 5.
source
SatelliteToolboxTransformations.r_mod_to_ers_iau2006Function
r_mod_to_ers_iau2006([T, ]jd_tt::Number, δΔϵ_2000::Number = 0, δΔΨ_2000::Number = 0) -> T

Compute the rotation that aligns the Mean of Date (MOD) reference frame with the Earth Reference System (ERS) at Julian day jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

Notice that one can provide corrections for the nutation in obliquity (δΔϵ_2000) and in longitude (δΔψ_2000) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop and compute_δΔϵ_δΔψ). This corrections are related to Free Core Nutation (FCN) that models the effect of a liquid Earth core.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the MOD frame with the ERS frame.
source
SatelliteToolboxTransformations.r_mod_to_gcrf_fk5Method
r_mod_to_gcrf_fk5([T, ]jd_tt::Number) -> T

Compute the rotation that aligns the Mean of Date (MOD) frame with the Geocentric Celestial Reference Frame (GCRF) at the Julian Day [Terrestrial Time] jd_tt. This algorithm uses the IAU-76/FK5 theory.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the MOD frame with the GCRF frame.

Remarks

The Mean of Date (MOD) frame is rotated into the Geocentric Celestial Reference Frame (GCRF) considering the IAU 1976 Precession model.

Notice that if the conversion TOD => MOD is performed without considering the EOP corrections, then the GCRF obtained by this rotation is what is usually called the J2000 reference frame.

source
SatelliteToolboxTransformations.r_mod_to_mj2000_iau2006Method
r_mod_to_mj2000_iau2006([T, ]jd_tt::Number) -> T

Compute the rotation that aligns the Mean of Date (MOD) reference frame with the J2000 mean equatorial frame at Julian day jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the MOD frame with the MJ2000 frame.

Remarks

The J2000 reference frame here is not equal to the previous definition in FK5 theory. It is the reason why it is internally called MJ2000. According to [1]:

The mean equinox of J2000.0 to be considered is not the “rotational dynamical mean equinox of J2000.0” as used in the past, but the “inertial dynamical mean equinox of J2000.0” to which the recent numerical or analytical solutions refer. The latter is associated with the ecliptic in the inertial sense, which is the plane perpendicular to the angular momentum vector of the orbital motion of the Earth-Moon barycenter as computed from the velocity of the barycenter relative to an inertial system. The rotational equinox is associated with the ecliptic in the rotational sense, which is perpendicular to the angular momentum vector computed from the velocity referred to the rotating orbital plane of the Earth-Moon barycenter. (The difference between the two angular momenta is the angular momentum associated with the rotation of the orbital plane.)

References

  • [1]: IERS (2010). Transformation between the International Terrestrial Reference System and the Geocentric Celestial Reference System. IERS Technical Note No. 36, Chapter 5.
source
SatelliteToolboxTransformations.r_mod_to_pef_fk5Function
r_mod_to_pef_fk5([T, ]jd_ut1::Number, jd_tt::Number[, δΔϵ_1980::Number, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the Mean of Date (MOD) reference frame with the Pseudo-Earth Fixed (PEF) frame at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory. Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The Julian Day in UT1 is used to compute the Greenwich Mean Sidereal Time (GMST) (see jd_to_gmst), whereas the Julian Day in Terrestrial Time is used to compute the nutation in the longitude. Notice that the Julian Day in UT1 and in Terrestrial Time must be equivalent, i.e. must be related to the same instant. This function does not check this.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the MOD frame with the PEF frame.
source
SatelliteToolboxTransformations.r_mod_to_temeFunction
r_mod_to_teme([T, ]jd_tt::Number[, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0]) -> T

Compute the rotation that aligns the Mean of Date (MOD) frame with the True Equator Mean Equinox (TEME) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233). Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop). .

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the MOD frame with the TEME frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_mod_to_tirs_iau2006Function
r_mod_to_tirs_iau2006([T, ]jd_ut1::Number, jd_tt::Number, δΔϵ_2000::Number = 0, δΔΨ_2000::Number = 0) -> T

Compute the rotation that aligns the Mean of Date (MOD) reference frame with the Terrestrial Intermediate Reference System (TIRS) at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

Notice that one can provide corrections for the nutation in obliquity (δΔϵ_2000) and in longitude (δΔψ_2000) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop and compute_δΔϵ_δΔψ). This corrections are related to Free Core Nutation (FCN) that models the effect of a liquid Earth core.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Info

This composed rotation TIRS <=> ERS <=> MOD is implemented as a new function because the single rotations TIRS <=> ERS and ERS <=> MOD call the function nutation_eo, which has a high computational burden. In this case, the composed algorithm is about 2x faster than calling those function separately.

Returns

  • T: The rotation that aligns the TIRS frame with the ERS frame.
source
SatelliteToolboxTransformations.r_mod_to_tod_fk5Function
r_mod_to_tod_fk5([T, ]jd_tt::Number[, δΔϵ_1980::Number, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the Mean of Date (MOD) frame with the True of Date (TOD) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory. Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the MOD frame with the TOD frame.

Remarks

The Mean of Date (MOD) frame is rotated into the True of Date (TOD) frame considering the 1980 IAU Theory of Nutation. The IERS EOP corrections must be added if one wants to make the rotation consistent with the Geocentric Celestial Reference Systems (GCRS).

source
SatelliteToolboxTransformations.r_pef_to_itrf_fk5Method
r_pef_to_itrf_fk5([T, ]x_p::Number, y_p::Number) -> T

Compute the rotation that aligns the Pseudo-Earth Fixed (PEF) with the International Terrestrial Reference Frame (ITRF) considering the polar motion represented by the angles x_p [rad] and y_p [rad] that are obtained from IERS EOP Data (see fetch_iers_eop).

x_p is the polar motion displacement about X-axis, which is the IERS Reference Meridian direction (positive south along the 0˚ longitude meridian). y_p is the polar motion displacement about Y-axis (90˚W or 270˚E meridian).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the PEF frame with the ITRF.

Remarks

The ITRF is defined based on the International Reference Pole (IRP), which is the location of the terrestrial pole agreed by international committees [1]. The Pseudo-Earth Fixed, on the other hand, is defined based on the Earth axis of rotation, or the Celestial Intermediate Pole (CIP). Hence, PEF XY-plane contains the True Equator. Furthermore, since the recovered latitude and longitude are sensitive to the CIP, then it should be computed considering the PEF frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_pef_to_mod_fk5Function
r_pef_to_mod_fk5([T, ]jd_ut1::Number, jd_tt::Number[, δΔϵ_1980::Number, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the Pseudo-Earth Fixed (PEF) frame with the Mean of Date (MOD) at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory. Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The Julian Day in UT1 is used to compute the Greenwich Mean Sidereal Time (GMST) (see jd_to_gmst), whereas the Julian Day in Terrestrial Time is used to compute the nutation in the longitude. Notice that the Julian Day in UT1 and in Terrestrial Time must be equivalent, i.e. must be related to the same instant. This function does not check this.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the PEF frame with the TOD frame.
source
SatelliteToolboxTransformations.r_pef_to_temeMethod
r_pef_to_teme([T, ]jd_tt::Number) -> T

Compute the rotation that aligns the Pseudo-Earth Fixed (PEF) frame with the True Equator Mean Equinox (TEME) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the PEF frame with the TEME frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_pef_to_tod_fk5Function
r_pef_to_tod_fk5([T, ]jd_ut1::Number, jd_tt::Number[, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the Pseudo-Earth Fixed (PEF) frame with the True of Date (TOD) frame at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory. Notice that one can provide correction for the nutation in longitude (δΔψ_1980) [rad] that is usually obtained from IERS EOP Data (see fetch_iers_eop).

The Julian Day in UT1 is used to compute the Greenwich Mean Sidereal Time (GMST) (see jd_to_gmst), whereas the Julian Day in Terrestrial Time is used to compute the nutation in the longitude. Notice that the Julian Day in UT1 and in Terrestrial Time must be equivalent, i.e. must be related to the same instant. This function does not check this.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the PEF frame with the TOD frame.

Remarks

The Pseudo-Earth Fixed (PEF) frame is rotated into the True of Date (TOD) frame considering the 1980 IAU Theory of Nutation. The IERS EOP corrections must be added if one wants to make the rotation consistent with the Geocentric Celestial Reference Systems (GCRS).

source
SatelliteToolboxTransformations.r_teme_to_gcrfFunction
r_teme_to_gcrf([T, ]jd_tt::Number[, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0]) -> T

Compute the rotation that aligns the True Equator Mean Equinox (TEME) frame with the Geocentric Celestial Reference Frame (GCRF) at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233). Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Info

The EOP data related to the nutation of the obliquity (δΔϵ_1980) and the nutation of the longitude (δΔψ_1980) can be omitted. In this case, the GCRF frame is what is usually called J2000 reference frame.

Returns

  • T: The rotation that aligns the TEME frame with the GCRF frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_teme_to_modFunction
r_teme_to_mod([T, ]jd_tt::Number[, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0]) -> T

Compute the rotation that aligns the True Equator Mean Equinox (TEME) frame with the Mean of Date (MOD) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233). Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TEME frame with the MOD frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_teme_to_pefMethod
r_teme_to_pef([T, ]jd_tt::Number) -> T

Compute the rotation that aligns the True Equator Mean Equinox (TEME) frame with the Pseudo-Earth Fixed (PEF) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TEME frame with the PEF frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_teme_to_todFunction
r_teme_to_tod([T, ]jd_tt::Number[, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0]) -> T

Compute the rotation that aligns the True Equator Mean Equinox (TEME) frame with the True of Date (TOD) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233). Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TEME frame with the TOD frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_tirs_to_cirs_iau2006Method
r_tirs_to_cirs_iau2006([T, ]jd_ut1::Number) -> T

Compute the rotation that aligns the Terrestrial Intermediate Reference System (TIRS) with the Celestial Intermediate Reference System (CIRS) at the Julian Day jd_ut1 [UT1]. This algorithm uses the IAU-2006 theory.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TIRS frame with the CIRS frame.

Remarks

The reference frames TIRS and CIRS are separated by a rotation about the Z-axis of the Earth Rotation Angle, which is the angle between the Conventional International Origin (CIO) and the Terrestrial Intermediate Origin (TIO) [1]. The latter is a reference meridian on Earth that is located about 100m away from Greenwich meridian along the equator of the Celestial Intermediate Pole (CIP) [1].

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_tirs_to_ers_iau2006Function
r_tirs_to_ers_iau2006([T, ]jd_ut1::Number, jd_tt::Number, δΔΨ_2000::Number = 0) -> T

Compute the rotation that aligns the Terrestrial Intermediate Reference System (TIRS) with the Earth Reference System (ERS) at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

Notice that one can provide corrections for the nutation in longitude (δΔψ_2000) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop and compute_δΔϵ_δΔψ). This corrections are related to Free Core Nutation (FCN) that models the effect of a liquid Earth core.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TIRS frame with the ERS frame.

Remarks

The reference frames TIRS and ERS are separated by a rotation about the Z-axis of the Greenwhich apparent sidereal angle (GAST). This angle is computed using the IAU-2006 theory, which consist of obtaining the Earth Rotation Angle (ERA) and subtracting the result of the Equation of Origins (EO).

source
SatelliteToolboxTransformations.r_tirs_to_itrf_iau2006Method
r_tirs_to_itrf_iau2006([T, ]jd_tt::Number, x_p::Number, y_p::Number) -> T

Compute the rotation that aligns the Terrestrial Intermediate Reference System (TIRS) with the International Terrestrial Reference Frame (ITRF) considering the polar motion represented by the angles x_p [rad] and y_p [rad] that are obtained from IERS EOP Data (see fetch_iers_eop).

x_p is the polar motion displacement about X-axis, which is the IERS Reference Meridian direction (positive south along the 0˚ longitude meridian). y_p is the polar motion displacement about Y-axis (90˚W or 270˚E meridian).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TIRS frame with the ITRF frame.

Remarks

The ITRF is defined based on the International Reference Pole (IRP), which is the location of the terrestrial pole agreed by international committees [1]. The Terrestrial Intermediate Reference Frame (TIRS), on the other hand, is defined based on the Earth axis of rotation, or the Celestial Intermediate Pole (CIP). Hence, TIRS XY-plane contains the True Equator. Furthermore, since the recovered latitude and longitude are sensitive to the CIP, then it should be computed considering the TIRS frame.

The TIRS and PEF (IAU-76/FK5) are virtually the same reference frame, but according to [1] it is convenient to separate the names as the exact formulae differ.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.r_tirs_to_mod_iau2006Function
r_tirs_to_mod_iau2006([T, ]jd_ut1::Number, jd_tt::Number, δΔϵ_2000::Number = 0, δΔΨ_2000::Number = 0) -> T

Compute the rotation that aligns the Terrestrial Intermediate Reference System (TIRS) with the Mean of Date (MOD) reference frame at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-2006 theory.

Notice that one can provide corrections for the nutation in obliquity (δΔϵ_2000) and in longitude (δΔψ_2000) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop and compute_δΔϵ_δΔψ). This corrections are related to Free Core Nutation (FCN) that models the effect of a liquid Earth core.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Info

This composed rotation TIRS <=> ERS <=> MOD is implemented as a new function because the single rotations TIRS <=> ERS and ERS <=> MOD call the function nutation_eo, which has a high computational burden. In this case, the composed algorithm is about 2x faster than calling those function separately.

Returns

  • T: The rotation that aligns the TIRS frame with the ERS frame.
source
SatelliteToolboxTransformations.r_tod_to_mod_fk5Function
r_tod_to_mod_fk5([T, ]jd_tt::Number[, δΔϵ_1980::Number, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the True of Date (TOD) frame with the Mean of Date (MOD) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory. Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TOD frame with the MOD frame.

Remarks

The True of Date (TOD) frame is rotated into the Mean of Date (MOD) frame considering the 1980 IAU Theory of Nutation. The IERS EOP corrections must be added if one wants to make the rotation consistent with the Geocentric Celestial Reference Systems (GCRS).

source
SatelliteToolboxTransformations.r_tod_to_pef_fk5Function
r_tod_to_pef_fk5([T, ]jd_ut1::Number, jd_tt::Number[, δΔψ_1980::Number]) -> T

Compute the rotation that aligns the True of Date (TOD) frame with the Pseudo-Earth Fixed (PEF) frame at the Julian Day jd_ut1 [UT1] and jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory. Notice that one can provide correction for the nutation in longitude (δΔψ_1980) [rad] that is usually obtained from IERS EOP Data (see fetch_iers_eop).

The Julian Day in UT1 is used to compute the Greenwich Mean Sidereal Time (GMST) (see jd_to_gmst), whereas the Julian Day in Terrestrial Time is used to compute the nutation in the longitude. Notice that the Julian Day in UT1 and in Terrestrial Time must be equivalent, i.e. must be related to the same instant. This function does not check this.

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TOD frame with the PEF frame.

Remarks

The True of Date (TOD) frame is rotated into the Pseudo-Earth Fixed (PEF) frame considering the 1980 IAU Theory of Nutation. The IERS EOP corrections must be added if one wants to make the rotation consistent with the Geocentric Celestial Reference Systems (GCRS).

source
SatelliteToolboxTransformations.r_tod_to_temeFunction
r_tod_to_teme([T, ]jd_tt::Number[, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0]) -> T

Compute the rotation that aligns the True of Date (TOD) frame with the True Equator Mean Equinox (TEME) frame at the Julian Day jd_tt [Terrestrial Time]. This algorithm uses the IAU-76/FK5 theory and TEME definition in [1](p. 233). Notice that one can provide corrections for the nutation in obliquity (δΔϵ_1980) [rad] and in longitude (δΔψ_1980) [rad] that are usually obtained from IERS EOP Data (see fetch_iers_eop).

The rotation type is described by the optional variable T. If it is DCM, then a DCM will be returned. Otherwise, if it is Quaternion, then a Quaternion will be returned. In case this parameter is omitted, then it falls back to DCM.

Returns

  • T: The rotation that aligns the TOD frame with the TEME frame.

References

  • [1]: Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications. Microcosm Press, Hawthorn, CA, USA.
source
SatelliteToolboxTransformations.read_iers_eopMethod
read_iers_eop(filename::String[, data_type]) -> EopIau1980 | EopIau2000A

Read IERS EOP data from the file filename. The user must specify if the data is related to the model IAU 1980 (data_type = Val(:IAU1980)), which is the default, or to the model IAU 2000A (data_type = Val(:IAU2000A)).

Note

The input file must be exactly the same as provided by IERS in CSV format. One can download it using the following commands:

  • IAU 1980

    curl -O https://datacenter.iers.org/data/csv/finals.all.csv wget https://datacenter.iers.org/data/csv/finals.all.csv

  • IAU 2000A

    curl -O https://datacenter.iers.org/data/csv/finals2000A.all.csv wget https://datacenter.iers.org/data/csv/finals2000A.all.csv

Returns

  • EopIau1980 or EopIau2000A, depending on data_type: An object with the interpolations of the EOP. Notice that the interpolation indexing is set to the Julian Day.
source
SatelliteToolboxTransformations.sv_ecef_to_eciMethod
sv_ecef_to_eci(sv::OrbitStateVector, ECEF, ECI, jd_utc[, eop]) -> OrbitStateVector

Convert the orbit state vector sv from the Earth-Centered, Earth-Fixed (ECEF) reference frame ECEF to the Earth-Centered Inertial (ECI) reference frame at the Julian day jd_utc [UTC]. The eop may be required depending on the selection of the input and output reference system. For more information, see the documentation of the function r_ecef_to_eci.

Info

It is assumed that the input velocity and acceleration in sv are obtained by an observer on the ECEF frame. Thus, the output will contain the velocity and acceleration as measured by an observer on the ECI frame.

source
SatelliteToolboxTransformations.sv_eci_to_ecefMethod
sv_eci_to_ecef(sv::OrbitStateVector, ECI, ECEF, jd_utc[, eop]) -> OrbitStateVector

Convert the orbit state vector sv from the Earth-Centered Inertial (ECI) reference frame ECI to the Earth-Centered, Earth-Fixed (ECEF) reference frame at the Julian day jd_utc [UTC]. The eop may be required depending on the selection of the input and output reference system. For more information, see the documentation of the function r_eci_to_ecef.

Info

It is assumed that the input velocity and acceleration in sv are obtained by an observer on the ECI frame. Thus, the output will contain the velocity and acceleration as measured by an observer on the ECEF frame.

source