FrameTransformations.Frames.ComputableAxesVectorType
ComputableAxesVector(from, to, order::Int)

Store the properties required to retrieve the i-th order components of a desired vector. Arguments from and to are the NAIFIDs or AbstractFramePoint instances that define the observer and target points.

Only orders between 1 and 3 are supported.

Example

julia> @point SSB 0 SolarSystemBarycenter

julia> @point Sun 10 

julia> ComputableAxesVector(SSB, Sun, 1)
ComputableAxesVector(0, 10, 1)

julia> ComputableAxesVector(0, 10, 1)
ComputableAxesVector(0, 10, 1)
FrameTransformations.Frames.FrameAxesNodeType
FrameAxesNode{O, T, N} <: AbstractJSMDGraphNode

Define a set of axes.

Fields

  • name – axes name
  • classSymbol representing the class of the axes
  • id – axes ID (equivalent of NAIFId for axes)
  • parentid – ID of the parent axes
  • comp – properties for computable axes
  • R – rotation matrix for fixed relative axes
  • fFrameAxesFunctions container
  • angles – vector storing the libration angles retrived from ephemerides
FrameTransformations.Frames.FramePointNodeType
FramePointNode{O, T, N} <: AbstractJSMDGraphNode

Define a frame system point.

Fields

  • name – point name
  • classSymbol representing the class of the point
  • axesid – ID of the axes in which the point coordinates are expressed
  • parentid – NAIF ID of the parent point
  • NAIFId – NAIF ID of the point
  • stv – vector storing the point state vectors
  • epochs – vector storing the epochs associated to stv
  • nzo – last order at which stv has been computed
  • fFramePointFunctions container
FrameTransformations.Frames.FrameSystemType
FrameSystem{O, T, S, E}

A FrameSystem instance manages a collection of user-defined FramePointNode and FrameAxesNode objects, enabling efficient computation of arbitrary transformations between them. It is created by specifying the maximum transformation order O, the outputs datatype T and an AbstractTimeScale instance S. Additionally, an AbstractEphemerisProvider instance E can be provided to compute transformations that require ephemeris data.

The following transformation orders are accepted:

  • 1: position
  • 2: position and velocity
  • 3: position, velocity and acceleration
  • 4: position, velocity, acceleration and jerk

By specifying the maximum transformation the FrameSystem memory usage and performance can be optimised and tailored to the user's needs.


FrameSystem{O, T}()

Create a FrameSystem object of order O and datatype T. The BarycentricDynamicalTime is automatically assigned as the default time scale. The resulting object is constructed with a NullEphemerisProvider, which does not allow the computation of transformation that involve ephemeris files.

Examples

julia> F = FrameSystem{2, Float64}();

julia> @axes ICRF 1 

julia> @axes ECLIPJ2000 17 

julia> add_axes_inertial!(F, ICRF)

julia> add_axes_eclipj2000!(F, ECLIPJ2000, ICRF)

julia> rotation6(F, ICRF, ECLIPJ2000, 0.0)
Rotation{2, Float64}
[...]

julia> rotation9(F, ICRF, ECLIPJ2000, 0.0)
ERROR: Insufficient frame system order: transformation requires at least order 3.

FrameSystem{O, T, S}()

Create a FrameSystem object of order O, datatype T and time scale S. The resulting object is constructed with a NullEphemerisProvider, which does not allow the computation of transformation that involve ephemeris files.

Examples

julia> F = FrameSystem{1, Float64, TerrestrialTime}();

julia> @axes ICRF 1 

julia> @axes ECLIPJ2000 17 

julia> add_axes_inertial!(F, ICRF)

julia> add_axes_eclipj2000!(F, ECLIPJ2000, ICRF)

julia> ep_tt = Epoch("2023-02-10T12:00:00 TT")
2023-02-10T12:00:00.000 TT

julia> rotation3(F, ICRF, ECLIPJ2000, ep_tt)
Rotation{1,Float64}([...])

julia> ep_tdb = Epoch("2023-02-10T12:00:00 TDB")
2023-02-10T12:00:00.000 TDB

julia> rotation3(F, ICRF, ECLIPJ2000, ep_tdb)
ERROR: ArgumentError: Incompatible epoch timescale: expected TerrestrialTime, found BarycentricDynamicalTime.
[...]

FrameSystem{O, T}(eph::AbstractEphemerisProvider)

Create a FrameSystem object of order O and datatype T by providing an instance of an AbstractEphemerisProvider subtype. The timescale is automatically set to the one associated to the ephemeris files loaded in eph. This constructor shall be used when the user desires to compute transformations that involve ephemeris data.

Note

All the kernels that will be used must be loaded within eph. Once the FrameSystem has been created, no additional kernel can be added nor removed.

Examples

julia> using Ephemerides 

julia> eph = EphemerisProvider(DE440_KERNEL_PATH);

julia> F = FrameSystem{2, Float64}(eph)
FrameSystem{2, Float64, BarycentricDynamicalTime, EphemerisProvider}(
  eph: 1-kernel EphemerisProvider,
  points: EMPTY
  axes: EMPTY
)

See also

See also add_axes_inertial!, add_point_root!, vector3 and rotation3

FrameTransformations.Frames.RotationType
Rotation{S, N}

A container to efficiently compute S-th order rotation matrices of type N between two set of axes. It stores the Direction Cosine Matrix (DCM) and its time derivatives up to the (S-1)-th order. Since this type is immutable, the data must be provided upon construction and cannot be mutated later.

The rotation of state vector between two set of axes is computed with an ad-hoc overload of the product operator. For example, a 3rd order Rotation object R, constructed from the DCM A and its time derivatives δA and δ²A rotates a vector v = [p, v, a] as:

̂v = [A*p, δA*p + A*v, δ²A*p + 2δA*v + A*a]

A Rotation object R call always be converted to a SMatrix or a MMatrix by invoking the proper constructor.

Examples

julia> A = angle_to_dcm(π/3, :Z)
DCM{Float64}:
  0.5       0.866025  0.0
 -0.866025  0.5       0.0
  0.0       0.0       1.0

julia> R = Rotation(A);

julia> SM = SMatrix(R)
3×3 SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
  0.5       0.866025  0.0
 -0.866025  0.5       0.0
  0.0       0.0       1.0

julia> MM = MMatrix(R)
3×3 MMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
  0.5       0.866025  0.0
 -0.866025  0.5       0.0
  0.0       0.0       1.0

Rotation(dcms::DCM...)

Create a Rotation object from a Direction Cosine Matrix (DCM) and any of its time derivatives. The rotation order is inferred from the number of inputs, while the rotation type is obtained by promoting the DCMs types.

Examples

julia> A = angle_to_dcm(π/3, :Z); 

julia> δA = DCM(0.0I); 

julia> δ²A = DCM(0.0I); 

julia> R = Rotation(A, δA, δ²A); 

julia> typeof(R) 
Rotation{3, Float64}

julia> R2 = Rotation(A, B, C, DCM(0.0im*I)); 

julia> typeof(R2)
Rotation{4, ComplexF64}

Rotation{S}(dcms::DCM...) where S

Create a Rotation object of order S. If the number of dcms is smaller than S, the remaining slots are filled with null DCMs, otherwise if the number of inputs is greater than S, only the first S DCMs are used.

Warning

Usage of this constructor is not recommended as it may yield unexpected results to unexperienced users.


Rotation{S1}(dcms::NTuple{S2, DCM{N}}) where {S1, S2, N}

Create a Rotation object from a tuple of Direction Cosine Matrix (DCM) and its time derivatives. If S1 < S2, only the first S1 DCMs are considered, otherwise the remaining orders are filled with null DCMs.

Examples

julia> A = angle_to_dcm(π/3, :Z);

julia> B = angle_to_dcm(π/4, π/6, :XY);

julia> R3 = Rotation{3}((A,B));

julia> R3[1] == A
true 

julia> R3[3] 
DCM{Float64}:
 0.0  0.0  0.0
 0.0  0.0  0.0
 0.0  0.0  0.0

Rotation{S}(u::UniformScaling{N}) where {S, N}
Rotation{S, N}(u::UniformScaling) where {S, N}

Create an S-order identity Rotation object of type N with identity position rotation and null time derivatives.

Examples

julia> Rotation{1}(1.0I) 
Rotation{1, Float64}(([1.0 0.0 0.0; 0.0 1.0 0.0; 0.0 0.0 1.0],))

julia> Rotation{1, Int64}(I)
Rotation{1, Int64}(([1 0 0; 0 1 0; 0 0 1],))

Rotation{S1}(rot::Rotation{S2, N}) where {S1, S2, N}
Rotation{S1, N}(R::Rotation{S2}) where {S1, S2, N}

Transform a Rotation object of order S2 to order S1 and type N. The behaviour of these functions depends on the values of S1 and S2:

  • S1 < S2: Only the first S1 components of rot are considered.
  • S1 > S2: The missing orders are filled with null DCMs.

Examples

julia> A = angle_to_dcm(π/3, :Z);

julia> B = angle_to_dcm(π/4, π/6, :XY);

julia> R1 = Rotation(A, B);

julia> order(R1)
2

julia> R2 = Rotation{1}(R1);

julia> order(R2)
1

julia> R2[1] == A 
true

julia> R3 = Rotation{3}(R1);

julia> R3[3]
DCM{Float64}:
 0.0  0.0  0.0
 0.0  0.0  0.0
 0.0  0.0  0.0

Rotation(m::DCM{N}, ω::AbstractVector) where N

Create a 2nd order Rotation object of type N to rotate between two set of axes a and b from a Direction Cosine Matrix (DCM) and the angular velocity vector ω of b with respect to a, expressed in b

See also

See also rotation3, rotation6 and rotation9.

Base.invMethod
inv(rot::Rotation)

Compute the invese of the rotation object rot. The operation is efficiently performed by taking the transpose of each rotation matrix within rot.

FrameTransformations.Frames._get_comp_axes_vector12Method
_get_comp_axes_vector12(frame, v::ComputableAxesVector, axesid, t)

Compute a 12-elements vector in the desired axes at the given time between two points of the frame system. The returned vector depends on the order in v as follows:

  • 1: position, velocity, acceleration, jerk
Warning

This function does not support orders of `v' higher than 1, because it would require the computation of vectors of order 5 and 6, which are currently not supported.

FrameTransformations.Frames._get_comp_axes_vector3Method
_get_comp_axes_vector3(frame, v::ComputableAxesVector, axesid, t)

Compute a 3-elements vector in the desired axes at the given time between two points of the frame system. The returned vector depends on the order in v as follows:

  • 1: position
  • 2: velocity
  • 3: acceleration
FrameTransformations.Frames._get_comp_axes_vector6Method
_get_comp_axes_vector6(frame, v::ComputableAxesVector, axesid, t)

Compute a 6-elements vector in the desired axes at the given time between two points of the frame system. The returned vector depends on the order in v as follows:

  • 1: position, velocity
  • 2: velocity, acceleration
  • 3: acceleration, jerk
FrameTransformations.Frames._get_comp_axes_vector9Method
_get_comp_axes_vector9(frame, v::ComputableAxesVector, axesid, t)

Compute a 9-elements vector in the desired axes at the given time between two points of the frame system. The returned vector depends on the order in v as follows:

  • 1: position, velocity, acceleration
  • 2: velocity, acceleration, jerk
Warning

This function does not support orders of `v' higher than 2, because it would require the computation of vectors of order 5, which is currently not supported.

FrameTransformations.Frames._two_vectors_basisMethod
_two_vectors_basis(a, b, seq::Symbol, fc::Function, fk::Function)

Generate a 3D right-handed orthogonal vector basis and/or its time derivatives from the vectors a and b, according to the directions specified in seq and the input cross function fc. fk is a function that filters a to guarantee type-stability.

The accepted sequence directions are: :XY, :YX, :XZ, :ZX, :YZ, :ZY

The standard basis, its 1st, 2nd-order and 3rd order time derivatives can be computed by passing cross, cross6, cross9 or cross12 to fc. The returned vectors will have a length of 3, 6 or 9, respectively.

FrameTransformations.Frames._two_vectors_to_rot6Method
_two_vectors_to_rot6(a, b, seq::Symbol)

Generate a direction cosine matrix and time-derivative, minimising the number of repeated computations.

See also

See twovectors_to_dcm and twovectors_to_δdcm for more information.

FrameTransformations.Frames._two_vectors_to_rot9Method
_two_vectors_to_rot9(a, b, seq::Symbol)

Generate a direction cosine matrix and its 1st and 2nd-order time-derivatives, minimising the number of repeated computations.

See also

See twovectors_to_dcm and twovectors_to_δ²dcm for more information.

FrameTransformations.Frames._twovectors_to_dcmMethod
_twovectors_to_dcm(a, b, seq::Symbol, fc::Function, fn::Function, fk::Function)

Generate a direction cosine matrix and/or its time derivatives from the vectors a and b, according to the directions specifeid in seq.

Notes

fc and fn are used to control the derivative order.

FrameTransformations.Frames.add_axes_bci2000!Method
add_axes_bci2000!(frames, axes::AbstractFrameAxes, center, data)

Add axes as a set of Body-Centered Inertial (BCI) axes at J2000 relative to the frames system. The center point (i.e., the reference body) is center and can either be the point ID or an AbstractFramePoint instance. data is a dictionary containing a parsed TPC file.

Warning

The parent axes are automatically set to the ICRF (ID = 1). If the ICRF is not defined in frames, an error is thrown.


add_axes_bci2000!(frames, name::Symbol, axesid::Int, cid::Int, data)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also add_axes_fixedoffset!, add_axes_bcrtod! and Orient.AXESID_ICRF.

FrameTransformations.Frames.add_axes_bcrtod!Method
add_axes_bcrtod!(frames, axes::AbstractFrameAxes, center::AbstractFramePoint, data)

Add axes as a set of Body-Centered Rotating (BCR), True-of-Date (TOD) axes to the frames system. The center point (i.e., the reference body) is center. data is a dictionary containing a parsed TPC file. These axes are the equivalent of SPICE's IAU_<BODY_NAME> frames.

Warning

The parent axes are automatically set to the ICRF (ID = 1). If the ICRF is not defined in frames, an error is thrown.


add_axes_bcrtod!(frames, name::Symbol, axesid::Int, cname::Symbol, cid::Int, data)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro and of an AbstractFramePoint via the @point macro.

See also

See also add_axes_rotating!, add_axes_bci2000! and Orient.AXESID_ICRF.

FrameTransformations.Frames.add_axes_computable!Method
add_axes_computable!(frames, axes, parent, v1, v2, seq::Symbol)

Add axes as a set of computable axes to frames. Computable axes differ from rotating axes because they are computed through two vectors that are defined within the frame system itself. Computable axes are the equivalent of SPICE's parameterized two-vector frames.

These axes are built such that the first vector, known as the primary vector, is parallel to one axis of the frame; the component of the secondary vector orthogonal to the first is parallel to another axis of the frame, and the cross product of the two vectors is parallel to the remaining axis.

The primary and secondary vectors, v1 and v2 are instances of ComputableAxesVector, which is used to define the NAIF IDs of the vector origin and target, and its order. Current accepted order values are: 1 (position), 2 (velocity) and 3 (acceleration).

For example, to define a primary vector that is parallel to the Sun's (NAIF ID = 10) velocity with respect to the Solary system barycenter (NAIF ID = 0), v1 must be set as: v1 = ComputableAxesVector(10, 0, 1).

seq is a combination of two letters that is used to identify the desired pointing directions of the primary and secondary vectors. Accepted sequences are: :XY, :YX, :XZ, :ZX, :YZ and :ZY.

Given a spacecraft registered as a point in the frame system, an example of a set of computable axes is the Local Vertical Local Horizon (LVLH), where the spacecraf's nadir direction and velocity direction define the axes orientation.

Note

Regardless of the original set of axes in which the primary and secondary vectors are defined, the axes orientation is automatically computed by rotating them to parent.


add_axes_computable!(frames, name::Symbol, axesid::Int, parentid::Int, v1, v2, seq)

Low-level function to add axes name with id axesid to frames as computable axes without requiring the creation of an AbstractFrameAxes type via the @axes macro.

Examples

julia> using Ephemerides 

julia> eph = EphemerisProvider(DE440_KERNEL_PATH);

julia> FRAMES = FrameSystem{4, Float64}(eph);

julia> @point SSB 0 SolarySystemBarycenter 

julia> @point Sun 10 

julia> @axes ICRF 1 InternationalCelestialReferenceFrame 

julia> add_axes_inertial!(FRAMES, ICRF)

julia> add_point_root!(FRAMES, SSB, ICRF)

julia> add_point_ephemeris!(FRAMES, Sun)

julia> @axes SunFrame 2

julia> v1 = ComputableAxesVector(10, 0, 1)
ComputableAxesVector(10, 0, 1)

julia> v2 = ComputableAxesVector(10, 0, 2)
ComputableAxesVector(10, 0, 2)

julia> add_axes_computable!(FRAMES, SunFrame, ICRF, v1, v2, :XY)

See also

See also ComputableAxesVector, add_axes_fixedoffset!, add_axes_inertial! and add_axes_computable!

FrameTransformations.Frames.add_axes_eclipj2000!Function
add_axes_eclipj2000!(frames, axes::AbstractFrameAxes, parent, iau_model::IAUModel=iau1980)

Add axes as a set of inertial axes representing the Ecliptic Equinox of J2000 (ECLIPJ2000) to frames. The obliquity of the ecliptic is computed using the IAU Model iau_model.

The admissed parent set of axes are the following:

  • ICRF: for the International Celestial Reference Frame, with ID = 1
  • MEME2000: the Mean Earth/Moon Ephemeris of J2000, with ID = 22

add_axes_eclipj2000!(frames, name::Symbol, parentid::Int, iau_model::IAUModel=iau1980, 
    axesid::Int = Orient.AXESID_ECLIPJ2000)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

Examples

julia> @axes ICRF 1 InternationalCelestialReferenceFrame

julia> @axes ECLIPJ2000 17 EclipticEquinoxJ2000 

julia> FRAMES = FrameSystem{1, Float64}();

julia> add_axes_inertial!(FRAMES, ICRF)

julia> add_axes_eclipj2000!(FRAMES, ECLIPJ2000, ICRF)

See also

See also add_axes_inertial! and Orient.DCM_ICRF_TO_MEME2000

FrameTransformations.Frames.add_axes_ephemeris!Method
add_axes_ephemeris!(frames, axes, rot_seq::Symbol)
add_axes_ephemeris!(frames, axes, fun, δfun=nothing, δ²fun=nothing, δ³fun=nothing)

Add axes as a set of ephemeris axes to frames. The orientation of these axes is computed with a series of 3 rotations specified by rot_seq. The euler angles and their derivatives needed by the rotation are extracted from the ephemeris kernels loaded in frames. The parent axes are automatically assigned to the axes with respect to which the orientation data has been written in the kernels.

This operation is only possible if the ephemeris kernels loaded within frames contain orientation data for the AXES ID associated to axes. An error is returned if the parent axes ID is yet to be added to frames.

The rotation sequence is defined by a Symbol specifing the rotation axes. This function assigns dcm = A3 * A2 * A1 in which Ai is the DCM related with the i-th rotation. The possible rot_seq values are: :XYX, XYZ, :XZX, :XZY, :YXY, :YXZ, :YZX, :YZY, :ZXY, :ZXZ, :ZYX, or :ZYZ.

Alternatively, one can provide custom functions that take a vector of euler angles and their derivatives and return the DCM and its derivatives.

The input functions must accept only time as argument and their outputs must be as follows:

  • fun: return a Direction Cosine Matrix (DCM).
  • δfun: return the DCM and its 1st order time derivative.
  • δ²fun: return the DCM and its 1st and 2nd order time derivatives.
  • δ³fun: return the DCM and its 1st, 2nd and 3rd order time derivatives.
Warning

It is expected that the axes ID assigned by the user are aligned with those used to generate the ephemeris kernels. No check are performed on whether these IDs represent the same physical axes that are intended in the kernels.


add_axes_ephemeris!(frames, name::Symbol, axesid::Int, rot_seq::Symbol)
add_axes_ephemeris!(frames, name::Symbol, axesid::Int, fun, δfun=nothing, δ²fun=nothing, δ³fun=nothing)

Low-level functions to add axes name with id axesid to frames as ephemeris axes without requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also ComputableAxesVector, add_axes_fixedoffset!, add_axes_inertial! and add_axes_computable!

FrameTransformations.Frames.add_axes_fixedoffset!Method
add_axes_fixedoffset!(frames::FrameSystem, axes, parent, dcm::DCM)

Add axes as a set of fixed offset axes to frames. Fixed offset axes have a constant orientation with respect to their parent axes, represented by dcm, a Direction Cosine Matrix (DCM).

Note

While inertial axes do not rotate with respect to the star background, fixed offset axes are only constant with respect to their parent axes, but might be rotating with respect to some other inertial axes.


add_axes_fixedoffset!(frames, name::Symbol, axesid::Int, parentid::Int, dcm:DCM)

Low-level function to add axes name with id axesid to frames with a fixed-offset from parentid without requiring the creation of an AbstractFrameAxes type via the @axes macro.

Examples

julia> FRAMES = FrameSystem{1, Float64}();

julia> @axes ICRF 1 InternationalCelestialReferenceFrame 

julia> add_axes_inertial!(FRAMES, ICRF)

julia> @axes ECLIPJ2000 17 

julia> add_axes_fixedoffset!(FRAMES, ECLIPJ2000, ICRF, angle_to_dcm(π/3, :Z))

See also

See also add_axes_rotating!, add_axes_inertial! and add_axes_computable!

FrameTransformations.Frames.add_axes_inertial!Method
add_axes_inertial!(frames, axes; parent=nothing, dcm=nothing)

Add axes as a set of inertial axes to frames. Only inertial axes can be used as root axes to initialise the axes graph. Only after the addition of a set of inertial axes, other axes classes may be added aswell. Once a set of root-axes has been added, parent and dcm become mandatory fields.

Note

The parent of a set of inertial axes must also be inertial. FixedOffset axes that are that only have inertial parents are also accepted.


add_axes_inertial!(frames, name::Symbol, axesid::Int; parentid=nothing, dcm=nothing)

Low-level function to add axes name with id axesid to frames without requiring the creation of an AbstractFrameAxes type via the @axes macro.

Examples

julia> FRAMES = FrameSystem{2, Float64}();

julia> @axes ICRF 1 InternationalCelestialReferenceFrame 

julia> add_axes_inertial!(FRAMES, ICRF)

julia> @axes ECLIPJ2000 17 

julia> add_axes_inertial!(FRAMES, ECLIPJ2000)
ERROR: A set of parent axes for ECLIPJ2000 is required
[...]

julia> add_axes_inertial!(FRAMES, ECLIPJ2000; parent=ICRF, dcm=angle_to_dcm(π/3, :Z))

See also

See also add_axes_rotating!, add_axes_fixedoffset! and add_axes_computable!

FrameTransformations.Frames.add_axes_itrf!Function
add_axes_itrf!(frames, axes::AbstractFrameAxes, parent, model::IAU2006Model=iau2006b)

Add axes as a set of axes representing the International Terrestrial Reference Frame (ITRF) to frames. Use the model argument to specify which IAU model model should be used for the computations. The default is set to iau2006b.

Warning

If the ID of the parent set of axes is neither the ICRF (ID = 1) nor the GCRF (ID = 23), an error is thrown.


add_axes_itrf!(frames, name::Symbol, parentid::Int, model::IAU2006Model=iau2006b, 
    axesid::Int = Orient.AXESID_ITRF)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also add_axes_rotating! and orient_rot3_itrf_to_gcrf.

FrameTransformations.Frames.add_axes_me421!Method
add_axes_me421!(frames, axes::AbstractFrameAxes, parent)

Add axes as fixed offset axes representing the DE421 Moon's Mean Earth/Mean Rotation (ME) to frames.

Warning

The parent set of axes must either the DE440 Principal Axes (PA440, ID =

  1. or the DE421 Principal Axes (PA421, ID =

31006), otherwise an error is thrown. Depending on that, the relative axes orientation will be automatically selected by this function.


add_axes_me421!(frames, name::Symbol, parentid::Int, axesid::Int = Orient.AXESID_MOONME_DE421)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also add_axes_pa440!, and add_axes_pa421!, Orient.DCM_MOON_PA421_TO_ME421 and Orient.DCM_MOON_PA421_TO_ME421,

FrameTransformations.Frames.add_axes_meme2000!Method
add_axes_meme2000!(frames, axes::AbstractFrameAxes, parent)

Add axes as a set of inertial axes representing the Mean Equator Mean Equinox of J2000 to frames.

Warning

The the axes ID of the parent set of axes must be 1 (ICRF) or 17 (ECLIPJ2000) otherwise and error is thrown.


add_axes_meme2000!(frames, name::Symbol, parentid::Int, axesid::Int = Orient.AXESID_MEME2000)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

Examples

julia> @axes ICRF 1 InternationalCelestialReferenceFrame

julia> @axes MEME2000 22 MeanEquatorMeanEquinoxJ2000 

julia> FRAMES = FrameSystem{1, Float64}();

julia> add_axes_inertial!(FRAMES, ICRF)

julia> add_axes_meme2000!(FRAMES, MEME2000, ICRF)

See also

See also add_axes_inertial! and Orient.DCM_ICRF_TO_MEME2000

FrameTransformations.Frames.add_axes_mod!Method
add_axes_mod!(frames, axes::AbstractFrameAxes, parent)

Add axes as a set of projected axes representing the Mean Equator and Equinox of Date (MOD) to frames.

Note

Despite this class of axes has a rotation matrix that depends on time, its derivatives are assumed null.

Warning

The ID of the parent set of axes must be 1 (ICRF), otherwise an error is thrown.


add_axes_mod!(frames, name::Symbol, axesid::Int, parentid::Int)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also add_axes_projected! and Orient.orient_rot3_icrf_to_mod

FrameTransformations.Frames.add_axes_pa421!Method
add_axes_pa421!(frames, axes::AbstractFrameAxes)

Add axes as a set of ephemeris axes representing the DE421 Moon's Principal Axes (PA) to frames. The libration angles are extracted from the ephemeris kernels loaded within frames, an error is thrown if such orientation data is not available.

Warning

The parent axes are automatically set to the ICRF (ID = 1). If such axes are not registered in the frame system, an error is thrown.

Warning

To properly read the ephemeris kernels, the ID associated to the input axes must match NAIF's FRAME ID for the Moon PA DE421 axes (31006).


add_axes_pa421!(frames, name::Symbol, axesid::Int = Orient.AXESID_MOONPA_DE421)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also Orient.AXESID_MOONPA_DE421, Orient.orient_rot3_icrf_to_pa421, add_axes_pa440!, and add_axes_me421!

FrameTransformations.Frames.add_axes_pa440!Method
add_axes_pa440!(frames, axes::AbstractFrameAxes)

Add axes as a set of ephemeris axes representing the DE440 Moon's Principal Axes (PA) to frames. The libration angles are extracted from the ephemeris kernels loaded within frames, an error is thrown if such orientation data is not available.

Warning

The parent axes are automatically set to the ICRF (ID = 1). If such axes are not registered in the frame system, an error is thrown.

Warning

To properly read the ephemeris kernels, the ID associated to the input axes must match NAIF's FRAME ID for the Moon PA DE440 axes (31008).


add_axes_pa440!(frames, name::Symbol, axesid::Int = Orient.AXESID_MOONPA_DE440)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also Orient.AXESID_MOONPA_DE440, Orient.orient_rot3_icrf_to_pa440 and add_axes_pa421!.

FrameTransformations.Frames.add_axes_pef!Method
add_axes_pef!(frames, axes::AbstractFrameAxes, parent)

Add axes as a set of rotating axes representing the Pseudo Earth Fixed (PEF) axes to frames.

Warning

The ID of the parent set of axes must be 1 (ICRF), otherwise an error is thrown.


add_axes_pef!(frames, name::Symbol, axesid::Int, parentid::Int)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also add_axes_rotating!, Orient.orient_rot3_icrf_to_pef and Orient.orient_rot6_icrf_to_pef.

FrameTransformations.Frames.add_axes_projected!Method
add_axes_projected!(frames, axes, parent, fun)

Add axes as a set of projected axes to frames. The orientation of these axes depends only on time and is computed through the custom functions provided by the user.

Projected axes are similar to rotating axis, except that all the positions, velocity, etc ... are rotated by the 0-order rotation (i.e. the derivatives of the rotation matrix are null, despite the rotation depends on time).

Warning

It is expected that the input function and their outputs have the correct signature. This function does not perform any checks on the output types.


add_axes_projected!(frames, name::Symbol, axesid::Int, parentid::Int, fun)

Low-level function to add axes name with id axesid to frames as projected axes without requiring the creation of an AbstractFrameAxes type via the @axes macro.

FrameTransformations.Frames.add_axes_rotating!Function
add_axes_rotating!(frames, axes, parent, fun, δfun=nothing, δ²fun=nothing, δ³fun=nothing)

Add axes as a set of rotating axes to frames. The orientation of these axes depends only on time and is computed through the custom functions provided by the user.

The input functions must accept only time as argument and their outputs must be as follows:

  • fun: return a Direction Cosine Matrix (DCM).
  • δfun: return the DCM and its 1st order time derivative.
  • δ²fun: return the DCM and its 1st and 2nd order time derivatives.
  • δ³fun: return the DCM and its 1st, 2nd and 3rd order time derivatives.

If δfun, δ²fun or δ³fun are not provided, they are computed via automatic differentiation.

Warning

It is expected that the input functions and their outputs have the correct signature. This function does not perform any checks on the output types.


add_axes_rotating!(frames, name::Symbol, axesid::Int, parentid::Int, fun, δfun=nothing, 
    δ²fun=nothing, δ³fun=nothing)

Low-level function to add axes name with id axesid to frames as rotating axes without requiring the creation of an AbstractFrameAxes type via the @axes macro.

Examples

julia> FRAMES = FrameSystem{3, Float64}();

julia> @axes Inertial 1

julia> add_axes_inertial!(FRAMES, Inertial)

julia> @axes Synodic 2 

julia> fun(t) = angle_to_dcm(t, :Z);

julia> add_axes_rotating!(FRAMES, Synodic, Inertial, fun)

julia> R = rotation6(FRAMES, Inertial, Synodic, π/6);

julia> R[1]
DCM{Float64}:
  0.866025  0.5       0.0
 -0.5       0.866025  0.0
  0.0       0.0       1.0

julia> R[2]
DCM{Float64}:
 -0.5        0.866025  0.0
 -0.866025  -0.5       0.0
  0.0        0.0       0.0

See also

See also add_axes_fixedoffset!, add_axes_inertial! and add_axes_computable!

FrameTransformations.Frames.add_axes_teme!Method
add_axes_teme!(frames, axes::AbstractFrameAxes, parent)

Add axes as a set of projected axes representing the True Equator, Mean Equinox (TEME) to frames.

Note

Despite this class of axes has a rotation matrix that depends on time, its derivatives are assumed null.

Warning

The ID of the parent set of axes must be 1 (ICRF), otherwise an error is thrown.


add_axes_teme!(frames, name::Symbol, axesid::Int, parentid::Int)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also add_axes_projected! and Orient.orient_rot3_icrf_to_teme

FrameTransformations.Frames.add_axes_tod!Method
add_axes_tod!(frames, axes::AbstractFrameAxes, parent)

Add axes as a set of projected axes representing the True Equator of Date (TOD) to frames.

Note

Despite this class of axes has a rotation matrix that depends on time, its derivatives are assumed null.

Warning

The ID of the parent set of axes must be 1 (ICRF), otherwise an error is thrown.


add_axes_tod!(frames, name::Symbol, axesid::Int, parentid::Int)

Low-level function to avoid requiring the creation of an AbstractFrameAxes type via the @axes macro.

See also

See also add_axes_projected! and Orient.orient_rot3_icrf_to_tod

FrameTransformations.Frames.add_axes_topocentric!Method
add_axes_topocentric!(frames, axes::AbstractFrameAxes, parent, λ::Number, ϕ::Number, type::Symbol)

Add axes as a set of fixed-offset topocentric axes to frames. The orientation relative to the parent axes parent is defined throuh the longitude λ, the geodetic latitude ϕ and the type type, which may be any of the following:

  • :NED (North, East, Down): the X-axis points North, the Y-axis is directed eastward and the Z-axis points inwards towards the nadir.
  • :SEZ (South, East, Zenith): the X-axis points South, the Y-axis is directed East, and the Z-axis points outwards towards the zenith.
  • :ENU (East, North, Up): the X-axis points East, the Y-axis is directed North and the Z-axis points outwards towards the zenith.
Warning

The parent axes must be a set of body-fixed reference axes. When this constraint is not satisfied, the results may be fundamentally wrong.


add_axes_topocentric!(frames, name::Symbol, axesid::Int, parentid::Int, λ, ϕ, type)

Low-level function to avoid requiring the creation of an AbstractFrameAxes via the @axes macro.

See also

See also add_axes_fixedoffset! and add_point_surface!.

FrameTransformations.Frames.add_point_dynamical!Function
add_point_dynamical!(frames, point::AbstractFramePoint, parent, axes, fun, δfun=nothing, 
    δ²fun=nothing, δ³fun=nothing)

Add point as a time point to frames. The state vector for these points depends only on time and is computed through the custom functions provided by the user. axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.

The input functions must accept only time as argument and their outputs must be as follows:

  • fun: return a 3-elements vector: position
  • δfun: return a 6-elements vector: position and velocity
  • δ²fun: return a 9-elements vector: position, velocity and acceleration
  • δ³fun: return a 12-elements vector: position, velocity, acceleration and jerk

If δfun, δ²fun or δ³fun are not provided, they are computed with automatic differentiation.

Warning

It is expected that the input functions and their ouputs have the correct signature. This function does not perform any checks on whether the returned vectors have the appropriate dimensions.


add_point_dynamical!(frames, name::Symbol, pointid::Int, parentid::Int, axes, fun, 
    δfun=nothing, δ²fun=nothing, δ³fun=nothing)

Low-level function to add point name with ID id to frames as a dynamical point without requiring the creation of an AbstractFramePoint type via the @point macro. axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.

Examples

julia> FRAMES = FrameSystem{2, Float64}()

julia> @axes ICRF 1 

julia> add_axes_inertial!(FRAMES, ICRF)

julia> @point Origin 0 

julia> add_point_root!(FRAMES, Origin, ICRF)

julia> @point Satellite 1 

julia> satellite_pos(t::T) where T = [cos(t), sin(t), 0]

julia> add_point_dynamical!(FRAMES, Satellite, Origin, ICRF, satellite_pos)

julia> vector6(FRAMES, Origin, Satellite, ICRF, π/6)
6-element SVector{6, Float64} with indices SOneTo(6):
  0.8660254037844387
  0.49999999999999994
  0.0
 -0.49999999999999994
  0.8660254037844387
  0.0

See also

See also add_point_root!, add_point_ephemeris!,add_point_fixed! and add_point_updatable!

FrameTransformations.Frames.add_point_ephemeris!Method
add_point_ephemeris!(frames::FrameSystem, point::AbstractFramePoint)

Add point as an ephemeris point to frames. This function is intended for points whose state-vector is read from ephemeris kernels (i.e., de440.bsp). The parent point is automatically assigned to the point with respect to which the ephemeris data is written in the kernels. If that point is not available as an ephemeris point in the frame system, an error is thrown.

Note

The axes in which the state-vector is expressed are taken from the ephemeris data: an error is returned if the axes ID is yet to be added to frames.

Warning

It is expected that the NAIF ID and the axes ID assigned by the user are aligned with those used to generate the ephemeris kernels. No check are performed on whether these IDs represent the same physical bodies and axes that are intended in the kernels.


add_point_ephemeris!(frames::FrameSystem, name::Symbol, naifid::Int)

Low-level function to add point name with ID id to frames as an ephemeris point without requiring the creation of an AbstractFramePoint type via the @point macro.

Examples

julia> using Ephemerides 

julia> eph = EphemerisProvider(DE440_KERNEL_PATH);

julia> FRAMES = FrameSystem{2, Float64}(eph);

julia> @axes ICRF 1 InternationalCelestialReferenceFrame

julia> add_axes_inertial!(FRAMES, ICRF)

julia> @point SSB 0 SolarSystemBarycenter

julia> @point Sun 10 

julia> add_point_root!(FRAMES, SSB, ICRF)

julia> add_point_ephemeris!(FRAMES, Sun)

julia> @point Jupiter 599

julia> add_point_ephemeris!(FRAMES, Jupiter)
ERROR: ArgumentError: Ephemeris data for NAIFID 599 is not available in the kernels loaded in the given FrameSystem.
[...]

See also

See also add_point_root!, add_point_fixed!, add_point_dynamical! and add_point_updatable!

FrameTransformations.Frames.add_point_fixed!Method
add_point_fixed!(frames, point::AbstractFramePoint, parent, axes, offset::AbstractVector)

Add point as a fixed point to frames. Fixed points are those whose positions have a constant offset with respect their parent points in the given set of axes. Thus, points eligible for this class must have null velocity and acceleration with respect to parent. axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.


add_point_fixed!(frames, name::Symbol, pointid::Int, parentid::Int, axes, offset::AbstractVector)

Low-level function to add point name with ID id to frames as a fixed-point with respect to point parentid without requiring the creation of an AbstractFramePoint type via the @point macro. axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.

Examples

julia> FRAMES = FrameSystem{2, Float64}();

julia> @axes SF -3000 SatelliteFrame

julia> add_axes_inertial!(FRAMES, SF)

julia> @point SC -10000 Spacecraft

julia> @point SolarArrayCenter -10001

julia> add_point_root!(FRAMES, SC, SF)

julia> sa_offset = [0.10, 0.15, 0.30];

julia> add_point_fixed!(FRAMES, SolarArrayCenter, SC, SF, sa_offset)

See also

See also add_point_root!, add_point_ephemeris!, add_point_dynamical! and add_point_updatable!

FrameTransformations.Frames.add_point_root!Method
add_point_root!(frames::FramesSystem, point::AbstractFramePoint, axes)

Add point as a root point to frames to initialize the points graph. Only after the addition of a root point, other points may be added aswell. This point is intended as the origin, i.e., its position will equal (0., 0., 0.). axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.

Note

This operation can be performed only once per FrameSystem object: multiple root points in the same graph are both inadmissible and meaningless.


add_point_root!(frames::FrameSystem, name::Symbol, pointid::Int, axes)

Low-level function to add point name with ID pointid to frames as a root-point without requiring the creation of an AbstractFramePoint type via the @point macro. axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.

Examples

julia> FRAMES = FrameSystem{2, Float64}();

julia> @axes ICRF 1 InternationalCelestialReferenceFrame

julia> add_axes_inertial!(FRAMES, ICRF)

julia> @point SSB 0 SolarSystemBarycenter 

julia> add_point_root!(FRAMES, SSB, ICRF)

julia> @point Sun 10

julia> add_point_root!(FRAMES, Sun, ICRF)
ERROR: ArgumentError: A root-point is already registed in the given FrameSystem.
[...]

See also

See also add_point_ephemeris!, add_point_fixed!, add_point_dynamical! and add_point_updatable!

FrameTransformations.Frames.add_point_surface!Function
add_point_surface!(frames, point::AbstractFramePoint, parent, axes, λ::Number, 
    ϕ::Number, R::Number, f::Number=0.0, h::Number=0.0)

Add point to frames as a fixed point on the surface of the parent point body. The relative position is specified by the longitude λ, the geodetic latitude ϕ, the reference radius of the ellipsoid R and its flattening f. The altitude over the reference surface of the ellipsoid h defaults to 0.

Warning

axes must be a set of body-fixed reference axes for the body represented by parent. When this constraint is not satisfied, the results may be fundamentally wrong.


add_point_surface!(frames, point::AbstractFramePoint, parent, axes, pck,  λ, ϕ, h::Number=0.0)

In this case, the ellipsoid parameters are extracted from the input TPC kernel pck using the NAIFId associated to the parent point.


add_point_surface!(frames, name::Symbol, pointid::Int, parentid::Int, axesid::Int, 
    λ::Number, ϕ::Number, R::Number, f::Number=0.0, h::Number=0.0,)

Low-level function to avoid requiring the creation of an AbstractFramePoint via the @point macro.

See also

See also add_point_fixed! and add_axes_topocentric!.

FrameTransformations.Frames.add_point_updatable!Method
add_point_updatable!(frames, point::AbstractFramePoint, parent, axes)

Add point as an updatable point to frames. Differently from all the other classes, the state vector for updatable points (expressed in the set of input axes) must be manually updated before being used for other computations. axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.

Note

This class of points becomes particularly useful if the state vector is not known a-priori, e.g., when it is the output of an optimisation process which exploits the frame system.


add_point_updatable!(frames, name::Symbol, pointid::Int, parentid::Int, axes)

Low-level function to add point name with ID id to frames as an updatable point without requiring the creation of an AbstractFramePoint type via the @point macro. axes is the ID or AbstractFrameAxes instance in which the point state-vector is expressed.

Examples

julia> FRAMES = FrameSystem{2, Float64}();

julia> @axes ICRF 1  

julia> add_axes_inertial!(FRAMES, ICRF)

julia> @point Origin 0

julia> @point Satellite 1 

julia> add_point_root!(FRAMES, Origin, ICRF)

julia> add_point_updatable!(FRAMES, Satellite, Origin, ICRF)

julia> y = [10000., 200., 300.]

julia> update_point!(FRAMES, Satellite, y, 0.1)

julia> vector3(FRAMES, Origin, Satellite, ICRF, 0.1)
3-element SVector{3, Float64} with indices SOneTo(3):
 10000.0
   200.0
   300.0

julia> vector3(FRAMES, Origin, Satellite, ICRF, 0.2)
ERROR: UpdatablePoint with NAIFId = 1 has not been updated at time 0.2 for order 1

julia> vector6(FRAMES, Origin, Satellite, ICRF, 0.1)
ERROR: UpdatablePoint with NAIFId = 1 has not been updated at time 0.1 for order 2
[...]

See also

See also update_point!, add_point_root!, add_point_ephemeris!, add_point_dynamical! and add_point_fixed!

FrameTransformations.Frames.build_axesMethod
build_axes(frames, name, id, class, funs; parentid, dcm, cax_prop)

Create and add a FrameAxesNode to frames based on the input parameters. Current supported classes are: :InertialAxes, :FixedOffsetAxes, :RotatingAxes, :ProjectedAxes, :EphemerisAxes and :ComputableAxes.

Inputs

  • frames – Target frame system
  • name – Axes name, must be unique within frames
  • id – Axes ID, must be unique within frames
  • class – Axes class.
  • funsFrameAxesFunctions object storing the functions to compute the DCM and, eventually, its time derivatives. It must match the type and order of frames.

Keywords

  • parentid – Axes ID of the parent axes. Not required only for the root axes.
  • dcm – DCM with respect to the parent axes. Required only for FixedOffsetAxes.
  • cax_propComputableAxesProperties, required only by ComputableAxes.
Warning

This is a low-level function and is NOT meant to be directly used. Instead, to add a set of axes to the frame system, see add_axes_inertial!, add_axes_rotating!, add_axes_computable!, add_axes_fixedoffset! and add_axes_projected!.

FrameTransformations.Frames.build_pointMethod
build_point(frames, name, NAIFId, class, axesid, funs; parentid, offset)

Create and add a FramePointNode to frames based on the input parameters. Current supported point classes are: :RootPoint, :TimePoint, :EphemerisPoint, :FixedPoint and :UpdatablePoint.

Inputs

  • frames – Target frame system
  • name – Point name, must be unique within frames
  • NAIFId – Point NAIF ID, must be unique within frames
  • class – Point class.
  • axesid – ID of the axes in which the state vector of the point is expressed.
  • funsFramePointFunctions object storing the functions to update the state vectors of the point. It must match the type and order of frames

Keywords

  • parentid – NAIF ID of the parent point. Not required only for the root point.
  • offset – Position offset with respect to a parent point. Required only for FixedPoints.
Warning

This is a low-level function and is NOT meant to be directly used. Instead, to add a point to the frame system, see add_point_ephemeris!, add_point_fixed!, etc...

FrameTransformations.Frames.is_inertialMethod
is_inertial(frame::FrameSystem, axes::AbstractFrameAxes)
is_inertial(frame::FrameSystem, axesid::Int)

Return true if the given axes are inertial, i.e., non rotating with respect to the root inertial axes.

Note

FixedOffsetAxes with respect to an inertial set of axes, are also consired inertial.

FrameTransformations.Frames.is_timefixedMethod
is_timefixed(frame::FrameSystem, axes::AbstractFrameAxes)
is_timefixed(frame::FrameSystem, axesid::Int)

Return true if the given axes are time-fixed, i.e., their orientation does not change in time with respect to the root inertial axes.

Note

Only :InertialAxes and :FixedOffsetAxes defined with respect to other inertial axes are here considered as time fixed.

FrameTransformations.Frames.rotation12Method
rotation12(frame::FrameSystem, from, to, ep::Epoch)

Compute the rotation that transforms a 12-elements state vector from one specified set of axes to another at a given epoch.

Requires a frame system of order ≥ 4.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the axes to transform from
  • to – ID or instance of the axes to transform to
  • epEpoch of the rotation. Its timescale must match that of the frame system.

Output

A Rotation object of order 4.

FrameTransformations.Frames.rotation12Method
rotation12(frame::FrameSystem, from, to, t::Number)

Compute the rotation that transforms a 12-elements state vector from one specified set of axes to another at a given time t, expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.rotation3Method
rotation3(frame::FrameSystem, from, to, ep::Epoch)

Compute the rotation that transforms a 3-elements state vector from one specified set of axes to another at a given epoch.

Requires a frame system of order ≥ 1.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the axes to transform from
  • to – ID or instance of the axes to transform to
  • epEpoch of the rotation. Its timescale must match that of the frame system.

Output

A Rotation object of order 1.

FrameTransformations.Frames.rotation3Method
rotation3(frame::FrameSystem, from, to, t::Number)

Compute the rotation that transforms a 3-elements state vector from one specified set of axes to another at a given time t, expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.rotation6Method
rotation6(frame::FrameSystem, from, to, ep::Epoch)

Compute the rotation that transforms a 6-elements state vector from one specified set of axes to another at a given epoch.

Requires a frame system of order ≥ 2.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the axes to transform from
  • to – ID or instance of the axes to transform to
  • epEpoch of the rotation. Its timescale must match that of the frame system.

Output

A Rotation object of order 2.

FrameTransformations.Frames.rotation6Method
rotation6(frame::FrameSystem, from, to, t::Number)

Compute the rotation that transforms a 6-elements state vector from one specified set of axes to another at a given time t, expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.rotation9Method
rotation9(frame::FrameSystem, from, to, ep::Epoch)

Compute the rotation that transforms a 9-elements state vector from one specified set of axes to another at a given epoch.

Requires a frame system of order ≥ 3.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the axes to transform from
  • to – ID or instance of the axes to transform to
  • epEpoch of the rotation. Its timescale must match that of the frame system.

Output

A Rotation object of order 3.

FrameTransformations.Frames.rotation9Method
rotation9(frame::FrameSystem, from, to, t::Number)

Compute the rotation that transforms a 9-elements state vector from one specified set of axes to another at a given time t, expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.twovectors_to_dcmMethod
twovectors_to_dcm(a, b, seq)

Generate a direction cosine matrix from two time-dependent vectors a and b, following the directions specified in seq.

Inputs

  • a – The primary vector that will be aligned with the first directions specified in seq.

  • b – The secondary vector. The component of this vector that is orthogonal to the primary vector is aligned with the second direction specified in the sequence seq.

  • seq – Accepted sequence directions are: :XY, :YX, :XZ, :ZX, :YZ, :ZY

Warning

The primary and secondary vectors do not have to be orthogonal. However, a great loss of precision happens when the two vectors are almost aligned. This function does not perform any check on the angular separation of the two vectors. The user should ensure that the primary and secondary vector differ of at least 1 milliradian.

FrameTransformations.Frames.twovectors_to_δdcmMethod
twovectors_to_δdcm(a, b, seq)

Compute the time derivative of a direction cosine matrix generated from two time-dependent state vectors a and b, following the directions specified in seq.

Inputs

  • a and b – 6-elements state vectors (position and velocity).
  • seq – Accepted sequence directions are: :XY, :YX, :XZ, :ZX, :YZ, :ZY
FrameTransformations.Frames.twovectors_to_δ²dcmMethod
twovectors_to_δ²dcm(a, b, seq)

Compute the 2nd-order time derivative of a direction cosine matrix generated from two time-dependent state vectors a and b, following the directions specified in seq.

Inputs

  • a and b – 9-elements state vectors (position velocity and acceleration).
  • seq – Accepted sequence directions are: :XY, :YX, :XZ, :ZX, :YZ, :ZY
FrameTransformations.Frames.twovectors_to_δ³dcmMethod
twovectors_to_δ³dcm(a, b, seq)

Compute the 3rd-order time derivative of a direction cosine matrix generated from two time-dependent state vectors a and b, following the directions specified in seq.

Inputs

  • a and b – 12-elements state vectors (position, velocity, acceleration and jerk).
  • seq – Accepted sequence directions are: :XY, :YX, :XZ, :ZX, :YZ, :ZY
FrameTransformations.Frames.update_point!Method
update_point!(frames::FrameSystem, point, stv::AbstractVector, time)

Update the state vector of point at the input time in frames. The only accepted length for the input vector stv are 3, 6, 9 or 12. The order is automatically inferred from the vector length.

Examples

julia> FRAMES = FrameSystem{2, Float64}();
  
julia> @axes ICRF 1  

julia> add_axes_inertial!(FRAMES, ICRF)

julia> @point Origin 0

julia> @point Satellite 1 

julia> add_point_root!(FRAMES, Origin, ICRF)

julia> add_point_updatable!(FRAMES, Satellite, Origin, ICRF)

julia> y = [10000., 200., 300.];

julia> update_point!(FRAMES, Satellite, y, 0.1)

julia> vector3(FRAMES, Origin, Satellite, ICRF, 0.1)
3-element SVector{3, StaticArraysCore.SVector{3, Float64} with indices SOneTo(3):
 10000.0
   200.0
   300.0

julia> vector3(FRAMES, Origin, Satellite, ICRF, 0.2)
ERROR: UpdatablePoint with NAIFId = 1 has not been updated at time 0.2 for order 1 
[...]

julia> vector6(FRAMES, Origin, Satellite, ICRF, 0.1)
ERROR: UpdatablePoint with NAIFId = 1 has not been updated at time 0.1 for order 2 
[...] 

See also

See also add_point_updatable!

FrameTransformations.Frames.vector12Method
vector12(frame::FrameSystem, from, to, axes, ep::Epoch)

Compute 12-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired epoch ep.

Requires a frame system of order ≥ 4.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the observing point
  • to – ID or instance of the target point
  • axes – ID or instance of the output state vector axes
  • epEpoch of the observer. Its timescale must match that of the frame system.
FrameTransformations.Frames.vector12Method
vector12(frame, from, to, axes, t::Number)

Compute 12-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired time t expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.vector3Method
vector3(frame::FrameSystem, from, to, axes, ep::Epoch)

Compute 3-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired epoch ep.

Requires a frame system of order ≥ 1.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the observing point
  • to – ID or instance of the target point
  • axes – ID or instance of the output state vector axes
  • epEpoch of the observer. Its timescale must match that of the frame system.
FrameTransformations.Frames.vector3Method
vector3(frame, from, to, axes, ep::Epoch, ltcorr, dir; <keyword arguments>)

Compute a light-time corrected 3-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired epoch ep, using the aberration flag ltcorr, which may be any of the following AbstractLightTimeCorrection:

  • LightTimeCorrection: it applies the one-way light time (planetary aberration) correction, using a Newtonian formulation.

  • PlanetaryAberrationCorrection: it applies the one-way light time and stellar aberration corrections using a Newtonian fromulation. It modifies the vector obtained with the LightTimeCorrection option to account for the observer velocity with respect to the Solar System Barycenter.

The integer argument dir is used to specify the correction direction, as follows:

  • -1: for Reception, in which photons depart from the target's location at the light-time corrected epoch ep-lt and arrive at the observer's location at ep.

  • +1: for Transmission, in which photons depart from the observer's location at ep and arrive at the target's location at the light-time corrected epoch ep+lt.

Keyword Arguments

  • iters::Int=1: the number of iterations used to find the solution to the light time equation. For the solar system bodies, a solution is usually found in 3 iterations.

  • axescenter: ID or instance of the center point for axes. This parameter is used only when axes have an orientation that depends on time. In these cases, the point is used to find the time at which the orientation of the axes should be computed. It defaults to from.

Note

If the PlanetaryAberrationCorrection is applied, the frame system must be at least one order higher than that of the requested transformation.

See also

See also LightTime, PlanetaryAberration and vector6.

References

FrameTransformations.Frames.vector3Method
vector3(frame, from, to, axes, t::Number, ltcorr, dir; <keyword arguments>)

Compute a light-time corrected 3-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired time t, expressed in seconds since J2000, using the aberration flag ltcorr and the direction dir.

FrameTransformations.Frames.vector3Method
vector3(frame, from, to, axes, t::Number)

Compute 3-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired time t expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.vector6Method
vector6(frame::FrameSystem, from, to, axes, ep::Epoch)

Compute 6-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired epoch ep.

Requires a frame system of order ≥ 2.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the observing point
  • to – ID or instance of the target point
  • axes – ID or instance of the output state vector axes
  • epEpoch of the observer. Its timescale must match that of the frame system.
FrameTransformations.Frames.vector6Method
vector6(frame, from, to, axes, ep::Epoch, ltcorr, dir; <keyword arguments>)

Compute a light-time corrected 6-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired epoch ep, using the aberration flag ltcorr, which may be any of the following AbstractLightTimeCorrection:

  • LightTimeCorrection: it applies the one-way light time (planetary aberration) correction, using a Newtonian formulation.

  • PlanetaryAberrationCorrection: it applies the one-way light time and stellar aberration corrections using a Newtonian fromulation. It modifies the vector obtained with the LightTimeCorrection option to account for the observer velocity with respect to the Solar System Barycenter.

The integer argument dir is used to specify the correction direction, as follows:

  • -1: for Reception, in which photons depart from the target's location at the light-time corrected epoch ep-lt and arrive at the observer's location at ep.

  • +1: for Transmission, in which photons depart from the observer's location at ep and arrive at the target's location at the light-time corrected epoch ep+lt.

Keyword Arguments

  • iters::Int=1: the number of iterations used to find the solution to the light time equation. For the solar system bodies, a solution is usually found in 3 iterations.

  • axescenter: ID or instance of the center point for axes. This parameter is used only when axes have an orientation that depends on time. In these cases, the point is used to find the time at which the orientation of the axes should be computed. It defaults to from.

Note

If the PlanetaryAberrationCorrection is applied, the frame system must be at least one order higher than that of the requested transformation.

See also

See also LightTime, PlanetaryAberration and vector6.

References

FrameTransformations.Frames.vector6Method
vector6(frame, from, to, axes, t::Number, ltcorr, dir; <keyword arguments>)

Compute a light-time corrected 6-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired time t, expressed in seconds since J2000, using the aberration flag ltcorr and the direction dir.

FrameTransformations.Frames.vector6Method
vector6(frame, from, to, axes, t::Number)

Compute 6-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired time t expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.vector9Method
vector9(frame::FrameSystem, from, to, axes, ep::Epoch)

Compute 9-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired epoch ep.

Requires a frame system of order ≥ 3.

Inputs

  • frame – The FrameSystem container object
  • from – ID or instance of the observing point
  • to – ID or instance of the target point
  • axes – ID or instance of the output state vector axes
  • epEpoch of the observer. Its timescale must match that of the frame system.
FrameTransformations.Frames.vector9Method
vector9(frame, from, to, axes, t::Number)

Compute 9-elements state vector of a target point relative to an observing point, in a given set of axes, at the desired time t expressed in seconds since J2000 if ephemerides are used.

FrameTransformations.Frames.@axesMacro
@axes(name, id, type=nothing)

Define a new axes instance to alias the given id. This macro creates an AbstractFrameAxes subtype and its singleton instance callen name. Its type name is obtained by appending Axes to either name or type (if provided).

Examples

julia> @axes ICRF 1 InternationalCelestialReferenceFrame

julia> typeof(ICRF)
InternationalCelestialReferenceFrameAxes

julia> axes_alias(ICRF) 
1

julia> @axes IAU_EARTH 10013

julia> typeof(IAU_EARTH)
IauEarthAxes

See also

See also @point and axes_alias.

FrameTransformations.Frames.@pointMacro
@point(name, id, type=nothing)

Define a new point instance to alias the given NAIFID id. This macro creates an AbstractFramePoint subtype and its singleton instance called name. Its type name is obtained by appending Point to either name or type (if provided).

Examples

julia> @point Venus 299

julia> typeof(Venus)
VenusPoint 

julia> point_alias(Venus)
299

julia> @point EMB 3 EarthMoonBarycenter

julia> typeof(EMB) 
EarthMoonBarycenterPoint

julia> point_alias(EMB) 
3 

See also

See also @axes and point_alias.

FrameTransformations.Orient.AXESID_GCRFConstant
AXESID_GCRF

Axes ID for the Geocentric Celestial Reference Frame (GCRFF)

Note

Although the ICRF and GCRF axes are identical, they are based upon a different timescale. A different ID is here assigned to provide a robust way of distinguishing between the two. 23 has been chosen because it is one the unassigned axes ID among the built-in SPICE frames.

FrameTransformations.Orient.AXESID_ITRFConstant
AXESID_ITRF

NAIF Axes ID for the International Terrestrial Reference Frame (ITRF)

Note

This ID is based upon the ID used to refer to the ITRF93 in NAIF's high-accuracy Earth rotation model PCK kernels.

FrameTransformations.Orient.AXESID_MEME2000Constant
AXESID_MEME2000

Axes ID for the Mean Dynamical Equator and Equinox of J2000.0.

Note

In SPICE the J2000 and ICRF axes are considered equal, thus there exist no specific NAIF ID for the MEME2000 axes. 22 has been chosen because it is the first unassigned axes ID among the built-in SPICE frames.

FrameTransformations.Orient.CPNcConstant
CPNc

The singleton instance of type CPNC, representing the concise CPNc from Capitaine & Wallace, Concise CIO based precession-nutation formulations, (2008). This model truncates the X, Y series to deliver an accuracy of few mas.

FrameTransformations.Orient.CPNdConstant
CPNd

The singleton instance of type CPND, representing the concise CPNd from Capitaine & Wallace, Concise CIO based precession-nutation formulations, (2008). This model truncates the X, Y series to deliver an accuracy of few arcseconds.

FrameTransformations.Orient.DCM_B1950_TO_ECLIPB1950Constant
DCM_B1950_TO_ECLIPB1950

DCM for the rotation from the Mean Equator and Dynamical Equinox of B1950 (MEMEB1950) to the Mean Ecliptic Equinox of B1950. This corresponds to the transformation B1950 -> ECLIPB1950 in the SPICE toolkit, and uses the mean obliquity of the ecliptic from the IAU 1976 theory.

References

FrameTransformations.Orient.DCM_B1950_TO_FK4Constant
DCM_B1950_TO_FK4

DCM for the rotation from the Mean Equator and Dynamical Equinox of B1950 (MEMEB1950) ot the FK4 reference frame.

Note

The FK4 reference frame is obtained from the B1950 frame by applying the equinox offset determined by Fricke.

References

FrameTransformations.Orient.DCM_ICRF_TO_MEME2000Constant
DCM_ICRF_TO_MEME2000

DCM for the rotation from the International Celestial Reference Frame (ICRF) and the Mean Equator and Equinox of J2000.0 (MEME2000). This corresponds to the J2000 frame in the SPICE toolkit.

Note

The frame bias is here computed using the IAU 2006 Precession model, similarly to ESA's GODOT. Some other software libraries, such as Orekit, use the frame bias of the IAU 2000 precession model. The two definitions differ of about 1 arcsecond.

Moreover, according to Hilton there are multiple possibilities to define the proper rotation between the ICRS and the MEME2000. The transformation implemented here correspond to Eq. 6 using the parameters in Table 3, line 1 (RIERS).

References

  • Hilton, James L., and Catherine Y. Hohenkerk. – Rotation matrix from the mean dynamical equator and equinox at J2000. 0 to the ICRS. – Astronomy & Astrophysics 513.2 (2004): 765-770. DOI: 10.1051/0004-6361:20031552
  • SOFA docs
FrameTransformations.Orient.DCM_MEME2000_TO_B1950Constant
DCM_MEME2000_TO_B1950

DCM for the rotation from the Mean Equator and Equinox of J2000.0 (MEME2000) to the Mean Equator and Dynamical Equinox of B1950 (B1950).

Note

This rotation is obtained by precessing the J2000 frame backwards from Julian year 2000 to Besselian year 1950, using the 1976 IAU precession model. The rotation values are taken from the SPICE toolkit.

References

FrameTransformations.Orient.DCM_MEME2000_TO_ECLIPJ2000Constant
DCM_MEME2000_TO_ECLIPJ2000

DCM for the rotation from the Mean Equator and Equinox of J2000 (MEME2000) to the Mean Ecliptic Equinox. This corresponds to the transformation J2000 -> ECLIPJ2000 in the SPICE toolkit, and uses the mean obliquity of the ecliptic from the IAU 1976 theory.

FrameTransformations.Orient.IERS_EOPConstant
IERS_EOP

Earth Orientation Parameters interpolators: x/y pole, UT1-UTC, LOD, dX, dY (smoothed values at 1-day intervals) with respect to IAU 2006/2000A precession-nutation model and consistent with ITRF2014.

See also: EOPInterpolator

FrameTransformations.Orient.IERS_EOP_DATAConstant
IERS_EOP_DATA

Earth Orientation Parameters Data: x/y pole, UT1-UTC, LOD, dX, dY (smoothed values at 1-day intervals) with respect to IAU 2006/2000A precession-nutation model and consistent with ITRF2014.

See also: EOPData

FrameTransformations.Orient.iau2006bConstant
iau2006b

The singleton instance of type IAU2006B, representing the IAU 2006B family of models.

Note

This is not an official IERS model.

FrameTransformations.Orient.EOPDataType
EOPData{T}

EOP Data for IAU 2000A.

Fields

  • filename : File where the EOP data are stored.
  • j2000 : Independent valiable - time - in UTC.
  • j2000_TT : Independent valiable - time - in TT.
  • x, y: Polar motion with respect to the crust (arcsec).
  • UT1_UTC: Irregularities of the rotation angle (s).
  • UT1_TT: Irregularities of the rotation angle (s) w.r.t. TT timescale.
  • LOD: Length of day offset (ms).
  • dX, dY: Celestial pole offsets referred to the model IAU2000A (milliarcsec).
FrameTransformations.Orient.EOPInterpolatorType
EOPInterpolator{T}

EOP Data for IAU 2000A.

Fields

  • init: A flag indicating whether the EOPInterpolator has been initialized.
  • j2000: Independent variable (time), in UTC.
  • j2000_TT: Independent variable (time), in TT.
  • x, y: Polar motion with respect to the crust (arcsec).
  • UT1_UTC: Irregularities of the rotation angle (s).
  • UT1_TT: Irregularities of the rotation angle (s) with respect to TT timescale.
  • LOD: Length of day offset (ms).
  • dX, dY: Celestial pole offsets referred to the model IAU2000A (milliarcsec).
  • x_TT, y_TT: Polar motion with respect to the crust (arcsec) parametrized by TT epoch.
  • UT1_TT: Irregularities of the rotation angle (s) parametrized by TT epoch.
  • LOD_TT: Length of day offset (ms) parametrized by TT epoch.
  • dX_TT, dY_TT: Celestial pole offsets referred to the model IAU2000A (milliarcsec) parametrized by TT epoch.
FrameTransformations.Orient.FundamentalArgumentsType
FundamentalArguments(t::Number, m::IAUModel=iau2006a)

Compute the fundamental luni-solar and planetary arguments, in radians, associated to the given IAU model m, at time t expressed in TDB Julian centuries since J2000.

The available IAU model options are: iau2006a and iau2000b.

Note

Only the planetary arguments are affected by the IAU Model choice.

FrameTransformations.Orient.FundamentalArgumentsType
FundamentalArguments{N <: Number}

Type storing the fundamental luni-solar and planetary arguments.

Fields

  • Mₐ – Mean anomaly of the Moon
  • Sₐ – Mean anomaly of the Sun
  • uₘ – Mean longitude of the Moon minus mean longitude of the ascending node F
  • Dₛ – Mean elongation of the Moon from the Sun
  • Ωₘ – Mean longitude of the Moon's ascending node
  • λ_Me – Mercury's mean heliocentric longitude
  • λ_Ve – Venus's mean heliocentric longitude
  • λ_Ea – Earth's mean heliocentric longitude
  • λ_Ma – Mars's mean heliocentric longitude
  • λ_Ju – Jupiter's mean heliocentric longitude
  • λ_Sa – Saturn's mean heliocentric longitude
  • λ_Ur – Uranus's mean heliocentric longitude
  • λ_Ne – Neptune's mean heliocentric longitude
  • pₐ – General precession in longitude
FrameTransformations.Orient.LuniSolarArgumentsMethod
LuniSolarArguments(t::Number, m::IAUModel)

Compute the fundamental (Delaunay) Luni-Solar arguments, in radians, associated to the desired IAU model m, at time t expressed in TDB Julian centuries since J2000.

The returned values depend on the input model as follows:

  • IAU2006A: the Delaunay expressions are taken from the IERS 2010 Conventions.

  • IAU2000B: the expressions are taken from Simon et al. (1994), following ERFA's implementation of nut00b.c

Outputs

  • Mₐ – Mean anomaly of the Moon
  • Sₐ – Mean anomaly of the Sun
  • uₘ – Mean longitude of the Moon minus mean longitude of the ascending node F
  • Dₛ – Mean elongation of the Moon from the Sun
  • Ωₘ – Mean longitude of the Moon's ascending node

References

  • Luzum, B. and Petit G. (2012). The IERS Conventions (2010), IERS Technical Note No. 36
  • Simon, J. et al., (1994), Numerical expressions for precession formulae and mean elements for the Moon and the planets.
  • ERFA software library
FrameTransformations.Orient.PlanetaryArgumentsMethod
PlanetaryArguments(t::Number)

Compute the fundamental planetary arguments and the general precession, in radians, at time t expressed in TDB Julian centuries since J2000.

Outputs

  • λ☿ – Mercury's mean heliocentric longitude.
  • λ♀ – Venus's mean heliocentric longitude.
  • λe – Earth's mean heliocentric longitude.
  • λ♂ – Mars's mean heliocentric longitude.
  • λ♃ – Jupiter's mean heliocentric longitude.
  • λ♄ – Saturn's mean heliocentric longitude.
  • λ⛢ – Uranus's mean heliocentric longitude.
  • λ♆ – Neptune's mean heliocentric longitude.
  • – General precession in longitude.

References

FrameTransformations.Orient.cio_locatorMethod
cio_locator(m::IAUModel, t::Number, x::Number, y::Number)

Compute the CIO Locator s in radians, according to the IAU Model m, given the CIP coordinates X and Y at time t expressed in TT Julian centuries since J2000

The function has been implemented for the IAU2000, IAU2006 and the CPN models.

References

  • Luzum, B. and Petit G. (2012), The IERS Conventions (2010), IERS Technical Note No. 36
  • Wallace P. T. and Capitaine N. (2006), Precession-nutation procedures consistent with IAU 2006 resolutions, DOI: 10.1051/0004-6361:20065897
  • Capitaine N. and Wallace P. T. (2008), Concise CIO based precession-nutation formulations
  • ERFA library
  • ERFA library
FrameTransformations.Orient.cipMethod
cip(m::IAUModel, t::Number)

Compute Celestial Intermediate Pole vector.

  • Wallace, P. T., & Capitaine, N. (2006). Precession-nutation procedures consistent with IAU

2006 resolutions. Astronomy & Astrophysics.

FrameTransformations.Orient.cip_coordsMethod
cip_coords(m::IAUModel, t::Number)

Computes the CIP X, Y coordinates, in radians, according to the IAU model m at time t expressed in TT Julian Centuries since J2000.

This function has been implemented for the IAU2000, IAU2006 and the CPN models.

References

  • Wallace P. T. and Capitaine N. (2006), Precession-nutation procedures consistent with IAU 2006 resolutions, DOI: 10.1051/0004-6361:20065897
  • Capitaine N. and Wallace P. T. (2008), Concise CIO based precession-nutation formulations
FrameTransformations.Orient.cip_motionFunction
cip_motion(m::IAUModel, t::Number, dx::Number=0.0, dy::Number=0.0)

Compute the CIRS to GCRS rotation matrix, according to the IAU Model m, at time t expressed in TT Julian centuries since J2000. Optional IERS corrections for free-core nutation and time depedented effects can be provided through dx and dy.

References

FrameTransformations.Orient.earth_rotation_angleMethod
earth_rotation_angle(t::Number)

Compute the Earth Rotation Angle (ERA) in radians, i.e., the angle between the Celestial Intermediate Origin (CIO) and the Terrestrial Intermediate Origin (TIO) at time t expressed as UT1 days since J2000.

Note

The function uses the fractional UT1 date to gain additional precision in the computations (0.002737.. instead of 1.002737..)

References

FrameTransformations.Orient.earth_rotation_rateMethod
earth_rotation_rate(LOD::Number)

Compute the true angular velocity of the Earth accounting for the Length of the Day, i.e., the instantaneous rate of change of UT1 with respect to a uniform time scale.

FrameTransformations.Orient.ecliptic_poleMethod
ecliptic_pole(m::IAU2006Model, t::Number)

Computes ecliptic pole of date in ICRF.

  • Wallace, P. T., & Capitaine, N. (2006). Precession-nutation procedures consistent with IAU

2006 resolutions. Astronomy & Astrophysics, Eq. 20.

FrameTransformations.Orient.equinoxes_equationMethod
equinoxes_equation(m::IAUModel, tt::Number)

Compute the Equation of the Equinoxes, in radians, according to the IAU Model m given time tt expressed in TT Julian centuries since J2000.

This function has been implemented for the IAU2006 and IAU2000 models.

Note

This function neglects the difference between TT and TDB.

References

FrameTransformations.Orient.frame_biasMethod
frame_bias(::IAU2000Model)

Compute the frame bias components of the IAU 2000 precession-nutation models, in radians.

Notes

  • The frame bias corrections in longitude and obliquity are required to correct for the offset between the GCRS pole and the mean J2000 pole. They define, with respect to the GCRS axes, a J2000 mean pole that is consistent with teh IAU 2000A precession-nutation model.

  • The function also returns an offset in right ascension taken from Chapront et al. (2002), necessary to completely describe the frame bias, but that is not part of the original IAU 2000A model.

References

  • ERFA software library
FrameTransformations.Orient.fw2xyMethod
fw2xy(ϵ::Number, ψ::Number, γ::Number, φ::Number)

Compute the CIP X and Y coordinates from Fukushima-Williams bias-precession-nutation angles, in radians.

Inputs

  • ϵ – F-W angle with IAU 2006A/B nutation corrections.
  • ψ – F-W angle with IAU 2006A/B nutation corrections.
  • γ – F-W angle
  • ϕ – F-W angle

References

FrameTransformations.Orient.fw_anglesMethod
fw_angles(m::IAU2006Model, t::Number)

Compute the precession angles in radians, following the IAU 2006 Fukushima-Williams 4-angle formulation at time t expressed in TT Julian centuries since J2000.

Outputs

  • γ – F-W 1st angle
  • ϕ – F-W 2nd angle
  • ψ – F-W 3rd angle
  • ε – F-W 4th angle

References

FrameTransformations.Orient.fw_matrixMethod
fw_matrix(γ, ϕ, ψ, ε)

Form the Nutation-Precession-Bias (NPB) rotation matrix given the Fukushima-Williams angles, expressed in radians.

The present function can construct three different matrices depending on which angles are supplied as arguments:

  • NPB: To obtain the Nutation-Precession-Bias matrix, generate the four standard FW precession angles (̄γ, ̄ϕ, ̄ψ, ϵₐ) then generate the nutation components Δψ and Δϵ and add them to ̄ψ, ϵₐ. Finally, call the present functions using those four angles as arguments.

  • PB: To obtain the precession-frame bias matrix, generate the four standard FW precession angles and call the present function.

  • B: To obtain the frame bias matrix, generate the four standard FW precession angles at date J2000.0 and call this function.

The remaining nutation-only and precession only matrices can be obtained by combining these three appropriately.

References

FrameTransformations.Orient.geoc2posMethod
geoc2pos(r::Number, λ::Number, ϕ::Number)
geoc2pos(geoc::AbstractArray)

Transform geocentric coordinates in a cartesian position vector, given the longitude λ, the geocentric latitude ϕ and the radius r.

FrameTransformations.Orient.geoc2pvMethod
geoc2pv(geoc::AbstractVector)

Transform a spherical geocentric 6-elements state vector (radius, longitude, geocentric latitude and their derivatives) into a cartesian 6-elements vector (position and velocity).

FrameTransformations.Orient.geod2posMethod
geod2pos(h::Number, λ::Number, ϕ::Number, R::Number, f::Number)

Transform longitude λ, geodetic latitude ϕ and altitude over the reference ellipsoid to a cartesian position vector, given the reference radius R and the flattening f.

References

  • Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications, Microcosm Press, Hawthorn, CA, USA.
FrameTransformations.Orient.init_eopMethod
init_eop(filename)

Initialize Earth Orientation Parameters (EOP) from file.

Warn

This function must be called to initialize the EOP data used by frames, in case Earth-associated frames are used.

Warn

This function accept only .eop.dat files. Please use prepare_eop to transform IERS EOP files in this format.

FrameTransformations.Orient.orient_bias_precessionMethod
orient_bias_precession(m::IAUModel, t::Number)

Form the precession-frame bias (PB) matrix that transforms vectors from the GCRS to the mean of date, following the IAU Model m at time t expressed as TT Julian centuries since J2000.

The function has been implemented for the IAU2000 and IAU2006 models.

References:

  • IAU: Trans. International Astronomical Union, Vol. XXIVB; Proc. 24th General Assembly, Manchester, UK. Resolutions B1.3, B1.6. (2000)
  • Wallace P. T. and Capitaine N. (2006), Precession-nutation procedures consistent with IAU 2006 resolutions, DOI: 10.1051/0004-6361:20065897
  • ERFA pmat06 function.
  • ERFA pmat00 function.
FrameTransformations.Orient.orient_bias_precession_nutationMethod
orient_bias_precession_nutation(m::IAUModel, t::Number)

Compute the equinox-based bias-precession-nutation matrix using the IAU Model m procedures at time t expressed in TT Julian centuries since J2000.

The function has been implemented for the IAU2000 and IAU2006 models.

Note

The computed matrix rotates a vector from the GCRS to the true equatorial triad of date.

References

FrameTransformations.Orient.orient_gastMethod
orient_gast(m::IAUModel, t::Number, θ::Number)

Compute the Greenwich Apparent Sidereal Time (GAST), in radians, given time t as TT Julian centuries since J2000 and the Earth Rotation Angle (ERA) θ, in radians, according to the IAU Model m.

The function has been implemented for the IAU2000 and IAU2006 models.

References

FrameTransformations.Orient.orient_gastMethod
orient_gast(m::IAUModel, t::Number)

Compute the Greenwich Apparent Sidereal Time (GAST), in radians, given time t as TT Julian centuries since J2000 according to the IAU Model m.

Note

For the IAU2000B model, as an approximation ERFA uses UT1 instead of TDB (or TT) to compute the precession component of GMST and the equation of the equinoxes. This approximation is not performed in this framework.

Note

This function computes the Earth Rotation Angle (ERA) by automatically converting TT to UT1. Lower-level interfaces are also available to minimise the number of repeated computations.

See also

See also orient_gmst.

FrameTransformations.Orient.orient_gmstMethod
orient_gmst(m::IAU1980Model, ut1::Number)

Compute the Greenwich Mean Sidereal Time (GMST), in radians, according to the IAU 1980 models, given time ut1 expressed in UT1 Julian centuries since J2000.

References

FrameTransformations.Orient.orient_gmstMethod
orient_gmst(m::IAUModel, tt::Number, θ::Number)

Compute the Greenwich Mean Sidereal Time (GMST), in radians, according to the IAU Model m, given the Earth Rotation Angle (ERA) θ, in radians and the time tt expressed in TT Julian centuries since J2000.

The function has been implemented for the IAU2000 and IAU2006 models.

References

FrameTransformations.Orient.orient_gmstMethod
orient_gmst(::IAUModel, t::Number)

Compute the Greenwich Mean Sidereal Time (GMST), in radians, according to the IAU Model m, given time tt expressed in TT Julian centuries since J2000.

Note

This function is valid for the IAU2000 and IAU2006 models only.

Note

This function computes the Earth Rotation Angle (ERA) by automatically converting TT to UT1. Lower-level interfaces are also available to minimise the number of repeated computations.

See also

See also orient_gast.

FrameTransformations.Orient.orient_nutationMethod
orient_nutation(m::IAUModel, t::Number)

Compute the nutation components in longitude and obliquity for the IAU Model m, in radians, at time t expressed in TT Julian Centuries since J2000.

Notes

  • Due to their theoretical basis, the original developments required t expressed as TDB. However, in practice, it is usually more convenient to use Terrestrial Time (TT) as it makes no significant differences (< 0.01 μas) in the final result.

  • For the IAU 2006A model, the function strictly follows the SOFA implementation. It first computes the IAU 2000A nutation, then applies adjustments for the consequences of the change in obliquity from the IAU 1980 ecliptic to the IAU 2006 ecliptic and (ii) for the secular variation in the Earth's dynamical form factor J2. These corrections ensure that the IAU 2000A nutation is consistent with the IAU 2006 precession model. Please note that the coefficients available on the IERS tables already include those corrections, and are retrieved by multiplying the amplitudes of the SOFA nutation in longitude coefficients by 1.00000047.

  • The computation of the free-core nutation and time dependent effects are excluded from this model. To achieve the < 1μas accuracy with the IAU 2006/2000 A precession-nutation models, such effects must be included a-posteriori (through dX and dY) using the IERS EOP data.

  • For the IAU 2000B model, the nutation series is truncated from nearly 1400 terms to only 77, yet it still delivers results of 1 mas accuracy at present epochs. In particular, it delivers a pole accurate to 1 mas from 1900 to 2100 (only very occasionally just outside 1 mas). The coefficients are taken from SOFA's implementation, which slighlty differ from those reported in McCarthy and Luzum (2003). Comparisons with IAU 2006A show that the SOFA version between 1995 and 2050 delivers 0.283 mas RMSE (0.994 mas in the worst case), whereas the IERS Conventions website version delivers 0.312 RMSE (1.125 mas in the worst case).

  • The IAU 2000B model includes constant planetary bias terms that compensate for long-period nutations. These amplitudes used in this implementation are optimised for a rigorous method, where frame bias, precession and nutation are applied separately and in that order (see SOFA's documentation for further insights).

  • A simplified version of the Fundamental Arguments, taken from Simon et al (1994) is exploited for IAU2000B as the error introduced is below the model accuracy ( ~0.1 mas).

References

  • Luzum, B. and Petit G. (2012), The IERS Conventions (2010), IERS Technical Note No. 36
  • Wallace P. T. and Capitaine N. (2006), Precession-nutation procedures consistent with IAU 2006 resolutions, DOI: 10.1051/0004-6361:20065897
  • Simon, J. et al., (1994), Numerical expressions for precession formulae and mean elements for the Moon and the planets.
  • ERFA nut06a and nut00b functions
FrameTransformations.Orient.orient_obliquityMethod
orient_obliquity(m::IAUModel, t::Number)

Compute the mean obliquity of the ecliptic at epoch, in radians, at time t expressed in TT Julian centuries since J2000.

Note

This function is implemented only for IAU1980 and IAU2006 models. IAU 2000 Models implement proper precession-rate corrections to the IAU1980 mean obliquity.

References

FrameTransformations.Orient.orient_rot12_itrf_to_gcrfFunction
orient_rot12_itrf_to_gcrf(m::IAUModel, tt, ut1, xₚ, yₚ, dX=0.0, dY=0.0, LOD=0.0)

Compute the rotation matrix from ITRF to GCRF and its derivatives up to order 3, according to the IAU Model m, at time tt and ut1 expressed in TT seconds and UT1 days since J2000, respectively.

This function has been implemented for IAU2000, IAU2006 and CPN models.

Note

All the input quantities xₚ, yₚ, dX and dY must be expressed in radians

FrameTransformations.Orient.orient_rot12_itrf_to_gcrfFunction
orient_rot12_itrf_to_gcrf(m::IAUModel, t::Number)

Compute the rotation matrix from ITRF to GCRF and its time derivatives up to order 3 at time t expressed as TT seconds since J2000, according to the IAU Model m.

FrameTransformations.Orient.orient_rot3_icrf_to_modMethod
orient_rot3_icrf_to_mod(tt::Number)

Compute the rotation matrix from the International Celestial Reference Frame (ICRF) to the Mean Equinox Mean Equator of Date at time tt, expressed in TT seconds since J2000.

Mean Equator Of Date is obtained applying frame bias and precession to the ICRF pole and origin. Fukushima-Williams parametrization for the equator and ecliptic precession is used. Consistent with the IAU2006 precession model.

FrameTransformations.Orient.orient_rot3_icrf_to_pa421Method
orient_rot3_icrf_to_pa421(eph::AbstractEphemerisProvider, t::Number)

Compute the rotation matrix from the ICRF to the DE421 Moon's Principal Axes at the given input time t, expressed in seconds since J2000.

Warning

This function is not optimised for performance (it allocates!). The user is suggested to retrieve the Principal axes orientation using the dedicated Frame System functions.

FrameTransformations.Orient.orient_rot3_icrf_to_pa440Method
orient_rot3_icrf_to_pa440(eph::AbstractEphemerisProvider, t::Number)

Compute the rotation matrix from the ICRF to the DE440 Moon's Principal Axes at the given input time t, expressed in seconds since J2000.

Warning

This function is not optimised for performance (it allocates!). The user is suggested to retrieve the Principal axes orientation using the dedicated Frame System functions.

FrameTransformations.Orient.orient_rot3_icrf_to_todMethod
orient_rot3_icrf_to_tod(tt::Number; [m]::IAU2006Model=iau2006a)

Compute the rotation matrix from the International Celestial Reference Frame (ICRF) to the True Equator of Date at time tt, expressed in TT seconds since J2000.

True Equator of Date is obtained applying frame bias, precession and nutation to the ICRF pole and origin.

FrameTransformations.Orient.orient_rot3_itrf_to_gcrfFunction
orient_rot3_itrf_to_gcrf(m::IAUModel, tt, ut1, xₚ, yₚ, dX=0.0, dY=0.0)

Compute the rotation matrix from ITRF to GCRF according to the IAU Model m, at time tt and ut1 expressed in TT seconds and UT1 days since J2000, respectively.

This function has been implemented for IAU2000, IAU2006 and CPN models.

Note

All the input quantities xₚ, yₚ, dX and dY must be expressed in radians

FrameTransformations.Orient.orient_rot3_itrf_to_gcrfMethod
orient_rot3_itrf_to_gcrf(m::IAUModel, t::Number)

Compute the rotation matrix from ITRF to GCRF at time t expressed as TT seconds since J2000, according to the IAU Model m, as follows:

  • IAU2000A: the pole coordinates (xₚ, yₚ) and the free-core nutation and time corrections to the CIP coordinates (dX, dY) are interpolated from the latest released IERS EOP data. The precession-nutation matrix is computed using the full IAU 2000A model.

  • IAU2000B: only the pole coordinates (xₚ, yₚ) are interpolated from the latest EOP data. The Free Core Nutation (FCN) corrections dX, dY are neglected. The precession-nutation matrix is computed following the IAU 2000 model but with truncated expressions for the nutation corrections.

  • IAU2006A: the pole coordinates (xₚ, yₚ) and the free-core nutation and time corrections to the CIP coordinates (dX, dY) are interpolated from the latest released IERS EOP data. The precession-nutation matrix is computed using the full IAU 2006/2000A model.

  • IAU2006B: only the pole coordinates (xₚ, yₚ) are interpolated from the latest EOP data. The Free Core Nutation (FCN) corrections dX, dY are neglected. The precession-nutation matrix is computed following the IAU 2006A model but with truncated expressions for the nutation corrections.

  • CPNc: a concise model with a cut-off at 2.5 mas of the X and Y series, delivering a worst-case accuracy of about 15 mas between 1995-2050. It does not take into account the Free Core Nutation (~0.2 mas).

  • CPNd: an extremely concise formulation with an accuracy of about 1 arcsec between 1995 and 2050. It neglects polar-motion (~0.25 arcsec), the FCN corrections and the CIO locator.

References

  • Luzum, B. and Petit G. (2012), The IERS Conventions (2010), IERS Technical Note No. 36
  • Wallace P. T. and Capitaine N. (2006), Precession-nutation procedures consistent with IAU 2006 resolutions, DOI: 10.1051/0004-6361:20065897
  • Capitaine N. and Wallace P. T. (2008), Concise CIO based precession-nutation formulations
FrameTransformations.Orient.orient_rot3_itrf_to_pefMethod
orient_rot3_itrf_to_pef(t::Number)

Compute the rotation matrix from the International Terrestrial Reference Frame (ITRF) to the Pseudo-Earth Fixed Frame at time t, expressed in TT seconds since J2000.

References

  • Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications, Microcosm Press, Hawthorn, CA, USA.

See also

See also IERS_EOP.

FrameTransformations.Orient.orient_rot3_mod_to_temeMethod
orient_rot3_mod_to_teme(tt::Number; [m]::IAU2006Model=iau2006a)

Compute the rotation matrix from the Mean Equator of Date (MOD) frame to the True Equator, Mean Equinox of date at time tt, expressed in TT seconds since J2000.

This is implemented with a small angle approx of Eq. 4 of Vallado, "Coordinate Frames of the US Space Object Catalogs."

FrameTransformations.Orient.orient_rot3_pef_to_todMethod
orient_rot3_pef_to_tod(tt::Number; [m]::IAUModel=iau2006a)

Compute the rotation matrix from the Pseudo-Earth Fixed (PEF) to the True Equator of Date (TOD) Frame at time tt, expressed in TT seconds since J2000.

This is using IAU-76/FK5 Reduction. This is a sidereal time only rotation. Eq. 3-80, Sec. 3.7.3 of Vallado (2013).

See also

See also orient_gast.

FrameTransformations.Orient.orient_rot3_tod_to_modMethod
orient_rot3_tod_to_mod(tt::Number; [m]::IAU2006Model=iau2006a)

Compute the rotation matrix from the True Equator of Date (TOD) to the Mean Equator of Date (MOD) Frame at time tt, expressed in TT seconds since J2000.

This is using IAU-76/FK5 Reduction. This is a nutation only rotation. Eq. 3-86, Sec. 3.7.3 of Vallado (2013).

FrameTransformations.Orient.orient_rot6_icrf_to_pefMethod
orient_rot6_icrf_to_pef(tt::Number; [m]::IAU2006Model=iau2006a)

Compute the rotation matrix the derivative of the transformation from the International Celestial Reference Frame (ICRF) to the Pseudo Earth Fixed (PEF) Frame at time tt, expressed in TT seconds since J2000.

See also

See also orient_rot3_icrf_to_pef.

FrameTransformations.Orient.orient_rot6_itrf_to_gcrfFunction
orient_rot6_itrf_to_gcrf(m::IAUModel, tt, ut1, xₚ, yₚ, dX=0.0, dY=0.0, LOD=0.0)

Compute the rotation matrix from ITRF to GCRF and its derivative, according to the IAU Model m, at time tt and ut1 expressed in TT seconds and UT1 days since J2000, respectively.

This function has been implemented for IAU2000, IAU2006 and CPN models.

Note

All the input quantities xₚ, yₚ, dX and dY must be expressed in radians

FrameTransformations.Orient.orient_rot6_pef_to_todMethod
orient_rot6_pef_to_tod(tt::Number; [m]::IAUModel=iau2006a)

Compute the rotation matrix and the derivative of the transformation from the Pseudo-Earth Fixed (PEF) to the True Equator of Date (TOD) Frame at time tt, expressed in TT seconds since J2000.

This is using IAU-76/FK5 Reduction. This is a sidereal time only rotation. Eq. 3-80, Sec. 3.7.3 of Vallado (2013).

See also

See also orient_rot3_pef_to_tod, orient_gast and IERS_EOP.

FrameTransformations.Orient.orient_rot9_itrf_to_gcrfFunction
orient_rot9_itrf_to_gcrf(m::IAUModel, tt, ut1, xₚ, yₚ, dX=0.0, dY=0.0, LOD=0.0)

Compute the rotation matrix from ITRF to GCRF and its derivatives up to order 2, according to the IAU Model m, at time tt and ut1 expressed in TT seconds and UT1 days since J2000, respectively.

This function has been implemented for IAU2000, IAU2006 and CPN models.

Note

All the input quantities xₚ, yₚ, dX and dY must be expressed in radians

FrameTransformations.Orient.orient_rot9_itrf_to_gcrfFunction
orient_rot9_itrf_to_gcrf(m::IAUModel, t::Number)

Compute the rotation matrix from ITRF to GCRF and its time derivatives up to order 2 at time t expressed as TT seconds since J2000, according to the IAU Model m.

FrameTransformations.Orient.origins_equationMethod
origins_equation(m::IAU2006Model, t::Number)

Compute the Equation of the Origins (EO), in radians, following the IAU2006 precession and IAU2000A nutation models, given time t expressed in TT Julian centuries since J2000.0.

References

FrameTransformations.Orient.polar_motionMethod
polar_motion(t::Number, xₚ::Number, yₚ::Number)

Compute the Polar Motion rotation matrix from ITRF to TIRS, according to the IERS 2010 Conventions, at time t expressed in TT Julian centuries since J2000. The function requires xp and yp, the Celestial Intermediate Pole (CIP) coordinates with respect to the International Celestial Reference Frame (ITFR).

References

FrameTransformations.Orient.pos2geocMethod
pos2geoc(pos::AbstractVector)

Transform a cartesian 3-elements position vector pos into radius, longitude and geocentric latitude, respectively.

FrameTransformations.Orient.pos2geodFunction
pos2geod(pos::AbstratVector, R::Number, f::Number, toll::Number=1e-12)

Transform a cartesian 3-elements position vector pos into longitude, geodetic latitude and altitude over the reference ellipsoid with radius R and flattening f.

References

  • Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications, Microcosm Press, Hawthorn, CA, USA.
FrameTransformations.Orient.precession_anglesMethod
precession_angles(m::IAU1980Model, t::Number)

Compute the precession angles from Lieske et al. 1977 model, in radians, at time t expressed in TT Julian centuries since J2000.

References

  • ERFA software library
FrameTransformations.Orient.precession_rateMethod
precession_rate(m::IAU2000Model, t::Number)

Compute the precession-rate part of the IAU 2000 precession-nutation models, in radians, at time t expressed as TT Julian centuries since J2000.

References

  • ERFA software library
FrameTransformations.Orient.prepare_eopFunction
prepare_eop(iers_file::AbstractString, output_filename::AbstractString="iau2000a")

Prepare Earth Orientation Parameters (EOP) data from IERS EOP C04 files to JSMD's eop.dat convenience format. The output_filename should not include the file extension, which is automatically added by this function.

# Save a new file called: test.eop.dat
prepare_eop("input.csv", "test")
FrameTransformations.Orient.pv2geocMethod
pv2geoc(pv::AbstractVector)

Transform a cartesian 6-elements state vector (position and velocity) into radius, longitude, geocentric latitude and their derivatives, respectively.

FrameTransformations.Orient.read_iers_eop_finalsMethod
read_iers_eop_finals(filename::AbstractString)

Read IERS EOP C04 files in csv format.

Returns

  • j2000_utc: Julian days since J2000 in UTC.
  • j2000_tt: Julian days since J2000 in TT (Terrestrial Time).
  • x_pole: Celestial pole offset in the X direction (arcsec).
  • y_pole: Celestial pole offset in the Y direction (arcsec).
  • ut1_utc: UT1-UTC time difference (s).
  • ut1_tt: UT1-TT time difference (s).
  • lod: Length of Day (s).
  • dX: Celestial pole offset rate in the X direction (milliarcsec).
  • dY: Celestial pole offset rate in the Y direction (milliarcsec).

The function reads IERS EOP C04 files in CSV format and extracts relevant Earth Orientation Parameters (EOP) data. It then updates predictions, filling missing values with zeros for LOD, dX, and dY. Finally, the function parametrizes EOP with respect to both UTC and TT time scales for convenience.

References

  • https://maia.usno.navy.mil/ser7/readme.finals2000A
  • http://hpiers.obspm.fr/eoppc/bul/bulb/explanatory.html
  • https://maia.usno.navy.mil
FrameTransformations.Orient.tio_locatorMethod
tio_locator(t::Number)

Compute the TIO locator s' at date, positioning the Terrestrial Intermediate Origin on the equator of the Celestial Intermediate Pole (CIP) at time t expressed as TT Julian centuries since J2000.

This function approximates the unpredictable motion of the TIO locator s' with its secular drift of ~0.47 μas/century.

References

FrameTransformations.Orient.xys2mMethod
xys2m(x::Number, y::Number, s::Number)

Compute the Intermediate-to-Celestial matrix given the CIP x, y' coordinates and the CIO locators`, all in radians.

References