Learn R Programming

asteRisk (version 1.0.0)

hpop: High-precision numerical orbital propagator

Description

Given the position and velocity of a satellite at a given time (in the GCRF system of coordinates), propagates its position by calculating its acceleration (based on a force model) and solving the resulting second-order ODE through numerical integration. This allows propagation of orbits with considerably higher accuracy than other propagators such as SGP4 and SDP4, but at the expense of a much higher computational cost. The forces and effects currently considered are gravitational attraction by the Earth (using a geopotential model based on spherical harmonics); effects of Earth ocean and solid tides; gravitational attraction by the Moon, Sun and planets (considered as point masses); solar radiation pressure; atmospheric drag, and relativistic effects. The force field is based on the forces described in Satellite Orbits: Models, Methods and Applications (Oliver Montenbruck and Eberhard Gill) and Fundamentals of Astrodynamics and Applications (David Vallado). Parts of this implementation are based on a previous MATLAB implementation by Meysam Mahooti. The NRLMSISE-00 model is used to calculate atmospheric density for the calculation of atmospheric drag. The high-precision numerical orbital propagator requires the asteRiskData package, which provides the data and coefficients required for calculation of the modeled forces. asteRiskData can be installed by running install.packages('asteRiskData', repos='https://rafael-ayala.github.io/drat/')

Usage

hpop(position, velocity, dateTime, times, satelliteMass, dragArea, 
     radiationArea, dragCoefficient, radiationCoefficient, ...)

Arguments

position

Initial position of the satellite in the GCRF system of coordinates. Should be provided as a numeric vector with 3 components that indicate the X, Y and Z components of the position in meters.

velocity

Initial velocity of the satellite in the GCRF system of coordinates. Should be provided as a numeric vector with 3 components that indicate the X, Y and Z components of the position in meters/second.

dateTime

Date time string in the YYYY-MM-DD HH:MM:SS format indicating the time corresponding to the initial position and velocity, in UTC time.

times

Vector with the times at which the position and velocity of the satellite should be calculated, in seconds since the initial time.

satelliteMass

Mass of the satellite in kilograms.

dragArea

Effective area of the satellite for atmospheric drag in squared meters.

radiationArea

Effective area of the satellite subject to the effect of radiation pressure in squared meters.

dragCoefficient

Drag coefficient (Cd) used for the calculation of atmospheric drag. For low Earth-orbiting satellites, a value of 2.2 is frequently employed if a better approximation is not available.

radiationCoefficient

Coefficient for the force resulting from radiation pressure. This parameter is usually referred to as reflectivity coefficient (Cr) and the value varies for different satellites and orbits. If unknown, a value of 1.2 is usually a decent approximation.

...

Additional parameters to be passed to ode to control how numerical integration is performed. By default, the RADAU5 solver is used.

Value

A matrix with the results of the numerical integration at the requested times. Each row contains the results for one of the requested times. The matrix contains seven columns: time (indicating the time for the corresponding row in seconds since the initial time), X, Y, Z (indicating the X, Y and Z components of the position for that time in meters), dX, dY and dZ (indicating the X, Y and Z components of the velocity for that time in meters/second). Positions and velocities are returned in the GCRF frame of reference.

References

Satellite Orbits: Models, Methods and Applications. Oliver Montenbruck and Eberhard Gill. Fundamentals of Astrodynamics and Applications. David Vallado. https://www.mathworks.com/matlabcentral/fileexchange/55167-high-precision-orbit-propagator https://ccmc.gsfc.nasa.gov/modelweb/models/nrlmsise00.php

Examples

Run this code
# NOT RUN {
if(requireNamespace("asteRiskData", quietly = TRUE)) {
# The following are the position and velocity in the GCRF frame of satellite
# MOLNIYA 1-83 the 25th of June, 2006 at 00:33:43 UTC.

initialPosition <-c(-14568679.5026116, -4366250.78287623, 9417.9289105405)
initialVelocity <- c(-3321.17428902497, -3205.49400830455, 4009.26862308806) 
initialTime <- "2006-06-25 00:33:43"

# Molniya satellites have a mass of approximately 1600 kg and a cross-section of
# 15 m2. Additionally, let<U+00B4>s use 2.2 and 1.2 as approximately values of the
# drag and reflectivity coefficients, respectively.

molniyaMass <- 1600
molniyaCrossSection <- 15
molniyaCr <- 1.2
molniyaCd <- 2.2

# Let<U+00B4>s calculate the position and velocity of the satellite for each minute of
# the following 10 minutes.

targetTimes <- seq(0, 600, by=60)
hpop_results <- hpop(initialPosition, initialVelocity, initialTime, targetTimes, 
                     molniyaMass, molniyaCrossSection, molniyaCrossSection,
                     molniyaCr, molniyaCd)
}
# }

Run the code above in your browser using DataLab