Module kiam

This Python module is a part of the KIAM Astrodynamics Toolbox developed in Keldysh Institute of Applied Mathematics (KIAM), Moscow, Russia.

The module serves as a safe and convenient interface to Fortran-compiled astrodynamical routines and provides instruments for performing translations between variables, coordinate systems, and time descriptions, propagating the trajectories in various models, and getting fast answers on typical questions about the two and n-body problems. It also contains some plotting routins and useful matrix linear algebra operations.

The toolbox is licensed under the MIT License.

The GitHub page of the project: https://github.com/shmaxg/KIAMToolbox.

Functions

def astro_const() ‑> tuple

Get astronomical constants.

Returns:

uni_const : dict

Universal constants containing the speed of light (SoL) in km/s, astronomical unit (AU) in km, constant of gravitation (G) in km^3/kg/s^2, standard acceleration due to gravity (g0) in m/s^2, degrees in 1 radian (RAD).

star : dict

Contains a dictionary that constants of the Sun: the gravitational parameter (GM) in km^3/s^2, the mean radius (MeanRadius) in km.

planet : dict

Contains constants of the planets (Mercury, Venus, Earth, Mars, Jupiter, Saturn, Uranus, Neptune). The keys of the dictionary are the names of the planets. Eack planet[planet_name] is also a dictionary that contains the gravitational parameter of the planet (GM) in km^3/s^2, the mean raidus (MeanRadius) in km, the equator radius (EquatorRadius) in km, the semi-major axis of the orbit around the Sun (SemimajorAxis) in km. For the Earth there are additionaly the obliquity of the ecliptic (Obliquity) in degrees and its time derivative (dObliquitydt) in arcsec/cy (cy = century years).

moon : dict

Contains constants of the moons (currently only of the Moon). The dictionary has a single key named Moon and moon['Moon'] is also a dictionary. That dictionary contains the gravitational parameter of the Moon (GM) in km^3/s^2, the mean raidus (MeanRadius) in km, the semi-major axis of the orbit around the Sun (SemimajorAxis) in km.

small_body : dict

Contains constants of the small celestial bodies (currently only of the Pluto). The dictionary has a single key named Pluto and small_body['Pluto'] is also a dictionary. That dictionary contains the gravitational parameter of the Pluto (GM) in km^3/s^2, the mean raidus (MeanRadius) in km, the equator radius (EquatorRadius) in km, the semi-major axis of the orbit around the Sun (SemimajorAxis) in km.

Examples:

uni_const, star, planet, moon, small_body = kiam.astro_const()  # If you need all the dicts

_, star, planet, _, _ = kiam.astro_const()  # If you need only star and planet dicts

print(star['Sun']['MeanRadius'])  # Mean radius of the Sun

print(planet['Earth']['GM'])  # Gravitational parameter of the Earth

print(planet['Mars']['SemimajorAxis'])  # Semi-major axis of the Mars's orbit.
def b1crs2b2crs(body1: str, body2: str, xb1crs: numpy.ndarray, jd: Union[float, numpy.ndarray], dist_unit: float, vel_unit: float) ‑> numpy.ndarray

Translate phase vectors from one CRS c/s to another CRS c/s.

Parameters:

body1 : str

The name of the first body.

Options: 'sun', 'mercury', 'venus', 'earth', 'moon', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune'.

body2 : str

The name of the second (target) body.

Options: 'sun', 'mercury', 'venus', 'earth', 'moon', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune'.

xb1crs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the CRS coordinate system of body1.

Vector structure: [x, y, z, vx, vy, vz]

jd : float, numpy.ndarray, shape (n,)

Julian dates corresponding to columns in xb1crs

dist_unit : float

The unit of distance in km

vel_unit : float

The unit of velocity in km/s

Returns:

xb2crs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the CRS coordinate system of body2.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

# Example 1 (6D -> 6D):

ku = kiam.units('sun', 'mars')

xb1crs = numpy.array([1, 0, 0, 0, 1, 0])  # wrt the Sun

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xb2crs = kiam.b1crs2b2crs('sun', 'mars', xb1crs, jd, ku['DistUnit'], ku['VelUnit'])  # wrt Mars

print(xb2crs)

# Example 2 (6x1 -> 6x1)

ku = kiam.units('sun', 'mars')

xb1crs = numpy.array([[1, 0, 0, 0, 1, 0]]).T  # wrt the Sun

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xb2crs = kiam.b1crs2b2crs('sun', 'mars', xb1crs, jd, ku['DistUnit'], ku['VelUnit'])  # wrt Mars

print(xb2crs)
def body_surface(body: str, radius: float = 1.0, quality: str = 'medium')

Return figure object for showing the surface of a celestial body (Earth, Moon).

Parameters:

body : str

The name of the celestial body.

Options: 'earth', 'moon'.

radius : float

The radius of the body.

Default: 1.0.

quality : str

The quality of the image.

Options: 'high', 'medium' (default), 'low'.

Returns:

fig : plotly.graph_objects.Figure

The Plotly figure object..

Examples:

fig = kiam.body_surface('earth')

fig.show()

fig = kiam.body_surface('moon')

fig.show()
def box_plot(*args)

Creates summary statistics with boxplots.

Parameters:

*args : Tuple[numpy.ndarray]

The 1D arrays. For each of them a boxplot is created.

Returns:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Examples:

y0 = numpy.random.randn(50) - 1

y1 = numpy.random.randn(50) + 1

fig = kiam.box_plot(y0, y1)

fig.show()
def cart2latlon(cart: numpy.ndarray) ‑> numpy.ndarray

Cartesian coordinates to latitude and longitude.

Parameters:

cart : numpy.ndarray, shape (3,), (3, n)

3D vector or array of column 3D vectors of Cartesian coordinates.

Vector structure: [x, y, z]

Returns:

latlon : numpy.ndarray, shape (2,), (2, n)

2D Vector or array of column 2D vectors of latitude and longitude pairs.

Vector structure: [lat, lon], where

lat in [-pi/2, pi/2],

lon in [-pi, pi].

Examples:

cart = numpy.array([1, 0, 0])

latlon = kiam.cart2latlon(cart)

print(latlon)

# [0. 0.]
def cart2sphere(cart: numpy.ndarray) ‑> numpy.ndarray

Cartesian coordinates to spherical coordinates.

Parameters:

cart : numpy.ndarray, shape (3,), (3, n)

3D vector or 3xn array of column 3D vectors of Cartesian coordinates

Vector structure: [x, y, z]

Returns:

sphere : numpy.ndarray, shape (3,), (3, n)

3D vector or 3xn array of column 3D vectors of spherical coordinates

Vector structure: [r, phi, theta], where


phi in [-pi, pi],

theta in [0, pi],

x = r*cos(theta)*cos(phi),

y = r*cos(theta)*sin(phi),

z = r*sin(theta).

Examples:

cart = numpy.array([1, 0, 0])

sphere = kiam.cart2sphere(cart)

print(sphere)

# [1.         0.         1.57079633]
def cosd(x: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

Cosine of a degree argument.

Parameters:

x : float, numpy.ndarray

Angle or an array of angles in degrees.

Returns:

s : float, numpy.ndarray

A cosine or array of cosines of angles in degrees.

Examples:

print(kiam.cosd(60))

# 0.5000000000000001
def cotand(x: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

Coangent of a degree argument.

Parameters:

x : float, numpy.ndarray

Angle or an array of angles in degrees.

Returns:

s : float, numpy.ndarray

A cotangent or array of cotangents of angles in degrees.

Examples:

print(kiam.cotand(45))

# 1.0000000000000002
def cr3bp_fb(t: float, s: numpy.ndarray, mu: float, stm_req: bool) ‑> numpy.ndarray

Right-hand side of the circular restricted three-body problem equations of motion wrt the first primary body.

Parameters:

t : float

Time

s : numpy.ndarray, shape (6,), (42,)

Phase state vector containing position and velocity and (if stm_req = True) vectorized state-transition matrix.

Vector structure:

[x, y, z, vx, vy, vz] if stm_req = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm_req = True

mu : float

Mass parameter of the three-body system

stm_req : bool

Flag to calculate the derivative of the state-transition matrix

Returns:

f : numpy.ndarray, shape (6,), (42,)

Gravitational acceleration according to the circular restricted three-body model of motion wrt the first primary body extended (if stm_req = True) by the derivative of the state-transition matrix.

Vector structure:

[fx, fy, fz, fvx, fvy, fvz] if stm_req = False

[fx, fy, fz, fvx, fvy, fvz, fm11, fm21, fm31, … ] if stm_req = True

Examples:

t0 = 0.0

s0 = numpy.array([0.5, 0, 0, 0, 1, 0])

mu = 1.2e-02

dsdt = kiam.cr3bp_fb(t0, s0, mu, False)

print(dsdt)

# [ 0.     1.     0.    -1.416 -0.    -0.   ]
def cr3bp_sb(t: float, s: numpy.ndarray, mu: float, stm_req: bool) ‑> numpy.ndarray

Right-hand side of the circular restricted three-body problem equations of motion wrt the secondary primary body.

Parameters:

t : float

Time

s : numpy.ndarray, shape (6,), (42,)

Phase state vector containing position and velocity and (if stm_req = True) vectorized state-transition matrix.

Vector structure:

[x, y, z, vx, vy, vz] if stm_req = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm_req = True

mu : float

Mass parameter of the three-body system

stm_req : bool

Flag to calculate the derivative of the state-transition matrix

Returns:

f : numpy.ndarray, shape (6,), (42,)

Gravitational acceleration according to the circular restricted three-body model of motion wrt the secondary primary body extended (if stm_req = True) by the derivative of the state-transition matrix.

Vector structure:

[fx, fy, fz, fvx, fvy, fvz] if stm_req = False

[fx, fy, fz, fvx, fvy, fvz, fm11, fm21, fm31, … ] if stm_req = True

Examples:

t0 = 0.0

s0 = numpy.array([0.5, 0, 0, 0, 1, 0])

mu = 1.2e-02

dsdt = kiam.cr3bp_sb(t0, s0, mu, False)

print(dsdt)

# [ 0.          1.          0.          3.00088889 -0.         -0.        ]
def deg2rad(deg: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

Degrees to radians conversion.

Parameters:

deg : float, numpy.ndarray

Angle or array of angles in degrees.

Returns:

rad : float, numpy.ndarray

Angle or array of angles in radians.

Examples:

print(kiam.deg2rad(180))

# 3.141592653589793
def dotainvb(a: numpy.ndarray, b: numpy.ndarray) ‑> numpy.ndarray

Calculate a*b^(-1) for matrices a and b.

Parameters:

a : numpy.ndarray, shape (n, n)

A square matrix.

b : numpe.ndarray, shape (n, n)

A square matrix.

Returns:

c : numpy.ndarray, shape (n, n)

The matrix that equals a*b^(-1)

Examples:

a = numpy.array([[1, 2], [3, 4]])

b = numpy.array([[1, 2], [3, 4]])

c = kiam.dotainvb(a, b)

print(c)

# [[1. 0.]
# [0. 1.]]
def ea2ta(ea: Union[float, numpy.ndarray], ecc: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

Eccentric anomaly to true anomaly.

Parameters:

ea : float, numpy.ndarray, shape (n,)

Scalar or array of eccentric anomalies.

ecc : float, numpy.ndarray, shape (n,)

Scalar or array of eccentricities. In case of array, the dimension should match the one of ea.

Returns:

ta : float, numpy.ndarray, shape (n,)

Scalar or array of true anomalies. Domain: (-pi, pi).

If ea and ecc are scalars, then ta is a scalar.

If ea is a scalar, ecc is a vector, then ta is a vector of the same size as ecc.

If ea is a vector, ecc is a scalar, then ta is a vector of the same size as ea.

If ea and ecc are vectors with the same size, then ta is a vector of the same size.

Examples:

ea = numpy.array([0.0, numpy.pi])

ecc = 0.1

ta = kiam.ea2ta(ea, ecc)
def ee2rv(ee: numpy.ndarray, mu: float, grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Equinoctial orbital elements to position and velocity.

Parameters:

ee : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column vectors of equinoctial orbital elements:

h = sqrt(p/mu),

ex = e*cos(Omega+omega),

ey = e*sin(Omega+omega),

ix = tan(i/2)*cos(Omega),

iy = tan(i/2)*sin(Omega),

L = theta + omega + Omega,

where

mu - gravitational parameter,

p - semi-latus rectum (focal parameter),

e - eccentricity,

Omega - right ascension of the ascending node,

omega - argument of pericenter,

i - inclination.

mu : float

Gravitational parameter.

grad_req : bool

Flag to calculate the derivatives of position and velocity wrt elemets.

Returns:

rv : numpy.ndarray, shape (6,), (6,n)

6D phase vector or array of 6D column vectors containing position and velocity.

Vector structure: [x, y, z, vx, dy, dz].

drv : numpy.ndarray, shape (6,6), (6,6,n)

6x6 matrix or 6x6xn array of partial derivatives of rv wrt ee (drv/dee).

Returns only if grad_req = True.

Examples:

ee = numpy.array([1, 0, 0, 0, 0, 0])

rv = kiam.ee2rv(ee, 1.0, False)

rv, drv = kiam.ee2rv(ee, 1.0, True)

print(rv)
def eye2vec(n: int) ‑> numpy.ndarray

Vector form of an identity matrix.

Parameters:

n : int

The number of rows and columns in the identity matrix.

Returns:

a : numpy.ndarray, shape (n**2,)

Vector form of the identity matrix.

Examples:

a = kiam.eye2vec(3)

print(a)

# [1. 0. 0. 0. 1. 0. 0. 0. 1.]
def gcrs2hcrs(xgcrs: numpy.ndarray, jd: Union[float, numpy.ndarray], dist_unit: float, vel_unit: float) ‑> numpy.ndarray

Translate phase vectors from GCRS c/s to HCRS c/s.

Parameters:

xgcrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the GCRS coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

jd : float, numpy.ndarray, shape (n,)

Julian dates corresponding to columns in xgcrs

dist_unit : float

The unit of distance in km

vel_unit : float

The unit of velocity in km/s

Returns:

xhcrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the HCRS coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

# Example 1 (6D -> 6D):

ku = kiam.units('sun', 'earth')

xgcrs = numpy.array([1, 0, 0, 0, 1, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xhcrs = kiam.gcrs2hcrs(xgcrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xhcrs)

# Example 2 (6x1 -> 6x1):

ku = kiam.units('sun', 'earth')

xgcrs = numpy.array([[1, 0, 0, 0, 1, 0]]).T

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xhcrs = kiam.gcrs2hcrs(xgcrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xhcrs)
def gcrs2itrs(xgcrs: numpy.ndarray, jd: Union[float, numpy.ndarray], grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Translate vector from GCRS c/s to ITRS c/s.

Parameters:

xgcrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the GCRS coordinate system

jd : float, numpy.ndarray, shape (n,)

Julian date(s) corresponding to column(s) in xgcrs

grad_req : bool

Flag to calculate the derivatives of the ITRS vector wrt the GCRS vector

Returns:

xitrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the ITRS coordinate system

dxitrs : numpy.ndarray, shape (3,3), (6,6), (3,3,n), (6,6,n)

3x3 or 6x6 matrix or 3x3xn or 6x6xn array of partial derivatives of xitrs wrt xgcrs (dxitrs/dxgcrs).

Returns only if grad_req = True.

Examples:

jd = kiam.juliandate(2022, 11, 22, 0, 0, 0)

xgcrs = numpy.array([1, 0, 0])

xitrs = kiam.gcrs2itrs(xgcrs, jd, False)

xitrs, dxitrs = kiam.gcrs2itrs(xgcrs, jd, True)

print(xitrs)
def gcrs2scrs(xgcrs: numpy.ndarray, jd: Union[float, numpy.ndarray], dist_unit: float, vel_unit: float) ‑> numpy.ndarray

Translate phase vectors from GCRS c/s to SCRS c/s.

Parameters:

xgcrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the GCRS coordinate system.

Vector structure: [x, y, z, vx, vy, vz]

jd : float, numpy.ndarray, shape (n,)

Julian dates corresponding to columns in xgcrs

dist_unit : float

The unit of distance in km

vel_unit : float

The unit of velocity in km/s

Returns:

xscrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the SCRS coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

# Example 1 (6D -> 6D):

ku = kiam.units('earth', 'moon')

xgcrs = numpy.array([1, 0, 0, 0, 1, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xscrs = kiam.gcrs2scrs(xgcrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xscrs)

# Example 2 (6x1 -> 6x1)

ku = kiam.units('earth', 'moon')

xgcrs = numpy.array([[1, 0, 0, 0, 1, 0]]).T

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xscrs = kiam.gcrs2scrs(xgcrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xscrs)
def get_altitude_km(period_hours: float, body: str) ‑> float

Calculate the altitude of a circular orbit with a given period around a specified celestial body.

Parameters:

period : float

The circular orbit period in hours.

body : str

The name of the celesial body.

Returns:

altitude_km : float

The altitude above the surface of the body in km.

Examples:

altitude_km = kiam.get_altitude_km(1.5, 'earth')

print(altitude_km)

# 281.5472668353086
def get_circular_velocity_km_s(altitude_km: float, body: str) ‑> float

Calculate the circular velocity at a given altitude around a specified celestial body.

Parameters:

altitude_km : float

The altitude above the surface of the body in km.

body : str

The name of the celesial body.

Returns:

velocity : float

The circular velocity at the given altitude.

Examples:

velocity_km_s = kiam.get_circular_velocity_km_s(200.0, 'earth')

print(velocity_km_s)

# 7.7884829462208724
def get_dv_hohmann(r1_nondim: float, r2_nondim: float) ‑> float

Calculate delta-v in a Hohmann transfer.

Parameters:

r1_nondim : float

Nondimensional distance to the center of mass of the central body at the start.

r2_nondim : float

Nondimensional distance to the center of mass of the central body at the end.

Returns:

dv : float

Nondimensional delta-v in the Hohmann transfer connecting r1_nondim and r2_nondim. It is assumed that the gravitational parameter equals 1.0.

Examples:

dv = kiam.get_dv_hohmann(1.0, 2.0)

print(dv)

# 0.024944026382329704
def get_order(altitude_thousands_km: float, approx_level: str = 'soft') ‑> int

The minimum order and degree of the complex lunar gravitational field at a given altitude according to the Trofimov–Shirobokov model.

Parameters:

altitude_thousands_km : float

The altitude above the lunar surface in km.

approx_level : str

The level of approximation, can be 'soft' or 'hard'.

Returns:

order : int

The order and degree of the complex lunar gravitational field.

Examples:

order = kiam.get_order(2.0, approx_level='soft')

print(order)

# 8.0
def get_period_hours(altitude_km: float, body: str) ‑> float

Calculate the circular orbit period with a given altitude around a specified celestial body.

Parameters:

altitude_km : float

The altitude above the surface of the body in km.

body : str

The name of the celesial body.

Returns:

period : float

The circular orbit period in hours.

Examples:

period_hours = kiam.get_period_hours(200.0, 'earth')

print(period_hours)

# 1.4725041134211172
def get_tof_hohmann(r1_nondim: float, r2_nondim: float) ‑> float

Calculate the time of flight in a Hohmann transfer.

Parameters:

r1_nondim : float

Nondimensional distance to the center of mass of the central body at the start.

r2_nondim : float

Nondimensional distance to the center of mass of the central body at the end.

Returns:

tof : float

Nondimensional time of flight in the Hohmann transfer connecting r1_nondim and r2_nondim. It is assumed that the gravitational parameter equals 1.0.

Examples:

tof = kiam.get_tof_hohmann(1.0, 2.0)

print(tof)

# 5.771474235728388
def grid_off(fig: plotly.graph_objs._figure.Figure)

Removes the grid.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Only Scatter and Scatter3d figure types are supported.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.grid_off(fig)

fig.show()

fig = kiam.grid_on(fig)

fig.show()
def grid_on(fig: plotly.graph_objs._figure.Figure)

Shows the grid.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Only Scatter and Scatter3d figure types are supported.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.grid_off(fig)

fig.show()

fig = kiam.grid_on(fig)

fig.show()
def hcrs2gcrs(xhcrs: numpy.ndarray, jd: Union[float, numpy.ndarray], dist_unit: float, vel_unit: float) ‑> numpy.ndarray

Translate phase vectors from HCRS c/s to GCRS c/s.

Parameters:

xhcrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the HCRS coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

jd : float, numpy.ndarray, shape (n,)

Julian dates corresponding to columns in xhcrs

dist_unit : float

The unit of distance in km

vel_unit : float

The unit of velocity in km/s

Returns:

xgcrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the GCRS coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

# Example 1 (6D -> 6D)

ku = kiam.units('sun', 'earth')

xhcrs = numpy.array([1, 0, 0, 0, 1, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xgcrs = kiam.hcrs2gcrs(xhcrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xgcrs)

# Example 2 (6x1 -> 6x1)

ku = kiam.units('sun', 'earth')

xhcrs = numpy.array([[1, 0, 0, 0, 1, 0]]).T

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xgcrs = kiam.hcrs2gcrs(xhcrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xgcrs)
def ine2rot(xine: numpy.ndarray, t: Union[float, numpy.ndarray], t0: Union[float, numpy.ndarray]) ‑> numpy.ndarray

Translate phase vectors from INE c/s to ROT c/s.

Parameters:

xine : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the INE coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

t : float, numpy.ndarray, shape (1,), (n,)

Time(s) corresponding to column(s) of xine

t0 : float, numpy.ndarray, shape (1,), (n,)

Time(s) of INE and ROT c/s coincidence for each column of xine.

Returns:

xrot : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the ROT coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

xine = numpy.array([1, 0, 0, 0, 1, 0])

t = 1.0

t0 = 0.0

xrot = kiam.ine2rot(xine, t, t0)

print(xrot)
def ine2rot_eph(xine: numpy.ndarray, jd: Union[float, numpy.ndarray], first_body: str, secondary_body: str, dist_unit: float, vel_unit: float) ‑> numpy.ndarray

Translate phase vectors from INEEPH c/s to ROTEPH c/s.

Parameters:

xine : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the INEEPH coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

jd : float, numpy.ndarray, shape (n,)

Julian date(s) corresponding to column(s) in xine

first_body : str

Name of the first primary body

Options: 'sun', 'mercury', 'venus', 'earth', 'moon', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune'

secondary_body : str

Name of the secondary primary body

Options: 'sun', 'mercury', 'venus', 'earth', 'moon', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune'

dist_unit : float

The unit of distance in km

vel_unit : float

The unit of velocity in km/s

Returns:

xrot : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the ROTEPH coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

xine = numpy.array([1, 0, 0, 0, 1, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

ku = kiam.units('earth', 'moon')

xrot = kiam.ine2rot_eph(xine, jd, 'earth', 'moon', ku['DistUnit'], ku['VelUnit'])

print(xrot)
def invadotb(a: numpy.ndarray, b: numpy.ndarray) ‑> numpy.ndarray

Calculate a^(-1)*b for matrices a and b.

Parameters:

a : numpy.ndarray, shape (n, n)

A square matrix.

b : numpe.ndarray, shape (n, n)

A square matrix.

Returns:

c : numpy.ndarray, shape (n, n)

The matrix that equals a^(-1)*b.

Examples:

a = numpy.array([[1, 2], [3, 4]])

b = numpy.array([[1, 2], [3, 4]])

c = kiam.invadotb(a, b)

print(c)

# [[1.00000000e+00 0.00000000e+00]
# [8.32667268e-17 1.00000000e+00]]
def is_visible(r_sat: numpy.ndarray, lat_deg: Union[int, float, numpy.ndarray], long_deg: Union[int, float, numpy.ndarray], body_radius: float, threshold_deg: Union[int, float, numpy.ndarray]) ‑> tuple

Get visibility statuses (0 or 1) of a vector from a point on a sphere surface.

Parameters:

r_sat : numpy.ndarray, shape (3,), (3,n)

Radius-vector(s) around a sphere surface

lat_deg : int, float, numpy.ndarray, shape (m,)

Latitude of a point(s) on a surface in degrees

long_deg : int, float, numpy.ndarray, shape (m,)

Longitude of a point(s) on a surface in degrees

body_radius : float

Body radius

threshold_deg : int, float, numpy.ndarray, shape (n,)

Minimum angle below which the vector is not visible.

Returns:

status : numpy.ndarray, shape (n,m)

Visibility statuses of the r_sat vectors from lat_deg/long_deg points.

n - number of vectors in r_sat

m - number of points on the surface

elev_deg : numpy.ndarray, shape (n,m)

Elevation angles in degrees

n - number of vectors in r_sat

m - number of points on the surface

azim_deg : numpy.ndarray, shape (n,m)

Azimuth angles in degrees

n - number of vectors in r_sat

m - number of points on the surface

Examples:

r_sat = numpy.array([2, 0, 0])

lat_deg = 0.0

long_deg = 0.0

body_radius = 1.0

threshold_deg = 5.0

status, elev_deg, azim_deg = kiam.is_visible(r_sat, lat_deg, long_deg, body_radius, threshold_deg)

print(status, elev_deg, azim_deg)
def itrs2gcrs(xitrs: numpy.ndarray, jd: Union[float, numpy.ndarray], grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Translate vector from ITRS c/s to GCRS c/s.

Parameters:

xitrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the ITRS coordinate system

jd : float, numpy.ndarray, shape (n,)

Julian date(s) corresponding to column(s) of xitrs

grad_req : bool

Flag to calculate the derivatives of the GCRS vector wrt the ITRS vector

Returns:

xgcrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the GCRS coordinate system

dxgcrs : numpy.ndarray, shape (3,3), (6,6), (3,3,n), (6,6,n)

3x3 or 6x6 matrix or 3x3xn or 6x6xn array of partial derivatives of xgcrs wrt xitrs (dxgcrs/dxitrs).

Returns only if grad_req = True.

Examples:

jd = kiam.juliandate(2022, 11, 22, 0, 0, 0)

xitrs = numpy.array([1, 0, 0])

xgcrs = kiam.itrs2gcrs(xitrs, jd, False)

xgcrs, dxgcrs = kiam.itrs2gcrs(xitrs, jd, True)

print(xgcrs)
def jd2time(jd: float) ‑> datetime.datetime

Julian date to usual date and time.

Parameters:

jd : float

Julian date

Returns:

time : datetime.datetime

Date and time object of type datetime.datetime

Examples:

print(kiam.jd2time(2459905.5))

# 2022-11-22 00:00:00
def juliandate(year: int, month: int, day: int, hour: int, minute: int, second: int) ‑> float

Usual date to Julian date.

Parameters:

year : int

Year

month : int

Month

day : int

Day

hour : int

Hour

minute : int

Minute

second : int

Second

Returns:

jd : float

Julian date

Examples:

jd = kiam.juliandate(2022, 11, 22, 0, 0, 0)

print(jd)

# 2459905.5
def kepler(mean_anomaly: float, ecc: float, atol: float = 1e-10, maxiter: int = 1000) ‑> float

Solve the Kepler equation.

Parameters:

mean_anomaly : float

The mean anomaly.

ecc : float

The eccentricity in range [0, 1).

atol : float

The absolute tolerance. Newton iterations will finish when absolute difference between two successive approximations to the solution of the Kepler equation will be lower than atol. Default if 1E-10.

maxiter : int

The maximum number of Newton iterations. Default is 1000.

Returns:

ea : float

Eccentric anolmaly, solution to the Kepler equation E - e*sin(E) = M.

Examples:

ea = kiam.kepler(2*numpy.pi, 0.1)
def latlon2cart(latlon: numpy.ndarray) ‑> numpy.ndarray

Latitude and longitude to Cartesian coordinates.

Parameters:

latlon : numpy.ndarray, shape (2,), (2, n)

2D Vector or array of column 2D vectors of latitude and longitude pairs.

Vector structure: [lat, lon], where

lat in [-pi/2, pi/2],

lon in [-pi, pi]

Returns:

cart : numpy.ndarray, shape (3,), (3, n)

3D vector or array of column 3D vectors of Cartesian coordinates.

Vector structure: [x, y, z].

Examples:

latlon = numpy.array([0, 0])

cart = kiam.latlon2cart(latlon)

print(cart)

# [1. 0. 0.]
def leap_second_count() ‑> int

Returns current leap seconds count.

Returns:

leap_second_count() : int

The current leap seconds count.

Examples:

lsc = kiam.leap_second_count()  # 37

print(lsc)
def legend_off(fig: plotly.graph_objs._figure.Figure)

Removes the legend.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.legend_off(fig)

fig.show()

fig = kiam.legend_on(fig)

fig.show()
def legend_on(fig: plotly.graph_objs._figure.Figure)

Shows the legend.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.legend_off(fig)

fig.show()

fig = kiam.legend_on(fig)

fig.show()
def load(filename: str) ‑> Any

Loads a variable from a specified file.

Parameters:

filename : str

A path to the file.

Returns:

var : Any

A variable contained in the file.

Examples:

a = kiam.load('variable_a')

print(a)
def lvlh2mer(xlvlh: numpy.ndarray, lat: float, lon: float) ‑> numpy.ndarray

Translate phase vector(s) from LVLH c/s to MER c/s.

Parameters:

xlvlh : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the LVLH coordinate system

lat : float

Latitude that specifies the LVLH c/s in radians

lon : float

Longitude that specifies the LVLH c/s in radians

Returns:

xmer : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the MER coordinate system

Examples:

xlvlh = numpy.array([1, 2, 3])

lat = 0.0

lon = 1.0

xmer = kiam.lvlh2mer(xlvlh, lat, lon)

print(xmer)
def mat2vec(a: numpy.ndarray) ‑> numpy.ndarray

Square matrix to vector form translation.

Parameters:

a : numpy.ndarray, shape (n, n)

A square matrix.

Returns:

v : numpy.ndarray, shape (n**2,)

Vector form of the matrix.

Vector structure (Fortran/MATLAB order): [a11, a21, a31, … ]

Examples:

a = numpy.array([[1, 2], [3, 4]])

v = kiam.mat2vec(a)

print(v)

# [1 3 2 4]
def mer2lvlh(xmer: numpy.ndarray, lat: float, lon: float) ‑> numpy.ndarray

Translate phase vector(s) from MER c/s to LVLH c/s.

Parameters:

xmer : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the MER coordinate system

lat : float

Latitude that specifies the LVLH c/s in radians

lon : float

Longitude that specifies the LVLH c/s in radians

Returns:

xlvlh : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the LVLH coordinate system

Examples:

xmer = numpy.array([1, 2, 3])

lat = 0.0

lon = 1.0

xlvlh = kiam.mer2lvlh(xmer, lat, lon)

print(xlvlh)
def mer2scrs(xmer: numpy.ndarray, jd: Union[float, numpy.ndarray], grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Translate vectors from MER c/s to SCRS c/s.

Parameters:

xmer : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the MER coordinate system

jd : float, numpy.ndarray, shape (n,)

Julian date(s) corresponding to vector or columns in xmer

grad_req : bool

Flag to calculate the derivatives of the SCRS vector wrt the MER vector

Returns:

xscrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the SCRS coordinate system

dxscrs : numpy.ndarray, shape (3,3), (6,6), (3,3,n), (6,6,n)

Array of matrices of partial derivatives of xscrs wrt xmer (dxscrs/dxmer).

Returns only if grad_req = True.

Examples:

xmer = numpy.array([1, 0, 0])

jd = kiam.juliandate(2022, 11, 22, 0, 0, 0)

xscrs = kiam.mer2scrs(xmer, jd, False)

dxscrs = kiam.mer2scrs(xmer, jd, True)

print(xscrs)
def nbp_ee_earth(t: float, s: numpy.ndarray, stm_req: bool, sources: dict, data: dict, units_data: dict) ‑> numpy.ndarray

Right-hand side of the n-body problem equations of motion wrt the Earth in terms of the equinoctial orbital elements.

Parameters:

t : float

Time

s : numpy.ndarray, shape (6,), (42,)

Phase state vector containing equinoctial elements and (if stm_req = True) vectorized state-transition matrix.

Vector structure:

[h, ex, ey, ix, iy, L] if stm_req = False,

[h, ex, ey, ix, iy, L, m11, m21, m31, …] if stm_req = True,

h = sqrt(p/mu),

ex = e*cos(Omega+omega),

ey = e*sin(Omega+omega),

ix = tan(i/2)*cos(Omega),

iy = tan(i/2)*sin(Omega),

L = theta + omega + Omega,

where

mu - gravitational parameter,

p - semi-latus rectum (focal parameter),

e - eccentricity,

Omega - right ascension of the ascending node,

omega - argument of pericenter,

i - inclination

stm_req : bool

Flag to calculate the derivative of the state-transition matrix

sources : dict

Dictionary that contains the perturbations that should be accounted.

The dictionary keys:

'atm' (Earth's atmosphere)

'j2' (Earth's J2)

'srp' (Solar radiation pressure)

'sun' (Gravitational acceleration of the Sun)

'mercury' (Gravitational acceleration of Mercury)

'venus' (Gravitational acceleration of Venus)

'earth' (Gravitational acceleration of the Earth)

'mars' (Gravitational acceleration of Mars)

'jupiter' (Gravitational acceleration of Jupiter)

'saturn' (Gravitational acceleration of Saturn)

'uranus' (Gravitational acceleration of Uranus)

'neptune' (Gravitational acceleration of Neptune)

'moon' (Gravitational acceleration of the Moon)

'cmplxmoon' (Complex gravitational acceleration of the Moon)

If sources[key] = True, the corresponding perturbation will be accounted.

If sources[key] = False, the corresponding perturbation will not be accounted.

For Earth's atmosphere, several levels are implemented.

If sources['atm'] == False, the atmosphere is not accounted.

If sources['atm'] == 'low', the low long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'mean', the mean long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'high', the high long term solar and geomagnetic activities are accounted.

The sources dictionary with all False values can be created by the kiam.prepare_sources_dict() function.

data : dict

A dictionary that contains auxilary data.

The dictionary keys:

'jd_zero' (Julian date that corresponds to t = 0)

'order' (Order of the lunar complex gravitational field)

'area' (Area of the spacecraft to account in atmospheric drag and SRP, m^2)

'mass' (Mass of the spacecraft to account in atmospheric drag and SRP, kg)

The data should be submitted even if the corresponding perturbations are not accounted.

units_data : dict

A dictionary that contains the units.

The dictionary keys:

'DistUnit' (The unit of distance in km)

'VelUnit' (The unit of velocity in km/s)

'TimeUnit' (The unit of time in days)

'AccUnit' (The unit of acceleration in m/s^2)

'RSun' (The radius of the Sun in the units of distance)

'REarth' (The radius of the Earth in the units of distance)

'RMoon' (The radius of the Moon in the units of distance)

'SunGM' (The nondimensional gravitational parameter of the Sun)

'MercuryGM' (The nondimensional gravitational parameter of Mercury)

'VenusGM' (The nondimensional gravitational parameter of Venus)

'EarthGM' (The nondimensional gravitational parameter of the Earth)

'MoonGM' (The nondimensional gravitational parameter of the Moon)

'MarsGM' (The nondimensional gravitational parameter of Mars)

'JupiterGM' (The nondimensional gravitational parameter of Jupiter)

'SaturnGM' (The nondimensional gravitational parameter of Saturn)

'UranusGM' (The nondimensional gravitational parameter of Uranus)

'NeptuneGM' (The nondimensional gravitational parameter of Neptune)

The units dictionary can be created by the kiam.prepare_units_dict() function.

The gravitational parameter in the specified units should be 1.0.

Returns:

f : numpy.ndarray, shape (6,), (42,)

Phase state time dericatives according to the specified n-body problem equations of motion extended (if stm_req = True) by the derivative of the state-transition matrix.

Vector structure:

[fh, fex, fey, fix, fiy, fL] if stm_req = False

[fh, fex, fey, fix, fiy, fL, fm11, fm21, fm31, … ] if stm_req = True

Examples:

t = 0.0

s = numpy.array([1, 0, 0, 0, 0, 0])

stm_req = False

sources = kiam.prepare_sources_dict()

data = kiam.prepare_data_dict()

data['jd_zero'] = kiam.juliandate(2022, 11, 1, 0, 0, 0)

data['area'] = 1.0

data['mass'] = 100.0

units_data = kiam.prepare_units_dict('earth')

dsdt = kiam.nbp_ee_earth(t, s, stm_req, sources, data, units_data)

print(dsdt)

# [0. 0. 0. 0. 0. 1.]
def nbp_ee_moon(t: float, s: numpy.ndarray, stm_req: bool, sources: dict, data: dict, units_data: dict) ‑> numpy.ndarray

Right-hand side of the n-body problem equations of motion wrt the Moon in terms of the equinoctial orbital elements.

Parameters:

t : float

Time

s : numpy.ndarray, shape (6,), (42,)

Phase state vector containing equinoctial elements and (if stm_req = True) vectorized state-transition matrix.

Vector structure:

[h, ex, ey, ix, iy, L] if stm_req = False,

[h, ex, ey, ix, iy, L, m11, m21, m31, …] if stm_req = True,

h = sqrt(p/mu),

ex = e*cos(Omega+omega),

ey = e*sin(Omega+omega),

ix = tan(i/2)*cos(Omega),

iy = tan(i/2)*sin(Omega),

L = theta + omega + Omega,

where

mu - gravitational parameter,

p - semi-latus rectum (focal parameter),

e - eccentricity,

Omega - right ascension of the ascending node,

omega - argument of pericenter,

i - inclination

stm_req : bool

Flag to calculate the derivative of the state-transition matrix

sources : dict

Dictionary that contains the perturbations that should be accounted. The dictionary keys:

'atm' (Earth's atmosphere)

'j2' (Earth's J2)

'srp' (Solar radiation pressure)

'sun' (Gravitational acceleration of the Sun)

'mercury' (Gravitational acceleration of Mercury)

'venus' (Gravitational acceleration of Venus)

'earth' (Gravitational acceleration of the Earth)

'mars' (Gravitational acceleration of Mars)

'jupiter' (Gravitational acceleration of Jupiter)

'saturn' (Gravitational acceleration of Saturn)

'uranus' (Gravitational acceleration of Uranus)

'neptune' (Gravitational acceleration of Neptune)

'moon' (Gravitational acceleration of the Moon)

'cmplxmoon' (Complex gravitational acceleration of the Moon)

If sources[key] = True, the corresponding perturbation will be accounted.

If sources[key] = False, the corresponding perturbation will not be accounted.

For Earth's atmosphere, several levels are implemented.

If sources['atm'] == False, the atmosphere is not accounted.

If sources['atm'] == 'low', the low long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'mean', the mean long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'high', the high long term solar and geomagnetic activities are accounted.

The sources dictionary with all False values can be created by the kiam.prepare_sources_dict() function.

data : dict

A dictionary that contains auxilary data.

The dictionary keys:

'jd_zero' (Julian date that corresponds to t = 0)

'order' (Order of the lunar complex gravitational field)

'area' (Area of the spacecraft to account in atmospheric drag and SRP, m^2)

'mass' (Mass of the spacecraft to account in atmospheric drag and SRP, kg)

The data should be submitted even if the corresponding perturbations are not accounted.

units_data : dict

A dictionary that contains the units.

The dictionary keys:

'DistUnit' (The unit of distance in km)

'VelUnit' (The unit of velocity in km/s)

'TimeUnit' (The unit of time in days)

'AccUnit' (The unit of acceleration in m/s^2)

'RSun' (The radius of the Sun in the units of distance)

'REarth' (The radius of the Earth in the units of distance)

'RMoon' (The radius of the Moon in the units of distance)

'SunGM' (The nondimensional gravitational parameter of the Sun)

'MercuryGM' (The nondimensional gravitational parameter of Mercury)

'VenusGM' (The nondimensional gravitational parameter of Venus)

'EarthGM' (The nondimensional gravitational parameter of the Earth)

'MoonGM' (The nondimensional gravitational parameter of the Moon)

'MarsGM' (The nondimensional gravitational parameter of Mars)

'JupiterGM' (The nondimensional gravitational parameter of Jupiter)

'SaturnGM' (The nondimensional gravitational parameter of Saturn)

'UranusGM' (The nondimensional gravitational parameter of Uranus)

'NeptuneGM' (The nondimensional gravitational parameter of Neptune)

The units dictionary can be created by the kiam.prepare_units_dict() function.

The gravitational parameter in the specified units should be 1.0.

Returns:

f : numpy.ndarray, shape (6,), (42,)

Phase state time dericatives according to the specified n-body problem equations of motion extended (if stm_req = True) by the derivative of the state-transition matrix.

Vector structure:

[fh, fex, fey, fix, fiy, fL] if stm_req = False

[fh, fex, fey, fix, fiy, fL, fm11, fm21, fm31, … ] if stm_req = True

Examples:

t = 0.0

s = numpy.array([1, 0, 0, 0, 0, 0])

stm_req = False

sources = kiam.prepare_sources_dict()

data = kiam.prepare_data_dict()

data['jd_zero'] = kiam.juliandate(2022, 11, 1, 0, 0, 0)

data['area'] = 1.0

data['mass'] = 100.0

units_data = kiam.prepare_units_dict('moon')

dsdt = kiam.nbp_ee_moon(t, s, stm_req, sources, data, units_data)

print(dsdt)

# [0. 0. 0. 0. 0. 1.]
def nbp_rv_body(body: str, t: float, s: numpy.ndarray, stm_req: bool, sources: dict, data: dict, units_data: dict) ‑> numpy.ndarray

Right-hand side of the n-body problem equations of motion wrt the specified body in terms of the position and velocity variables.

Parameters:

body : str

The body wrt that the right-hand side of the equations of motion is calculated.

Options: Sun, Mercury, Venus, Mars, Jupiter, Saturn, Uranus, Neptune

t : float

Time

s : numpy.ndarray, shape (6,), (42,)

Phase state vector containing position and velocity and (if stm_req = True) vectorized state-transition matrix.

Vector structure:

[x, y, z, vx, vy, vz] if stm_req = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm_req = True

stm_req : bool

Flag to calculate the derivative of the state-transition matrix

sources : dict

Dictionary that contains the perturbations that should be accounted.

The dictionary keys:

'srp' (Solar radiation pressure)

'sun' (Gravitational acceleration of the Sun)

'mercury' (Gravitational acceleration of Mercury)

'venus' (Gravitational acceleration of Venus)

'earth' (Gravitational acceleration of the Earth)

'mars' (Gravitational acceleration of Mars)

'jupiter' (Gravitational acceleration of Jupiter)

'saturn' (Gravitational acceleration of Saturn)

'uranus' (Gravitational acceleration of Uranus)

'neptune' (Gravitational acceleration of Neptune)

If sources[key] = True, the corresponding perturbation will be accounted.

If sources[key] = False, the corresponding perturbation will not be accounted.

The sources dictionary with all False values can be created by the kiam.prepare_sources_dict() function.

data : dict

A dictionary that contains auxilary data.

The dictionary keys:

'jd_zero' (Julian date that corresponds to t = 0)

'order' (Order of the lunar complex gravitational field)

'area' (Area of the spacecraft to account in atmospheric drag and SRP, m^2)

'mass' (Mass of the spacecraft to account in atmospheric drag and SRP, kg)

The data should be submitted even if the corresponding perturbations are not accounted.

units_data : dict

A dictionary that contains the units.

The dictionary keys:

'DistUnit' (The unit of distance in km)

'VelUnit' (The unit of velocity in km/s)

'TimeUnit' (The unit of time in days)

'AccUnit' (The unit of acceleration in m/s^2)

'RSun' (The radius of the Sun in the units of distance)

'REarth' (The radius of the Earth in the units of distance)

'RMoon' (The radius of the Moon in the units of distance)

'SunGM' (The nondimensional gravitational parameter of the Sun)

'MercuryGM' (The nondimensional gravitational parameter of Mercury)

'VenusGM' (The nondimensional gravitational parameter of Venus)

'EarthGM' (The nondimensional gravitational parameter of the Earth)

'MoonGM' (The nondimensional gravitational parameter of the Moon)

'MarsGM' (The nondimensional gravitational parameter of Mars)

'JupiterGM' (The nondimensional gravitational parameter of Jupiter)

'SaturnGM' (The nondimensional gravitational parameter of Saturn)

'UranusGM' (The nondimensional gravitational parameter of Uranus)

'NeptuneGM' (The nondimensional gravitational parameter of Neptune)

The units dictionary can be created by the kiam.prepare_units_dict() function.

The gravitational parameter in the specified units should be 1.0.

Returns:

f : numpy.ndarray, shape (6,), (42,)

Gravitational acceleration according to the specified n-body problem equations of motion extended (if stm_req = True) by the derivative of the state-transition matrix.

Vector structure:

[fx, fy, fz, fvx, fvy, fvz] if stm_req = False

[fx, fy, fz, fvx, fvy, fvz, fm11, fm21, fm31, … ] if stm_req = True

Examples:

t = 0.0

s = numpy.array([1, 0, 0, 0, 1, 0])

stm_req = False

sources = kiam.prepare_sources_dict()

sources['jupiter'] = True

data = kiam.prepare_data_dict()

data['jd_zero'] = kiam.juliandate(2022, 11, 1, 0, 0, 0)

data['area'] = 1.0

data['mass'] = 100.0

units_data = kiam.prepare_units_dict('sun')

dsdt = kiam.nbp_rv_body('sun', t, s, stm_req, sources, data, units_data)

print(dsdt)
def nbp_rv_earth(t: float, s: numpy.ndarray, stm_req: bool, sources: dict, data: dict, units_data: dict) ‑> numpy.ndarray

Right-hand side of the n-body problem equations of motion wrt the Earth in terms of the position and velocity variables.

Parameters:

t : float

Time

s : numpy.ndarray, shape (6,), (42,)

Phase state vector containing position and velocity and (if stm_req = True) vectorized state-transition matrix.

Vector structure:

[x, y, z, vx, vy, vz] if stm_req = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm_req = True

stm_req : bool

Flag to calculate the derivative of the state-transition matrix

sources : dict

Dictionary that contains the perturbations that should be accounted.

The dictionary keys:

'atm' (Earth's atmosphere)

'j2' (Earth's J2)

'srp' (Solar radiation pressure)

'sun' (Gravitational acceleration of the Sun)

'mercury' (Gravitational acceleration of Mercury)

'venus' (Gravitational acceleration of Venus)

'earth' (Gravitational acceleration of the Earth)

'mars' (Gravitational acceleration of Mars)

'jupiter' (Gravitational acceleration of Jupiter)

'saturn' (Gravitational acceleration of Saturn)

'uranus' (Gravitational acceleration of Uranus)

'neptune' (Gravitational acceleration of Neptune)

'moon' (Gravitational acceleration of the Moon)

'cmplxmoon' (Complex gravitational acceleration of the Moon)

If sources[key] = True, the corresponding perturbation will be accounted.

If sources[key] = False, the corresponding perturbation will not be accounted.

For Earth's atmosphere, several levels are implemented.

If sources['atm'] == False, the atmosphere is not accounted.

If sources['atm'] == 'low', the low long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'mean', the mean long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'high', the high long term solar and geomagnetic activities are accounted.

The sources dictionary with all False values can be created by the kiam.prepare_sources_dict() function.

data : dict

A dictionary that contains auxilary data.

The dictionary keys:

'jd_zero' (Julian date that corresponds to t = 0)

'order' (Order of the lunar complex gravitational field)

'area' (Area of the spacecraft to account in atmospheric drag and SRP, m^2)

'mass' (Mass of the spacecraft to account in atmospheric drag and SRP, kg)

The data should be submitted even if the corresponding perturbations are not accounted.

units_data : dict

A dictionary that contains the units.

The dictionary keys:

'DistUnit' (The unit of distance in km)

'VelUnit' (The unit of velocity in km/s)

'TimeUnit' (The unit of time in days)

'AccUnit' (The unit of acceleration in m/s^2)

'RSun' (The radius of the Sun in the units of distance)

'REarth' (The radius of the Earth in the units of distance)

'RMoon' (The radius of the Moon in the units of distance)

'SunGM' (The nondimensional gravitational parameter of the Sun)

'MercuryGM' (The nondimensional gravitational parameter of Mercury)

'VenusGM' (The nondimensional gravitational parameter of Venus)

'EarthGM' (The nondimensional gravitational parameter of the Earth)

'MoonGM' (The nondimensional gravitational parameter of the Moon)

'MarsGM' (The nondimensional gravitational parameter of Mars)

'JupiterGM' (The nondimensional gravitational parameter of Jupiter)

'SaturnGM' (The nondimensional gravitational parameter of Saturn)

'UranusGM' (The nondimensional gravitational parameter of Uranus)

'NeptuneGM' (The nondimensional gravitational parameter of Neptune)

The units dictionary can be created by the kiam.prepare_units_dict() function.

The gravitational parameter in the specified units should be 1.0.

Returns:

f : numpy.ndarray, shape (6,), (42,)

Gravitational acceleration according to the specified n-body problem equations of motion extended (if stm_req = True) by the derivative of the state-transition matrix.

Vector structure:

[fx, fy, fz, fvx, fvy, fvz] if stm_req = False

[fx, fy, fz, fvx, fvy, fvz, fm11, fm21, fm31, … ] if stm_req = True

Examples:

t = 0.0

s = numpy.array([1, 0, 0, 0, 1, 0])

stm_req = False

sources = kiam.prepare_sources_dict()

data = kiam.prepare_data_dict()

data['jd_zero'] = kiam.juliandate(2022, 11, 1, 0, 0, 0)

data['area'] = 1.0

data['mass'] = 100.0

units_data = kiam.prepare_units_dict('earth')

dsdt = kiam.nbp_rv_earth(t, s, stm_req, sources, data, units_data)

print(dsdt)

# [ 0.  1.  0. -1. -0. -0.]
def nbp_rv_moon(t: float, s: numpy.ndarray, stm_req: bool, sources: dict, data: dict, units_data: dict) ‑> numpy.ndarray

Right-hand side of the n-body problem equations of motion wrt the Moon in terms of the position and velocity variables.

Parameters:

t : float

Time

s : numpy.ndarray, shape (6,), (42,)

Phase state vector containing position and velocity and (if stm_req = True) vectorized state-transition matrix.

Vector structure:

[x, y, z, vx, vy, vz] if stm_req = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm_req = True

stm_req : bool

Flag to calculate the derivative of the state-transition matrix

sources : dict

Dictionary that contains the perturbations that should be accounted.

The dictionary keys:

'atm' (Earth's atmosphere)

'j2' (Earth's J2)

'srp' (Solar radiation pressure)

'sun' (Gravitational acceleration of the Sun)

'mercury' (Gravitational acceleration of Mercury)

'venus' (Gravitational acceleration of Venus)

'earth' (Gravitational acceleration of the Earth)

'mars' (Gravitational acceleration of Mars)

'jupiter' (Gravitational acceleration of Jupiter)

'saturn' (Gravitational acceleration of Saturn)

'uranus' (Gravitational acceleration of Uranus)

'neptune' (Gravitational acceleration of Neptune)

'moon' (Gravitational acceleration of the Moon)

'cmplxmoon' (Complex gravitational acceleration of the Moon)

If sources[key] = True, the corresponding perturbation will be accounted.

If sources[key] = False, the corresponding perturbation will not be accounted.

For Earth's atmosphere, several levels are implemented.

If sources['atm'] == False, the atmosphere is not accounted.

If sources['atm'] == 'low', the low long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'mean', the mean long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'high', the high long term solar and geomagnetic activities are accounted.

The sources dictionary with all False values can be created by the kiam.prepare_sources_dict() function.

data : dict

A dictionary that contains auxilary data.

The dictionary keys:

'jd_zero' (Julian date that corresponds to t = 0)

'order' (Order of the lunar complex gravitational field)

'area' (Area of the spacecraft to account in atmospheric drag and SRP, m^2)

'mass' (Mass of the spacecraft to account in atmospheric drag and SRP, kg)

The data should be submitted even if the corresponding perturbations are not accounted.

units_data : dict

A dictionary that contains the units.

The dictionary keys:

'DistUnit' (The unit of distance in km)

'VelUnit' (The unit of velocity in km/s)

'TimeUnit' (The unit of time in days)

'AccUnit' (The unit of acceleration in m/s^2)

'RSun' (The radius of the Sun in the units of distance)

'REarth' (The radius of the Earth in the units of distance)

'RMoon' (The radius of the Moon in the units of distance)

'SunGM' (The nondimensional gravitational parameter of the Sun)

'MercuryGM' (The nondimensional gravitational parameter of Mercury)

'VenusGM' (The nondimensional gravitational parameter of Venus)

'EarthGM' (The nondimensional gravitational parameter of the Earth)

'MoonGM' (The nondimensional gravitational parameter of the Moon)

'MarsGM' (The nondimensional gravitational parameter of Mars)

'JupiterGM' (The nondimensional gravitational parameter of Jupiter)

'SaturnGM' (The nondimensional gravitational parameter of Saturn)

'UranusGM' (The nondimensional gravitational parameter of Uranus)

'NeptuneGM' (The nondimensional gravitational parameter of Neptune)

The units dictionary can be created by the kiam.prepare_units_dict() function.

The gravitational parameter in the specified units should be 1.0.

Returns:

f : numpy.ndarray, shape (6,), (42,)

Gravitational acceleration according to the specified n-body problem equations of motion extended (if stm_req = True) by the derivative of the state-transition matrix.

Vector structure:

[fx, fy, fz, fvx, fvy, fvz] if stm_req = False

[fx, fy, fz, fvx, fvy, fvz, fm11, fm21, fm31, … ] if stm_req = True

Examples:

t = 0.0

s = numpy.array([1, 0, 0, 0, 1, 0])

stm_req = False

sources = kiam.prepare_sources_dict()

data = kiam.prepare_data_dict()

data['jd_zero'] = kiam.juliandate(2022, 11, 1, 0, 0, 0)

data['area'] = 1.0

data['mass'] = 100.0

units_data = kiam.prepare_units_dict('moon')

dsdt = kiam.nbp_rv_moon(t, s, stm_req, sources, data, units_data)

print(dsdt)

# [ 0.  1.  0. -1. -0. -0.]
def oe2rv(oe: numpy.ndarray, mu: float, grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Classical orbital elements to position and velocity.

Parameters:

oe : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column vectors of classical orbital elements:

a (semi-major axis),

e (eccentricity),

i (inclination),

Omega (right ascension of the ascending node),

omega (argument of pericenter),

theta (true anomaly).

mu : float

Gravitational parameter.

grad_req : bool

Flag to calculate the derivatives of position and velocity wrt elements.

Returns:

rv : numpy.ndarray, shape (6,), (6,n)

6D phase vector or array of 6D column vectors containing position and velocity.

Vector structure: [x, y, z, vx, dy, dz].

drv : numpy.ndarray, shape (6,6), (6,6,n)

6x6 matrix or 6x6xn array of partial derivatives of rv wrt oe (drv/doe).

Returns only if grad_req = True.

Examples:

oe = numpy.array([1, 0.1, 1.0, 0.0, 0.0, 0.0])

rv = kiam.oe2rv(oe, 1.0, False)

rv, drv = kiam.oe2rv(oe, 1.0, True)

print(rv)
def planet_state(jd: float, center: str, target: str) ‑> numpy.ndarray

Gives position and velocity of the planet at specified julian date wrt to the specified center (planet).

Parameters:

jd : float

Julian date

center : str

Name of the center planet

target : str

Name of the target planet

Returns:

state : numpy.ndarray, shape(6,)

State of the target planet wrt the center planet.

Position in km, velocity in km/s.

Examples:

s = kiam.planet_state(kiam.juliandate(2022, 12, 3, 0, 0, 0), 'Earth', 'Moon')

print(s)

# [ 3.76623766e+05  7.07472988e+04  1.01213236e+04
#  -1.36269070e-01 8.97864551e-01  4.72492325e-01 ]
def plot(x: numpy.ndarray, y: numpy.ndarray, fig: plotly.graph_objs._figure.Figure = None, xlabel: str = 'x', ylabel: str = 'y', name: str = '', axis_equal: bool = False, grid: str = 'on')

Creates a 2D line plot.

Parameters:

x : numpy.ndarray, shape (n,)

The x-axis nodes.

y : numpt.ndarray, shape (n,)

The y-axis data.

fig : plotly.graph_objects.Figure

The Plotly figure object. If provided (not None), then the line plot will be added to the existing figure in fig.

xlabel : str

The x-axis label.

ylabel : str

The y-axis label

name : str

The name of the plot to be indicated in the legend.

'axis_equal' : bool

Sets axis to be equal. False by default.

'grid' : str

Grid option setting.

Options:

None, 'plotly' – grid used by default by plotly (not used in papers normally)

'on' – white background, black dashed grid

'off' – disables the grid

Returns:

fig : plotly.graph_objects.Figure

The (updated) Plotly figure object.

Examples:

# Example 1 (minimal):

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig.show()

# Example 2:

x = numpy.array([1, 2, 3, 4, 5])

y1 = numpy.array([2, 3, 0, 1, 2])

y2 = numpy.array([3, 4, 1, 2, 3])

fig = kiam.plot(x, y1, name='blue')

fig = kiam.plot(x, y2, fig, name='red')  # add to the existing figure

fig.show()
def plot3(x: numpy.ndarray, y: numpy.ndarray, z: numpy.ndarray, fig: plotly.graph_objs._figure.Figure = None, xlabel: str = 'x', ylabel: str = 'y', zlabel: str = 'z', name: str = '', axis_equal: bool = False, grid: str = 'on')

Creates a 3D line plot.

Parameters:

x : numpy.ndarray, shape (n,)

The x-axis data.

y : numpt.ndarray, shape (n,)

The y-axis data.

z : numpt.ndarray, shape (n,)

The z-axis data.

fig : plotly.graph_objects.Figure

The Plotly figure object. If provided (not None), then the line plot will be added to the existing figure in fig.

xlabel : str

The x-axis label.

ylabel : str

The y-axis label

zlabel : str

The z-axis label

name : str

The name of the plot to be indicated in the legend.

'axis_equal' : bool

Sets axis to be equal. False by default.

'grid' : str

Grid option setting.

Options:

None, 'plotly' – grid used by default by plotly (not used in papers normally)

'on' – white background, black dashed grid

'off' – disables the grid

Returns:

fig : plotly.graph_objects.Figure

The (updated) Plotly figure object.

Examples:

# Example 1 (minimal):

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

z = numpy.array([3, 4, 1, 2, 3])

fig = kiam.plot3(x, y, z)

fig.show()

# Example 2:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

z1 = numpy.array([3, 4, 1, 2, 3])

z2 = numpy.array([4, 5, 2, 3, 4])

fig = kiam.plot3(x, y, z1, name='blue')

fig = kiam.plot3(x, y, z2, fig, name='red')

fig.show()
def polar_plot(r: numpy.ndarray, theta_deg: numpy.ndarray, mode: str = 'lines')

Plots a line in polar coordinates.

Parameters:

r : numpy.ndarray, shape (n,)

The radiuses.

theta_deg : numpy.ndarray, shape (n,)

The angles in degrees.

mode : str

The line display mode.

Options: 'lines' (default), 'markers', 'lines+markers'

Returns:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Examples:

r = numpy.array([0.5, 1, 2, 2.5, 3, 4])

theta = numpy.array([35, 70, 120, 155, 205, 240])

fig = kiam.polar_plot(r, theta)

fig.show()
def prepare_data_dict() ‑> dict

Auxilary function that returns a dictionary of data.

Returns:

data : dict

The data structure used in right-hand side of the equations of motion.

The dictionary keys:

'jd_zero' (Julian date that corresponds to t = 0)

'order' (Order of the lunar complex gravitational field)

'area' (Area of the spacecraft to account in atmospheric drag and SRP, m^2)

'mass' (Mass of the spacecraft to account in atmospheric drag and SRP, kg)

All the values are assigned to zero.

Examples:

print(kiam.prepare_data_dict())

# {'jd_zero': 0.0, 'order': 0, 'area': 0.0, 'mass': 0.0}
def prepare_sources_dict() ‑> dict

Auxilary function that returns a dictionary of perturbations.

Returns:

sources : dict

The dictionary keys:

'atm' (Earth's atmosphere)

'j2' (Earth's J2)

'srp' (Solar radiation pressure)

'sun' (Gravitational acceleration of the Sun)

'mercury' (Gravitational acceleration of Mercury)

'venus' (Gravitational acceleration of Venus)

'earth' (Gravitational acceleration of the Earth)

'mars' (Gravitational acceleration of Mars)

'jupiter' (Gravitational acceleration of Jupiter)

'saturn' (Gravitational acceleration of Saturn)

'uranus' (Gravitational acceleration of Uranus)

'neptune' (Gravitational acceleration of Neptune)

'moon' (Gravitational acceleration of the Moon)

'cmplxmoon' (Complex gravitational acceleration of the Moon)

Examples:

print(kiam.prepare_sources_dict())

# {'sun': False, 'mercury': False, 'venus': False, 'earth': False, 'moon': False, 'mars': False, 'jupiter': False, 'saturn': False, 'uranus': False, 'neptune': False, 'srp': False, 'cmplxmoon': False, 'atm': False, 'j2': False}
def prepare_units_dict(units_name: str) ‑> dict

Auxilary function that returns a dictionary of units.

Parameters:

units_name : str

A name of the units that should be used.

Options:

'dim': a dictionary with unity values will be returned

'earth': a dictionary of the earth units will be returned

'moon': a dictionary of the moon units will be returned

Returns:

units_dict : dict

A dictionary that containes the following keys:

'DistUnit': the unit of distance, km

'VelUnit': the unit of velocity, km/s

'TimeUnit': the unit of time, days

'AccUnit': the units of acceleration, m/s^2

'RSun': the mean radius of the Sun in DistUnit units

'REarth': the mean radius of the Earth in DistUnit units

'RMoon': the mean radius of the Moon in DistUnit units

Examples:

print(kiam.prepare_units_dict('earth'))

# {'DistUnit': 6371.0084, 'VelUnit': 7.909787126714006, 'TimeUnit': 0.009322440916154166, 'AccUnit': 9.820224438870717, 'RSun': 109.19778413728038, 'REarth': 1.0, 'RMoon': 0.27270408244949107}
def propagate_br4bp(central_body: str, tspan: numpy.ndarray, x0: numpy.ndarray, mu: float, gm4b, a4b: float, theta0: float, stm: bool, atol: float = 1e-10, rtol: float = 1e-10) ‑> tuple

Propagate trajectory in the bi-circular restricted four-body model of motion.

Parameters:

central_body : str

First or secondary primary body or baricenter as the origin of the coordinate system

Options: 'first', 'secondary', 'center'

tspan : numpy.ndarray, shape (n,)

Time nodes at which the solution is required

x0 : numpy.ndarray, shape (6,), (42,)

Initial state containing:

position and velocity (if stm = False),

position and velocoty extended by vectorized state-transition matrix (if stm = True),

Vectory structure:

[x, y, z, vx, vy, vz] if stm = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm = True

mu : float

Mass parameter of the three-body system

gm4b : float

Scaled gravitational parameter of the fourth (perturbing) body

a4b : float

Distance from the center of mass of the primary bodies to the fourth body in units where the distance between the primaries equals 1.

theta0 : float

Initial value of the synodic phase - the angle between the direction to the fourth body from the center of mass of the primaries and the line connecting the primaties.

stm : bool

Flag to calculate the derivative of the state-transition matrix.

atol : float

Absolute tolerance when integrating the equations. Default is 1e-10.

rtol : float

Relative tolerance when integrating the equations. Default is 1e-10.

Returns:

t : numpy.ndarray, shape(n,)

Times (nodes) in tspan at which the solution is obtained.

y : numpy.ndarray, shape(6, n)

Array of column trajectory phase states.

Vector structure:

[x, y, z, vx, vy, vz] if stm = False.

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm = True.

Examples:

central_body = 'first'

tspan = numpy.linspace(0, 10, 1000)

x0 = numpy.array([0.5, 0, 0, 0, 0.5, 0])

mu = 1.215e-02

stm = False

gm4b =  3.289005596e+05

a4b = 389.170375544352

theta0 = 0.0

t, y = kiam.propagate_br4bp(central_body, tspan, x0, mu, gm4b, a4b, theta0, stm)

print(t[-1], y[:, -1])
def propagate_cr3bp(central_body: str, tspan: numpy.ndarray, x0: numpy.ndarray, mu: float, stm: bool, atol: float = 1e-10, rtol: float = 1e-10) ‑> tuple

Propagate trajectory in the circular restricted three-body model of motion.

Parameters:

central_body : str

First or secondary primary body or barycenter as the origin of the coordinate system

Options: 'first', 'secondary', 'center'

tspan : numpy.ndarray, shape (n,)

Time nodes at which the solution is required

x0 : numpy.ndarray, shape (6,), (42,)

Initial state containing:

position and velocity (if stm = False),

position and velocoty extended by vectorized state-transition matrix (if stm = True),

Vectory structure:

[x, y, z, vx, vy, vz] if stm = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm = True

mu : float

Mass parameter of the three-body system

stm : bool

Flag to calculate the derivative of the state-transition matrix

atol : float

Absolute tolerance when integrating the equations. Default is 1e-10.

rtol : float

Relative tolerance when integrating the equations. Default is 1e-10.

Returns:

t : numpy.ndarray, shape(n,)

Times (nodes) in tspan at which the solution is obtained

y : numpy.ndarray, shape(6, n)

Array of column trajectory phase states.

Vector structure:

[x, y, z, vx, vy, vz] if stm = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if stm = True

Examples:

central_body = 'first'

tspan = numpy.linspace(0, 10, 1000)

x0 = numpy.array([0.5, 0, 0, 0, 0.5, 0])

mu = 1e-02

stm = False

t, y = kiam.propagate_cr3bp(central_body, tspan, x0, mu, stm)

print(t[-1], y[:, -1])
def propagate_hill(tspan: numpy.ndarray, x0: numpy.ndarray, atol: float = 1e-10, rtol: float = 1e-10) ‑> tuple

Propagate trajectory in Hill's model of motion.

Parameters:

tspan : numpy.ndarray, shape (n,)

Time nodes at which the solution is required

x0 : numpy.ndarray, shape (6,)

Initial state containing position and velocity.

Vector structure: [x, y, z, vx, vy, vz]

Returns:

t : numpy.ndarray, shape(n,)

Times (nodes) in tspan at which the solution is obtained

y : numpy.ndarray, shape(6, n)

Array of column trajectory phase states.

Vector structure: [x, y, z, vx, vy, vz].

atol : float

Absolute tolerance when integrating the equations. Default is 1e-10.

rtol : float

Relative tolerance when integrating the equations. Default is 1e-10.

Examples:

tspan = numpy.linspace(0, 2*numpy.pi, 1000)

x0 = numpy.array([-0.5, 0, 0, 0, 2.0, 0])

t, y = kiam.propagate_hill(tspan, x0)

print(t[-1], y[:, -1])
def propagate_nbp(central_body: str, tspan: numpy.ndarray, x0: numpy.ndarray, sources_dict: dict, dat_dict: dict, units_dict: dict, stm: bool, variables: str, atol: float = 1e-10, rtol: float = 1e-10, control_function: Callable = None) ‑> tuple

Propagate trajectory in the n-body model of motion.

Parameters:

central_body : str

Name of the central body

tspan : numpy.ndarray, shape (n,)

Time nodes at which the solution is required

x0 : numpy.ndarray, shape (6,), (42,)

Initial state containing:

position and velocity (if variables = 'rv', stm = False),

position and velocoty extended by vectorized state-transition matrix (if variables = 'rv_stm', stm = True),

equinoctial orbital elements (if variables = 'ee', stm = False),

equinoctial orbital elements extended by vectorized state-transition matrix (if variables = 'ee_stm', stm = True),

Vector structure:

[x, y, z, vx, vy, vz] if variables = 'rv' and stm = False

[x, y, z, vx, vy, vz, m11, m21, m31, …] if variables = 'rv_stm' and stm = True

[h, ex, ey, ix, iy, L] if variables = 'ee' and stm = False

[h, ex, ey, ix, iy, L, m11, m21, m31, …] if variables = 'ee_stm' and stm = True

h = sqrt(p/mu),

ex = e*cos(Omega+omega),

ey = e*sin(Omega+omega),

ix = tan(i/2)*cos(Omega),

iy = tan(i/2)*sin(Omega),

L = theta + omega + Omega,

where

mu - gravitational parameter,

p - semi-latus rectum (focal parameter),

e - eccentricity,

Omega - right ascension of the ascending node,

omega - argument of pericenter,

i - inclination

sources_dict : dict

Dictionary that contains the perturbations that should be accounted.

The dictionary keys:

'atm' (Earth's atmosphere)

'j2' (Earth's J2)

'srp' (Solar radiation pressure)

'sun' (Gravitational acceleration of the Sun)

'mercury' (Gravitational acceleration of Mercury)

'venus' (Gravitational acceleration of Venus)

'earth' (Gravitational acceleration of the Earth)

'mars' (Gravitational acceleration of Mars)

'jupiter' (Gravitational acceleration of Jupiter)

'saturn' (Gravitational acceleration of Saturn)

'uranus' (Gravitational acceleration of Uranus)

'neptune' (Gravitational acceleration of Neptune)

'moon' (Gravitational acceleration of the Moon)

'cmplxmoon' (Complex gravitational acceleration of the Moon)

If sources[key] = True, the corresponding perturbation will be accounted.

If sources[key] = False, the corresponding perturbation will not be accounted.

For Earth's atmosphere, several levels are implemented.

If sources['atm'] == False, the atmosphere is not accounted.

If sources['atm'] == 'low', the low long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'mean', the mean long term solar and geomagnetic activities are accounted.

If sources['atm'] == 'high', the high long term solar and geomagnetic activities are accounted.

The sources dictionary with all False values can be created by the kiam.prepare_sources_dict() function.

dat_dict : dict

A dictionary that contains auxilary data.

The dictionary keys:

'jd_zero' (Julian date that corresponds to t = 0)

'order' (Order of the lunar complex gravitational field)

'area' (Area of the spacecraft to account in atmospheric drag and SRP, m^2)

'mass' (Mass of the spacecraft to account in atmospheric drag and SRP, kg)

The data should be submitted even if the corresponding perturbations are not accounted.

'units_dict' : dict

A dictionary that contains the units of distance, velocity, time, acceleration, and the gravitational parameters of the bodies.

This variable can be generated by kiam.prepare_units_dict function.

stm : bool

Flag to calculate the derivative of the state-transition matrix

variables : str

Type of variables used to propagate the trajectory.

If stm = False, then variables should be 'rv' or 'ee'.

If stm = True, then variables should be 'rv_stm' or 'ee_stm'.

control_function : Callable

The control function that returns force vector, specific impulse and (if stm is True) force vector time derivative,

force vector state derivative, specific impulce time derivative, specific impulse state derivative. The control function

should take two arguments: the time and the phase state that corresponds to variables.

None by default.

atol : float

Absolute tolerance when integrating the equations. Default is 1e-10.

rtol : float

Relative tolerance when integrating the equations. Default is 1e-10.

Returns:

t : numpy.ndarray, shape(n,)

Times (nodes) in tspan at which the solution is obtained

y : numpy.ndarray, shape(6, n), shape(42, n)

Array of column trajectory phase states extended (if stm = True) by vectorized state-transition matrices.

Vector structure:

[x, y, z, vx, vy, vz] if stm = False and variables = 'rv'

[x, y, z, vx, vy, vz, m11, m21, m31, … ] if stm = True and variables = 'rv_stm'

[h, ex, ey, ix, iy, L] if stm = False and variables = 'ee'

[h, ex, ey, ix, iy, L, m11, m21, m31, … ] if stm = True and variables = 'ee_stm'

h = sqrt(p/mu),

ex = e*cos(Omega+omega),

ey = e*sin(Omega+omega),

ix = tan(i/2)*cos(Omega),

iy = tan(i/2)*sin(Omega),

L = theta + omega + Omega,

where

mu - gravitational parameter,

p - semi-latus rectum (focal parameter),

e - eccentricity,

Omega - right ascension of the ascending node,

omega - argument of pericenter,

i - inclination

Examples:

central_body = 'earth'

tspan = numpy.linspace(0, 100, 10000)

x0 = numpy.array([1, 0, 0, 0, 1, 0])

sources_dict = kiam.prepare_sources_dict()

dat_dict = kiam.prepare_data_dict()

units_dict = kiam.prepare_units_dict('earth')

stm = False

variables = 'rv'

t, y = kiam.propagate_nbp(central_body, tspan, x0, sources_dict, dat_dict, units_dict, stm, variables)

print(t[-1], y[:, -1])

Examples with using the control function can be found on GitHub: https://github.com/shmaxg/KIAMToolbox/tree/master/examples

def propagate_r2bp(tspan: numpy.ndarray, x0: numpy.ndarray, atol: float = 1e-10, rtol: float = 1e-10) ‑> tuple

Propagate trajectory in the two-body model of motion.

Parameters:

tspan : numpy.ndarray, shape (n,)

Time nodes at which the solution is required

x0 : numpy.ndarray, shape (6,)

Initial state containing position and velocity.

Vector structure: [x, y, z, vx, vy, vz]

Returns:

t : numpy.ndarray, shape(n,)

Times (nodes) in tspan at which the solution is obtained

y : numpy.ndarray, shape(6, n)

Array of column trajectory phase states.

Vector structure: [x, y, z, vx, vy, vz].

atol : float

Absolute tolerance when integrating the equations. Default is 1e-10.

rtol : float

Relative tolerance when integrating the equations. Default is 1e-10.

Examples:

tspan = numpy.linspace(0, 100, 10000)

x0 = numpy.array([1, 0, 0, 0, 1, 0])

t, y = kiam.propagate_r2bp(tspan, x0)

print(t[-1], y[:, -1])
def r2bp(t: float, s: numpy.ndarray) ‑> numpy.ndarray

Right-hand side of the restricted two-body problem equations of motion.

Parameters:

t : float

Time

s : numpy.ndarray, shape (6,)

Phase state vector containing position and velocity.

Vector structure [x, y, z, vx, vy, vz].

Returns:

f : numpy.ndarray, shape (6,)

Gravity acceleration according to the two-body model.

Vector structure [fx, fy, fz, fvx, fvy, fvz].

Examples:

t0 = 0.0

s0 = numpy.array([1, 0, 0, 0, 1, 0])

print(kiam.r2bp(t0, s0))

# [ 0.  1.  0. -1. -0. -0.]
def rad2deg(rad: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

Radians to degrees conversion.

Parameters:

rad : float, numpy.ndarray

Angle or array of angles in radians.

Returns:

deg : float, numpy.ndarray

Angle or array of angles in degrees.

Examples:

print(kiam.rad2deg(3.141592))

# 179.99996255206332
def rot2ine(xrot: numpy.ndarray, t: Union[float, numpy.ndarray], t0: Union[float, numpy.ndarray]) ‑> numpy.ndarray

Translate phase vectors from ROT c/s to INE c/s.

Parameters:

xrot : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the ROT coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

t : float, numpy.ndarray, shape (n,)

Time(s) corresponding to column(s) of xrot

t0 : float, numpy.ndarray, shape (1,), (n,)

Time(s) of of INE and ROT c/s coincidence for each column of xrot.

Returns:

xine : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the INE coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

xrot = numpy.array([1, 0, 0, 0, 0, 0])

t = 1.0

t0 = 0.0

xine = kiam.rot2ine(xrot, t, t0)

print(xine)
def rot2ine_eph(xrot: numpy.ndarray, jd: Union[float, numpy.ndarray], first_body: str, secondary_body: str, dist_unit: float, vel_unit: float) ‑> numpy.ndarray

Translate phase vectors from ROTEPH c/s to INEEPH c/s.

Parameters:

xrot : numpy.ndarray, shape (6,), (6,n)

6D vector array or 6D column phase vectors in the ROTEPH coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

jd : float, numpy.ndarray, shape (n,)

Julian date(s) corresponding to column(s) in xrot

first_body : str

Name of the first primary body

Options: 'sun', 'mercury', 'venus', 'earth', 'moon', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune'

secondary_body : str

Name of the secondary primary body

Options: 'sun', 'mercury', 'venus', 'earth', 'moon', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune'

dist_unit : float

The unit of distance in km

vel_unit : float

The unit of velocity in km/s

Returns:

xine : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the INEEPH coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

xrot = numpy.array([1, 0, 0, 0, 1, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

ku = kiam.units('earth', 'moon')

xine = kiam.ine2rot_eph(xrot, jd, 'earth', 'moon', ku['DistUnit'], ku['VelUnit'])

print(xine)
def rv2ee(rv: numpy.ndarray, mu: float, grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Position and velocity to equinoctial orbital elements.

Parameters:

rv : numpy.ndarray, shape (6,), (6,)

6D phase vector or array of 6D column vectors containing position and velocity.

Vector structure: [x, y, z, vx, dy, dz]

mu : float

Gravitational parameter

grad_req : bool

Flag to calculate the derivatives of elements wrt position and velocity

Returns:

ee : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column vectors of equinoctial orbital elements:

h = sqrt(p/mu),

ex = e*cos(Omega+omega),

ey = e*sin(Omega+omega),

ix = tan(i/2)*cos(Omega),

iy = tan(i/2)*sin(Omega),

L = theta + omega + Omega,

where

mu - gravitational parameter,

p - semi-latus rectum (focal parameter),

e - eccentricity,

Omega - right ascension of the ascending node,

omega - argument of pericenter,

i - inclination.

dee : numpy.ndarray, shape (6,6), (6,6,n)

6x6 matrix or 6x6xn array of partial derivatives of ee wrt rv (dee/drv).

Returns only if grad_req = True.

Examples:

rv = numpy.array([1, 0, 0, 0, 1, 0])

ee = kiam.rv2ee(rv, 1.0, False)

ee, dee = kiam.rv2ee(rv, 1.0, True)

print(ee)
def rv2oe(rv: numpy.ndarray, mu: float, grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Position and velocity to classical orbital elements.

Parameters:

rv : numpy.ndarray, shape (6,), (6,n)

6D phase vector or array of column 6D phase vectors containing position and velocity.

Vector structure: [x, y, z, vx, dy, dz].

mu : float

Gravitational parameter.

grad_req : bool

Flag to calculate the derivatives of elements wrt position and velocity.

Returns:

oe : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column vectors of classical orbital elements:

a (semi-major axis),

e (eccentricity),

i (inclination),

Omega (right ascension of the ascending node),

omega (argument of pericenter),

theta (true anomaly)

doe : numpy.ndarray, shape (6,6), (6,6,n)

6x6 matrix or 6x6xn array of partial derivatives of oe wrt rv (doe/drv).

Returns only if grad_req = True.

Examples:

rv = numpy.array([1, 0, 0, 0.1, 1, 0.1])

oe = kiam.rv2oe(rv, 1.0, False)

oe, doe = kiam.rv2oe(rv, 1.0, True)

print(oe)
def save(variable: Any, filename: str) ‑> None

Saves a variable into a specified file.

Parameters:

variable : Any

Variable to be saved.

For limitations on variables see the pickle package https://docs.python.org/3/library/pickle.html

filename : str

A path to the file.

Examples:

a = numpy.random.rand(10, 10)

kiam.save(a, 'variable_a')
def save_figure(fig: plotly.graph_objs._figure.Figure, filename: str)

Saves the figure as an interactive html file.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

filename : str

The file name (or file path) with extension to which the figure should be saved.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

kiam.save_figure(fig, 'myfig.html')
def save_image(fig: plotly.graph_objs._figure.Figure, filename: str, scale: int = 2)

Saves the figure as a static image (PNG, PDF, etc).

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

filename : str

The file name (or file path) with extension to which the figure should be saved.

scale : int

The scale parameter controls dpi. Default is 2.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

kiam.save_image(fig, 'myfig.png')
def scrs2gcrs(xscrs: numpy.ndarray, jd: Union[float, numpy.ndarray], dist_unit: float, vel_unit: float) ‑> numpy.ndarray

Translate phase vectors from SCRS c/s to GCRS c/s.

Parameters:

xscrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the SCRS coordinate system

Vector structure: [x, y, z, vx, vy, vz]

jd : float, numpy.ndarray, shape (n,)

Julian dates corresponding to columns in xscrs

dist_unit : float

The unit of distance in km

vel_unit : float

The unit of velocity in km/s

Returns:

xgcrs : numpy.ndarray, shape (6,), (6,n)

6D vector or array of 6D column phase vectors in the GCRS coordinate system.

Vector structure: [x, y, z, vx, vy, vz].

Examples:

# Example 1 (6D -> 6D):

ku = kiam.units('earth', 'moon')

xscrs = numpy.array([1, 0, 0, 0, 1, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xgcrs = kiam.scrs2gcrs(xscrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xgcrs)

# Example 2 (6x1 -> 6x1):

ku = kiam.units('earth', 'moon')

xscrs = numpy.array([[1, 0, 0, 0, 1, 0]]).T

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xgcrs = kiam.scrs2gcrs(xscrs, jd, ku['DistUnit'], ku['VelUnit'])

print(xgcrs)
def scrs2mer(xscrs: numpy.ndarray, jd: Union[float, numpy.ndarray], grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Translate vectors from SCRS c/s to MER c/s.

Parameters:

xscrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the SCRS coordinate system

jd : float, numpy.ndarray, shape (n,)

Julian date(s) corresponding to vector or columns in xscrs

grad_req : bool

Flag to calculate the derivatives of the MER vector wrt the SCRS vector

Returns:

xmer : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the MER coordinate system.

dxmer : numpy.ndarray, shape (3,3), (6,6), (3,3,n), (6,6,n)

Matrix or array of matrices of partial derivatives of xmer wrt xscrs (dxmer/dxscrs).

Returns only if grad_req = True.

Examples:

xscrs = numpy.array([1, 0, 0])

jd = kiam.juliandate(2022, 11, 22, 0, 0, 0)

xmer = kiam.scrs2mer(xscrs, jd, False)

dxmer = kiam.scrs2mer(xscrs, jd, True)

print(xmer)
def scrs2pa(xscrs: numpy.ndarray, jd: Union[float, numpy.ndarray], grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Translate vector from SCRS c/s to PA c/s.

Parameters:

xscrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the SCRS coordinate system

jd : float, numpy.ndarray, shape (n,)

Julian date(s) corresponding to column(s) in xscrs

grad_req : bool

Flag to calculate the derivatives of the PA vector wrt the SCRS vector

Returns:

xpa : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the PA coordinate system

dxpa : numpy.ndarray, shape (3,3), (6,6), (3,3,n), (6,6,n)

3x3 or 6x6 matrix or 3x3xn or 6x6xn array of partial derivatives of xpa wrt xscrs (dxpa/dxscrs).

Returns only if grad_req = True.

Examples:

jd = kiam.juliandate(2022, 11, 22, 0, 0, 0)

xscrs = numpy.array([1, 0, 0])

xpa = kiam.scrs2pa(xscrs, jd, False)

xpa, dxpa = kiam.scrs2pa(xscrs, jd, True)

print(xpa)
def scrs2sors(xscrs: numpy.ndarray, jd: Union[float, numpy.ndarray], grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Translate vectors from SCRS c/s to SORS c/s.

Parameters:

xscrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the SCRS coordinate system

jd : float, numpy.ndarray, shape (n,)

Julian dates corresponding to columns in xscrs

grad_req : bool

Flag to calculate the derivatives of the SORS vector wrt the SCRS vector

Returns:

xsors : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the SORS coordinate system

dxsors : numpy.ndarray, shape (3,3), (6,6), (3,3,n), (6,6,n)

Matrix or array of matrices of partial derivatives of xsors wrt xscrs (dxsors/dxscrs).

Returns only if grad_req = True.

Examples:

xscrs = numpy.array([1, 0, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xsors = kiam.scrs2sors(xscrs, jd, False)

xsors, dxsors = kiam.scrs2sors(xscrs, jd, True)

print(xsors)
def set_axis_equal(fig: plotly.graph_objs._figure.Figure)

Sets axis to be equal.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Only Scatter and Scatter3d figure types are supported.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.set_axis_equal(fig)

fig.show()
def set_default_grid(fig: plotly.graph_objs._figure.Figure)

Sets the default grid.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Only Scatter and Scatter3d figure types are supported.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.grid_off(fig)

fig.show()

fig = kiam.set_default_grid(fig)

fig.show()
def set_xlabel(fig: plotly.graph_objs._figure.Figure, xlabel: str)

Set a custom x-axis label.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Only Scatter and Scatter3d figure types are supported.

xlabel : str

The new x-axis label.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Only Scatter and Scatter3d figure datatypes supported.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.set_xlabel(fig, 'x variable')

fig.show()
def set_ylabel(fig: plotly.graph_objs._figure.Figure, ylabel: str)

Set a custom y-axis label.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Only Scatter and Scatter3d figure types are supported.

ylabel : str

The new y-axis label.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot(x, y)

fig = kiam.set_ylabel(fig, 'y variable')

fig.show()
def set_zlabel(fig: plotly.graph_objs._figure.Figure, zlabel: str)

Set a custom y-axis label.

Parameters:

fig : plotly.graph_objects.Figure

The Plotly figure object.

Only Scatter3d figure types are supported.

zlabel : str

The new z-axis label.

Returns:

fig : plotly.graph_objects.Figure

The updated Plotly figure object.

Examples:

x = numpy.array([1, 2, 3, 4, 5])

y = numpy.array([2, 3, 0, 1, 2])

z = numpy.array([2, 3, 0, 1, 2])

fig = kiam.plot3(x, y, z)

fig = kiam.set_zlabel(fig, 'z variable')

fig.show()
def sind(x: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

Sine of a degree argument.

Parameters:

x : float, numpy.ndarray

Angle or an array of angles in degrees.

Returns:

s : float, numpy.ndarray

A sine or array of sines of angles in degrees.

Examples:

print(kiam.sind(30))

# 0.49999999999999994
def sors2scrs(xsors: numpy.ndarray, jd: Union[float, numpy.ndarray], grad_req: bool = False) ‑> Union[numpy.ndarray, tuple[numpy.ndarray, numpy.ndarray]]

Translate vectors from SORS c/s to SCRS c/s.

Parameters:

xsors : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the SORS coordinate system

jd : float, numpy.ndarray, shape (n,)

Julian dates corresponding to columns in xsors

grad_req : bool

Flag to calculate the derivatives of the SCRS vector wrt the SORS vector

Returns:

xscrs : numpy.ndarray, shape (3,), (6,), (3,n), (6,n)

3D vector, 6D vector or array of 3D or 6D column vectors in the SCRS coordinate system

dxscrs : numpy.ndarray, shape (3,3), (6,6), (3,3,n), (6,6,n)

Matrix or array of matrices of partial derivatives of xscrs wrt xsors (dxscrs/dxsors).

Returns only if grad_req = True.

Examples:

xsors = numpy.array([1, 0, 0])

jd = kiam.juliandate(2022, 12, 6, 0, 0, 0)

xscrs = kiam.sors2scrs(xsors, jd, False)

xscrs, dxscrs = kiam.sors2scrs(xsors, jd, True)

print(xscrs)
def sphere2cart(sphere: numpy.ndarray) ‑> numpy.ndarray

Spherical coordinates to Cartesian coordinates.

Parameters:

sphere : numpy.ndarray, shape (3,), (3, n)

3D vector or 3xn array of column 3D vectors of spherical coordinates

Vector structure: [r, phi, theta], where


phi in [-pi, pi],

theta in [0, pi],

x = r*cos(theta)*cos(phi),

y = r*cos(theta)*sin(phi),

z = r*sin(theta)

Returns:

cart : numpy.ndarray, shape (3,), (3, n)

3D vector or 3xn array of column 3D vectors of Cartesian coordinates.

Vector structure: [x, y, z].

Examples:

sphere = numpy.array([1, 0, 0])

cart = kiam.sphere2cart(sphere)

print(cart)

# [0. 0. 1.]
def sphere_coordinates(radius: float, nlat: int, nlon: int)

Get x, y, z coordinates on a sphere.

Parameters:

radius : float

The radius of the sphere.

nlat : int

The number of latitude angles in a grid.

nlon : int

The number of longitude angles in a grid.

Returns:

x : numpy.ndarray, shape(nlat, nlon)

The x-coordinates.

y : numpy.ndarray, shape(nlat, nlon)

The y-coordinates.

z : numpy.ndarray, shape(nlat, nlon)

The z-coordinates.

Examples:

x, y, z = kiam.sphere_coordinates(1.0, 100, 100)
def ta2ea(ta: Union[float, numpy.ndarray], ecc: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

True anomaly to eccentric anomaly.

Parameters:

ta : float, numpy.ndarray, shape (n,)

Scalar or array of true anomalies.

ecc : float, numpy.ndarray, shape (n,)

Scalar or array of eccentricities. In case of array, the dimension should match the one of ta.

Returns:

ea : float, numpy.ndarray, shape (n,)

Scalar or array of eccentric anomalies. Domain: (-pi, pi).

If ta and ecc are scalars, then ea is a scalar.

If ta is a scalar, ecc is a vector, then ea is a vector of the same size as ecc.

If ta is a vector, ecc is a scalar, then ea is a vector of the same size as ta.

If ta and ecc are vectors with the same size, then ea is a vector of the same size.

Examples:

ta = numpy.array([0.0, numpy.pi])

ecc = 0.1

ea = kiam.ta2ea(ta, ecc)
def tai2utc(tai: datetime.datetime) ‑> datetime.datetime

TAI to UTC conversion.

Parameters:

tai : datetime.datetime

The International Atomic Time (TAI).

Returns:

utc : datetime.datetime

The Coordinated Universal Time (UTC).

Examples:

tai = datetime.datetime(2023, 3, 8, 12, 0, 0)  # 2023-03-18 12:00:00

utc = kiam.tai2utc(tai)  # 2023-03-08 11:59:27.816000

print(utc)
def tand(x: Union[float, numpy.ndarray]) ‑> Union[float, numpy.ndarray]

Tangent of a degree argument.

Parameters:

x : float, numpy.ndarray

Angle or an array of angles in degrees.

Returns:

s : float, numpy.ndarray

A tangent or array of tangents of angles in degrees.

Examples:

print(kiam.tand(45))

0.9999999999999999
def time2jd(time: datetime.datetime) ‑> float

Usual date and time to Julian date.

Parameters:

time : datetime.datetime

Date and time object of type datetime.datetime

Returns:

jd : float

Julian date

Examples:

jd = kiam.time2jd(datetime.datetime(2022, 11, 22, 0, 0, 0, 0))

print(jd)

# 2459905.5
def to_float(*args: Any) ‑> tuple

Convert all arguments to the float64 type.

Parameters:

*args

Arguments separated by comma to convert to float64.

Returns:

float_args : tuple

Tuple of numpy arrays with components converted to float64 arguments.

Examples:

f = kiam.to_float([1, 2], 3, [4, 5, 6])

print(f)

# (array([1., 2.]), array(3.), array([4., 5., 6.]))
def tt2utc(tt: datetime.datetime) ‑> datetime.datetime

TT to UTC conversion.

Parameters:

tt : datetime.datetime

The Terrestrial Time (TT).

Returns:

utc : datetime.datetime

The Coordinated Universal Time (UTC).

Examples:

tt = datetime.datetime(2023, 3, 8, 12, 0, 0)  # 2023-03-18 12:00:00

utc = kiam.tt2utc(tt)  # 2023-03-08 11:58:50.816000

print(utc)
def units(*args: str) ‑> dict

Get units of distance, velocity, time, and gravitational parameters.

Parameters:

*args

Name or a pair of names of a celestial bodies

Options for a single argument: 'earth', 'moon', 'sun', 'mercury', 'venus', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune', 'pluto'

Options for two arguments: ('earth', 'moon'), ('sun', 'earth')

Returns:

units_dict : dict

A dictionary containing the units of distance, velocity, and time.

'DistUnit' – the unit of distance, km

'VelUnit' – the unit of velocity, km/s

'TimeUnit' – the unit of time, days

'AccUnit' – the unit of acceleration, m/s^2

'SunGM' – the nondimensional gravitational parameter of the Sun

'MercuryGM' – the nondimensional gravitational parameter of Mercury

'VenusGM' – the nondimensional gravitational parameter of Venus

'EarthGM' – the nondimensional gravitational parameter of the Earth

'MoonGM' – the nondimensional gravitational parameter of the Moon

'EarthMoonGM' – the nondimensional gravitational parameter of the Earth+Moon system

'MarsGM' – the nondimensional gravitational parameter of Mars

'JupiterGM' – the nondimensional gravitational parameter of Jupiter

'SaturnGM' – the nondimensional gravitational parameter of Saturn

'UranusGM' – the nondimensional gravitational parameter of Uranus

'NeptuneGM' – the nondimensional gravitational parameter of Neptune

Examples:

un = kiam.units('earth')

DU = un['DistUnit']  # Unit of distance for the earth system of units

print(DU)

un = kiam.units('earth', 'moon')

VU = un['VelUnit']  # Unit of velocity for the Earth-Moon system of units

print(VU)
def utc2tai(utc: datetime.datetime) ‑> datetime.datetime

UTC to TAI conversion.

Parameters:

utc : datetime.datetime

The Coordinated Universal Time (UTC).

Returns:

tai : datetime.datetime

The International Atomic Time (TAI).

Examples:

utc = datetime.datetime(2023, 3, 8, 12, 0, 0)  # 2023-03-18 12:00:00

tai = kiam.utc2tai(utc)  # 2023-03-08 12:00:32.184000

print(tai)
def utc2tt(utc: datetime.datetime) ‑> datetime.datetime

UTC to TT conversion.

Parameters:

utc : datetime.datetime

The Coordinated Universal Time (UTC).

Returns:

tt : datetime.datetime

The Terrestrial Time (TT).

Examples:

utc = datetime.datetime(2023, 3, 8, 12, 0, 0)  # 2023-03-18 12:00:00

tt = kiam.utc2tt(utc)  # 2023-03-08 12:01:09.184000

print(tt)
def vec2mat(v: numpy.ndarray) ‑> numpy.ndarray

Vector to square matrix translation.

Parameters:

v : numpy.ndarray, shape (n**2,)

A vector.

Returns:

a : numpy.ndarray, shape (n, n)

A square matrix.

Matrix structure (Fortran/MATLAB order): [[v1, v2, ..., vn], [v_(n+1), ...]].T

Examples:

v = numpy.array([1, 2, 3, 4])

m = kiam.vec2mat(v)

print(m)

# [[1 3]
# [2 4]]