celestial.OrbitalEl - EranOfek/AstroPack GitHub Wiki

Description

celestial.OrbitalEl is a container class for orbital elements, including manipulations and ephemeris generation. The class can be used to store orbital elements or to upload orbital elements from the JPL website, and use them to generate ephemeris.

Related Functionality

Additional functionality related to orbital elements, and ephemeris ia available in the following packages, classes, and functions:

  • celestial.INPOP - A class that store and use the INPOP high-accuracy Solar System ephemeris.
  • celestial.SolarSys - General tools for calculating low-accuracy and medium-accuracy Solar System ephemeris.
  • celestial.Kepler - Tools for solving the Kepler equation.
  • celestial.SolarSys.jpl_horizons - A function for downloading ephemeris tables from the JPL horizons website.

Uploading JPL orbital elements.

For downloading the JPL orbital elements you can use the Installer class:

I = Installer;
I.install('MinorPlanets');
I.install('MinorPlanetsCT');

By default, the orbital elements files will be stored in the ~/matlab/data/SolarSys/MinorPlanets/ directory.

Properties

The class contains the following properties:

Each property contains a vector or array in which different rows refer to different objects.

  • Number - Minor planet number.
  • Designation - Minor planet designation.
  • Node - Longitude of the ascending node.
  • W - Argument of the periastron.
  • Incl - Inclination
  • A - Semi Major Axis.
  • PeriDist - Periastron distance.
  • Eccen - Eccentricity.
  • Tp - Time of periastron.
  • Epoch - Epoch of orbital elements.
  • Mepoch - Mean anomaly at epoch.
  • MagPar - Magnitude parameters. For example, for 'MagType' HG, these are the [H G] parameters.
  • MagType - Magnitude type (e.g., 'HG').
  • Equinox - Equinox.
  • AngUnits - Angular units (e.g., 'deg')
  • LenUnits - Length units (e.g., 'au')
  • TimeUnits - Time units (e.g., 'day')
  • K - Gaussian gravitational constant
  • Ref - JPL reference

Methods

In the following, we cover some of the available methods and show some examples of how to use them:

Loading the JPL orbital elements of minor planets and comets.

In order to load all the known minor planets and comets.

E=celestial.OrbitalEl.loadSolarSystem;

In this specific case, the object E will contain 3 elements:

  1. Orbital elements of numbered minor planets.
  2. Orbital elements of unnumbered minor planets.
  3. Orbital elements of comets.

However, in these lists the objects does not have a common epoch. For the function that uses orbital integration over multiple object simoltanously, a common epoch is needed. In order to load such a list, you can use:

E=celestial.OrbitalEl.loadSolarSystem('merge');

Utilities

% Return the number of orbital elements in each object element.
E.numEl

% select specific orbital elements by object index, or flag of logicals
% The following example will return the first object
Res   = selectFlag(E(1), 1, true);

% You can search and select objects using:
F = strcmp(E(1).Designation, 'Ceres');
Res   = selectFlag(E(1), F, true);

Additional functionality

  • merge - Merge the orbital elements in several elements of the OrbitalEl object.

Derived parameters

% Calculate the mean motion of all objects in OrbitalEl object:
E(1).meanMotion

% Calculate the orbital periods
E(1).period

% Calculate the semi-latus rectum
E(1).semiLatusRectum

Ephemeris

OrbElA = celestial.OrbitalEl.loadSolarSystem;
% Solve the Kepler equation to all minor planets for a single time (JD=2451545)
[Nu, R, E, Vel, M] = keplerSolve(OrbElA(1), 2451545)

% load a single asteroid
OrbEl = celestial.OrbitalEl.loadSolarSystem([],9804);
% Solve the Kepler equation for a single asteroid at multiple times:
[Nu, R, E, Vel, M] = keplerSolve(OrbEl(1), 2451545+(1:1:10)')

% Calculate rectangular ecliptic position from Nu and R:
[x,y,z] = trueAnom2rectPos(OrbElA(1), 1, 1)

% Rectangular position and velocity:
[V,X] = trueAnom2rectVel(OrbElA(1), 1)

To calculate the ephemeris for a Geocentric or Topocentric observer:

% load a single asteroid:
OrbEl = celestial.OrbitalEl.loadSolarSystem([],9804);
JD = celestial.time.julday([9 9 2021]);
% Generate ephemeris for a Geocentric observer, for a single asteroid at multiple times:
% using Kepler equation
Cat = celestial.ephem.ephemKepler(OrbEl, JD +(1:1:100)');

% load all numbered asteroids:
OrbEl = celestial.OrbitalEl.loadSolarSystem('num');
RAD = 180./pi;
% Calculate ephemeris for all asteroids at a given time for a topocentric observer:
% using Kepler equation:
CatE = celestial.ephem.ephemKepler(OrbEl, JD, 'GeoPos',[35./RAD 30./RAD 415],'MaxIterLT',2,'IncludeMag',false);

% There are several function for direct integration (much more accurate).
% For example, for multiple object at a single time:
OrbEl = celestial.OrbitalEl.loadSolarSystem('merge');
JD = 2460000;
Result = celestial.ephem.ephemKeplerMultiObj(OrbEl, JD)

Additional functionality:

  • nu2posOrbitalFrame - Convert true anomaly to position in the orbital frame.
  • magnitude - Calculate the magnitude for an OrbitalEl object

Search minor planets near sky position

The following example shows how to search for all known minor planets and comets within some sky coordinate at some time.

OrbEl= celestial.OrbitalEl.loadSolarSystem;               
% Search for minor planets (numbered, unnumbered, comets) within 1000 arcsec from RA=180, Dec=10 at JD=2451545                   
[Result, Names] = searchMinorPlanetsNearPosition(OrbEl, 2451545, 180, 10, 1000);
% The output Result in this case is a 3-element AstroCatalog object, with the ephemeris of minor planets in the search region

Basic parameters for ephemeris calculations

% Eccentric anomaly to radius vector
Result = eccAnom2radius(E(1), 1.1);

% True anomaly to radius vector
Result = trueAnom2radius(E(1), 1.1);

% Eccentric anomaly to True anaomaly
Result = eccAnom2trueAnom(E(1), 1.1);

% True anomaly to eccentric anomaly
Result = trueAnom2eccAnom(E(1), 1.1);

% calculate dNu/dt (time derivative of true anomaly) at nu
E(1).nuDot(1.1)

% calculate dr/dt (time derivative of radius vector) at nu
E(1).rDot(1.1)

% Radius vector to velocity
E(1).r2vel(1.1)

Conversions and utilities

To convert orbital elements to Thiele-Innes elements:

TI=thiele_innes(OrbEl(1));

To generate a table of orbital elements:

E=celestial.OrbitalEl.loadSolarSystem;
table(E)

To generate random orbital elements:

% see help for details:
R = celestial.OrbitalEl.randomElements;

Additional examples

Plot the apparent magnitude/positions of all numbered minor planets at a specific epoch

OrbElA = celestial.OrbitalEl.loadSolarSystem;
C=OrbElA(1).ephem(2451545);                  
hist(C.Catalog.Mag)
plot(C.Catalog.RA, C.Catalog.Dec,'.')