Library

Union of all Earth-Centered Earth-Fixed (ECEF) frames supported.

Union of all Earth-Centered Inertial (ECI) frames supported.

Union of all of date Earth-Centered Inertial (ECI) frames supported.

Union of all supported rotation descriptions.

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)]

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)]

Structure to store the information about a gravity model.

Structure to store the information contained in ICGEM files.

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.

Low level J2 orbit propagator structure.

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].

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].

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

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.

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

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].

function 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.

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

Fields

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

Structure that holds the information related to the SGP4 propagator.

Fields

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

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

Fields

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

Gravitational constants for SGP4.

Fields

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

Low level SGP4 structure.

This structure contains the same elements of the TLE with the same units.

Fields

• name: Name of the satellite.

First line

• sat_num: Satellite number.
• classification: Classification ('U', 'C', or 'S').
• int_designator: International designator.
• epoch_year: Epoch year (two digits).
• epoch_day: Epoch day (day + fraction of the day).
• epoch: The epoch represented in Julian Day.
• dn_o2: 1st time derivative of mean motion / 2 [rev/day²].
• ddn_o6: 2nd time derivative of mean motion / 6 [rev/day³].
• bstar: B* drag term.
• elem_set_number: Element set number.
• checksum_l1: Checksum of the line 1 (modulo 10).

Second line

• i: Inclination [deg].
• Ω: Right ascension of the ascending node [deg].
• e: Eccentricity.
• ω: Argument of perigee [deg].
• M: Mean anomaly [deg].
• n: Mean motion [rev/day].
• rev_num: Revolution number at epoch [rev].
• checksum_l2: Checksum of the line 2 (modulo 10).

Low level Two Body orbit propagator structure.

function DatetoJD(dateTime::DateTime)

Convert the date and time dateTime to Julian Day.

function DatetoJD(date::Date)

Convert the date date to Julian Day.

function 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).

function 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].

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

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

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

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

function 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].

function 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].

function 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.

function 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.

function 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.

function 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.

function 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.

function 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.

function 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.

function 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.

function 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].

function 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.

function 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.

function 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.

function 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.

function 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.

function 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.

function angvel(a::Number, e::Number, i::Number, pert::Symbol = :J2)
function 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.

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

function change_oe_frame(a::Number, e::Number, i::Number, Ω::Number, ω::Number, f::Number, conv_args...)
function 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 ˚

function 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.

function compute_U(gm_coefs::GravityModel_Coefs{T}, r::AbstractVector, n_max::Number = 0) 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. If n_max is less or equal 0, then the maximum available degree will be used. If n_max is omitted, then it defaults to 0.

function compute_dU(gm_coefs::GravityModel_Coefs{T}, r::AbstractVector, n_max::Number = 0) 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. If n_max is less or equal 0, then the maximum available degree will be used. If n_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.

function compute_g(gm_coefs::GravityModel_Coefs{T}, r::AbstractVector [, n_max::Number]) 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. If n_max is less or equal 0, then the maximum available degree will be used. If n_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.

function 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.

function 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].

function 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.

function 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].

function 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.

function 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).

function dArgPer(a::Number, e::Number, i::Number, pert::Symbol = :J2)
function 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.

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

function dRAAN(a::Number, e::Number, i::Number, pert::Symbol = :J2)
function 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.

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

function dlegendre([N,] ϕ::Number, n_max::Number, 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.

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(ϕ)].

function 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. Notice, however, that dP must be a square 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.

function dlegendre([N,] ϕ::Number, P::Matrix, ph_term::Bool)

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

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

This algorithm needs the matrix P with the associated Legendre function. This can be computed using the function legendre. The maximum order of the computed derivatives will be selected according to the dimensions of the matrix P. Notice that this matrix must be computed using the same normalization (see T) 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.

Returns

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

function dlegendre_conventional(ϕ::Number, n_max::Number, 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.

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(ϕ)].

function dlegendre_conventional(ϕ::Number, P::AbstractMatrix, ph_term = 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ϕ

This algorithm needs the matrix P with the conventional associated Legendre function. This can be computed using the function legendre_schmidt_quasi_normalized. The maximum order of the computed derivatives will be selected according to the dimensions of the 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.

Returns

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

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.

function 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. Notice, however, that dP must be a square 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.

function dlegendre_fully_normalized(ϕ::Number, n_max::Number, ph_term::Bool = false)

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.

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(ϕ)].

function dlegendre_fully_normalized(ϕ::Number, P::AbstractMatrix, ph_term = false)

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ϕ

This algorithm needs the matrix P with the Schmidt fully normalized associated Legendre function. This can be computed using the function legendre_fully_normalized. The maximum order of the computed derivatives will be selected according to the dimensions of the 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.

Returns

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

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.

function 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. Notice, however, that dP must be a square 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.

function dlegendre_schmidt_quasi_normalized(ϕ::Number, n_max::Number, 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 maximum degree that will be computed is 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(ϕ)].

function dlegendre_schmidt_quasi_normalized(ϕ::Number, P::AbstractMatrix, ph_term = 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ϕ

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. The maximum order of the computed derivatives will be selected according to the dimensions of the 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.

Returns

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

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.

function 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. Notice, however, that dP must be a square 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.

function dsinit(epoch::T, nll_0::T, all_0::T, e_0::T, i_0::T, Ω_0::T, ω_0::T, M_0::T, dotM::T, dotω::T, dotΩ::T) where T<:Number

Initialize the deep space structure. This function performs the initial computations and save the values at an instance of the structure SGP4_DeepSpace. Those will be used when calling the functions dsper! and dpsec!.

Args

• epoch: Epoch of the initial orbit [Julian Day].
• nll_0: Initial mean motion [rad/min].
• all_0: Initial semi-major axis [ER].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascencion of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• M_0: Initial mean motion [rad].
• dotM: Time-derivative of the mean motion [rad/min].
• dotω: Time-derivative of the argument of perigee [rad/min].
• dotΩ: Time-derivative of the RAAN [rad/min].

Returns

An instance of the structure SGP4_DeepSpace with the initalized values.

function dsper!(sgp4_ds::SGP4_DeepSpace{T}, e_k::T, i_k::T, Ω_k::T, ω_k::T, M_k::T, Δt:Number) where T<:Number

Compute the effects caused by Lunar-Solar periodics.

Notice that the values in the structure sgp4_ds will be modified.

Args

• sgp4_ds: Deep space structure (see SGP4_DeepSpace).
• e_k: Current eccentricity.
• i_k: Current inclination [rad].
• Ω_k: Current right ascension of the ascending node [rad].
• ω_k: Current argument of perigee [rad].
• M_k: Current mean anomaly [rad].
• Δt: Time interval since the epoch [min].

Returns

The following elements perturbed by lunar-solar periodics.

• Eccentricity.
• Inclination [rad].
• Right ascension of the ascending node [rad].
• Argument of perigee [rad].
• Mean anomaly [rad].

function dssec!(sgp4_ds::SGP4_DeepSpace{T}, nll_0::T, e_0::T, i_0::T, ω_0::T, Ω_k::T, ω_k::T, M_k::T, dotω::T, Δt::Number) where T<:Number

Compute the secular effects.

Notice that the values in the structure sgp4_ds will be modified.

Args

• sgp4_ds: Deep space structure (see SGP4_DeepSpace).
• nll_0: Initial mean motion [rad/min].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• ω_0: Initial argument of perigee [rad].
• Ω_k: Current right ascension of the ascending node [rad].
• ω_k: Current argument of perigee [rad].
• M_k: Current mean anomaly [rad].
• dotω: Time-derivative of the argument of perigee [rad/min].
• Δt: Time interval since the epoch [min].

Returns

The following elements perturbed by the secular effects:

• Mean motion [rad/min].
• Eccentricity.
• Inclination [rad].
• Right ascension of the ascending node [rad].
• Argument of perigee [rad].
• Mean anomaly [rad].

function 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].

function 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.

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

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

function f_to_E(orb::Orbit)

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

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

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

function f_to_M(orb::Orbit)

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

function 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().

function 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().

function 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.

function 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).

function 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).

function 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).

function 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.

function 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.

function 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.

function 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.

function 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.

function init_orbit_propagator(orbp_type::Type{Val{:J2}}, orb_0::Orbit, dn_o2::Number = 0, ddn_o6::Number = 0, j2_gc::J2_GravCte = j2_gc_wgs84)
function init_orbit_propagator(orbp_type::Type{Val{:sgp4}}, orb_0::Orbit, bstar::Number = 0, sgp4_gc::SGP4_GravCte = sgp4_gc_wgs84)
function init_orbit_propagator(orbp_type::Type{Val{:twobody}}, orb_0::Orbit, μ::Number = m0)

Initialize the orbit propagator orbp_type, which can be:

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

Args

• orb_0: Initial orbital elements (see Orbit).

J2 orbit propagator

• 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_wgs84).

SGP4 orbit propagator

• bstar: (OPTIONAL) B* parameter of the SGP4 (Default = 0).
• sgp4_gc: (OPTIONAL) Gravitational constants (Default = sgp4_gc_wgs84).

Two-body orbit propagator

• μ: (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

SGP4 Orbit Propagator

Notice that the orbit elements specified in orb_0 must be represented in TEME frame.

This 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.

function init_orbit_propagator(orbp_type::Type{Val{:J2}}, tle::TLE, j2_gc::J2_GravCte = j2_gc_wgs84)
function init_orbit_propagator(orbp_type::Type{Val{:sgp4}}, tle::TLE, sgp4_gc::SGP4_Structure = sgp4_gc_wgs84)
function init_orbit_propagator(orbp_type::Type{Val{:twobody}}, tle::TLE, μ::Number = m0)

Initialize the orbit propagator orbp_type, which can be:

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

Args

• tle: TLE that will be used to initialize the propagator.

J2 orbit propagator

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

SGP4 orbit propagator

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

Two-body orbit propagator

• μ: (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

SGP4 Orbit Propagator

This 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.

function init_orbit_propagator(orbp_type::Type{Val{:J2}}, epoch::Number, n_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, M_0::Number, dn_o2::Number = 0, ddn_o6::Number = 0, j2_gc::J2_GravCte{T} = j2_gc_wgs84) where T
function init_orbit_propagator(orbp_type::Type{Val{:sgp4}}, epoch::Number, n_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, M_0::Number, bstar::Number = 0, sgp4_gc::SGP4_GravCte{T} = sgp4_gc_wgs84) where T
function init_orbit_propagator(orbp_type::Type{Val{:twobody}}, epoch::Number, n_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, M_0::Number, μ::T = m0) where T

Initialize the orbit propagator orbp_type, which can be:

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

Args

• epoch: Initial orbit epoch [Julian Day].
• n_0: Initial angular velocity [rad/s].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascension of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• M_0: Initial mean anomaly [rad].

J2 orbit propagator

• 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_wgs84).

SGP4 orbit propagator

• bstar: (OPTIONAL) Initial B* parameter of the SGP4 (Default = 0).
• sgp4_gc: (OPTIONAL) Gravitational constants (Default = sgp4_gc_wgs84).

Two-body orbit propagator

• μ: (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

SGP4 Orbit Propagator

Notice that the orbit elements must be represented in TEME frame.

This 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.

function 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).

function 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].

function 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.

function j2_init(j2_gc::J2_GravCte{T}, epoch::Number, n_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, M_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].
• n_0: Mean motion at epoch [rad/s].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascension of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• M_0: Initial mean 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.

function jb2008(JD::Number, glat::Number, glon::Number, h::Number)
function 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.

function 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.

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

Convert the Keplerian elements (a, e, i, Ω, ω, and f) to a Cartesian representation (position vector r and velocity vector v).

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].

function legendre([N,] ϕ::Number, n_max::Number, ph_term::Bool = false)

Compute the associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is 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 square matrix with the Legendre associated functions P_n,m[cos(ϕ)].

function 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. Notice, however, that P must be a square matrix.

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.

function legendre_conventional(ϕ::Number, n_max::Number, ph_term::Bool = false)

Compute the conventional associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is 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 square matrix with the Legendre associated functions P_n,m[cos(ϕ)].

function 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. Notice, however, that P must be a square matrix.

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.

function legendre_fully_normalized(ϕ::Number, n_max::Number, ph_term::Bool = false)

Compute the fully normalized associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is 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 square 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.

function 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. Notice, however, that P must be a square matrix.

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.

function legendre_schmidt_quasi_normalized(ϕ::Number, n_max::Number, ph_term::Bool = false)

Compute the fully normalized associated Legendre function P_n,m[cos(ϕ)]. The maximum degree that will be computed is 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 square 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.

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. Notice, however, that P must be a square matrix.

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.

function 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] | 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.

function 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.

function 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].

function 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].

function 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].

function 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].

function 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.

function 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].

function parse_icgem(filename::AbstractString)

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

function period(a::Number, e::Number, i::Number, pert::Symbol = :J2)
function 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.

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

function 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].

function propagate!(orbp, t::Number) where T
function 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.

function propagate_to_epoch!(orbp, JD::Number) where T
function 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.

function 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, only IAU-76/FK5 is supported.

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. Currently, only the IAU-76/FK5 model is supported.

ECEF Frame

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

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

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.

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

function 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, only IAU-76/FK5 is supported.

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. Currently, only the IAU-76/FK5 model is supported.

ECEF Frame

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

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

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).

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. 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.

MOD and TOD

In this function, MOD and TOD frames are always defined with IERS EOP corrections. Hence, if one wants to obtain the MOD and TOD frames according to the original IAU-76/FK5 theory, it is necessary to use the low-level functions in file ./src/transformations/fk5/fk5.jl.

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

function 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, only IAU-76/FK5 is supported.

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. Currently, only the IAU-76/FK5 model is supported.

ECI Frame

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

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

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.

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. 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.

MOD and TOD

In this function, MOD and TOD frames are always defined with IERS EOP corrections. Hence, if one wants to obtain the MOD and TOD frames according to the original IAU-76/FK5 theory, it is necessary to use the low-level functions in file ./src/transformations/fk5/fk5.jl.

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

function rECEFtoECI([T,] ECIo, ECIf, JD_UTC::Number [, eop_data])
function 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, only IAU-76/FK5 is supported.

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. Currently, only the IAU-76/FK5 model is supported.

ECI Frame

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

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

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. The following table shows the requirements for EOP data given the selected frames.

MOD and TOD

In this function, MOD and TOD frames are always defined with IERS EOP corrections. Hence, if one wants to obtain the MOD and TOD frames according to the original IAU-76/FK5 theory, it is necessary to use the low-level functions in file ./src/transformations/fk5/fk5.jl.

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

function 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.

function 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.

function 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.

function 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.

function 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.

function 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.

function 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.

function 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.

function 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).

function 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.

function 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.

function 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.

function 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).

function 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.

function 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.

function 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.

function 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.

function 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).

function 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).

function 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.

function 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

function read_tle(tle_filename::String, verify_checksum::Bool = true)

Read the TLEs in the file tle_filename and return an array of TLE with the parsed TLEs.

If verify_checksum if true, then the checksum of both TLE lines will be verified. Otherwise, the checksum will not be checked. If verify_checksum is omitted, then it defaults to true.

function read_tle_from_string(tles::String, verify_checksum::Bool = true)
function read_tle_from_string(tle_l1::String, tle_l2::String, verify_checksum::Bool = false)

Parse a set of TLEs in the string tles or one TLE with first line tle_l1 and second line tle_l2. This function returns an array of TLE with the parsed TLEs.

If verify_checksum if true, then the checksum of both TLE lines will be verified. Otherwise, the checksum will not be checked. If verify_checksum is omitted, then it defaults to true.

function 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].

function rv_to_kepler(r::Vector, v::Vector)

Convert a Cartesian representation (position vector r [m] and velocity vector v [m/s²]) to the Keplerian elements.

Returns

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

References

The algorithm was adapted from [1].

function satellite_beta_angle(JD0::Number, a::Number, e::Number, i::Number, RAAN::Number, numDays::Integer)

Compute the beta angle of a satellite.

Args

• JD0: Initial instant for the analysis [Julian day].
• a: Semi-major axis of the orbit [m].
• e: Orbit eccentricity.
• i: Orbit inclination [rad].
• RAAN: Right ascension of the ascending node at JD0 [rad].
• numDays: Number of days of the analysis.

Returns

The beta angle [deg] computed for each day.

function 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.

function satellite_check_station(r_e::Vector, rs_e::Vector, minElev::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 minElev [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.

function satellite_check_station(r_e::Vector, lat_s::Number, lon_s::Number, h_s::Number, minElev::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 minElev [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.

function satellite_eclipse_time(JD0::Number, a::Number, e::Number, i::Number, w::Number, RAAN::Number, numDays::Integer, relative::Bool = false)

Compute the eclipse time of a satellite.

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 of the analysis.
• relative: Compute the eclipse time relative to the nodal period.

Returns

The following table:

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

function 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.

function 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.

function 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 [-π,π].

function 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 [-π,π].

function 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 [-π,π].

function 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 [-π,π].

function 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.

function 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.

function sgp4!(sgp4d::SGP4_Structure{T}, t::Number) where T

Propagate the orbit defined in sgp4d (see SGP4_Structure) until the time t [min]. Notice that the values in sgp4d will be modified.

Returns

• The position vector represented in TEME frame at time t [km].
• The velocity vector represented in TEME frame at time t [km/s].

function sgp4_init(spg4_gc::SGP4_GravCte{T}, epoch::Number, n_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, M_0::Number, bstar::Number) where T

Initialize the data structure of SGP4 algorithm.

Args

• spg4_gc: SPG4 gravitational constants (see SGP4_GravCte).
• epoch: Epoch of the orbital elements [Julian Day].
• n_0: SGP type "mean" mean motion at epoch [rad/min].
• e_0: "Mean" eccentricity at epoch.
• i_0: "Mean" inclination at epoch [rad].
• Ω_0: "Mean" longitude of the ascending node at epoch [rad].
• ω_0: "Mean" argument of perigee at epoch [rad].
• M_0: "Mean" mean anomaly at epoch [rad].
• bstar: Drag parameter (B*).

Returns

The structure SGP4_Structure with the initialized parameters.

function 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] |

sortlistssorbitsbyheight(ssorbits::Matrix)

Sort the list of Sun-synchronous orbits ss_orbits (see list_ss_orbits_by_rep_period) by height and return a new matrix.

function 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.

function 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].

function 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].

function 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].

function 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.

function twobody_init(epoch::Number, n_0::Number, e_0::Number, i_0::Number, Ω_0::Number, ω_0::Number, M_0::Number, μ::T) where T

Initialize the data structure of Two Body orbit propagator algorithm.

Args

• epoch: Epoch of the orbital elements [s].
• n_0: Mean motion at epoch [rad/s].
• e_0: Initial eccentricity.
• i_0: Initial inclination [rad].
• Ω_0: Initial right ascension of the ascending node [rad].
• ω_0: Initial argument of perigee [rad].
• M_0: Initial mean anomaly.
• μ: Standard gravitational parameter of the central body [m^3/s^2].

Returns

The structure TwoBody_Structure with the initialized parameters.

macro 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.

macro tle_str(str)

Parse a set of TLEs in the string str and return as an array of TLE. This version verifies the checksum of the TLE. If the checksum verification is not desired, see @tlenc_str.

Example

julia> tles = tle"""
       CBERS 4
       1 40336U 14079A   18166.15595376 -.00000014  00000-0  10174-4 0  9993
       2 40336  98.4141 237.7928 0001694  75.7582 284.3804 14.35485112184485
       SCD 1
       1 22490U 93009B   18165.62596833  .00000225  00000-0  11410-4 0  9991
       2 22490  24.9690 231.7852 0042844 200.7311 292.7198 14.44524498338066
       SCD 2
       1 25504U 98060A   18165.15074951  .00000201  00000-0  55356-5 0  9994
       2 25504  24.9961  80.1303 0017060 224.4822 286.6438 14.44043397 37312
       """

macro tlenc_str(str)

Parse a set of TLEs in the string str and return as an array of TLE. This version does not verify the checksum of the TLE. If the checksum verification is required, see @tle_str.

Example

julia> tles = tlenc"""
       CBERS 4
       1 40336U 14079A   18166.15595376 -.00000014  00000-0  10174-4 0  9993
       2 40336  98.4141 237.7928 0001694  75.7582 284.3804 14.35485112184485
       SCD 1
       1 22490U 93009B   18165.62596833  .00000225  00000-0  11410-4 0  9991
       2 22490  24.9690 231.7852 0042844 200.7311 292.7198 14.44524498338066
       SCD 2
       1 25504U 98060A   18165.15074951  .00000201  00000-0  55356-5 0  9994
       2 25504  24.9961  80.1303 0017060 224.4822 286.6438 14.44043397 37312
       """

Scale height for the exponential atmospheric model [km].

Base altitude for the exponential atmospheric model [km].

Nominal density for the exponential atmospheric model [kg/m³].

Constants for the Jacchia-Roberts 1971 Atmospheric Model.

Index of the species for the Jacchia-Roberts 1971 Atmospheric Model.

Structure with the constants for the Jacchia-Roberts 1971 Atmospheric Model.

Structure to store the interpolations of the data in DTCFILE.TXT file.

Fields

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

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.

Structure to store the interpolations of the data in WDC files.

Fields

• Kp: Kp index.
• Ap: Ap index.

function _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.

function _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.

function _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 where N2 where N3

Compute the temperature and density profiles for lower atmosphere.

Returns

• The density.
• The temperature.

function _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 where N

Compute the temperature and density profiles for MSIS models.

This algorithm uses new lower thermo polynomial.

Returns

• The density.
• The temperature.

function _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.

function _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 (?).

function _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 (?).

function _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.

function _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.

function _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.

function _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.

function _jb2008_M(z::R) where R

Compute the mean molecular mass at altitude z [km] using the empirical profile in eq. 1 [3].

function _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.

function _jb2008_grav(z::R) where R

Compute the gravity [m/s] at altitude z [km] according to the model Jacchia 1971 [3].

function _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].

function _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.

function _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₁₀ρ.

function _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].

function _jb2008_δf1(z, Tx, T∞)

Auxiliary function to compute the integrand in _jb2008_int.

function _jb2008_δf2(z, Tx, T∞)

Auxiliary function to compute the integrand in _jb2008_int.

function _jr1971_M(z::R) where R

Compute the mean molecular mass at altitude z [km] using the empirical profile in eq. 1 [3,4].

function _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.

function _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.

function _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.

function _parse_fluxtable(path::AbstractString)

Parse the fluxtable.txt file in path and return an instance of the structure _fluxtable_Structure with the initialize interpolations.

function _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

function _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!

function _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.

function _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.

function _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.

function _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.

function compute_checksum(str::AbstractString)

Compute the checksum of the line str modulo 10.

The algorithm is simple: add all the numbers in the line, ignoring letters, spaces, periods, and plus signs, but assigning +1 to the minus signs. The checksum is the remainder of the division by 10.

function 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.

function print_tle(io::IO, tle::TLE, color::Bool = true)

Print the TLE tle in the IO io.

If color is true, then the text will be printed using colors. If color is omitted, then it defaults to true.

macro _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.

macro _parse_float(input)

Parse the input to float substituting all Ds and ds to e, so that we can convert numbers in FORTRAN format.

macro parsevalue(T, str, linenum)

Parse the string str using the type T. If it is not succeeded, then throw an error indicating the line line_num with the problem.