Library

T_ECEFs

Union of all Earth-Centered Earth-Fixed (ECEF) frames supported by the IAU-76/FK5 theory.

T_ECEFs_IAU_2006

Union of all Earth-Centered Earth-Fixed (ECEF) frames supported by IAU-2006/2010 theory.

T_ECIs

Union of all Earth-Centered Inertial (ECI) frames supported by the IAU-76/FK5 theory.

T_ECIs_IAU_2006

Union of all Earth-Centered Inertial (ECI) frames supported by IAU-2006/2010 theory.

T_ECIs_of_date

Union of all of date Earth-Centered Inertial (ECI) frames supported by the IAU-76/FK5 theory.

T_ROT

Union of all supported rotation descriptions.

EOPData_IAU1980{T}

EOP Data for IAU 1980.

Fields

• x, y: Polar motion with respect to the crust [arcsec].
• UT1_UTC: Irregularities of the rotation angle [s].
• LOD: Length of day offset [s].
• dPsi, dEps: Celestial pole offsets referred to the model IAU1980 [arcsec].
• *_err: Errors in the components [same unit as the component].

Remarks

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

x[DatestoJD(2018,19,06,0,0,0)]

EOPData_IAU2000A{T}

EOP Data for IAU 2000A.

Fields

• x, y: Polar motion with respect to the crust [arcsec].
• UT1_UTC: Irregularities of the rotation angle [s].
• LOD: Length of day offset [s].
• dX, dY: Celestial pole offsets referred to the model IAU2000A [arcsec].
• *_err: Errors in the components [same unit as the component].

Remarks

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

x[DatestoJD(2018,19,06,0,0,0)]

GravityModel_Coefs{T}

Structure to store the information about a gravity model.

ICGEM

Structure to store the information contained in ICGEM files.

J2_GravCte{T}

Gravitational constants for J2 orbit propagator.

Fields

• R0: Earth equatorial radius [m].
• μm: √GM [er/s]^(3/2).
• J2: The second gravitational zonal harmonic of the Earth.

J2_Structure{T}

Low level J2 orbit propagator structure.

J4_GravCte{T}

Gravitational constants for J4 orbit propagator.

Fields

• R0: Earth equatorial radius [m].
• μm: √GM [er/s]^(3/2).
• J2: The second gravitational zonal harmonic of the Earth.
• J4: The fourth gravitational zonal harmonic of the Earth.

J4_Structure{T}

Low level J4 orbit propagator structure.

JB2008_Output

Output structure of the Jacchia-Bowman 2008.

Fields

• nN2: Number density of N₂ [1/m³].
• nO2: Number density of O₂ [1/m³].
• nO: Number density of O [1/m³].
• nAr: Number density of Ar [1/m³].
• nHe: Number density of He [1/m³].
• nH: Number density of H [1/m³].
• rho: Total density [kg/m³].
• T_exo: Exospheric temperature [K].
• Tz: Temperature at the selected altitude [K].

JR1971_Output

Output structure of the Jacchia-Roberts 1971 model.

Fields

• nN2: Number density of N₂ [1/m³].
• nO2: Number density of O₂ [1/m³].
• nO: Number density of O [1/m³].
• nAr: Number density of Ar [1/m³].
• nHe: Number density of He [1/m³].
• nH: Number density of H [1/m³].
• rho: Total density [kg/m³].
• T_exo: Exospheric temperature [K].
• Tz: Temperature at the selected altitude [K].

NRLMSISE00_Flags

Flags to configure NRLMSISE-00.

Fields

• output_m_kg
• F107_Mean
• time_independent
• sym_annual
• sym_semiannual
• asym_annual
• asyn_semiannual
• diurnal
• semidiurnal
• daily_ap
• all_ut_long_effects
• longitudinal
• ut_mixed_ut_long
• mixed_ap_ut_long
• terdiurnal
• departures_from_eq
• all_tinf_var
• all_tlb_var
• all_tn1_var
• all_s_var
• all_tn2_var
• all_nlb_var
• all_tn3_var
• turbo_scale_height
• use_ap_array

NRLMSISE00_Output

Output structure for NRLMSISE00 model.

Fields

• den_N: Nitrogen number density [U].
• den_N2: N₂ number density [U].
• den_O: Oxygen number density [U].
• den_aO: Anomalous Oxygen number density [U].
• den_O2: O₂ number density [U].
• den_H: Hydrogen number density [U].
• den_He: Helium number density [U].
• den_Ar: Argon number density [U].
• den_Total: Total mass density [T/U] (this value has different meanings for routines gtd7 and gtd7d).
• T_exo: Exospheric temperature [K].
• T_alt: Temperature at the selected altitude [K].
• flags: Flags used to compute NRLMSISE-00 model.

Notice that:

• If flags.output_m_kg is false, then [U] is [cm⁻³] and [T] is [g/cm⁻³].
• If flags.output_m_kg is true, then [U] is [m⁻³] and [T] is [kg/m⁻³].

Remarks

Anomalous oxygen is defined as hot atomic oxygen or ionized oxygen that can become appreciable at high altitudes (> 500 km) for some ranges of inputs, thereby affection drag on satellites and debris. We group these species under the term Anomalous Oxygen, since their individual variations are not presently separable with the drag data used to define this model component.

NRLMSISE00_Structure{T}

Structure with the configuration parameters for NRLMSISE-00 model. It can be created using the function conf_nrlmsise00.

Orbit{T1,T2}

This structure defines the orbit in terms of the Keplerian elements.

Fields

• t: Orbit epoch.
• a: Semi-major axis [m].
• e: Eccentricity.
• i: Inclination [rad].
• Ω: Right ascension of the ascending node [rad].
• ω: Argument of perigee [rad].
• f: True anomaly [rad].

Orbit(a::Number, e::Number, i::Number, Ω::Number, ω::Number, f::Number)

Create an orbit with semi-major axis a [m], eccentricity e, inclination i [rad], right ascension of the ascending node Ω [rad], argument of perigee ω [rad], and true anomaly f [rad].

Returns

An object of type Orbit with the specified orbit. The orbit epoch is defined as 0.0.

OrbitPropagator{T}

Abstract type of the orbit propagator. Every propagator structure must be a subtype of this type and must implement the following API functions:

propagate!(orbp, t::Number)
propagate!(orbp, t::AbstractVector)
propagate_to_epoch!(orbp, JD::Number)
propagate_to_epoch!(orbp, JD::AbstractVector)
step!(orbp, Δt::Number)

OrbitPropagatorJ2{T} <: OrbitPropagator{T}

Structure that holds the information related to the J2 orbit propagator.

Fields

• orb: Mean orbital elements (see Orbit).
• j2d: Structure that stores the J2 orbit propagator data (see J2_Structure).

OrbitPropagatorJ4{T} <: OrbitPropagator{T}

Structure that holds the information related to the J4 orbit propagator.

Fields

• orb: Mean orbital elements (see Orbit).
• j4d: Structure that stores the J4 orbit propagator data (see J4_Structure).

OrbitPropagatorSGP4{T} <: OrbitPropagator{T}

Structure that holds the information related to the SGP4 propagator.

Fields

• orb: Mean orbital elements (see Orbit).
• sgp4_gc: Gravitational contents of the SGP4 algorithm (see SGP4_GravCte).
• sgp4d: Structure that stores the SGP4 data (see SGP4_Structure).

OrbitPropagatorTwoBody{T} <: OrbitPropagator{T}

Structure that holds the information related to the Two Body orbit propagator.

Fields

• orb: Mean orbital elements (see Orbit).
• tbd: Structure that stores the Two Body orbit propagator data (see TwoBody_Structure).

TwoBody_Structure{T}

Low level Two Body orbit propagator structure.

DatetoJD(dateTime::DateTime)

Convert the date and time dateTime to Julian Day.

DatetoJD(date::Date)

Convert the date date to Julian Day.

DatetoJD(Y::Int, M::Int, D::Int, h::Int, m::Int, s::Number)

Convert a date represented using the Gregorian Calendar (Year = y, Month = M (1-12), Day = D, Hour = h (0-24), minute = m, and second = s) to Julian Day.

Remarks

The algorithm was obtained from [2] (Accessed on 2018-04-11).

ECEFtoGeodetic(r_e::AbstractVector)

Convert the vector r_e [m] represented in the Earth-Centered, Earth-Fixed (ECEF) reference frame into Geodetic coordinates (WGS-84).

Returns

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

Remarks

Based on algorithm in [3].

E_to_M(e::Number, E::Number)

Compute the mean anomaly (0,2π) [rad] given the eccentricity e and the eccentric anomaly E [rad].

E_to_f(e::Number, E::Number)

Compute the true anomaly (0,2π) [rad] given the eccentricity e and the eccentric anomaly E [rad].

GeodetictoECEF(lat::Number, lon::Number, h::Number)

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

Remarks

Based on algorithm in [3].

GeodetictoGeocentric(ϕ_gd::Number, h::Number)

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

Returns

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

Remarks

Based on algorithm in [4, p. 3].

J2000toGMST(J2000_UT1::Number)

Compute the Greenwich Mean Sideral Time (GMST) [rad] given the instant J2000_UT1 in J2000.0 reference [UT1].

Remarks

Based on algorithm in 2, accessed at 2015-12-01.

JD_TTtoUTC(JD_TT::Number, ΔAT::Number = 37)

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

JD_UT1toUTC(JD_UT1::Number, ΔUT1::Number)

Convert the Julian Day in UT1 JD_UT1 to the Julian Day in UTC using the accumulated difference ΔUT1, which is provided by IERS EOP Data.

JD_UTCtoUT1(JD_UTC::Number, eop::Union{EOPData_IAU1980,EOPData_IAU2000A})

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

JD_UTCtoTT(JD_UTC::Number [, ΔAT::Number])

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

JD_UTCtoUT1(JD_UTC::Number, ΔUT1::Number)

Convert the Julian Day in UTC JD_UTC to the Julian Day in UT1 using the accumulated difference ΔUT1, which is provided by IERS EOP Data.

JD_UTCtoUT1(JD_UTC::Number, eop::Union{EOPData_IAU1980,EOPData_IAU2000A})

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

JDtoDate([T,] JD::Number)

Convert a date represented in Julian Day JD to Gregorian Calendar. The optional parameter T defines the return type. If T is omitted, then it defaults to Int.

Returns

If T is omitted or Int, then a tuple with the following data will be returned:

• Year.
• Month (1 => January, 2 => February, ...).
• Day.
• Hour (0 - 24).
• Minute (0 - 59).
• Second (0 - 59).

Notice that if T is Int, then the seconds field will be Integer. Otherwise, it will be floating point.

If T is Date, then it will return the Julia structure Date. Notice that the hours, minutes, and seconds will be neglected because the structure Date does not handle them.

If T is DateTime, then it will return the Julia structure DateTime.

Remarks

The algorithm was obtained from [2] (Accessed on 2018-04-11). In [2], there is the following warning:

Note: This method will not give dates accurately on the Gregorian Proleptic Calendar, i.e., the calendar you get by extending the Gregorian calendar backwards to years earlier than 1582. using the Gregorian leap year rules. In particular, the method fails if Y<400.

JDtoGMST(JD_UT1::Number)

Compute the Greenwich Mean Sideral Time (GMST) [rad] for the Julian Day JD_UT1 [UT1].

Remarks

Based on algorithm in [1, pp. 188].

M_to_E(e::Number, M::Number, tol::Number = 1e-10)

Compute the eccentric anomaly (0,2π) [rad] given the eccentricity e and the mean anomaly M [rad]. This function uses the Newton-Raphson algorithm and the tolerance to accept the solution is tol.

M_to_f(e::Number, M::Number, tol::Number = 1e-10)

Compute the true anomaly (0,2π) [rad] given the eccentricity e and the mean anomaly M [rad]. This function uses the Newton-Raphson algorithm and the tolerance to accept the solution is tol.

adjacent_track_angle_grss(h::Number, T::Number, i::Number, To::Int, lat::Number)

Compute the angle between two adjacent ground tracks [rad] in a given latitude lat [rad] measured from the satellite position for a ground repeating, Sun-synchronous orbit with altitude in the Equator h [m], period T [s], inclination i [rad], and orbit cycle To [days].

Remarks

The functions does not check if the orbit is a GRSS orbit.

adjacent_track_angle_grss(h::Number, a::Number, e::Number, i::Number, To::Int, lat::Number)

Compute the angle between two adjacent ground tracks [rad] in a given latitude lat [rad] measured from the satellite position for a ground repeating, Sun-synchronous orbit with altitude in the Equator h [m], semi-major axis a [m], eccentricity e, inclination i [rad], and orbit cycle To [days].

Remarks

The functions does not check if the orbit is a GRSS orbit.

adjacent_track_distance_grss(T::Number, i::Number, To::Int, lat::Number)

Compute the distance between adjacent ground tracks [m] at a given latitude lat [rad] for a ground repeating, Sun-synchronous orbit with period T [s], inclination i [rad], and orbit cycle To [days].

Remarks

The functions does not check if the orbit is a GRSS orbit.

adjacent_track_distance_grss(a::Number, e::Number, i::Number, To::Int, lat::Number)

Compute the distance between adjacent ground tracks [m] at a given latitude lat [rad] for a ground repeating, Sun-synchronous orbit with semi-major axis a [m], eccentricity e, inclination i [rad], and orbit cycle To [days].

Remarks

The functions does not check if the orbit is a GRSS orbit.

angvel(a::Number, e::Number, i::Number, pert::Symbol = :J2)
angvel(orb::Orbit, pert::Symbol = :J2)

Compute the angular velocity [rad/s] of an object in an orbit with semi-major axis a [m], eccentricity e, and inclination i [rad], using the perturbation terms specified by the symbol pert. The orbit can also be specified by orb, which is an instance of the structure Orbit.

pert can be:

• :J0: Consider a Keplerian orbit.
• :J2: Consider the perturbation terms up to J2.
• :J4: Consider the perturbation terms J2, J4, and J2².

If pert is omitted, then it defaults to :J2.

angvel_to_a(n::Number, e::Number, i::Number, pert::Symbol = :J2; μ::Number = m0, max_iter::Int = 20, tol::Number = 1e-10)

Compute the semi-major axis that will provide an angular velocity n [rad/s] in an orbit with eccentricity e and inclination i [rad], using the perturbation terms specified by the symbol pert.

Notice that the angular velocity n is related to the nodal period, i.e. the time between two consecutive passages by the ascending node.

pert can be:

• :J0: Consider a Keplerian orbit.
• :J2: Consider the perturbation terms up to J2.
• :J4: Consider the perturbation terms J2, J4, and J2².

If pert is omitted, then it defaults to :J2.

Keyword

• μ: Standard gravitational parameter for Earth [m^3/s^2]. (Default = m0)
• max_iter: Maximum number of iterations allowed in the Newton-Raphson algorithm. (Default = 20)
• tol: Tolerance to stop the Newton-Raphson algorithm. (Default = 1e-10)

beta_angle(JD₀::Number, a::Number, e::Number, i::Number, RAAN::Number, Δt::Integer, pert::Symbol = :J2)

Compute the beta angle of an orbit with semi-major axis a [m], eccentricity e, inclination i [rad], and initial right ascension of the ascending node RAAN [rad]. The orbit epoch, which is also the day in which the analysis will begin, is JD₀ [Julian Day]. The analysis will be performed for each day during Δt days.

The argument pert can be used to select the perturbation terms that must be used when propagating the right ascencion of the ascending node. The possible values are:

• :J0: Consider a Keplerian orbit.
• :J2: Consider the perturbation terms up to J2.
• :J4: Consider the perturbation terms J2, J4, and J2².

If pert is omitted, then it defaults to :J2.

Returns

An array with two columns. The first one contains the days of the analysis and the second one contains the beta angle [rad] for each day.

change_oe_frame(a::Number, e::Number, i::Number, Ω::Number, ω::Number, f::Number, conv_args...)
change_oe_frame(oe::Orbit, conv_args...)

Change the reference frame of orbit elements. The orbit elements can be specified by a, e, i, Ω, ω, and f, or the structure oe (see Orbit).

The conversion arguments conv_args are the same arguments that one should pass to the function rECItoECI to convert between the desired frames. For more information, see the documentation of the function rECItoECI.

Args

• a: Semi-major axis [m].

• e: Excentricity.

• i: Inclination [rad].

• Ω: Right-ascension of the ascending node [rad].

• ω: Argument of perigee [rad].

• f: True anomaly [rad].

• conv_args...: Conversion arguments, which are the same arguments that one would pass to the function rECItoECI to convert between the desired frames.

• oe: An instance of the structure Orbit with the orbit elements that will be converted [SI units].

Returns

An instance of the structure Orbit with the Keplerian elements [SI units] converted to the new frame.

Examples

julia> eop = get_iers_eop(:IAU1980);

julia> teme_epoch = DatetoJD(2016,6,1,11,0,0);

julia> tod_epoch  = DatetoJD(2016,1,1,0,0,0);

julia> oe_teme    = Orbit(0,
                          7130.982e3,
                          0.001111,
                          98.405*pi/180,
                          227.336*pi/180,
                          90*pi/180,
                          320*pi/180)

                 Orbit
  =====================================
                  t =        0.0
    Semi-major axis =     7130.9820 km
       Eccentricity =        0.001111
        Inclination =       98.4050 ˚
               RAAN =      227.3360 ˚
    Arg. of Perigee =       90.0000 ˚
       True Anomaly =      320.0000 ˚

julia> oe_j2000 = change_oe_frame(oe_teme, TEME(), J2000(), teme_epoch, eop)

                 Orbit
  ======================================
                  t =        0.0
    Semi-major axis =     7130.9820 km
       Eccentricity =        0.001111
        Inclination =       98.3365 ˚
               RAAN =      227.1345 ˚
    Arg. of Perigee =       90.0604 ˚
       True Anomaly =      320.0000 ˚

julia> oe_tod   = change_oe_frame(oe_teme, TEME(), teme_epoch, TOD(), tod_epoch, eop)

                 Orbit
  ======================================
                  t =        0.0
    Semi-major axis =     7130.9820 km
       Eccentricity =        0.001111
        Inclination =       98.4037 ˚
               RAAN =      227.3306 ˚
    Arg. of Perigee =       90.0014 ˚
       True Anomaly =      320.0000 ˚

compute_RAAN_lt(JD::Number, asc_node_lt::Number)

Compute the RAAN (0,2π) [rad] so that the orbit plane local time is asc_node_lt [hour] at the Julian day JD.

compute_U(gm_coefs::GravityModel_Coefs{T}, r::AbstractVector, n_max::Number = -1, m_max::Number = -1) where T<:Number

Compute the gravitational potential [J/kg] at r (ITRF) [m] using the coefficients gm_coefs (see GravityModel_Coefs). The maximum degree that will be used while computing the spherical harmonics will be n_max and the maximum order is m_max.

If n_max is negative, then the maximum available degree will be used. If n_max is omitted, then it defaults to 0.

If m_max is negative or if it is greater than n_max, then it will be set to n_max. If m_max is omitted, then it defaults to 0.

compute_dU(gm_coefs::GravityModel_Coefs{T}, r::AbstractVector, n_max::Number = -1, m_max::Number = -1) where T<:Number

Compute the derivatives w.r.t. the spherical coordinates of the gravitational field (∂U/∂r, ∂U/∂ϕ, ∂U/∂λ) defined by the coefficients gm_coefs (see GravityModel_Coefs) at the position r [m] in ITRF. The maximum degree that will be used while computing the spherical harmonics will be n_max and the maximum order is m_max.

If n_max is negative, then the maximum available degree will be used. If n_max is omitted, then it defaults to 0.

If m_max is negative or if it is greater than n_max, then it will be set to n_max. If m_max is omitted, then it defaults to 0.

Returns

• The derivative of the gravitational field w.r.t. the radius (∂U/∂r).
• The derivative of the gravitational field w.r.t. the latitude (∂U/∂ϕ).
• The derivative of the gravitational field w.r.t. the longitude (∂U/∂λ).

Remarks

In this case, ϕ is the geocentric latitude and λ is the geocentric longitude.

compute_g(gm_coefs::GravityModel_Coefs{T}, r::AbstractVector, n_max::Number = -1, m_max::Number = -1) where T<:Number

Compute the gravitational acceleration (ITRF) [m/s²] at position r [m] (ITRF) using the coefficients gm_coefs (see GravityModel_Coefs). The maximum degree that will be used while computing the spherical harmonics will be n_max and the maximum order it m_max.

If n_max is negative, then the maximum available degree will be used. If n_max is omitted, then it defaults to 0.

If m_max is negative or if it is greater than n_max, then it will be set to n_max. If m_max is omitted, then it defaults to 0.

Remarks

Notice that this function computes the gravitational acceleration. Hence, the acceleration due to Earth rotation rate is not included.

compute_ss_orbit_by_ang_vel(n::Number, e::Number)

Compute the Sun-synchronous orbit given the angular velocity n [rad/s] and the eccentricity e.

Returns

• The semi-major axis [m].
• The inclination [rad].
• The residues of the two functions.
• A boolean variable that indicates if the numerical algorithm converged.

compute_ss_orbit_by_inclination(i::Number, e::Number)

Compute the Sun-synchronous orbit given the inclination i [rad] and the eccentricity e.

Returns

The semi-major axis of the Sun-synchronous orbit [m].

compute_ss_orbit_by_num_rev_per_day(numRevPD::Number, e::Number)

Compute the Sun-synchronous orbit given the number of revolutions per day numRevPD and the eccentricity e.

Returns

• The semi-major axis [m].
• The inclination [rad].
• The residues of the two functions.
• A boolean variable that indicates if the numerical algorithm converged.

compute_ss_orbit_by_semi_major_axis(a::Number, e::Number)

Compute the Sun-synchronous orbit given the semi-major axis a [m] and the eccentricity e.

Returns

The inclination of the Sun-synchronous orbit [rad].

conf_nrlmsise00(year::Int, doy::Int, sec::Number, alt::Number, g_lat::Number, g_long::Number, lst::Number, f107A::Number, f107::Number, ap::[Number, AbstractVector], flags::NRLMSISE00_Flags = NRLMSISE00_Flags())

Create the structure with the proper configuration to call the NRLMSISE-00 model.

Notice that the input variables have the same units of the original model.

Args

• year: Year (currently ignored).
• doy: Day of year.
• sec: Seconds in day [UT].
• alt: Altitude [km].
• g_lat: Geodetic latitude [deg].
• g_long: Geodetic longitude [deg].
• lst: Local apparent solar time (hours).
• f107A: 81 day average of F10.7 flux (centered on day of year doy).
• f107: Daily F10.7 flux for previous day.
• ap: Magnetic index (daily) if it is a number. If it is an array, then see Remarks.
• flags: (OPTIONAL) An instance of the structure NRLMSISE00_Flags with the configuration flags for NRLMSISE00. If omitted, then the default configurations will be used.

Returns

An instance of the structure NRLMSISE00_Structure.

Remarks

If ap is a Vector, then it must be a vector with 7 dimensions as described below:

Notes on input variables

UT, Local Time, and Longitude are used independently in the model and are not of equal importance for every situation. For the most physically realistic calculation these three variables should be consistent (lst=sec/3600 + g_long/15). The Equation of Time departures from the above formula for apparent local time can be included if available but are of minor importance.

f107 and f107A values used to generate the model correspond to the 10.7 cm radio flux at the actual distance of the Earth from the Sun rather than the radio flux at 1 AU. The following site provides both classes of values:

ftp://ftp.ngdc.noaa.gov/STP/SOLAR_DATA/SOLAR_RADIO/FLUX/

f107, f107A, and ap effects are neither large nor well established below 80 km and these parameters should be set to 150, 150, and 4 respectively.

create_gravity_model_coefs(icgem::ICGEM)

Return an instance of the structure GravityModel_Coefs based on the information obtained from an ICGEM file in icgem (see parse_icgem).

dArgPer(a::Number, e::Number, i::Number, pert::Symbol = :J2)
dArgPer(orb::Orbit, pert::Symbol = :J2)

Compute the time-derivative of the argument of perigee [rad/s] of an orbit with semi-major axis a [m], eccentricity e, and inclination i [rad], using the perturbation terms specified by the symbol pert. The orbit can also be specified by orb, which is an instance of the structure Orbit.

pert can be:

• :J0: Consider a Keplerian orbit.
• :J2: Consider the perturbation terms up to J2.
• :J4: Consider the perturbation terms J2, J4, and J2².

If pert is omitted, then it defaults to :J2.

dRAAN(a::Number, e::Number, i::Number, pert::Symbol = :J2)
dRAAN(orb::Orbit, pert::Symbol = :J2)

Compute the time-derivative of the right ascension of the ascending node [rad/s] of an orbit with semi-major axis a [m], eccentricity e, and inclination i [rad], using the perturbation terms specified by the symbol pert. The orbit can also be specified by orb, which is an instance of the structure Orbit.

pert can be:

• :J0: Consider a Keplerian orbit.
• :J2: Consider the perturbation terms up to J2.
• :J4: Consider the perturbation terms J2, J4, and J2².

If pert is omitted, then it defaults to :J2.

dlegendre([N,] ϕ::Number, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false)

Compute the first-order derivative of the associated Legendre function P_n,m[cos(ϕ)] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

The optional parameter N can be used to select the normalization. The following values are valid:

• Val{:full}: Compute the fully normalized associated Legendre function (see legendre_fully_normalized).
• Val{:schmidt}: Compute the Schmidt quasi-normalized associated Legendre function (see legendre_schmidt_quasi_normalized).
• Val{:conv}: Compute the conventional associated Legendre function (see dlegendre_conventional!).

If N is omitted, then the full normalization will be used (Val{:full}).

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the first-order derivative of the Legendre associated functions P_n,m[cos(ϕ)].

dlegendre!([N,] dP::AbstractMatrix, ϕ::Number, P::AbstractMatrix, ph_term::Bool = false)

Compute the first-order derivative of the associated Legendre function P_n,m[x] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The derivatives will be stored in the matrix dP. Hence, the maximum degree and order that will be computed are given by the dimensions of this matrix.

This algorithm needs the matrix P with the associated Legendre function. This can be computed using the function legendre. Notice that this matrix must be computed using the same normalization (see N) as the one selected here.

The optional parameter N can be used to select the normalization. The following values are valid:

• Val{:full}: Compute the fully normalized associated Legendre function (see dlegendre_fully_normalized!).
• Val{:schmidt}: Compute the Schmidt quasi-normalized associated Legendre function (see dlegendre_schmidt_quasi_normalized!).
• Val{:conv}: Compute the conventional associated Legendre function (see dlegendre_conventional!).

If N is omitted, then the full normalization will be used (Val{:full}).

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

dlegendre_conventional!(dP::AbstractMatrix, ϕ::Number, P::AbstractMatrix, ph_term::Bool = false)

Compute the first-order derivative of the conventional associated Legendre function P_n,m[cos(ϕ)] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The derivatives will be stored in the matrix dP. Hence, the maximum degree and order that will be computed are given by the dimensions of this matrix.

This algorithm needs the matrix P with the conventional associated Legendre function. This can be computed using the function legendre_conventional.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Remarks

The user is responsible to pass a matrix P with the correct values. For example, if ph_term is true, then P must also be computed with ph_term set to true.

dlegendre_conventional(ϕ::Number, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false)

Compute the first-order derivative of the conventional associated Legendre function P_n,m[cos(ϕ)] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the first-order derivative of the Legendre associated functions P_n,m[cos(ϕ)].

dlegendre_fully_normalized!(dP::AbstractMatrix, ϕ::Number, P::AbstractMatrix, ph_term::Bool = false)

Compute the first-order derivative of the fully normalized associated Legendre function P_n,m[cos(ϕ)] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The derivatives will be stored in the matrix dP. Hence, the maximum degree and order that will be computed are given by the dimensions of this matrix.

This algorithm needs the matrix P with the fully normalized associated Legendre function. This can be computed using the function legendre_fully_normalized.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Remarks

The user is responsible to pass a matrix P with the correct values. For example, if ph_term is true, then P must also be computed with ph_term set to true.

dlegendre_fully_normalized(ϕ::T, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false) where T<:AbstractFloat

Compute the first-order derivative of the Schmidt fully normalized associated Legendre function P_n,m[cos(ϕ)] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the first-order derivative of the Legendre associated functions P_n,m[cos(ϕ)].

dlegendre_schmidt_quasi_normalized!(dP::AbstractMatrix, ϕ::Number, P::AbstractMatrix, ph_term::Bool = false)

Compute the first-order derivative of the Schmidt quasi-normalized associated Legendre function P_n,m[cos(ϕ)] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The derivatives will be stored in the matrix dP. Hence, the maximum degree and order that will be computed are given by the dimensions of this matrix.

This algorithm needs the matrix P with the Schmidt quasi-normalized associated Legendre function. This can be computed using the function legendre_schmidt_quasi_normalized.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Remarks

The user is responsible to pass a matrix P with the correct values. For example, if ph_term is true, then P must also be computed with ph_term set to true.

dlegendre_schmidt_quasi_normalized(ϕ::T, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false) where T<:AbstractFloat

Compute the first-order derivative of the Schmidt quasi-normalized associated Legendre function P_n,m[cos(ϕ)] w.r.t. ϕ [rad]:

dP_n,m[cos(ϕ)]
--------------
      dϕ

The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the first-order derivative of the Legendre associated functions P_n,m[cos(ϕ)].

eclipse_time_summary(JD₀::Number, a::Number, e::Number, i::Number, RAAN::Number, w::Number, Δd::Integer, relative::Bool = false, Δt₀::AbstractFloat = -1.0)

Compute the eclipse time of an orbit with semi-major axis a [m], eccentricity e, inclination i [rad], initial right ascension of the ascending node RAAN [rad], and initial argument of perigee w [rad]. The orbit epoch, which is also the day in which the analysis will begin, is JD₀ [Julian Day]. The analysis will be performed for each day during Δd days.

This function will compute the eclipse time of one orbit per day.

If the argument relative is true, then the computed times will be relative to the nodal period [%]. Otherwise, they will be computed in seconds. By default, relative = false.

The argument Δt₀ can be used to select the time step in which the orbit will be propagated. Notice that this algorithm performs a numerical search to find the beginning of each section (sunlight, penumbra, and umbra) with millisecond precision. Thus, selecting a high number for Δt₀ will make the analysis faster, but the accuracy is lost if a region time span is smalled than Δt₀. If this parameter is omitted or if it is negative, then the time step will be selected automatically to match a mean anomaly step of 5°.

All the analysis is performed using a J2 orbit propagator.

Returns

The following table:

    day | Sunlight Time | Penumbra Time | Umbra Time
   -----+---------------+---------------+------------

epoch(orbp)

Return the epoch of the propagator orbp [JD].

equation_of_time(JD::Number)

Compute the difference between the Sun apparent local time and the Sun mean local time [rad], which is called Equation of Time, at the Julian Day JD. The algorithm was adapted from [1, p. 178, 277-279].

expatmosphere(h::Number)

Compute the atmospheric density [kg/m³] at the altitude h [m] (above the ellipsoid) using the exponential atmospheric model:

                ┌            ┐
                │    h - h₀  │
ρ(h) = ρ₀ ⋅ exp │ - ──────── │ ,
                │      H     │
                └            ┘

in which ρ₀, h₀, and H are parameters obtained from tables that depend only on h.

f_to_E(e::Number,f::Number)

Compute the eccentric anomaly (0,2π) [rad] given the eccentricity e and the true anomaly f [rad].

f_to_E(orb::Orbit)

Compute the eccentric anomaly (0,2π) [rad] given the orbit orb (see Orbit).

f_to_M(e::Number, f::Number)

Compute the mean anomaly (0,2π) [rad] given the eccentricity e and the true anomaly f [rad].

f_to_M(orb::Orbit)

Compute the mean anomaly (0,2π) [rad] given the orbit orb (see Orbit).

geomag_dipole(r_e::AbstractVector, pole_lat::Number, pole_lon::Number, m::Number)

Compute the geomagnetic field [nT] using the simplified dipole model at position r_e (ECEF reference frame). This function considers that the latitude of the South magnetic pole (which lies in the North hemisphere) is pole_lat [rad] and the longitude is pole_lon [rad]. Furthermore, the dipole moment is considered to be m [A.m²].

geomag_dipole(r_e::AbstractVector, year::Number = 2019)

Compute the geomagnetic field [nT] using the simplified dipole model at position r_e (ECEF reference frame). This function uses the year year to obtain the position of the South magnetic pole (which lies in the North hemisphere) and the dipole moment. If year is omitted, then it will be considered as 2019.

Remarks

In both functions, the output vector will be represented in the ECEF reference frame.

get_Ap(JD::Number; mean::Tuple{Int} = (), daily = false)

Return the Ap index.

If mean is a tuple of two integers (hi, hf), then the average between hi and hf previous hours will be computed.

If mean is empty and daily is true, then the day average will be computed.

If mean keyword is empty, and daily keyword is false, then the Ap at Julian day JD will be computed.

By default, mean is empty and daily is false.

get_DstΔTc(JD::Number)

Get the value of the index DstΔTc at Julian Day JD.

This function requires the initialization of the variable _dtcfile_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_F10(JD::Number)

Get the value of the index F10 at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_F81a(JD::Number)

Get the value of the index F81a at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_Kp(JD::Number)

Return the Kp index at Julian Day JD.

get_M10(JD::Number)

Get the value of the index M10 at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_M81a(JD::Number)

Get the value of the index M81a at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_S10(JD::Number)

Get the value of the index S10 at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_S81a(JD::Number)

Get the value of the index S81a at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_Y10(JD::Number)

Get the value of the index Y10 at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_Y81a(JD::Number)

Get the value of the index Y81a at Julian Day JD.

This function requires the initialization of the variable _solfsmy_data. Otherwise, an exception will be raised. To initialize it, run init_space_indices().

get_iers_eop(data_type::Symbol = :IAU1980; force_download = false)

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

• IAU1980: Get IERS EOP C04 IAU1980 data.
• IAU2000A: Get IERS EOP C04 IAU2000A data.

If data_type is omitted, then it defaults to IAU1980.

The files are downloaded using the RemoteFile package with daily updates. Hence, if one desires to force a download before the scheduled time, then set the keyword force_download to true.

Returns

A structure (EOPData_IAU1980 or EOPData_IAU2000A, depending on data_type) with the interpolations of the EOP parameters. Notice that the interpolation indexing is set to the Julian Day.

get_iers_eop_iau_1980(url::String = "https://datacenter.iers.org/data/latestVersion/223_EOP_C04_14.62-NOW.IAU1980223.txt")

Get the IERS EOP C04 IAU1980 data from the URL url. If url is omitted, then it defaults to https://datacenter.iers.org/data/latestVersion/223EOPC04_14.62-NOW.IAU1980223.txt

The file is downloaded using the RemoteFile package with daily updates. Hence, if one desires to force a download before the scheduled time, then set the keyword force_download to true.

Returns

The structure EOPData_IAU1980 with the interpolations of the EOP parameters. Notice that the interpolation indexing is set to the Julian Day.

Remarks

For every field in EOPData_IAU1980 to interpolation between two points in the grid is linear. If extrapolation is needed, then if will use the nearest value (flat extrapolation).

get_iers_eop_iau_2000A(url::String = "https://datacenter.iers.org/data/latestVersion/224_EOP_C04_14.62-NOW.IAU2000A224.txt"; force_download = false)

Get the IERS EOP C04 IAU2000A data from the URL url. If url is omitted, then it defaults to https://datacenter.iers.org/data/latestVersion/224EOPC04_14.62-NOW.IAU2000A224.txt

The file is downloaded using the RemoteFile package with daily updates. Hence, if one desires to force a download before the scheduled time, then set the keyword force_download to true.

Returns

The structure EOPData_IAU2000A with the interpolations of the EOP parameters. Notice that the interpolation indexing is set to the Julian Day.

Remarks

For every field in EOPData_IAU2000A to interpolation between two points in the grid is linear. If extrapolation is needed, then if will use the nearest value (flat extrapolation).

get_space_index(T, JD::Number; ...)

Return the space index T at the day JD [Julian Day]. T can be:

Daily 10.7-cm solar flux

The daily 10.7-cm solar flux can be obtained using:

• F10(): 10.7-cm adjusted solar flux [10⁻²² W/(M² Hz)].
• F10adj(): 10.7-cm adjusted solar flux [10⁻²² W/(M² Hz)].
• F10obs(): 10.7-cm observed solar flux [10⁻²² W/(M² Hz)].

These indices require fluxtable (see init_space_indices).

Daily average 10.7-cm solar flux

The daily average 10.7-cm solar flux, centered at JD, can be obtained using:

• F10M(): 10.7-cm adjusted solar flux [10⁻²² W/(M² Hz)].
• F10Madj(): 10.7-cm adjusted solar flux [10⁻²² W/(M² Hz)].
• F10Mobs(): 10.7-cm observed solar flux [10⁻²² W/(M² Hz)].

In this case, the keyword window::Int can be passed to select the size of the window. By default, it is selected as 81.

These indices require fluxtable (see init_space_indices).

Daily Kp and Ap

• Kp(): Kp index (daily mean).
• Kp_vect(): A vector containing the Kp index for the following hours of the day: 0-3h, 3-6h, 6-9h, 9-12h, 12-15h, 15-18h, 18-20h, 20-23h.
• Ap(): Ap index (daily mean).
• Ap_vect(): A vector containing the Ap index for the following hours of the day: 0-3h, 3-6h, 6-9h, 9-12h, 12-15h, 15-18h, 18-20h, 20-23h.

These indices require wdcfiles (see init_space_indices).

Daily S10, M10, and Y10

• S10(): EUV index (26-34 nm) scaled to F10.7.
• M10(): MG2 index scaled to F10.7.
• Y10(): Solar X-ray & Lya index scaled to F10.7.

These indices require solfsmy (see init_space_indices).

81-day centered average of S10, M10, and Y10.

• S81a: EUV 81-day averaged centered index.
• M81a: MG2 81-day averaged centered index.
• Y81a: Solar X-ray & Lya 81-day averaged centered index.

These indices require solfsmy (see init_space_indices).

Exospheric temperature variation due to Dst

• DstΔTc: Exospheric temperature variation due to Dst [K].

This index requires dtcfile (see init_space_indices).

get_ΔAT(JD::Number)

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

Remarks

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

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

ground_station_accesses(orbp, vrs_e,     Δt, ECI, ECEF, vargs...; kwargs...)
ground_station_accesses(orbp, [(WGS84)], Δt, ECI, ECEF, vargs...; kwargs...)

Compute the accesses of a satellite with orbit propagator orbp (see init_orbit_propagator) to the ground stations defined in the vector vrs_e. The analysis interval begins in the propagator epoch and lasts Δt [s].

The ground stations can be specified by an array of 3×1 vectors describing the ground stations position in an ECEF frame vrs_e or by an array of tuples containing the WGS84 position of each ground station [(WGS84)]:

(latitude [rad], longitude [rad], altitude [m])

Args

• ECI: Earth-Centered Inertial frame in which the state vector of the propagator is represented.
• ECEF: Earth-Centered, Earth-fixed frame to be used for the analysis. It must be the same frame used to compute the ground station position vector.
• vargs...: list of additional arguments to be passed to the function rECItoECEF when converting the ECI frame to the ECEF.

Keywords

• θ: Minimum elevation angle for communication between the satellite and the ground stations [rad]. (Default = 10ᵒ)
• reduction: A function that receives a boolean vector with the visibility between the satellite and each ground station. It must return a boolean value indicating if the access must be computed or not. This is useful to merge access time between two or more stations. (Default = v->|(v...) i.e. compute the access if at least one ground station is visible)

ground_station_gaps(args...; kwargs...)

Compute the gaps between the accesses of ground stations. The arguments and keywords are the same as the ones used in the function ground_station_accesses.

Notice that the gap analysis starts in the orbit propagator epoch and ends in the instant defined by the argument Δt.

ground_station_visible(r_e::AbstractVector, rs_e::AbstractVector, θ::Number)

Check if the satellite with position vector r_e (ECEF) is inside the visibility circle of a ground station with position vector rs_e (ECEF) and a minimum elevation angle of θ [rad].

Notice that r_e and rs_e must be represented in the same ECEF frame, and must have the same unit.

Returns true if the satellite is inside the visibility circle, or false otherwise.

ground_station_visible(r_e::AbstractVector, lat_s::Number, lon_s::Number, h_s::Number, θ::Number)

Check if the satellite with position vector r_e (ECEF) is inside the visibility circle of a ground station with latitude lat_s [rad], longitude lon_s [rad], altitude h_s (WGS-84), and a minimum elevation angle of θ [rad].

Notice that the units of r_e and h_s must be the same.

Returns true if the satellite is inside the visibility circle, or false otherwise.

ground_trace(orbp::OrbitPropagator{N}, eop_data::Union{Nothing, EOPData_IAU1980, EOPData_IAU2000A} = nothing; ECI = TEME(), ECEF = PEF(), span = 1.0) where N

Compute the ground trace of the object with orbit defined by orbp.

By default, it considers that the orbit elements on the propagator are represented in the True Equator, Mean Equinox (TEME) reference frame and the ground trace will be computed in the Pseudo-Earth Fixed (PEF) reference frame. Hence, no EOP data is needed. However, this can be changed by the keywords presented as follows.

Keywords

• eop_data: EOP data that will be used to convert the ECI reference frame to the ECEF reference frame. If nothing, then it will not be used (see rECItoECEF). (Default = nothing)
• ECI: ECI frame in which the orbit elements in orbp are represented. (Default = TEME())
• ECEF: ECEF frame that will be used to compute the ground trace. (Default = PEF())
• span: Defines for how much time the ground trace will be computed. The unit is the orbit period. (Default = 1.0)
• dt: Time interval between two samples [s]. (Default = 10.0)

Returns

A vector of tuples with the pairs (latitude,longitude) of the ground trace.

gtd7(nrlmsise00d::NRLMSISE00_Structure{T}) where T<:Number

NRLMSISE-00

Neutral Atmosphere Empirical Model from the surface to lower exosphere.

This routine computes the NRLMSISE-00 outputs (see NRLMSISE00_Output) using the configurations in the structure nrlmsise00 (see NRLMSISE00_Structure).

Args

• nrlmsise00d: An instance of NRLMSISE00_Structure.

Returns

An instance of structure NRLMSISE00_Output with the outputs.

In this case, the total mass den_Total (see NRLMSISE00_Output) is the sum of the mass densities of the species He, O, N₂, O₂, Ar, H, and N, but does not include anomalous oxygen.

Remarks

1. The densities of O, H, and N are set to 0 below 72.5 km.
2. The exospheric temperature T_exo is set to global average for altitudes below 120 km. The 120 km gradient is left at global average value for altitudes below 72.5 km.
3. Anomalous oxygen is defined as hot atomic oxygen or ionized oxygen that can become appreciable at high altitudes (> 500 km) for some ranges of inputs, thereby affection drag on satellites and debris. We group these species under the term Anomalous Oxygen, since their individual variations are not presently separable with the drag data used to define this model component.

gtd7d(nrlmsise00d::NRLMSISE00_Structure{T}) where T<:Number

NRLMSISE-00

Neutral Atmosphere Empirical Model from the surface to lower exosphere.

This routine computes the NRLMSISE-00 outputs (see NRLMSISE00_Output) using the configurations in the structure nrlmsise00 (see NRLMSISE00_Structure).

Args

• nrlmsise00d: An instance of NRLMSISE00_Structure.

Returns

An instance of structure NRLMSISE00_Output with the outputs.

In this case, the total mass den_Total (see NRLMSISE00_Output) is the effective total mass density for drag and is the sum of the mass densities of all species in this model including the anomalous oxygen.

Remarks

1. The densities of O, H, and N are set to 0 below 72.5 km.
2. The exospheric temperature T_exo is set to global average for altitudes below 120 km. The 120 km gradient is left at global average value for altitudes below 72.5 km.
3. Anomalous oxygen is defined as hot atomic oxygen or ionized oxygen that can become appreciable at high altitudes (> 500 km) for some ranges of inputs, thereby affection drag on satellites and debris. We group these species under the term Anomalous Oxygen, since their individual variations are not presently separable with the drag data used to define this model component.

igrf12(date::Number, r::Number, λ::Number, Ω::Number, T; show_warns = true)

IGRF v12 Model

Compute the geomagnetic field vector [nT] at the date date [Year A.D.] and position (r, λ, Ω).

The position representation is defined by T. If T is Val{:geocentric}, then the input must be geocentric coordinates:

1. Distance from the Earth center r [m];
2. Geocentric latitude λ (-π/2, +π/2) [rad]; and
3. Geocentric longitude Ω (-π, +π) [rad].

If T is Val{:geodetic}, then the input must be geodetic coordinates:

1. Altitude above the reference ellipsoid h (WGS-84) [m];
2. Geodetic latitude λ (-π/2, +π/2) [rad]; and
3. Geodetic longitude Ω (-π, +π) [rad].

If T is omitted, then it defaults to Val{:geocentric}.

Notice that the output vector will be represented in the same reference system selected by the parameter T (geocentric or geodetic). The Y-axis of the output reference system always points East. In case of geocentric coordinates, the Z-axis points toward the center of Earth and the X-axis completes a right-handed coordinate system. In case of geodetic coordinates, the X-axis is tangent to the ellipsoid at the selected location and points toward North, whereas the Z-axis completes a right-hand coordinate system.

Keywords

• show_warns: Show warnings about the data (Default = true).

Remarks

The date must be greater or equal to 1900 and less than or equal 2025. Notice that a warning message is printed for dates greater than 2020.

Disclaimer

This function is an independent implementation of the IGRF model. It contains a more readable code than the original one in FORTRAN, because it uses features available in Julia language.

igrf12syn(isv::Int, date::Number, itype::Int, alt::Number, colat::Number, elong::Number; show_warns = true)

This is a Julia implementation of the official IGRF source code, which was written in Fortran [2]. The input and output variables are exactly the same as the ones described in the function igrf12syn in [2].

Args

• isv: 0 if main-field values are required, 1 if secular variation values are required.
• date: Year A.D.
• itype: 1 if geodetic (spheroid), 2 if geocentric (sphere).
• alt: Height above sea level [km] if itype = 1, or distance from the center of Earth [km] if itype = 2 (must be > 3485 km).
• colat: Colatitude (0 - 180) [˚].
• elong: East-Longitude (0 - 360) [˚].

Keywords

• show_warns: Show warnings about the data (Default = true).

Returns

• The north component [nT] if isv = 0, or [nT/year] if isv = 1.
• The east component [nT] if isv = 0, or [nT/year] if isv = 1.
• The vertical component [nT] if isv = 0, or [nT/year] if isv = 1.
• The total intensity if isv = 0, or rubbish if isv = 1.

Remarks

• The date must be greater or equal to 1900 and less than or equal 2025.

Notice that a warning message is printed for dates grated than 2020.

init_orbit_propagator(T, tle::TLE, ...)

Initialize the orbit propagator T using the TLE tle. The propagator type T can be:

• Val{:J2}: J2 orbit propagator;
• Val{:J4}: J4 orbit propagator;
• Val{:twobody}: Two-body orbit propagator; or
• Val{:sgp4}: SGP4 orbit propagator.

Additional optional arguments for the J2 orbit propagator

The initialization function of the J2 orbit propagator can receive the following optional parameter:

• j2_gc: (OPTIONAL) J2 orbit propagator gravitational constants (Default = j2_gc_egm08).

Additional optional arguments for the J4 orbit propagator

The initialization function of the J4 orbit propagator can receive the following optional parameter:

• j4_gc: (OPTIONAL) J4 orbit propagator gravitational constants (Default = j4_gc_egm08).

Additional optional arguments for the two body orbit propagator

The initialization function of the two body orbit propagator can receive the following optional parameter:

• μ: (OPTIONAL) Standard gravitational parameter of the central body [m^3/s^2] (Default = m0).

Additional optional arguments for the SGP4 orbit propagator

The initialization function of the SGP4 orbit propagator can receive the following optional parameter:

• sgp4_gc: (OPTIONAL) Gravitational constants (Default = sgp4_gc_wgs84).

Returns

A new instance of the orbit propagator structure that stores the information of the orbit propagator.

Remarks

The SGP4 implementation includes also the deep space perturbations, which was originally called SDP4 algorithm. Modern approaches, such as [2] and [3], identifies if the selected orbit must be propagated using the deep space perturbations and automatically applied them. This is sometimes called SGDP4 algorithm.

init_orbit_propagator(T, epoch::Number, a_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, f_0::Number, ...)
init_orbit_propagator(T, orb_0::Orbit, ...)

Initialize the orbit propagator T using the initial mean orbital elements. The propagator type T can be:

• Val{:J2}: J2 orbit propagator;
• Val{:J4}: J4 orbit propagator; or
• Val{:twobody}: Two-body orbit propagator.

The mean orbital elements can be passed individually of using an instance of the structure Orbit.

Args

• epoch: Initial orbit epoch [Julian Day].
• a_0: Initial semi-major axis [m].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascension of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• f_0: Initial true anomaly [rad].
• n_0: Initial angular velocity [rad/s].
• M_0: Initial mean anomaly [rad].
• orb_0: Instance of the structure Orbit with the initial mean orbital elements [SI].

Additional optional arguments for the J2 orbit propagator

The initialization function of the J2 orbit propagator can receive the following optional parameters:

• dn_o2: (OPTIONAL) First time derivative of mean motion divided by 2 [rad/s²] (Default = 0).
• ddn_o6: (OPTIONAL) Second time derivative of mean motion divided by 6 [rad/s³] (Default = 0).
• j2_gc: (OPTIONAL) J2 orbit propagator gravitational constants (Default = j2_gc_egm08).

Additional optional arguments for the J4 orbit propagator

The initialization function of the J4 orbit propagator can receive the following optional parameters:

• dn_o2: (OPTIONAL) First time derivative of mean motion divided by 2 [rad/s²] (Default = 0).
• ddn_o6: (OPTIONAL) Second time derivative of mean motion divided by 6 [rad/s³] (Default = 0).
• j4_gc: (OPTIONAL) J4 orbit propagator gravitational constants (Default = j4_gc_egm08).

Additional optional arguments for the two body orbit propagator

The initialization function of the two body orbit propagator can receive the following optional parameter:

• μ: (OPTIONAL) Standard gravitational parameter of the central body [m^3/s^2] (Default = m0).

Returns

A new instance of the orbit propagator structure that stores the information of the orbit propagator.

Remarks

If the orbit is defined in terms of the angular velocity (mean motion) instead of the semi-major axis, then it is possible to use the function angvel_to_a to convert.

init_space_indices(...)

Initialize all space indices. The files that will be initialized must be indicated by the array of symbols passed to the keyword argument enabled_files. If this is nothing, which is the default, then all files will be initialized. The symbol related to each file is described next.

Notice that the initialization process can be changed by a set of keywords as described next.

DTCFILE

Symbol: :dtcfile

This file contains the exospheric temperature variation caused by the Dst index. This is used for the JB2008 atmospheric model.

Keywords

• dtcfile_path: Path for the file DTCFILE.TXT. If nothing, then it will be downloaded. (Default = nothing)
• dtcfile_force_download: If true, then the file will always be downloaded if the path is not specified. (Default = false).

fluxtable

Symbol: :fluxtable

This file contains the F10.7 flux data in different formats.

Keywords

• fluxtable_path: Path for the file fluxtable.txt. If nothing, then it will be downloaded. (Default = nothing)
• fluxtable_force_download: If true, then the file will always be downloaded if the path is not specified. (Default = false).

SOLFSMY

Symbol: :solfsmy

This files contains the indices necessary for the JB2008 atmospheric model.

Keywords

• solfsmy_path: Path for the file SOLFSMY.TXT. If nothing, then it will be downloaded. (Default = nothing)
• solfsmy_force_download: If true, then the file will always be downloaded if the path is not specified. (Default = false).

WDC Files

Symbol: :wdcfiles

This set of files contain the Kp and Ap indices.

Keywords

• wdcfiles_path: Path for the directory with the WDC files. If nothing, then they will be downloaded. (Default = nothing)
• wdcfiles_force_download: If true, then the files will always be downloaded if the path is not specified. (Default = false).
• wdcfiles_oldest_year: Oldest year in which the WDC file will be obtained. (Default = past 3 years).
• wdcfiles_newest_year: Newest year in which the WDC file will be obtained. If it is nothing, then it defaults to the current year. (Default = nothing).

is_leap_year(year::Int)

Check if the year year is a leap year. It returns true if year is a leap year, or false otherwise.

Remarks

This algorithm was based on [3].

j2!(j2d::J2_Structure{T}, t::Number) where T

Propagate the orbit defined in j2d (see J2_Structure) until the time t [s]. Notice that the values in j2d will be modified.

Returns

• The position vector represented in the inertial frame at time t [m].
• The velocity vector represented in the inertial frame at time t [m/s]

Remarks

The inertial frame in which the output is represented depends on which frame it was used to generate the orbit parameters. If the orbit parameters are obtained from a TLE, then the inertial frame will be TEME. Notice, however, that the perturbation theory requires an inertial frame with true equator.

j2_init(j2_gc::J2_GravCte{T}, epoch::Number, a_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, f_0::Number, dn_o2::Number, ddn_o6::Number) where T

Initialize the data structure of J2 orbit propagator algorithm.

Args

• j2_gc: J2 orbit propagator gravitational constants (see J2_GravCte).
• epoch: Epoch of the orbital elements [Julian Day].
• a_0: Initial semi-major axis [m].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascension of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• f_0: Initial true anomaly [rad].
• dn_o2: First time derivative of the mean motion divided by two [rad/s^2].
• ddn_o6: Second time derivative of the mean motion divided by six [rad/s^3].

Returns

The structure J2_Structure with the initialized parameters.

Remarks

The inputs are the mean orbital elements.

j4!(j4d::J4_Structure{T}, t::Number) where T

Propagate the orbit defined in j4d (see J4_Structure) until the time t [s]. Notice that the values in j4d will be modified.

Returns

• The position vector represented in the inertial frame at time t [m].
• The velocity vector represented in the inertial frame at time t [m/s]

Remarks

The inertial frame in which the output is represented depends on which frame it was used to generate the orbit parameters. If the orbit parameters are obtained from a TLE, then the inertial frame will be TEME. Notice, however, that the perturbation theory requires an inertial frame with true equator.

j4_init(j4_gc::J4_GravCte{T}, epoch::Number, n_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, M_0::Number) where T

Initialize the data structure of J4 orbit propagator algorithm.

Args

• j4_gc: J4 orbit propagator gravitational constants (see J4_GravCte).
• epoch: Epoch of the orbital elements [Julian Day].
• a_0: Initial semi-major axis [m].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascension of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• f_0: Initial true anomaly [rad].
• dn_o2: First time derivative of the mean motion divided by two [rad/s^2].
• ddn_o6: Second time derivative of the mean motion divided by six [rad/s^3].

Returns

The structure J4_Structure with the initialized parameters.

Remarks

The inputs are the mean orbital elements.

jb2008(JD::Number, glat::Number, glon::Number, h::Number)
jb2008(JD::Number, glat::Number, glon::Number, h::Number, F10::Number, F10ₐ::Number, S10::Number, S10ₐ::Number, M10::Number, M10ₐ::Number, Y10::Number, Y10ₐ::Number, DstΔTc::Number)

Compute the atmospheric density using the Jacchia-Bowman 2008 (JB2008) model.

If the space indices are not provided (first call), then they will be obtained from the online database. In this case, the function init_space_indices() must be called first and the function will throw an exception if the selected JD is outside of the available data.

This model is a product of the Space Environment Technologies, more information can be seen in the websites:

http://sol.spacenvironment.net/jb2006/

http://sol.spacenvironment.net/jb2008/

Args

• JD: Julian day.

• glat: Geocentric latitude [rad].

• glon: Geocentric longitude [rad].

• h: Altitude [m].

• F10: 10.7-cm solar flux [10⁻²² W/(M² Hz)] (Tabular time 1 day earlier).

• F10ₐ: 10.7-cm averaged solar flux, 81-day centered on input time (Tabular time 1 day earlier).

• S10: EUV index (26-34 nm) scaled to F10.7 (Tabular time 1 day earlier).

• S10ₐ: EUV 81-day averaged centered index (Tabular time 1 day earlier).

• M10: MG2 index scaled to F10.7 (Tabular time 2 days earlier).

• M10ₐ: MG2 81-day averaged centered index (Tabular time 2 days earlier).

• Y10: Solar X-ray & Lya index scaled to F10.7 (Tabular time 5 days earlier).

• Y10ₐ: Solar X-ray & Lya 81-day averaged centered index (Tabular time 5 days earlier).

• DstΔTc: Temperature variation related to the Dst.

Returns

An instance of the structure JB2008_Output with the computed values.

jr1971(JD::Number, glat::Number, glon::Number, h::Number, F10::Number, F10ₐ::Number, Kp::Number)

Compute the atmospheric density using the Jacchia-Roberts 1971 model.

Args

• JD: Julian day.
• glat: Geodetic latitude [rad].
• glon: Geodetic longitude [rad].
• h: Altitude [m].
• F10: 10.7-cm solar flux [10⁻²² W/(M² Hz)].
• F10ₐ: 10.7-cm averaged solar flux, 81-day centered on input time.
• Kp: Kp geomagnetic index (with a delay of 3 hours).

Returns

An instance of the structure JR1971_Output with the computed values.

kepler_to_rv(a::Number, e::Number, i::Number, Ω::Number, ω::Number, f::Number)
kepler_to_rv(o::Orbit)

Convert the Keplerian elements (a, e, i, Ω, ω, and f) to a Cartesian representation (position vector r and velocity vector v). The Keplerian elements can also be passed inside an instance of the Orbit structure.

Args

• a: Semi-major axis [m].
• e: Excentricity.
• i: Inclination [rad].
• Ω: Right ascension of the ascending node [rad].
• ω: Argument of perigee [rad].
• f: True anomaly [rad].

Returns

• The position vector represented in the inertial reference frame [m].
• The velocity vector represented in the inertial reference frame [m].

References

This algorithm was adapted from [1] and [3, p. 37-38].

kepler_to_sv(orb::Orbit)

Convert the Keplerian elements in the structure orb to a state vector.

legendre([N,] ϕ::Number, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false)

Compute the associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

The optional parameter N can be used to select the normalization. The following values are valid:

• Val{:full}: Compute the fully normalized associated Legendre function (see legendre_fully_normalized).
• Val{:schmidt}: Compute the Schmidt quasi-normalized associated Legendre function (see legendre_schmidt_quasi_normalized).
• Val{:conv}: Compute the conventional associated Legendre function (see legendre_conventional).

If N is omitted, then the full normalization will be used (Val{:full}).

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the Legendre associated functions P_n,m[cos(ϕ)].

legendre!([N,] P::AbstractMatrix, ϕ::Number, ph_term::Bool = false)

Compute the associated Legendre function P_n,m[cos(ϕ)]. The maximum degree and order that will be computed are given by the dimensions of matrix P.

The result will be stored at matrix P.

The optional parameter N can be used to select the normalization. The following values are valid:

• Val{:full}: Compute the fully normalized associated Legendre function (see legendre_fully_normalized!).
• Val{:schmidt}: Compute the Schmidt quasi-normalized associated Legendre function (see legendre_schmidt_quasi_normalized!).
• Val{:conv}: Compute the conventional associated Legendre function (see legendre_conventional!).

If N is omitted, then the full normalization will be used.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

legendre_conventional!(P::AbstractMatrix, ϕ::Number, ph_term::Bool = false)

Compute the conventional associated Legendre function P_n,m[cos(ϕ)]. The maximum degree and order that will be computed are given by the dimensions of matrix P:

maximum degree -> number of rows
maximum order  -> number of columns

The result will be stored at matrix P.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

legendre_conventional(ϕ::T, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false) where T<:AbstractFloat

Compute the conventional associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the Legendre associated functions P_n,m[cos(ϕ)].

legendre_fully_normalized!(P::AbstractMatrix, ϕ::Number, ph_term::Bool = false)

Compute the fully normalized associated Legendre function P_n,m[cos(ϕ)]. The maximum degree and order that will be computed are given by the dimensions of matrix P:

maximum degree -> number of rows
maximum order  -> number of columns

The result will be stored at matrix P.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Remarks

This algorithm was based on [1]. Our definition of fully normalized associated Legendre function can be seen in [2, p. 546]. The conversion is obtained by:

             _                     -
            |  (n-m)! . k . (2n+1)  |      k = 1 if m  = 0
K_n,m = sqrt| --------------------- |,     k = 2 if m != 0
            |         (n+m)!        |
             -                     -
_
P_n,m = P_n,m * K_n,m,

      _
where P_n,m is the fully normalized Legendre associated function.

legendre_fully_normalized(ϕ::T, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false) where T<:AbstractFloat

Compute the fully normalized associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the Legendre associated functions P_n,m[cos(ϕ)].

Remarks

This algorithm was based on [1]. Our definition of fully normalized associated Legendre function can be seen in [2, p. 546]. The conversion is obtained by:

             _                     -
            |  (n-m)! . k . (2n+1)  |      k = 1 if m  = 0
K_n,m = sqrt| --------------------- |,     k = 2 if m != 0
            |         (n+m)!        |
             -                     -
_
P_n,m = P_n,m * K_n,m,

      _
where P_n,m is the fully normalized Legendre associated function.

legendre_schmidt_quasi_normalized!(P::AbstractMatrix, ϕ::Number, ph_term::Bool = false)

Compute the Schmidt quasi-normalized associated Legendre function P_n,m[cos(ϕ)] [3,4]. The maximum degree and order that will be computed are given by the dimensions of matrix P:

maximum degree -> number of rows
maximum order  -> number of columns

The result will be stored at matrix P.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Remarks

This algorithm was based on [3,4]. The conversion is obtained by:

             _           -
            |     (n-m)!  |    k = 1 if m  = 0
K_n,m = sqrt| k. -------- |,   k = 2 if m != 0
            |     (n+m)!  |
             -           -

=
P_n,m = P_n,m * K_n,m,

      =
where P_n,m is the quasi-normalized normalized Legendre associated function.

legendre_schmidt_quasi_normalized(ϕ::T, n_max::Integer, m_max::Integer = -1, ph_term::Bool = false) where T<:AbstractFloat

Compute the Schmidt quasi-normalized associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is n_max and the maximum order is m_max. Notice that if m_max is higher than n_max or negative, than it is set to n_max.

If ph_term is set to true, then the Condon-Shortley phase term (-1)ᵐ will be included. If ph_term is not present, then it defaults to false.

Returns

A matrix with the Legendre associated functions P_n,m[cos(ϕ)].

Remarks

This algorithm was based on [3,4]. The conversion is obtained by:

             _           -
            |     (n-m)!  |    k = 1 if m  = 0
K_n,m = sqrt| k. -------- |,   k = 2 if m != 0
            |     (n+m)!  |
             -           -

=
P_n,m = P_n,m * K_n,m,

      =
where P_n,m is the quasi-normalized normalized Legendre associated function.

list_ground_station_accesses(io, vargs...; kwargs...)

Print the ground station accesses to the io io. The arguments vargs... and keywords kwargs... are those of the function ground_station_accesses.

Additionally, the following keywords can be used to modify the behavior of this function:

• format: If :pretty, then a formatted table will be printed. If :csv, then the access data will be printed using the CSV format. (Default = :pretty)
• time_scale: Select the time scale of the access duration (:s for seconds, :m for minutes, and :h for hours). (Default = :m)

list_ground_station_gaps(io, vargs...; kwargs...)

Print the ground station gaps to the io io. The arguments vargs... and keywords kwargs... are those of the function ground_station_gaps.

Additionally, the following keywords can be used to modify the behavior of this function:

• format: If :pretty, then a formatted table will be printed. If :csv, then the access data will be printed using the CSV format. (Default = :pretty)
• time_scale: Select the time scale of the access duration (:s for seconds, :m for minutes, and :h for hours). (Default = :m)

list_ss_orbits_by_rep_period(minRep::Int, maxRep::Int, minAlt::Number=-1.0, maxAlt::Number=-1.0, e::Number=0.0)

Compute a list of repeating Sun-synchronous orbits.

Args

• minRep: Minimum repetition time of the orbit [days].
• maxRep: Maximum repetition time of the orbit [days].
• minAlt: Minimum altitude of the orbits on the list [m].
• maxAlt: Minimum altitude of the orbits on the list [m].
• e: Eccentricity.

Returns

A matrix containing the orbits found with the following format:

Semi-major axis [m] | Altitude [m] | Inclination [rad] | Period [s] | Int | Num | Den
--------------------|--------------|-------------------|------------|-----|-----|-----

in which the period is Int + Num/Den.

Remarks

If minAlt or maxAlt is < 0.0, then the altitude will not be checked when a orbit is added to the list.

load_gravity_model(T)

Load an embedded gravity model coefficients T and return an instance of the structure GravityModel_Coefs with the parsed values.

The current supported values for T are:

For other models, you can downlad the gfc file at

http://icgem.gfz-potsdam.de/home

and load it using the functions parse_icgem and create_gravity_model_coefs.

minimum_half_FOV_grss(h::Real, T::Real, i::Real, To::Integer)

Compute the minimum half FOV of a ground repeating Sun-synchronous (GRSS) orbit to cover the entire Equator within the revisit interval.

Args

• h: Orbit altitude in the Equator [m].
• T: Orbit period [s].
• i: Inclination [rad].
• To: Orbit cycle [days].

Returns

The minimum half FOV [rad].

minimum_half_FOV_grss(h::Real, a::Real, e::Real, i::Real, To::Integer)

Compute the minimum half FOV of a ground repeating Sun-synchronous (GRSS) orbit to cover the entire Equator within the revisit interval.

Args

• h: Orbit altitude in the Equator [m].
• a: Semi-major axis [m].
• e: Eccentricity.
• i: Inclination [rad].
• To: Orbit cycle [days].

Returns

The minimum half FOV [rad].

minimum_swath_grss(T::Real, i::Real, To::Integer)

Compute the minimum swath of a ground repeating Sun-synchronous (GRSS) orbit to cover the entire Equator within the revisit interval.

Args

• T: Orbit period [s].
• i: Inclination [rad].
• To: Orbit cycle [days].

Returns

The minimum swath [m].

minimum_swath_grss(a::Real, e::Real, i::Real, To::Integer)

Compute the minimum swath of a ground repeating Sun-synchronous (GRSS) orbit to cover the entire Equator within the revisit interval.

Args

• a: Semi-major axis [m].
• e: Eccentricity.
• i: Inclination [rad].
• To: Orbit cycle [days].

Returns

The minimum swath [m].

moon_position_i(JD_TDB::Number)

Compute the Moon position represented in the IAU-76/FK5 (mean-equator, mean-equinox), also called as J2000, at the Julian Day JD. The algorithm was adapted from [1, p. 288].

nrlmsise00(JD::Number, alt::Number, g_lat::Number, g_long::Number [, f107A::Number, f107::Number, ap::Union{Number,AbstractVector}]; output_si::Bool = true, dversion::Bool = true)

NRLMSISE-00

Neutral Atmosphere Empirical Model from the surface to lower exosphere.

This routine computes the NRLMSISE-00 outputs (see NRLMSISE00_Output) using the configurations in the structure nrlmsise00 (see NRLMSISE00_Structure).

Notice that the NRLMSISE-00 will be run using the default flags (see NRLMSISE00_DEFAULT_FLAGS). The user can only change the value of flags[:output_m_kg] using the keyword output_si to select whether the output must be converted to SI units. If more control is needed, then the user must manually call the function conf_nrlmsise00 and then call gtd7 or gtd7d with the desired flags.

If the space indices f107A, f107, and ap are missing, then they will be obtained from the online databases (see init_space_indices()).

Args

• JD: Julian Day [UTC].
• alt: Altitude [m].
• g_lat: Geodetic latitude [rad].
• g_long: Geodetic longitude [rad].
• f107A: 81 day average of F10.7 flux (centered on day of year JD).
• f107: Daily F10.7 flux for previous day.
• ap: Magnetic index (daily) if it is a number. If it is an array, then see Remarks.

Keywords

• output_si: (OPTIONAL) If true, then the output units will be [m⁻³] for species number density and [kg/m⁻³] for the total density. Otherwise, the units will be [cm⁻³] and [g/cm⁻³], respectively.
• dversion: (OPTIONAL) If true, run gtd7d. Otherwise, run gtd7 (see Remarks).

Returns

An instance of the structure NRLMSISE00_Output. The result in variable den_Total depends on the value of dversion (see Remarks, Notes on input variables).

Remarks

1. The densities of O, H, and N are set to 0 below 72.5 km.
2. The exospheric temperature T_exo is set to global average for altitudes below 120 km. The 120 km gradient is left at global average value for altitudes below 72.5 km.
3. Anomalous oxygen is defined as hot atomic oxygen or ionized oxygen that can become appreciable at high altitudes (> 500 km) for some ranges of inputs, thereby affection drag on satellites and debris. We group these species under the term Anomalous Oxygen, since their individual variations are not presently separable with the drag data used to define this model component.

AP

If ap is a Vector, then it must be a vector with 7 dimensions as described below:

Notes on input variables

f107 and f107A values used to generate the model correspond to the 10.7 cm radio flux at the actual distance of the Earth from the Sun rather than the radio flux at 1 AU. The following site provides both classes of values:

ftp://ftp.ngdc.noaa.gov/STP/SOLAR_DATA/SOLAR_RADIO/FLUX/

f107, f107A, and ap effects are neither large nor well established below 80 km and these parameters should be set to 150, 150, and 4 respectively.

If dversion is true, then the total mass den_Total (see NRLMSISE00_Output) is the sum of the mass densities of the species He, O, N₂, O₂, Ar, H, and N, but does not include anomalous oxygen.

If dversion is false, then total mass den_Total (see NRLMSISE00_Output) is the effective total mass density for drag and is the sum of the mass densities of all species in this model including the anomalous oxygen.

nutation_fk5(JD_TT::Number, n_max::Number = 106, nut_coefs_1980::Matrix = nut_coefs_1980)

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

an1  an2  an3  an4  an5  Ai  Bi  Ci  Di

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

Returns

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

parse_icgem(filename::AbstractString)

Parse the ICGEM file filename and return an instance of the structure ICGEM with the parsed data.

period(a::Number, e::Number, i::Number, pert::Symbol = :J2)
period(orb::Orbit, pert::Symbol = :J2)

Compute the period [s] of an object in an orbit with semi-major axis a [m], eccentricity e, and inclination i [rad], using the perturbation terms specified by the symbol pert. The orbit can also be specified by orb, which is an instance of the structure Orbit.

pert can be:

• :J0: Consider a Keplerian orbit.
• :J2: Consider the perturbation terms up to J2.
• :J4: Consider the perturbation terms J2, J4, and J2².

If pert is omitted, then it defaults to :J2.

precession_fk5(JD_TT::Number)

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

Returns

The angles (ζ, Θ, z) as described in [1, p. 226-228].

precession_nutation_iau2006(JD_TT::Number)

Compute the coordinates X, Y, and s related to the Celestial Intermediate Pole (CIP) with respect to the Geocentric Celestial Reference Frame (GCRF). This accounts for the effects of both precession and nutation of the CIP.

Returns

• The coordinate X of the CIP w.r.t. the GCRF;
• The coordinate Y of the CIP w.r.t. the GCRF;
• The CIO locator s that provides the position of the CIO on the Equator of the CIP corresponding to the kinematical definition of the non-rotation origin in the GCRS when the CIP is moving with respect to the GCRS between the reference epoch and the epoch due to precession and nutation [1, p. 214].

propagate!(orbp, t::Number) where T
propagate!(orbp, t::AbstractVector) where T

If t is a number, then propagate orbp by t [s] from the orbit epoch. Otherwise, if t is an array, then propagate the orbit in orbp using the time instants defined in the vector t [s].

In both cases, the orbit propagator algorithm is the one related to the structure orbp.

The structure orbp will contain the elements at the last propagation instant.

Returns

• The mean Keplerian elements represented in inertial frame in each time instant (see Orbit) [SI units].
• The position vector represented in inertial frame in each time instant [m].
• The velocity vector represented in inertial frame in each time instant [m].

If t is an array, then those values will be an array containing the information related to each epoch in t.

Remarks

The inertial frame in which the output is represented depends on which frame it was used to generate the orbit parameters. If the orbit parameters are obtained from a TLE, then the inertial frame will be TEME. Notice, however, that the perturbation theory requires an inertial frame with true equator.

propagate_to_epoch!(orbp, JD::Number) where T
propagate_to_epoch!(orbp, JD::AbstractVector) where T

If t is a number, then propagate orbp until the epoch JD [Julian Day]. Otherwise, if JD is an array, then propagate the orbit in orbp using the epochs defined in the vector t [Julian Day].

In both cases, the orbit propagator algorithm is the one related to the structure orbp.

The structure orbp will contain the elements at the last propagation instant.

Returns

• The mean Keplerian elements represented in inertial frame in each time instant (see Orbit) [SI units].
• The position vector represented in inertial frame in each time instant [m].
• The velocity vector represented in inertial frame in each time instant [m].

If JD is an array, then those values will be an array containing the information related to each epoch in JD.

Remarks

The inertial frame in which the output is represented depends on which frame it was used to generate the orbit parameters. If the orbit parameters are obtained from a TLE, then the inertial frame will be TEME. Notice, however, that the perturbation theory requires an inertial frame with true equator.

rCIRStoGCRF_iau2006([T::Type,] JD_TT::Number, dX::Number = 0, dY::Number = 0)

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

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

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

Returns

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

rCIRStoTIRS_iau2006([T::Type,] JD_UT1::Number)

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

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

Returns

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

Remarks

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

rECEFtoECEF([T,] ECEFo, ECEFf, JD_UTC::Number, eop_data)

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

Rotation description

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

The possible values are:

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

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

Conversion model

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

ECEF Frame

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

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

EOP Data

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

Returns

The rotation description represented by T that rotates the ECEF reference frame into alignment with the ECI reference frame.

Examples

julia> eop_IAU1980 = get_iers_eop(:IAU1980);

julia> rECEFtoECEF(PEF(), ITRF(), DatetoJD(1986,6,19,21,35,0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
  1.0          0.0         4.35684e-7
  0.0          1.0         1.44762e-6
 -4.35684e-7  -1.44762e-6  1.0

julia> rECEFtoECEF(Quaternion, PEF(), ITRF(), DatetoJD(1986,6,19,21,35,0), eop_IAU1980)
Quaternion{Float64}:
  + 0.9999999999997147 - 7.236343481310813e-7.i + 2.1765518308012794e-7.j + 0.0.k

julia> eop_IAU2000A = get_iers_eop(:IAU2000A);

julia> rECEFtoECEF(TIRS(), ITRF(), DatetoJD(1986,6,19,21,35,0), eop_IAU2000A)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
  1.0          3.08408e-11  -4.3531e-7
 -3.14708e-11  1.0          -1.44727e-6
  4.3531e-7    1.44727e-6    1.0

julia> rECEFtoECEF(Quaternion, TIRS(), ITRF(), DatetoJD(1986,6,19,21,35,0), eop_IAU2000A)
Quaternion{Float64}:
  + 0.9999999999997146 - 7.236343481345639e-7.i + 2.176551830689726e-7.j + 1.5577911634233308e-11.k

rECEFtoECI([T,] ECEF, ECI, JD_UTC::Number [, eop_data])

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

Rotation description

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

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

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

Conversion model

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

ECEF Frame

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

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

ECI Frame

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

• TEME(): ECI will be selected as the True Equator Mean Equinox (TEME) reference frame.
• TOD(): ECI will be selected as the True of Date (TOD).
• MOD(): ECI will be selected as the Mean of Date (MOD).
• J2000(): ECI will be selected as the J2000 reference frame.
• GCRF(): ECI will be selected as the Geocentric Celestial Reference Frame (GCRF).
• CIRS(): ECEF will be selected as the Celestial Intermediate Reference System (CIRS).

EOP Data

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

| Model | ECEF | ECI | EOP Data | |:–––––––|:–––-|:––––|:–––––––-=| | IAU-76/FK5 | ITRF | GCRF | EOP IAU1980 | | IAU-76/FK5 | ITRF | J2000 | EOP IAU1980 | | IAU-76/FK5 | ITRF | MOD | EOP IAU1980 | | IAU-76/FK5 | ITRF | TOD | EOP IAU1980 | | IAU-76/FK5 | ITRF | TEME | EOP IAU1980 | | IAU-76/FK5 | PEF | GCRF | EOP IAU1980 | | IAU-76/FK5 | PEF | J2000 | Not required¹ | | IAU-76/FK5 | PEF | MOD | Not required¹ | | IAU-76/FK5 | PEF | TOD | Not required¹ | | IAU-76/FK5 | PEF | TEME | Not required¹ | | IAU-2006/2010 | ITRF | CIRS | EOP IAU2000A | | IAU-2006/2010 | ITRF | GCRF | EOP IAU2000A | | IAU-2006/2010 | TIRS | CIRS | Not required¹ | | IAU-2006/2010 | TIRS | GCRF | Not required¹ ² |

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

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

MOD and TOD

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

Returns

The rotation description represented by T that rotates the ECEF reference frame into alignment with the ECI reference frame.

Examples

julia> eop_IAU1980 = get_iers_eop(:IAU1980);

julia> rECEFtoECI(DCM, ITRF(), GCRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267      0.78518     -0.00132979
 -0.78518      -0.619267     3.33492e-5
 -0.000797313   0.00106478   0.999999

julia> rECEFtoECI(ITRF(), GCRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267      0.78518     -0.00132979
 -0.78518      -0.619267     3.33492e-5
 -0.000797313   0.00106478   0.999999

julia> rECEFtoECI(PEF(), J2000(), DatetoJD(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619271      0.785176    -0.00133066
 -0.785177     -0.619272     3.45854e-5
 -0.000796885   0.00106622   0.999999

julia> rECEFtoECI(PEF(), J2000(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267      0.78518     -0.00133066
 -0.78518      -0.619267     3.45854e-5
 -0.000796879   0.00106623   0.999999

julia> rECEFtoECI(Quaternion, ITRF(), GCRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
Quaternion{Float64}:
  + 0.4363098936462618 - 0.0005909969666939257.i + 0.00030510511316206974.j + 0.8997962182293519.k

julia> eop_IAU2000A = get_iers_eop(:IAU2000A);

julia> rECEFtoECI(ITRF(), GCRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU2000A)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267      0.78518     -0.00132979
 -0.78518      -0.619267     3.33502e-5
 -0.000797312   0.00106478   0.999999

julia> rECEFtoECI(TIRS(), GCRF(), DatetoJD(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619271      0.785176    -0.00133066
 -0.785177     -0.619272     3.45884e-5
 -0.000796885   0.00106623   0.999999

julia> rECEFtoECI(Quaternion, ITRF(), GCRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU2000A)
Quaternion{Float64}:
  + 0.4363098936309669 - 0.000590996988144556.i + 0.0003051056555230158.j + 0.8997962182365703.k

rECItoECEF([T,] ECI, ECEF, JD_UTC::Number [, eop_data])

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

Rotation description

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

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

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

Conversion model

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

ECI Frame

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

• TEME(): ECI will be selected as the True Equator Mean Equinox (TEME) reference frame.
• TOD(): ECI will be selected as the True of Date (TOD).
• MOD(): ECI will be selected as the Mean of Date (MOD).
• J2000(): ECI will be selected as the J2000 reference frame.
• GCRF(): ECI will be selected as the Geocentric Celestial Reference Frame (GCRF).
• CIRS(): ECEF will be selected as the Celestial Intermediate Reference System (CIRS).

ECEF Frame

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

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

EOP Data

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

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

²: In this case, the terms that account for the free-core nutation and time dependent effects of the Celestial Intermediate Pole (CIP) position with respect to the GCRF will not be available, reducing the precision. The conversion between the frames depends on EOP Data (see get_iers_eop and read_iers_eop). If IAU-76/FK5 model is used, then the type of eop_data must be EOPData_IAU1980. The following table shows the requirements for EOP data given the selected frames.

MOD and TOD

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

Returns

The rotation description represented by T that rotates the ECI reference frame into alignment with the ECEF reference frame.

Examples

julia> eop_IAU1980 = get_iers_eop(:IAU1980);

julia> rECItoECEF(DCM, GCRF(), ITRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267    -0.78518     -0.000797313
  0.78518     -0.619267     0.00106478
 -0.00132979   3.33492e-5   0.999999

julia> rECItoECEF(GCRF(), ITRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267    -0.78518     -0.000797313
  0.78518     -0.619267     0.00106478
 -0.00132979   3.33492e-5   0.999999

julia> rECItoECEF(J2000(), PEF(), DatetoJD(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619271    -0.785177    -0.000796885
  0.785176    -0.619272     0.00106622
 -0.00133066   3.45854e-5   0.999999

julia> rECItoECEF(J2000(), PEF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267    -0.78518     -0.000796879
  0.78518     -0.619267     0.00106623
 -0.00133066   3.45854e-5   0.999999

julia> rECItoECEF(Quaternion, GCRF(), ITRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU1980)
Quaternion{Float64}:
  + 0.4363098936462618 + 0.0005909969666939257.i - 0.00030510511316206974.j - 0.8997962182293519.k

julia> eop_IAU2000A = get_iers_eop(:IAU2000A);

julia> rECItoECEF(GCRF(), ITRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU2000A)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619267    -0.78518     -0.000797312
  0.78518     -0.619267     0.00106478
 -0.00132979   3.33502e-5   0.999999

julia> rECItoECEF(GCRF(), TIRS(), DatetoJD(1986, 06, 19, 21, 35, 0))
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 -0.619271    -0.785177    -0.000796885
  0.785176    -0.619272     0.00106623
 -0.00133066   3.45884e-5   0.999999

julia> rECItoECEF(Quaternion, GCRF(), ITRF(), DatetoJD(1986, 06, 19, 21, 35, 0), eop_IAU2000A)
Quaternion{Float64}:
  + 0.4363098936309669 + 0.000590996988144556.i - 0.0003051056555230158.j - 0.8997962182365703.k

rECEFtoECI([T,] ECIo, ECIf, JD_UTC::Number [, eop_data])
rECEFtoECI([T,] ECIo, JD_UTCo::Number, ECIf, JD_UTCf::Number [, eop_data])

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

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

Rotation description

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

The possible values are:

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

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

Conversion model

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

ECI Frame

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

• TEME(): ECI will be selected as the True Equator Mean Equinox (TEME) reference frame.
• TOD(): ECI will be selected as the True of Date (TOD).
• MOD(): ECI will be selected as the Mean of Date (MOD).
• J2000(): ECI will be selected as the J2000 reference frame.
• GCRF(): ECI will be selected as the Geocentric Celestial Reference Frame (GCRF).
• CIRS(): ECEF will be selected as the Celestial Intermediate Reference System (CIRS).

EOP Data

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

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

MOD and TOD

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

Returns

The rotation description represented by T that rotates the origin ECI reference frame into alignment with the destination ECI reference frame.

Examples

julia> eop_IAU1980 = get_iers_eop(:IAU1980);

julia> rECItoECI(DCM, GCRF(), J2000(), DatetoJD(1986, 6, 19, 21, 35, 0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
  1.0          -2.45469e-12   4.56602e-10
  2.45466e-12   1.0          -1.84455e-9
 -4.56602e-10   1.84455e-9    1.0

julia> rECItoECI(Quaternion, TEME(), GCRF(), DatetoJD(1986, 6, 19, 21, 35, 0), eop_IAU1980)
Quaternion{Float64}:
  + 0.9999986335698654 + 1.8300414020900853e-5.i + 0.0006653038276169474.j - 0.0015132396749411375.k

julia> rECItoECI(TOD(), DatetoJD(1986,6,19,21,35,0), TOD(), DatetoJD(1987,5,19,3,0,0), eop_IAU1980)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 1.0          -0.000224087  -9.73784e-5
 0.000224086   1.0          -5.79859e-6
 9.73797e-5    5.77677e-6    1.0

julia> rECItoECI(Quaternion, TOD(), JD_J2000, MOD(), JD_J2000, eop_IAU1980)
Quaternion{Float64}:
  + 0.9999999993282687 - 1.400220690336851e-5.i + 1.3473593746216003e-5.j - 3.107834312843103e-5.k

julia> rECItoECI(J2000(), TEME(), DatetoJD(1986,6,19,21,35,0))
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
  0.999995    0.0030265    0.00133055
 -0.00302645  0.999995    -3.86125e-5
 -0.00133066  3.45854e-5   0.999999

julia> eop_IAU2000A = get_iers_eop(:IAU2000A);

julia> rECItoECI(CIRS(), GCRF(), DatetoJD(1986,6,19,21,35,0), eop_IAU2000A)
3×3 StaticArrays.SArray{Tuple{3,3},Float64,2,9}:
 0.999999     3.88379e-8  -0.00133066
 7.18735e-9   1.0          3.45882e-5
 0.00133066  -3.45882e-5   0.999999

julia> rECItoECI(Quaternion, CIRS(), GCRF(), DatetoJD(1986,6,19,21,35,0), eop_IAU2000A)
Quaternion{Float64}:
  + 0.9999997785177528 + 1.7294102099105917e-5.i + 0.0006653310148723835.j + 7.912627369563795e-9.k

rGCRFtoCIRS_iau2006([T::Type,] JD_TT::Number, dX::Number = 0, dY::Number = 0)

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

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

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

Returns

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

rGCRFtoITRF_fk5([T,] JD_UT1::Number, JD_TT::Number, x_p::Number, y_p::Number [, δΔϵ_1980::Number, δΔψ_1980::Number])

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

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

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

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

Returns

The rotation that aligns the GCRF frame with the ITRF frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rGCRFtoMOD_fk5([T,] JD_TT::Number)

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

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

Returns

The rotation that aligns the GCRF frame with the MOD frame. The rotation representation is selected by the optional parameter T.

Remarks

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

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

rGCRFtoTEME([T,] JD_TT::Number [, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0])

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

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

Returns

The rotation that aligns the GCRF frame with the TEME frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rITRFtoGCRF_fk5([T,] JD_UT1::Number, JD_TT::Number, x_p::Number, y_p::Number [, δΔϵ_1980::Number, δΔψ_1980::Number])

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

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

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

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

Returns

The rotation that aligns the ITRF frame with the GCRF frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rITRFtoPEF_fk5([T,] x_p::Number, y_p::Number)

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

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

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

Returns

The rotation that aligns the ITRF frame with the PEF frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rITRFtoTIRS_iau2006([T::Type,] JD_TT::Number, x_p::Number, y_p::Number)

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

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

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

Returns

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

Remarks

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

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

rMODtoGCRF_fk5([T,] JD_TT::Number)

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

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

Returns

The rotation that aligns the MOD frame with the GCRF frame. The rotation representation is selected by the optional parameter T.

Remarks

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

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

rMODtoPEF_fk5([T,] JD_UT1::Number, JD_TT::Number [, δΔϵ_1980::Number, δΔψ_1980::Number])

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

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

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

Returns

The rotation that aligns the MOD frame with the PEF frame. The rotation representation is selected by the optional parameter T.

rMODtoTEME([T,] JD_TT::Number [, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0])

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

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

Returns

The rotation that aligns the MOD frame with the TEME frame. The rotation representation is selected by the optional parameter T.

rMODtoTOD_fk5([T,] JD_TT::Number [, δΔϵ_1980::Number, δΔψ_1980::Number])

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

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

Returns

The rotation that aligns the MOD frame with the TOD frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rPEFtoITRF_fk5([T,] x_p::Number, y_p::Number)

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

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

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

Returns

The rotation that aligns the PEF frame with the ITRF. The rotation representation is selected by the optional parameter T.

Remarks

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

rPEFtoMOD_fk5([T,] JD_UT1::Number, JD_TT::Number [, δΔϵ_1980::Number, δΔψ_1980::Number])

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

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

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

Returns

The rotation that aligns the PEF frame with the TOD frame. The rotation representation is selected by the optional parameter T.

rPEFtoTEME([T,] JD_TT::Number)

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

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

Returns

The rotation that aligns the PEF frame with the TEME frame. The rotation representation is selected by the optional parameter T.

rPEFtoTOD_fk5([T,] JD_UT1::Number, JD_TT::Number [, δΔψ_1980::Number])

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

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

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

Returns

The rotation that aligns the PEF frame with the TOD frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rTEMEtoGCRF([T,] JD_TT::Number [, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0])

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

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

Returns

The rotation that aligns the TEME frame with the GCRF frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rTEMEtoMOD([T,] JD_TT::Number [, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0])

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

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

Returns

The rotation that aligns the TEME frame with the MOD frame. The rotation representation is selected by the optional parameter T.

rTEMEtoPEF([T,] JD_TT::Number)

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

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

Returns

The rotation that aligns the TEME frame with the PEF frame. The rotation representation is selected by the optional parameter T.

rTEMEtoTOD([T,] JD_TT::Number [, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0])

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

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

Returns

The rotation that aligns the TEME frame with the TOD frame. The rotation representation is selected by the optional parameter T.

rTIRStoCIRS_iau2006([T::Type,] JD_UT1::Number)

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

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

Returns

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

Remarks

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

rTIRStoITRF_iau2006([T::Type,] JD_TT::Number, x_p::Number, y_p::Number)

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

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

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

Returns

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

Remarks

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

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

rTODtoMOD_fk5([T,] JD_TT::Number [, δΔϵ_1980::Number, δΔψ_1980::Number])

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

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

Returns

The rotation that aligns the TOD frame with the MOD frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rTODtoPEF_fk5([T,] JD_UT1::Number, JD_TT::Number [, δΔψ_1980::Number])

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

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

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

Returns

The rotation that aligns the TOD frame with the PEF frame. The rotation representation is selected by the optional parameter T.

Remarks

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

rTODtoTEME([T,] JD_TT::Number [, δΔϵ_1980::Number = 0, δΔψ_1980::Number = 0])

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

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

Returns

The rotation that aligns the TOD frame with the TEME frame. The rotation representation is selected by the optional parameter T.

read_iers_eop(filename::String, data_type::Symbol = :IAU1980)

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

Returns

A structure (EOPData_IAU1980 or EOPData_IAU2000A, depending on data_type) with the interpolations of the EOP parameters. Notice that the interpolation indexing is set to the Julian Day.

Remarks

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

• IAU 1980

    curl -O https://datacenter.iers.org/data/latestVersion/223_EOP_C04_14.62-NOW.IAU1980223.txt
    wget https://datacenter.iers.org/data/latestVersion/223_EOP_C04_14.62-NOW.IAU1980223.txt

• IAU 2000A

    curl -O https://datacenter.iers.org/data/latestVersion/224_EOP_C04_14.62-NOW.IAU2000A224.txt
    wget https://datacenter.iers.org/data/latestVersion/224_EOP_C04_14.62-NOW.IAU2000A224.txt

rv_to_kepler(r_i::AbstractVector, v_i::AbstractVector, t::Number = 0)

Convert a Cartesian representation (position vector r [m] and velocity vector v [m/s²]) to the Keplerian elements. Optionally, the user can specify the epoch of the returned elements using the parameter t. It it is omitted, then it default to 0.

Returns

An instance of the structure Orbit with the Keplerian elements [SI units].

Remarks

The special cases are treated as follows:

• Circular and equatorial: the right ascension of the ascending node and the argument of perigee are set to 0. Hence, the true anomaly is equal to the true longitude.
• Elliptical and equatorial: the right ascension of the ascending node is set to 0. Hence, the argument of perigee is equal to the longitude of periapsis.
• Circular and inclined: the argument of perigee is set to 0. Hence, the true anomaly is equal to the argument of latitude.

References

The algorithm was adapted from [1].

rv_to_kepler(x::Number, y::Number, z::Number, vx::Number, vy::Number, vz::Number)

Convert a Cartesian representation (position vector [x;y;z] [m] and velocity vector [vx;vy;vz] [m/s²]) to the Keplerian elements.

Returns

An instance of the structure Orbit with the Keplerian elements [SI units].

rv_to_mean_elements_sgp4(vJD::AbstractVector{T}, vr::AbstractVector{Tv}, vv::AbstractVector{Tv}, W = I; max_it = 10000, sgp4_gc = sgp4_gc_wgs84) where {T,Tv<:AbstractVector}

Compute the mean elements for SGP4 based on the position vr and velocity vectors vr represented in TEME reference frame. The epoch of those measurements [Julian Day] must be in vJD.

The matrix W defined the weights for the least-square algorithm.

The variable max_it defines the maximum allowed number of iterations.

The variable sgp4_gc defines which constants should be used when running SGP4.

Returns

• The epoch of the elements [Julian Day].
• The mean elements for SGP4 algorithm:
  • Mean motion [rad/s];
  • Eccentricity [];
  • Inclination [rad];
  • Right ascension of the ascending node [rad];
  • Argument of perigee [rad];
  • Mean anomaly [rad];
  • BSTAR.
• The covariance matrix of the mean elements estimation.

rv_to_tle(args...; name::String = "UNDEFINED", sat_num::Int = 9999, classification::Char = 'U', int_designator = "999999", elem_set_number::Int = 0, rev_num, kwargs...)

Convert a set of position and velocity vectors represented in TEME reference frame to a TLE. The arguments args and keywords kwargs are the same as those described in the function rv_to_mean_elements_sgp4.

Additionally, the user can specify some parameters of the generated TLE.

This function prints the TLE to stdout using the function print_tle and also returns the TLE string.

satellite_check_Brazil(lat::Number, lon::Number)

Verify if a point described by latitude lat [rad] and longitude lon [rad] is inside Brazil. Returns true if the point is inside Brazil, of false otherwise.

Remarks

This function was based on the algorithm sent by Renato Branco to Ronan Arraes by e-mail at 2016-02-16.

satellite_lighting_condition(r_i::AbstractVector, s_i::AbstractVector)

Compute the satellite lighting condition given the Sun unitary vector s_i [m] and the satellite position vector r_i [m].

Returns

• SAT_LIGHTING_SUNLIGHT: Satellite is under sunlight.
• SAT_LIGHTING_PENUMBRA: Satellite is at penumbra region.
• SAT_LIGHTING_UMBRA: Satellite is at umbra region.

satellite_position_i(a::Number, e::Number, i::Number, RAAN::Number, w::Number, f::Number)

Compute the satellite position in the Earth-Centered Inertial (ECI) reference frame given the orbital elements a, e, i, RAAN, w, and f.

Notice that the ECI frame used will be the same as the frame of the orbital elements.

Args

• a: Semi-major axis.
• e: Eccentricity.
• i: Inclination [rad].
• RAAN: Right ascension of the ascending node [rad].
• w: Argument of perigee [rad].
• f: True anomaly [rad].

Returns

• The satellite position vector represented in the ECI reference frame.
• The unit vector perpendicular to the satellite position vector that lies on the orbit plane represented in the ECI reference frame.

Remarks

The satellite position vector will have the same unit of the semi-major axis.

satellite_sun_angle_earth_pointing(JD0::Number, a::Number, e::Number, i::Number, RAAN::Number, w::Number, numDays::Integer, fN_k::Function, meanAnomaly::Bool = false, step::Number = 0.1*pi/180.0)

Compute the Sun angle on a satellite surface for an Earth-pointing mission.

Args

• JD0: Initial instant for the analysis [Julian day].
• a: Semi-major axis of the orbit [m].
• e: Orbit eccentricity.
• i: Orbit inclination [rad].
• w: Argument of perigee [rad].
• RAAN: Right ascension of the ascending node at JD0 [rad].
• numDays: Number of days for the analysis.
• fN_k: Function f(s_b) that describes the solar panel normal at each k-th sampling step. Notice that s_b is the Sun vector represented in the body coordinate frame.
• meanAnomaly: (OPTIONAL) If true, compute using angular steps in the mean anomaly instead of in the orbit latitude (Default: false).
• step: (OPTIONAL) Mean anomaly step (Default: 0.1 deg).

Returns

A matrix containing the sun angle [rad] for each position in orbit for each day.

NOTE: if the Sun angle is larger than 90 deg or if the satellite is in eclipse, then NaN is returned in the matrix.

Remarks

The body reference frame is defined as:

• Z axis points towards the center of Earth;
• Y axis points towards the negative direction of orbit normal;
• X axis completes the right-hand reference frame.

If the mean anomaly is used, then the average value of the output is the average sun radiation received by the satellite surface, because every angular steps have a fixed time interval.

If the mean anomaly is used, then the angle interval is [0, 2π]. Otherwise, the angle interval is [-π,π].

satellite_sun_angle_earth_pointing(JD0::Number, a::Number, e::Number, i::Number, RAAN::Number, w::Number, numDays::Integer, N::AbstractVector, step::Number = 0.1*pi/180.0)

Compute the Sun angle on a satellite surface for an Earth-pointing mission.

Args

• JD0: Initial instant for the analysis [Julian day].
• a: Semi-major axis of the orbit [m].
• e: Orbit eccentricity.
• i: Orbit inclination [rad].
• w: Argument of perigee [rad].
• RAAN: Right ascension of the ascending node at JD0 [rad].
• numDays: Number of days for the analysis.
• N: Vector normal to the surface represented in the body reference frame.
• meanAnomaly: (OPTIONAL) If true, compute using angular steps in the mean anomaly instead of in the orbit latitude (Default: false).
• step: (OPTIONAL) Mean anomaly step (Default: 0.1 deg).

Returns

A matrix containing the Sun angle for each position in orbit for each day.

NOTE: if the Sun angle is larger than 90 deg or if the satellite is in eclipse, then NaN is returned in the matrix.

Remarks

The body reference frame is defined as:

• Z axis points towards the center of Earth;
• Y axis points towards the negative direction of orbit normal;
• X axis completes the right-hand reference frame.

If the mean anomaly is used, then the average value of the output is the average sun radiation received by the satellite surface, because every angular steps have a fixed time interval.

If the mean anomaly is used, then the angle interval is [0, 2π]. Otherwise, the angle interval is [-π,π].

satellite_sun_radiation_earth_pointing(JD0::Number, a::Number, e::Number, i::Number, RAAN::Number, w::Number, numDays::Integer, N::Vector, meanAnomaly::Bool = false, step::Number = 0.1*pi/180.0)

Compute the Sun radiation on a satellite surface for an Earth-pointing mission.

Args

• JD0: Initial instant for the analysis [Julian day].
• a: Semi-major axis of the orbit [m].
• e: Orbit eccentricity.
• i: Orbit inclination [rad].
• w: Argument of perigee [rad].
• RAAN: Right ascension of the ascending node at JD0 [rad].
• numDays: Number of days for the analysis.
• N: Vector normal to the surface represented in the body reference frame.
• meanAnomaly: (OPTIONAL) If true, compute using angular steps in the mean anomaly instead of in the orbit latitude (Default: false).
• step: (OPTIONAL) Mean anomaly step (Default: 0.1 deg).

Returns

A matrix containing the Sun radiation [W/m²] for each position in orbit for each day.

NOTE: if the Sun angle is larger than 90 deg or if the satellite is in eclipse, then NaN is returned in the matrix.

Remarks

The body reference frame is defined as:

• Z axis points towards the center of Earth;
• Y axis points towards the negative direction of orbit normal;
• X axis completes the right-hand reference frame.

If the mean anomaly is used, then the average value of the output is the average sun radiation received by the satellite surface, because every angular steps have a fixed time interval.

If the mean anomaly is used, then the angle interval is [0, 2π]. Otherwise, the angle interval is [-π,π].

satellite_sun_radiation_earth_pointing(JD0::Number, a::Number, e::Number, i::Number, RAAN::Number, w::Number, numDays::Integer, fN_k::Function, meanAnomaly::Bool = false, step::Number = 0.1*pi/180.0)

Compute the Sun radiation on a satellite surface for an Earth-pointing mission.

Args

• JD0: Initial instant for the analysis [Julian day].
• a: Semi-major axis of the orbit [m].
• e: Orbit eccentricity.
• i: Orbit inclination [rad].
• w: Argument of perigee [rad].
• RAAN: Right ascension of the ascending node at JD0 [rad].
• numDays: Number of days for the analysis.
• fN_k: Function f(s_b) that describes the solar panel normal at each k-th sampling step. Notice that s_b is the Sun vector represented in the body coordinate frame.
• meanAnomaly: (OPTIONAL) If true, compute using angular steps in the mean anomaly instead of in the orbit latitude (Default: false).
• step: (OPTIONAL) Mean anomaly step (Default: 0.1 deg).

Returns

A matrix containing the Sun radiation [W/m²] for each position in orbit for each day.

NOTE: if the Sun angle is larger than 90 deg or if the satellite is in eclipse, then NaN is returned in the matrix.

Remarks

The body reference frame is defined as:

• Z axis points towards the center of Earth;
• Y axis points towards the negative direction of orbit normal;
• X axis completes the right-hand reference frame.

If the mean anomaly is used, then the average value of the output is the average sun radiation received by the satellite surface, because every angular steps have a fixed time interval.

If the mean anomaly is used, then the angle interval is [0, 2π]. Otherwise, the angle interval is [-π,π].

satellite_sun_radiation_earth_pointing_mean(JD0::Number, a::Number, e::Number, i::Number, RAAN::Number, w::Number, numDays::Integer, N::AbstractVector, step::Number = 0.1*pi/180.0)

Compute the mean Sun radiation on a satellite surface for an Earth-pointing mission.

Args

• JD0: Initial instant for the analysis [Julian day].
• a: Semi-major axis of the orbit [m].
• e: Orbit eccentricity.
• i: Orbit inclination [rad].
• w: Argument of perigee [rad].
• RAAN: Right ascension of the ascending node at JD0 [rad].
• numDays: Number of days for the analysis.
• N: Vector normal to the surface represented in the body reference frame.
• meanAnomaly: (OPTIONAL) If true, compute using angular steps in the mean anomaly instead of in the orbit latitude (Default: false).
• step: (OPTIONAL) Mean anomaly step (Default: 0.1 deg).

Returns

The mean Sun radiation on a surface [W/m²].

Remarks

For more details, see satellitesunradiationearthpointing.

satellite_sun_radiation_earth_pointing_mean(JD0::Number, a::Number, e::Number, i::Number, RAAN::Number, w::Number, numDays::Integer, fN_k::Function, step::Number = 0.1*pi/180.0)

Compute the mean Sun radiation on a satellite surface for an Earth-pointing mission.

Args

• JD0: Initial instant for the analysis [Julian day].
• a: Semi-major axis of the orbit [m].
• e: Orbit eccentricity.
• i: Orbit inclination [rad].
• w: Argument of perigee [rad].
• RAAN: Right ascension of the ascending node at JD0 [rad].
• numDays: Number of days for the analysis.
• fN_k: Function f(s_b) that describes the solar panel normal at each k-th sampling step. Notice that s_b is the Sun vector represented in the body coordinate frame.
• meanAnomaly: (OPTIONAL) If true, compute using angular steps in the mean anomaly instead of in the orbit latitude (Default: false).
• step: (OPTIONAL) Mean anomaly step (Default: 0.1 deg).

Returns

The mean Sun radiation on a surface [W/m²].

Remarks

For more details, see satellitesunradiationearthpointing.

satsv(t::T1, r::AbstractVector{T2}, v::AbstractVector{T3} = [0,0,0], a::AbstractVector{T4} = [0,0,0]) where {T1<:Number, T2<:Number, T3<:Number, T4<:Number}
satsv(t::T1, vec::AbstractVector{T2}) where {T1<:Number, T2<:Number}

Create a new satellite state vector (see SatelliteStateVector) using the position r, velocity v, and acceleration a. It is also possible to pass a vector vec with the information concatenated.

sim_RAAN_J2(a::Number, e::Number, i::Number, RAAN_0::Number, numDays::Integer)

Simulate the RAAN of an orbit with semi-major axis a [m], eccentricity e, inclination i [rad], and initial RAAN RAAN_0 [rad] considering J2 perturbations. The analysis is performed for numDays days.

Returns

A numDays × 2 matrix in which the i-th line is:

| day | RAAN (0,2π) [rad] |

sort_list_ss_orbits_by_height(ss_orbits::Matrix)

Sort the list of Sun-synchronous orbits ss_orbits (see list_ss_orbits_by_rep_period) by height and return a new matrix.

step!(orbp, Δt::Number)

Propagate the orbit in orbp by Δt [s] using the algorithm of orbp. The new parameters will be written in orbp.

Returns

• The Keplerian elements represented in the inertial frame after the step (see Orbit) [SI units].
• The position vector represented in the inertial frame after the step [m].
• The velocity vector represented in the inertial frame after the step [m].

Remarks

The inertial frame in which the output is represented depends on which frame it was used to generate the orbit parameters. If the orbit parameters are obtained from a TLE, then the inertial frame will be TEME. Notice, however, that the perturbation theory requires an inertial frame with true equator.

sun_position_i(JD::Number)

Compute the Sun position represented in the Mean Equinox of Date (MOD) at the Julian Day JD. The algorithm was adapted from [3, p. 277-279].

sun_velocity_i(JD::Number)

Compute the Sun velocity represented in the Mean Equinox of Date (MOD) at the Julian Day JD. The algorithm was obtained by computing the time derivative of the Sun position in [3, p. 277-279].

svECEFtoECEF(sv::SatelliteStateVector, args...)

Convert the satellite state vector sv from an ECEF frame to another ECEF frame. The arguments args... must match those of the function rECEFtoECEF wihtout the rotation representation.

svECEFtoECI(sv::SatelliteStateVector, ECEF, ECI, JD_UTC [, eop_data])

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

svECItoECEF(sv::SatelliteStateVector, ECI, ECEF, JD_UTC [, eop_data])

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

svECItoECI(sv::SatelliteStateVector, args...)

Convert the satellite state vector sv from an ECI frame to another ECI frame. The arguments args... must match those of the function rECItoECI wihtout the rotation representation.

sv_to_kepler(sv::SatelliteStateVector)

Convert the state vector sv to Keplerian elements represented by an instance of the structure Orbit.

swath_width(h::real, HalfFOV::real)

Compute the swath width given the orbit altitude and the half FOV.

Args

• h: Orbit altitude [m].
• HalfFOV: Half field of view [rad].

Returns

The swath width [m].

twobody!(tbd::TwoBody_Structure, t::Number)

Propagate the orbit defined in tbd (see TwoBody_Structure) until the time t [s]. Notice that the values in tbd will be modified.

Returns

• The position vector represented in the inertial frame at time t [m].
• The velocity vector represented in the inertial frame at time t [m/s]

Remarks

The inertial frame in which the output is represented depends on which frame it was used to generate the orbit parameters. If the orbit parameters are obtained from a TLE, then the inertial frame will be TEME.

twobody_init(epoch::Number, a_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, f_0::Number, μ::T) where T

Initialize the data structure of Two Body orbit propagator algorithm.

Args

• epoch: Epoch of the orbital elements [s].
• a_0: Initial semi-major axis [m].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascension of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• f_0: Initial true anomaly.
• μ: Standard gravitational parameter of the central body [m^3/s^2].

Returns

The structure TwoBody_Structure with the initialized parameters.

@check_orbit(a, e)

Verify if the orbit with semi-major axis a [m] and eccentricity e is valid. This macro throws an exception if the orbit is not valid.

Return true is the orbit is valid, and false otherwise.

_expatmosphere_H

Scale height for the exponential atmospheric model [km].

_expatmosphere_h₀

Base altitude for the exponential atmospheric model [km].

_expatmosphere_ρ₀

Nominal density for the exponential atmospheric model [kg/m³].

_jr1971_constants

Constants for the Jacchia-Roberts 1971 Atmospheric Model.

_jr1971_id

Index of the species for the Jacchia-Roberts 1971 Atmospheric Model.

JR1971_CONSTANTS{T}

Structure with the constants for the Jacchia-Roberts 1971 Atmospheric Model.

SatelliteStateVector{T}

Store the state vector of the satellite.

Fields

• t: Epoch [Julian Day].
• r: Position vector [m].
• v: Velocity vector [m/s].
• a: Acceleration vector [m/s²].

_DTCFILE_Structure

Structure to store the interpolations of the data in DTCFILE.TXT file.

Fields

• DstΔTc: Temperature variation due to Dst [K].

_SOLFSMY_Structure

Structure to store the interpolations of the data in SOLFSMY.TXT file.

Fields

• F10: 10.7-cm solar flux [10⁻²² W/(m² Hz)].
• F81a: 10.7-cm averaged solar flux, 81-day centered on input time.
• S10: EUV index.
• S81a: EUV 81-day averaged centered index.
• M10: MG2 index scaled to F10.
• M81a: MG2 81-day averaged centered index.
• Y81a: Solar X-ray & Lya 81-day averaged centered index.
• Y81a: Solar X-ray & Lya 81-day averaged centered index.

_WDC_Structure

Structure to store the interpolations of the data in WDC files.

Fields

• Kp: Kp index.
• Ap: Ap index.

_ccor(alt::T, r::T, h1::T, zh::T) where T<:Number

Chemistry / Dissociation correction for MSIS models.

Args

• alt: Altitude.
• r: Target ratio.
• h1: Transition scale length.
• zh: Altitude of 1/2 r.

Returns

The chemistry / dissociation correction.

_ccor2(alt::T, r::T, h1::T, zh::T, h2::T) where T<:Number

Chemistry / Dissociation correction for MSIS models.

Args

• alt: Altitude.
• r: Target ration.
• h1: Transition scale length.
• zh: Altitude of 1/2 r.
• h2: Transition scale length 2.

Returns

The chemistry / dissociation correction.

_densm(re::T, gsurf::T, alt::T, d0::T, xm::T, tz::T, zn3::StaticVector{N3,T}, tn3::AbstractVector{T}, tgn3::AbstractVector{T}, zn2::StaticVector{N2,T}, tn2::AbstractVector{T}, tgn2::AbstractVector{T}) where {T<:Number,N2,N3}

Compute the temperature and density profiles for lower atmosphere.

Returns

• The density.
• The temperature.

_densu(re::T, gsurf::T, alt::T, dlb::T, tinf::T, tlb::T, xm::T, alpha::T, zlb::T, s2::T, zn1::StaticVector{N,T}, tn1::AbstractVector{T}, tgn1::AbstractVector{T}) where {T<:Number,N}

Compute the temperature and density profiles for MSIS models.

This algorithm uses new lower thermo polynomial.

Returns

• The density.
• The temperature.

_dnet(dd::T, dm::T, zhm::T, xmm::T, xm::T) where T<:Number

Turbopause correction for MSIS models.

Args

• dd: Diffusive density.
• dm: Full mixed density.
• zhm: Transition scale length.
• xmm: Full mixed molecular weight.
• xm: Species molecular weight.

Returns

The combined density.

_glob7s(p::AbstractVector{T}, nrlmsise00d::NRLMSISE00_Structure{T}) where T<:Number

Version of Globe for lower atmosphere (1999-10-26).

Args

• p: Vector with the coefficients.
• nrlmsise00d: NRLMSISE-00 structure (see NRLMSISE00_Structure).

Returns

The temperature (?).

_globe7!(p::AbstractVector{T}, nrlmsise00d::NRLMSISE00_Structure{T}) where T<:Number

Compute G(L) function.

Notice that the parameters apt and apdf of structure nrlmsise00d are modified.

Args

• p: Vector with the coefficients.
• nrlmsise00d: NRLMSISE-00 structure (see NRLMSISE00_Structure).

Returns

The temperature (?).

_init_dctfile(;force_download = false, local_path = nothing)

Initialize the data in the file DTCFILE.TXT by creating _dtcfile_data. The initialization process is composed of:

1. Download the file, if it is necessary;
2. Parse the file;
3. Create the interpolations and the structures.

If the keyword force_download is true, then the file will always be downloaded.

The user can also specify a location for the file using the keyword local_path. If it is nothing, which is the default, then the file will be downloaded.

_init_fluxtable(;force_download = false, local_path = nothing)

Initialize the data in the file fluxtable.txt by creating _fluxtable_data. The initialization process is composed of:

1. Download the file, if it is necessary;
2. Parse the file;
3. Create the interpolations and the structures.

If the keyword force_download is true, then the file will always be downloaded.

The user can also specify a location for the file using the keyword local_path. If it is nothing, which is the default, then the file will be downloaded.

_init_solfsmy(;force_download = false, local_path = nothing)

Initialize the data in the file SOLFSMY.TXT by creating _solfsmy_data. The initialization process is composed of:

1. Download the file, if it is necessary;
2. Parse the file;
3. Create the interpolations and the structures.

If the keyword force_download is true, then the file will always be downloaded.

The user can also specify a location for the file using the keyword local_path. If it is nothing, which is the default, then the file will be downloaded.

_init_wdcfiles(;force_download = false, local_dir = nothing, wdcfiles_oldest_year = year(now())-3)

Initialize the data in the WDC files by creating _wdcfiles_data. The initialization process is composed of:

1. Download the files, if it is necessary;
2. Parse the files;
3. Create the interpolations and the structures.

If the keyword force_download is true, then the files will always be downloaded.

The user can also specify a location for the directory with the WDC files using the keyword local_dir. If it is nothing, which is the default, then the file will be downloaded.

The user can select what is the oldest year in which the data will be downloaded by the keyword wdcfiles_oldest_year. By default, it will download the data from 3 previous years.

The user can select what is the newest year in which the data will be downloaded by the keyword wdcfiles_newest_year. It it is nothing, which is the default, then it is set to the current year.

_jb2008_M(z::R) where R

Compute the mean molecular mass at altitude z [km] using the empirical profile in eq. 1 [3].

_jb2008_T(z::R, Tx::R, T∞::R) where R<:Number

Compute the temperature [K] at height z [km] given the temperature Tx [K] at the inflection point, and the exospheric temperature T∞ [K] according to the theory of the model Jacchia 1971 [3].

The inflection point is considered to by z = 125 km.

_jb2008_grav(z::R) where R

Compute the gravity [m/s] at altitude z [km] according to the model Jacchia 1971 [3].

_jb2008_highaltitude(h::Number, F10ₐ::Number)

Compute the high altitude exospheric density correction factor in altitude h [km] and the averaged 10.7-cm solar flux (81-day centered on input time) [10⁻²² W/(M² Hz)].

This function uses the model in Section 6.2 of [2].

_jb2008_int(z₀::Number, z₁::Number, R::Number, Tx::Number, T∞::Number, δf::Function)

Compute the integral of the function δf between z₀ and z₁ using the Newton-Cotes 4th degree method. R is a number that defines the step size, Tx is the temperature at the inflection point, and T∞ is the exospheric temperature.

The signature of the function δf is:

δf(z, Tx, T∞)

and it must be _jb2008_δf1 or _jb2008_δf2.

This function returns a tuple containing the integral and last value of z used in the numerical algorithm.

_jb2008_semiannual(doy::Number, h::Number, F10ₐ::Number, S10ₐ::Number, M10ₐ::Number)

Compute the semiannual variation of the density considering the JB2008 model [1].

Args

• doy: Day of the year + fraction of the day.
• h: Height [km].
• F10ₐ: Averaged 10.7-cm flux (81-day centered on input-time) [10⁻²² W/(M² Hz)].
• S10ₐ: EUV 81-day averaged centered index.
• M10ₐ: MG2 81-day averaged centered index.

Returns

• Semiannual F(z) heigh function.
• Semiannual G(t) yearly periodic function.
• Semiannual variation of the density Δsalog₁₀ρ.

_jb2008_ΔTc(F10::Number, lst::Number, glat::Number, h::Number)

Compute the correction in the Tc for Jacchia-Bowman model.

This correction is mention in [2]. However, the equations do not seem to match those in the source-code. The ones implemented here are exactly the same as in the source-code.

Args

• F10: F10.7 flux.
• lst: Local solar time (0 - 24) [hr].
• glat: Geocentric latitude [rad].
• h: Altitude [km].

Returns

The correction ΔTc [K].

_jb2008_δf1(z, Tx, T∞)

Auxiliary function to compute the integrand in _jb2008_int.

_jb2008_δf2(z, Tx, T∞)

Auxiliary function to compute the integrand in _jb2008_int.

_jr1971_M(z::R) where R

Compute the mean molecular mass at altitude z [km] using the empirical profile in eq. 1 [3,4].

_jr1971_T(z::R, Tx::R, T∞::R) where R<:Number

Compute the temperature [K] at height z [km] given the temperature Tx [K] at the inflection point, and the exospheric temperature T∞ [K] according to the theory of the model Jacchia-Roberts 1971 [1,3,4].

The inflection point is considered to by z = 125 km.

_jr1971_roots(p::Polynomial{R}) where R

Compute the roots of the polynomial p necessary to compute the density below 125 km. It returns the value r₁, r₂, x, and y.

_parse_dtcfile(path::AbstractString)

Parse the DTCFILE.TXT file in path and return an instance of the structure _DTCFILE_Structure with the initialized interpolations.

The format of the file DTCFILE.TXT must be:

DTC YYYY DOY DTC_0h DTC_1h DTC_2h ... DTC_22h DTC_23h

in which DOY is the day of the year and DTC_Xh is the ΔTc at hour X.

_parse_fluxtable(path::AbstractString)

Parse the fluxtable.txt file in path and return an instance of the structure _fluxtable_Structure with the initialize interpolations.

_parse_solfsmy(path::AbstractString)

Parse the SOLFSMY.TXT file in path and retur an instance of the structure _SOLFSMY_Structure with the initialized interpolations.

The format of the file SOLFSMY.TXT must be:

YYYY DDD   JulianDay  F10   F81c  S10   S81c  M10   M81c  Y10   Y81c  Ssrc

_parse_wdcfiles(filepaths::Vector{String}, years::Vector{Int})

Parse the WDC files with paths in filepaths related to the years in years.

Notice that the files must be sorted by the year!

_prepare_wdc_remote_files(oldest_year::Number, newest_year::Number)

Configure all the WDC remote files between newest_year and oldest_year. Notice that previous years will never be updated whereas the current year will be updated daily.

If oldest_year is greater than current year, then only the files from the current year will be downloaded.

If newest_year is smaller than oldest_year, then only the files from the oldest_year will be downloaded.

This function modifies the global variable _wdcfiles.

_spline(x::StaticVector{N,T}, y::StaticVector{N,T}, yp1::T, ypn::T) where {T<:Number,N}

Compute the 2nd derivatives of cubic spline interpolation function tabulated by x and y given the 2nd derivatives values at x[1] (yp1) and at x[N] (ypn).

This function was adapted from Numerical Recipes.

Args

• x: X components of the tabulated function in ascending order.
• y: Y components of the tabulated function evaluated at x.
• yp1: 2nd derivative value at x[1].
• ypn: 2nd derivative value at x[N].

Returns

The 2nd derivative of cubic spline interpolation function evaluated at x.

Remarks

Values higher than 1e30 in the 2nd derivatives at the borders (yp1 and ypn) are interpreted as 0.

_splini(xa::StaticVector{N,T}, ya::StaticVector{N,T}, y2a::StaticVector{N,T}, x::T) where {T<:Number,N}

Compute the integral of the cubic spline function from xa[1] to x.

Args

• xa: X components of the tabulated function in ascending order.
• ya: Y components of the tabulated function evaluated at xa.
• y2a: Second derivatives.
• x: Abscissa endpoint for integration.

Returns

The integral of cubic spline function from xa[1] to x.

_splint(xa::StaticVector{N,T}, ya::StaticVector{N,T}, y2a::StaticVector{N,T}, x::T) where {T<:Number,N}

Compute the cubic spline interpolation value at x.

This function was adapted from Numerical Recipes.

Args

• xa: X components of the tabulated function in ascending order.
• ya: Y components of the tabulated function evaluated at xa.
• y2a: Second derivatives.
• x: Abscissa endpoint for interpolation.

Returns

The cubic spline interpolation value at x.

find_crossing(f::Function, t₀::Number, t₁::Number, s₀, s₁; Δ = 1e-3, max = 100)

Return the crossing time tc in which the function f(t) goes from the state s₀ to the state s₁. It is assumed that f(t₀) = s₀ and f(t₁) = s₁.

If the computed interval is smalled than Δ or if the number of iterations is higher than max, then the algorithm stops.

gts7(nrlmsise00d::NRLMSISE00_Structure{T}) where T<:Number

Thermospheric portion of NRLMSISE-00. This function should not be called to compute NRLMSISE-00. Use gtd7 or gtd7d instead.

Args

• nrlmsise00d: An instance of NRLMSISE00_Structure.

Returns

An instance of structure NRLMSISE00_Structure with the outputs.

@_keyword_found(keyword, keywords_found, current_line)

Check if the keyword exists in the list keywords_found. If true, then throw an error indicating that the problem occurred on the current_line.

@_parse_float(input)

Parse the input to float substituting all Ds and ds to e, so that we can convert numbers in FORTRAN format.