environment

Functionalities of environment objects.

This module provides functionalities for environment objects. Specifically, it contains classes and functions that perform computations related to environment models of natural and artificial bodies. Much of the functionality in this module concerns classes stored inside Body objects, a list of which is in turn stored in a SystemOfBodies object. Note that the classes in this module are rarely created manually, but are instead created by the functionality in the environment_setup module.

Functions

transform_to_inertial_orientation(...)

Function to convert a Cartesian state vector from a body-fixed to an inertial frame

save_vehicle_mesh_to_file(...[, ...])

Function to save the mesh used for a hypersonic local inclination analysis to a file.

save_vehicle_mesh_to_file(local_inclination_analysis_object: tudatpy.kernel.numerical_simulation.environment.HypersonicLocalInclinationAnalysis, output_directory: str, output_file_prefix: str = '') None

Function to save the mesh used for a hypersonic local inclination analysis to a file.

Function to save the mesh used for a hypersonic local inclination analysis to a file. This function saves two files to the specified directory, with filenames: “ShapeFile.dat” and “SurfaceNormalFile.dat”, where these files names may be prefixed by an optional string (see below). The first of these files contains four columns defining the surface points that define mesh, with Column 0: point index; Column 1: x-position of point; Column 1: y-position of point; Column 2: z-position of point. The second file contains four columns with Column 0: point index; Column 1: x-component of surface normal; Column 1: y-position of surface normal; Column 2: z-position of surface normal.

Parameters:
  • local_inclination_analysis_object (HypersonicLocalInclinationAnalysis) – Object used to calculate the aerodynamics of the vehicle

  • output_directory (str) – Directory to which the files are to be saved

  • output_file_prefix (str, default='') – Optional prefix of output file names

transform_to_inertial_orientation(state_in_body_fixed_frame: numpy.ndarray[numpy.float64[6, 1]], current_time: float, rotational_ephemeris: tudatpy.kernel.numerical_simulation.environment.RotationalEphemeris) numpy.ndarray[numpy.float64[6, 1]]

Function to convert a Cartesian state vector from a body-fixed to an inertial frame

Function to convert a Cartesian state vector from a body-fixed to an inertial frame, using a RotationalEphemeris object as a model for the rotation. The body-fixed frame from which the conversion takes place is the body_fixed_frame_name frame, the (assumedly) inertial frame to which the conversion is done is inertial_frame_name.

This function calls rotate_state_to_frame() (with frame \(A\) the inertial frame, and frame \(B\) the body-fixed frame). The present function computes the required rotation matrix and its time derivative from the rotational_ephemeris input given here.

Parameters:
  • state_in_body_fixed_frame (numpy.ndarray[numpy.float64[6, 1]]) – Cartesian state (position and velocity) in the body-fixed frame

  • current_time (float) – Time at which the transformation is to be computed

  • rotational_ephemeris (RotationalEphemeris) – Boy rotation model that is to be used to convert the body-fixed state to inertial state

Returns:

Cartesian state transformed to inertial frame, using rotational_ephemeris model, from body-fixed state_in_body_fixed_frame

Return type:

numpy.ndarray[numpy.float64[6, 1]]

Enumerations

AerodynamicsReferenceFrames

Enumeration of reference frame identifiers typical for aerodynamic calculations.

AerodynamicCoefficientFrames

Enumeration of reference frames used for definition of aerodynamic coefficients.

AerodynamicCoefficientsIndependentVariables

Enumeration of the independent variables that can be used to compute aerodynamic coefficients.

AtmosphericCompositionSpecies

No documentation found.

class AerodynamicsReferenceFrameAngles

No documentation found.

Members:

latitude_angle

longitude_angle

heading_angle

flight_path_angle

angle_of_attack

angle_of_sideslip

bank_angle

property name
class AerodynamicsReferenceFrames

Enumeration of reference frame identifiers typical for aerodynamic calculations.

Enumeration of reference frame identifiers typical for aerodynamic calculations. Note that the frames are also defined in the absence of any aerodynamic forces and/or atmosphere. They define frames of a body w.r.t. a central body, with the details given by Mooij (1994). The chain of frames starts from the inertial frame, to the frame fixed to the central body (corotating), to the vertical frame (defined by the body’s relative position), the trajectory and aerodynamic frames (defined by the body’s relative velocity) and finally the body’s own body-fixed frame.

Members:

inertial_frame :

The global orientation (which is by definition inertial).

corotating_frame :

The body-fixed frame of the central body.

vertical_frame :

Frame with z-axis pointing towards origin of central body, the x-axis lies in the meridian plane and points towards the central-body-fixed z-axis (the y-axis completes the frame).

trajectory_frame :

The (airspeed-based) trajectory frame has the x-axis in the direction of the velocity vector relative to the atmosphere (airspeed-based velocity vector), z-axis lies in the vertical plane and points downwards (the y-axis completes the frame).

aerodynamic_frame :

The (airspeed-based) aerodynamic frame has the x-axis in the direction of the velocity vector relative to the atmosphere (airspeed-based velocity vector), z-axis co-linear with the aerodynamic lift vector, pointing in the opposite direction (the y-axis completes the frame)..

body_frame :

The body-fixed frame of the body itself.

property name
class AerodynamicCoefficientFrames

Enumeration of reference frames used for definition of aerodynamic coefficients.

Enumeration of reference frames used for definition of aerodynamic coefficients. There is a partial overlap between this enum and the AerodynamicsReferenceFrames. This enum combines a subset of those frames (which are typically used for aerodynamic coefficient definition), and a swap in sign. For instance, aerodynamic force coefficients are often defined positive along negative axes of the aerodynamic frame (drag, side force and lift coefficients)

Members:

positive_body_fixed_frame_coefficients :

The coefficients are defined in the body-fixed frame, with the directions the same as the body-fixed axes. For aerodynamic forces and moments, this results in the typical \(C_{x}, C_{y}, C_{y}\) (force) and \(C_{l}, C_{m}, C_{n}\) (moment) coefficients

negative_body_fixed_frame_coefficients :

Same as positive_body_fixed_frame_coefficients, but opposite in direction (so axes along negative body-fixed frame axes)

positive_aerodynamic_frame_coefficients :

Same as negative_aerodynamic_frame_coefficients, but opposite in direction (so axes along positive aerodynamic frame axes)

negative_aerodynamic_frame_coefficients :

The coefficients are defined in aerodynamic frame, with the directions the same as the negative axes. For aerodynamic forces, this results in the typical \(C_{D}, C_{S}, C_{D}\) force coefficients

property name
class AerodynamicCoefficientsIndependentVariables

Enumeration of the independent variables that can be used to compute aerodynamic coefficients.

Members:

mach_number_dependent :

angle_of_attack_dependent :

sideslip_angle_dependent :

altitude_dependent :

time_dependent :

temperature_dependent :

velocity_dependent :

he_number_density_dependent :

o_number_density_dependent :

n2_number_density_dependent :

o2_number_density_dependent :

ar_number_density_dependent :

h_number_density_dependent :

n_number_density_dependent :

anomalous_o_number_density_dependent :

control_surface_deflection_dependent : No documentation found.

undefined_independent_variable :

Can be used for a custom coefficient interface with other variables, at the expense of being able to use the FlightConditions class to automatically updates the aerodynamic coefficients during propagation.

property name
class AtmosphericCompositionSpecies

No documentation found.

Members:

o_species : No documentation found.

o2_species : No documentation found.

n2_species : No documentation found.

he_species : No documentation found.

h_species : No documentation found.

ar_species : No documentation found.

n_species : No documentation found.

anomalous_o_species : No documentation found.

property name

Classes

Ephemeris

Object that computes the state of a body as a function of time

RotationalEphemeris

Object that stores the rotational state of the bodies.

GcrsToItrsRotationModel

Object for high-accuracy GCRS<->ITRS rotation.

EarthOrientationAnglesCalculator

Object for computing high-accuracy Earth orientation angles

GravityFieldModel

Object that provides the gravity field of a body

SphericalHarmonicsGravityField

Object that provides a spherical harmonic gravity field of a body.

BodyShapeModel

Object that provides a shape model for a natural body.

RigidBodyProperties

Object that defines the mass, center of mass, and inertia tensor as a function of time.

AtmosphereModel

Object that provides the atmospheric properties of the body.

AerodynamicCoefficientInterface

Base class for computing the current aerodynamic coefficients of the body

HypersonicLocalInclinationAnalysis

GroundStation

Object used to define and store properties of a ground station.

GroundStationState

Object that performs computations of the current (body-fixed) position and frame conversions of the ground station.

FlightConditions

Object that calculates various state-derived quantities typically relevant for flight dynamics.

AtmosphericFlightConditions

Object that calculates various state-derived quantities typically relevant for flight dynamics, for flight in an atmosphere.

AerodynamicAngleCalculator

Object to calculate (aerodynamic) orientation angles, and frame transformations, from current vehicle state.

VehicleSystems

Object used to store physical (hardware) properties of a vehicle.

EngineModel

Body

Object that stores the environment properties and current state of a single body.

SystemOfBodies

Object that contains a set of Body objects and associated frame information.

class Ephemeris

Object that computes the state of a body as a function of time

Object (typically stored inside a Body object) that computes the state of a body as a function of time, both outside of a propagation, and during a propagation if the given body’s translational state is not propagated. Note that this object computes the state w.r.t. its own origin (defined by frame_origin), which need not be the same as the global frame origin of the environment.

cartesian_position(self: tudatpy.kernel.numerical_simulation.environment.Ephemeris, current_time: float) numpy.ndarray[numpy.float64[3, 1]]

As cartesian_state, but only the three position components

Parameters:

current_time (float) – Time (in seconds since J2000 in TDB time scale) at which the state is to be computed.

Returns:

Requested Cartesian position

Return type:

numpy.ndarray

cartesian_state(self: tudatpy.kernel.numerical_simulation.environment.Ephemeris, current_time: float) numpy.ndarray[numpy.float64[6, 1]]

This function returns the Cartesian state (position and velocity) at the given time, w.r.t. the frame_origin.

Parameters:

current_time (float) – Time (in seconds since J2000 in TDB time scale) at which the state is to be computed.

Returns:

Requested Cartesian state

Return type:

numpy.ndarray

cartesian_velocity(self: tudatpy.kernel.numerical_simulation.environment.Ephemeris, current_time: float) numpy.ndarray[numpy.float64[3, 1]]

As cartesian_state, but only the three velocity components

Parameters:

current_time (float) – Time (in seconds since J2000 in TDB time scale) at which the state is to be computed.

Returns:

Requested Cartesian velocity

Return type:

numpy.ndarray

property frame_orientation

read-only

Name of the frame orientation w.r.t which this object provides its states

Type:

str

property frame_origin

read-only

Name of the reference body/point w.r.t. which this object provides its states

Type:

str

class RotationalEphemeris

Object that stores the rotational state of the bodies.

Object that stores the rotational state of the bodies. This object can be used to calculate rotation matrices, which are used to transform coordinates between reference frames.

angular_velocity_in_body_fixed_frame(self: tudatpy.kernel.numerical_simulation.environment.RotationalEphemeris, time: float) numpy.ndarray[numpy.float64[3, 1]]

Function to get the body’s angular velocity vector, expressed in the body-fixed frame.

Function to get the body’s angular velocity vector \(\boldsymbol{\omega}^{(B)}\), expressed in the body-fixed frame \(B\). The calculation of the angular velocity depends on the specific rotation model that has been defined, either from an a priori definition (see rotation_model submodule) or from processing the results of propagation of the rotational equations of motion. Note that when numerically propagating rotational dynamics, this angular velocity vector is typically directly defined in the last three entries of the state vector.

Parameters:

current_time (float) – The time (in seconds since epoch J2000, TDB time scale) at which the angular velocity vector is evaluated

Returns:

Angular velocity vector of the body \(\boldsymbol{\omega}^{(B)}\) expressed in the body-fixed frame \(B\)

Return type:

numpy.ndarray

angular_velocity_in_inertial_frame(self: tudatpy.kernel.numerical_simulation.environment.RotationalEphemeris, time: float) numpy.ndarray[numpy.float64[3, 1]]

Function to get the body’s angular velocity vector, expressed in the inertial frame.

Function to get the body’s angular velocity vector \(\boldsymbol{\omega}^{(I)}\), expressed in the body-fixed frame \(I\). This quantity is computed from \(\mathbf{R}^{I/B}\boldsymbol{\omega}^{(B)}\), see the angular_velocity_in_body_fixed_frame and body_fixed_to_inertial_rotation functions.

Parameters:

current_time (float) – The time (in seconds since epoch J2000, TDB time scale) at which the angular velocity vector is evaluated

Returns:

Angular velocity vector of the body \(\boldsymbol{\omega}^{(B)}\) expressed in the body-fixed frame \(B\)

Return type:

numpy.ndarray

body_fixed_to_inertial_rotation(self: tudatpy.kernel.numerical_simulation.environment.RotationalEphemeris, time: float) numpy.ndarray[numpy.float64[3, 3]]

Function to get rotation matrix from body-fixed frame to inertial frame over time.

Function to get rotation matrix from body-fixed (target) frame to inertial (base) frame over time. The calculation of this rotation matrix depends on the specific rotation model that has been defined, either from an a priori definition (see rotation_model submodule) or from processing the results of propagation of the rotational equations of motion.

Parameters:

current_time (float) – The time (in seconds since epoch J2000, TDB time scale) at which the rotation matrix is evaluated

Returns:

Rotation matrix \(\mathbf{R}^{I/B}\) from body-fixed frame \(B\) to inertial frame I

Return type:

numpy.ndarray

inertial_to_body_fixed_rotation(self: tudatpy.kernel.numerical_simulation.environment.RotationalEphemeris, time: float) numpy.ndarray[numpy.float64[3, 3]]

Function to get rotation matrix from inertial frame to body-fixed frame over time.

Function computes the inverse (equal to transpose) rotation of the body_fixed_to_inertial_rotation function.

Parameters:

current_time (float) – The time (in seconds since epoch J2000, TDB time scale) at which the rotation matrix is evaluated

time_derivative_body_fixed_to_inertial_rotation(self: tudatpy.kernel.numerical_simulation.environment.RotationalEphemeris, time: float) numpy.ndarray[numpy.float64[3, 3]]

Function to get time derivative of rotation matrix from body-fixed frame to inertial frame over time.

Function to get time derivative of rotation matrix from body-fixed frame to inertial frame over time (see body_fixed_to_inertial_rotation), denoted \(\dot{\mathbf{R}}^{(I/B)}\),

Parameters:

current_time (float) – The time (in seconds since epoch J2000, TDB time scale) at which the rotation matrix derivative is evaluated

Returns:

Rotation matrix \(\dot{\mathbf{R}}^{I/B}\) from body-fixed frame \(B\) to inertial frame I

Return type:

numpy.ndarray

time_derivative_inertial_to_body_fixed_rotation(self: tudatpy.kernel.numerical_simulation.environment.RotationalEphemeris, time: float) numpy.ndarray[numpy.float64[3, 3]]

Function to get time derivative of rotation matrix from inertial frame to body-fixed frame over time.

Function to get time derivative of rotation matrix from inertial frame to body-fixed frame over time (see inertial_to_body_fixed_rotation), denoted \(\dot{\mathbf{R}}^{(B/I)}\),

Parameters:

current_time (float) – The time (in seconds since epoch J2000, TDB time scale) at which the rotation matrix derivative is evaluated

Returns:

Rotation matrix \(\dot{\mathbf{R}}^{B/I}\) from inertial frame I to body-fixed frame \(B\)

Return type:

numpy.ndarray

property body_fixed_frame_name

read-only

The identifier of the body-fixed frame, used in other parts of the simulation to identify it.

Type:

str

property inertial_frame_name

read-only

The identifier of the inertial frame, used in other parts of the simulation to identify it.

Type:

str

class GcrsToItrsRotationModel

Object for high-accuracy GCRS<->ITRS rotation.

Object derived from RotationalEphemeris that implements the high-accuracy GCRS<->ITRS rotation as per the IERS 2010 Conventions. The details of the model are described in gcrs_to_itrs() With the exception of \(s'\), the list of angles used to compute the full rotation are computed by an object of type EarthOrientationAnglesCalculator (which can be retrieved from this rotation model through angles_calculator.

property angles_calculator

read-only

Object that computes the Earth rotation angles \(X,Y,s,\theta_{E},x_{p},y_{p}\)

Type:

EarthOrientationAnglesCalculator

class EarthOrientationAnglesCalculator

Object for computing high-accuracy Earth orientation angles

get_gcrs_to_itrs_rotation_angles(self: tudatpy.kernel.numerical_simulation.environment.EarthOrientationAnglesCalculator, epoch: float, time_scale: tudatpy.kernel.astro.time_conversion.TimeScales = <TimeScales.tdb_scale: 2>) tuple[numpy.ndarray[numpy.float64[5, 1]], float]

Function to compute high-accuracy Earth orientation angles

Function to compute high-accuracy Earth orientation angle quantities \(X,Y,s,x_{p},y_{p}\) and UT1 (from which \(\theta_{E}\) is computed) as described in gcrs_to_itrs()

Parameters:
  • epoch (float) – Epoch at which the Earth orientation angles are to be compute

  • time_scale (TimeScales) – Time scale in which the input epoch is given

Returns:

Pair (tuple of size two) with the first entry a list of orientation angles \(X,Y,s,x_{p},y_{p}\) (in that order) and the second entry the current UT1.

Return type:

tuple[list[float],float]

class GravityFieldModel

Object that provides the gravity field of a body

Object (typically stored inside a Body object) that provides the gravity field of a body, typically (but not exclusively) for use in gravitational acceleration and torque models. This base class allows access to the gravitational parameter of the body. Specific derived classes are implemented to provide models for more detailed gravity field models (e.g. spherical harmonics, polyhedron).

get_gravitational_parameter(self: tudatpy.kernel.numerical_simulation.environment.GravityFieldModel) float
property gravitational_parameter

Value of the gravity field’s gravitational parameters \(\mu\)

Type:

float

class SphericalHarmonicsGravityField

Object that provides a spherical harmonic gravity field of a body.

Object (typically stored inside a Body object) that provides a spherical harmonic gravity field of a body, typically (but not exclusively) for use in gravitational acceleration and torque models. This class is derived from GravityFieldModel. This object is typically created using the spherical_harmonic() settings function. If any time variations of the gravity field are provided, an object of the derived class TimeVariableSphericalHarmonicsGravityField is created.

property cosine_coefficients

Matrix with cosine spherical harmonic coefficients \(\bar{C}_{lm}\) (geodesy normalized). Entry \((i,j)\) denotes coefficient at degree \(i\) and order \(j\).

Type:

numpy.ndarray[numpy.float64[l, m]]

property maximum_degree

read-only

Maximum spherical harmonic degree \(l_{max}\) for which coefficients are defined

Type:

int

property maximum_order

read-only

Maximum spherical harmonic order \(m_{max}\) for which coefficients are defined

Type:

int

property reference_radius

read-only

Reference radius \(R\) of the gravity field

Type:

float

property sine_coefficients

Matrix with sine spherical harmonic coefficients \(\bar{S}_{lm}\) (geodesy normalized). Entry \((i,j)\) denotes coefficient at degree \(i\) and order \(j\).

Type:

numpy.ndarray[numpy.float64[l, m]]

class BodyShapeModel

Object that provides a shape model for a natural body.

Object (typically stored inside a Body object) that provides a shape model for a body, for instance to compute the altitude from a body-centered state, or w.r.t. which to place ground stations. This shape model is typically only associated with natural bodies. Shape models for spacecraft (for non-conservative force models) use properties stored inside the VehicleSystems object.

property average_radius

read-only

Average radius of the body, for use in computations that assume a spherical body shape.

Type:

float

class RigidBodyProperties

Object that defines the mass, center of mass, and inertia tensor as a function of time.

Object that defines the mass, center of mass, and inertia tensor as a function of time, typically used for evaluation of torques and non-conservative forces in numerical state propagation. Note that this object does not define properties of a gravity field (it defines the inertial mass rather than the gravitational mass)

update(self: tudatpy.kernel.numerical_simulation.environment.RigidBodyProperties, time: float) None

Function to update the body properties to the current time. This function is called automatically during a propagation loop. In case these properties are not time-dependent (e.g. when using the constant_rigid_body_properties() settings) this function does nothing (since no update is needed).

Parameters:

current_time (float) – Time (in seconds since J2000 TDB) to which this object is to be updated

property current_center_of_mass

Position of the center of mass of the object (in the body-centered, body-fixed frame), as set by the latest call to the update function of this object.

property current_inertia_tensor

Inertia tensor of the object (with axes along those of the body-fixed frame), as set by the latest call to the update function of this object.

property current_mass

Mass of the object, as set by the latest call to the update function of this object.

class AtmosphereModel

Object that provides the atmospheric properties of the body.

Object that provides the atmospheric properties of the body, as a function of altitude, latitude, longitude and time. Depending on the implementation, the this dependence may be limited to altitude-only (e.g. standard atmosphere models). This object can be accessed directly by the user to compute atmospheric properties outside the loop of the propagation by calling one of its member functions. During the propagation, each body undergoing aerodynamic forces has a AtmosphericFlightConditions object associated with it (accessed from a boyd through flight_condition) that links the atmosphere model to the aerodynamic model.

get_density(self: tudatpy.kernel.numerical_simulation.environment.AtmosphereModel, altitude: float, longitude: float, latitude: float, time: float) float

Function to compute the atmospheric freestream density at a given location.

Parameters:
  • altitude (float) – Local altitude above the body surface at which the property is to be computed

  • latitude (float) – Geographic latitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • longitude (float) – Geographic longitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • time (float) – Time (in seconds since J2000 TDB) at which the property is to be computed.

Returns:

Freestream density at the given time and location

Return type:

float

get_number_density(self: tudatpy.kernel.numerical_simulation.environment.AtmosphereModel, species: tudatpy.kernel.numerical_simulation.environment.AtmosphericCompositionSpecies, altitude: float, longitude: float, latitude: float, time: float) float

Function to compute the atmospheric freestream number density of a given specie at a given location.

Parameters:
  • species (AtmosphericCompositionSpecies) – Atmospheric species for which the number density is to be computed

  • altitude (float) – Local altitude above the body surface at which the property is to be computed

  • latitude (float) – Geographic latitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • longitude (float) – Geographic longitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • time (float) – Time (in seconds since J2000 TDB) at which the property is to be computed.

Returns:

Freestream number density of the requested specie at the given time and location

Return type:

float

get_pressure(self: tudatpy.kernel.numerical_simulation.environment.AtmosphereModel, altitude: float, longitude: float, latitude: float, time: float) float

Function to compute the atmospheric freestream static pressure at a given location.

Parameters:
  • altitude (float) – Local altitude above the body surface at which the property is to be computed

  • latitude (float) – Geographic latitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • longitude (float) – Geographic longitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • time (float) – Time (in seconds since J2000 TDB) at which the property is to be computed.

Returns:

Freestream static pressure at the given time and location

Return type:

float

get_speed_of_sound(self: tudatpy.kernel.numerical_simulation.environment.AtmosphereModel, altitude: float, longitude: float, latitude: float, time: float) float

Function to compute the atmospheric freestream speed of sound at a given location.

Parameters:
  • altitude (float) – Local altitude above the body surface at which the property is to be computed

  • latitude (float) – Geographic latitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • longitude (float) – Geographic longitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • time (float) – Time (in seconds since J2000 TDB) at which the property is to be computed.

Returns:

Freestream speed of sound at the given time and location

Return type:

float

get_temperature(self: tudatpy.kernel.numerical_simulation.environment.AtmosphereModel, altitude: float, longitude: float, latitude: float, time: float) float

Function to compute the atmospheric freestream temperature at a given location.

Parameters:
  • altitude (float) – Local altitude above the body surface at which the property is to be computed

  • latitude (float) – Geographic latitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • longitude (float) – Geographic longitude (in the body-fixed frame of the body with the atmosphere) at which the property is to be computed

  • time (float) – Time (in seconds since J2000 TDB) at which the property is to be computed.

Returns:

Freestream temperature at the given time and location

Return type:

float

class AerodynamicCoefficientInterface

Base class for computing the current aerodynamic coefficients of the body

Base class for computing the current aerodynamic coefficients of the body. The implementation of the computation depends on the choice of aerodynamic coefficient model (see aerodynamic_coefficients for available options). During the propagation, this object is automatically updated to the current state by the AtmosphericFlightConditions object. The user may override the current aerodynamic coefficients when using, for instance, a custom aerodynamic guidance model (see here for an example). using the member functions of this class.

current_control_surface_force_coefficient_increment(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicCoefficientInterface, control_surface_name: str) numpy.ndarray[numpy.float64[3, 1]]

Function to get the contribution from a single control surface to the aerodynamic force coefficient, as compute by last call to update_full_coefficients()

Parameters:

control_surface_name (str) – The name of the control surface for which the contribution is to be retrieved

Returns:

Contribution from the requested control surface to the aerodynamic force coefficient

Return type:

numpy.ndarray

current_control_surface_moment_coefficient_increment(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicCoefficientInterface, control_surface_name: str) numpy.ndarray[numpy.float64[3, 1]]

Function to get the contribution from a single control surface to the aerodynamic moment coefficients, as compute by last call to update_full_coefficients()

Parameters:

control_surface_name (str) – The name of the control surface for which the contribution is to be retrieved

Returns:

Contribution from the requested control surface to the aerodynamic moment coefficients

Return type:

numpy.ndarray

set_control_surface_increments(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicCoefficientInterface, control_surface_list: dict[str, tudat::aerodynamics::ControlSurfaceIncrementAerodynamicInterface]) None

No documentation found.

update_coefficients(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicCoefficientInterface, independent_variables: list[float], time: float) None

Function to update the aerodynamic coefficients of the body only

Function to update the aerodynamic coefficients of the body only (without the control surface contribution), based on the current state. This function may be called by the user, but will set only the current_force_coefficients and current_moment_coefficients (while leaving the current_control_surface_free_force_coefficients and current_control_surface_free_moment_coefficients unchanged)

Parameters:
  • independent_variables (list[float]) – List of inputs from which the aerodynamic coefficients are to be computed, with each entry corresponding to the value of the physical variable defined by the independent_variable_names attribute.

  • time (float) – Current time (in seconds since J2000)

Returns:

Contribution from the requested control surface to the aerodynamic moment coefficients

Return type:

numpy.ndarray

update_full_coefficients(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicCoefficientInterface, independent_variables: list[float], control_surface_independent_variables: dict[str, list[float]], time: float, check_force_contribution: bool = True) None

Function to update the aerodynamic coefficients, from both the body and its control surfaces

Function to update the aerodynamic coefficients of both the body and its control surfaces, based on the current state. This function will call the update_coefficients() function to update the body coefficients. This function may be called by the user, and will set the following attributes: current_force_coefficients, current_moment_coefficients , current_control_surface_free_force_coefficients and current_control_surface_free_moment_coefficients. In addition, it will modify the coefficients returned by the current_control_surface_force_coefficient_increment() and current_control_surface_moment_coefficient_increment() functions

Parameters:
  • independent_variables (list[float]) – List of inputs from which the aerodynamic coefficients of the body are to be computed, with each entry corresponding to the value of the physical variable defined by the independent_variable_names attribute.

  • control_surface_independent_variables (dict[str,list[float]]) – List of inputs from which the control surface aerodynamic coefficients are to be computed (with dictionary key the control surface name), with each entry corresponding to the value of the physical variable defined by the control_surface_independent_variable_names attribute.

  • time (float) – Current time (in seconds since J2000)

  • check_force_contribution (bool, default = True) – Boolean that determines if the force contribution to the aerodynamic moments should be added. Note that this input is only used if the add_force_contribution_to_moments attribute is set to True.

property control_surface_independent_variable_names

read-only

List of independent variables from which the aerodynamic coefficients of each control surface are computed, with dictionary key being the control surface name (e.g. required input to update_full_coefficients() function).

Type:

dict[str,list[AerodynamicCoefficientsIndependentVariables]]

property current_coefficients

read-only

Concatenation of current_force_coefficients and current_moment_coefficients

Type:

numpy.ndarray

property current_control_surface_free_force_coefficients

read-only

Same as current_force_coefficients, but without contribution (if any) from control surfaces

Type:

numpy.ndarray

property current_control_surface_free_moment_coefficients

read-only

Same as current_moment_coefficients, but without contribution (if any) from control surfaces

Type:

numpy.ndarray

property current_force_coefficients

read-only

The current aerodynamic force coefficients, in the frame defined by the force_coefficient_frame attribute, as computed by the last call to the update_coefficients() function.

Type:

numpy.ndarray

property current_moment_coefficients

read-only

The current aerodynamic moment coefficients, in the frame defined by the moment_coefficient_frame attribute, as computed by the last call to the update_coefficients() function.

Type:

numpy.ndarray

property force_coefficient_frame

read-only

Reference frame in which the current_force_coefficients are defined

Type:

AerodynamicCoefficientFrames

property independent_variable_names

read-only

List of independent variables from which the aerodynamic coefficients are computed (e.g. required input to update_coefficients() function).

Type:

list[AerodynamicCoefficientsIndependentVariables]

property moment_coefficient_frame

read-only

Reference frame in which the current_moment_coefficients are defined

Type:

AerodynamicCoefficientFrames

property reference_area

read-only

The aerodynamic reference area \(A\) of the coefficients

Type:

float

class HypersonicLocalInclinationAnalysis
__init__(self: tudatpy.kernel.numerical_simulation.environment.HypersonicLocalInclinationAnalysis, independent_variable_points: list[list[float]], body_shape: tudatpy.kernel.math.geometry.SurfaceGeometry, number_of_lines: list[int], number_of_points: list[int], invert_orders: list[bool], selected_methods: list[list[int]], reference_area: float, reference_length: float, moment_reference_point: numpy.ndarray[numpy.float64[3, 1]], save_pressure_coefficients: bool = False) None

Class constructor, taking the shape of the vehicle, and various analysis options as input.

Parameters:
  • independent_variable_points (list[list[float]]) – List containing three lists, with each sublist containing the data points of each of the independent variables for the coefficient generation. The physical meaning of each of the three independent variables is: 0 = mach number, 1 = angle of attack, 2 = angle of sideslip. Each of the subvectors must be sorted in ascending order.

  • body_shape (SurfaceGeometry) – Class that defines the shape of the vehicle as a continuous surface. The local inclination analysis discretizes the surface of the vehicle into quadrilateral panels, defined by the other inputs to this constructor. In case the tudat.geometry.SurfaceGeometry object is made up of multiple sub-shapes, different settings may be used for each

  • number_of_lines (List[ float ]) – Number of discretization points in the first independent surface variable of each of the subparts of body_shape. The size of this list should match the number of parts of which the body_shape is composed. The first independent variable of a subpart typically runs along the longitudinal vehicle direction

  • number_of_points (List[ float ]) – Number of discretization points in the second independent surface variable of each of the subparts of body_shape. The size of this list should match the number of parts of which the body_shape is composed. The first independent variable of a subpart typically runs along the lateral vehicle direction

  • invert_orders (List[ bool ]) – Booleans to denote whether the surface normals of the panels of each discretized body_shape subpart are to be inverted (i.e. inward-facing->outward facing or vice versa). The size of this list should match the number of parts of which the body_shape is composed.

  • selected_methods (List[ List[ int ] ]) – Double list of selected local inclination methods, the first index (outer list) represents compression or expansion (0 and 1), the second index (inner list) denotes the vehicle part index. The size of this inner list should match the number of parts of which the body_shape is composed. The int defining the method type is interpreted as follows. For the compression methods, the following are available: * 0: Newtonian Method. * 1: Modified Newtonian. * 2 and 3: not available at this moment. * 4: Tangent-wedge method. * 5: Tangent-cone method. * 6: Modified Dahlem-Buck method. * 7: VanDyke unified pressure method. * 8: Smyth Delta Wing method. * 9: Hankey flat surface method The expansion method has the following options: * 0: Vacuum Pressure coefficient method. * 1: Zero Pressure function. * 4: High Mach base pressure method. * 3 or 5: Prandtl-Meyer method. * 6: ACM empirical pressure coefficient.

  • reference_area (float) – Reference area used to non-dimensionalize aerodynamic forces and moments.

  • moment_reference_point (numpy.ndarray) – Reference point wrt which aerodynamic moments are calculated.

  • save_pressure_coefficients (bool) – Boolean denoting whether to save the pressure coefficients that are computed to files

clear_data(self: tudatpy.kernel.numerical_simulation.environment.HypersonicLocalInclinationAnalysis) None
class GroundStation

Object used to define and store properties of a ground station.

Object (typically stored inside a Body object) used to define and store properties of a ground station, typically used in modelling tracking observations to/from a ground station.

set_pressure_function(self: tudatpy.kernel.numerical_simulation.environment.GroundStation, pressure_function: Callable[[float], float]) None
set_relative_humidity_function(self: tudatpy.kernel.numerical_simulation.environment.GroundStation, relative_humidity_function: Callable[[float], float]) None
set_temperature_function(self: tudatpy.kernel.numerical_simulation.environment.GroundStation, temperature_function: Callable[[float], float]) None
set_timing_system(self: tudatpy.kernel.numerical_simulation.environment.GroundStation, timing_system: tudatpy.kernel.numerical_simulation.environment.TimingSystem) None
set_transmitting_frequency_calculator(self: tudatpy.kernel.numerical_simulation.environment.GroundStation, transmitting_frequency_calculator: tudat::ground_stations::StationFrequencyInterpolator) None
set_water_vapor_partial_pressure_function(self: tudatpy.kernel.numerical_simulation.environment.GroundStation, water_vapor_partial_pressure_function: Callable[[float], float]) None
property pointing_angles_calculator

read-only

Object that performs computations of the azimuth and elevation of an arbitrary target as observed by the ground station

Type:

PointingAnglesCalculator

property pressure_function

Function that provides the local pressure at the ground station (typically use for media corrections) as a function of time

Type:
type:

Callable[[float], float]

property relative_humidity_function

Function that provides the local relative humidity at the ground station (typically use for media corrections) as a function of time

Type:
type:

Callable[[float], float]

property station_state

read-only

Object that performs computations of the current (body-fixed) position and frame conversions of the ground station.

Type:

GroundStationState

property temperature_function

Function that provides the local temperature at the ground station (typically use for media corrections) as a function of time

Type:
type:

Callable[[float], float]

property transmitting_frequency_calculator

Object that provides the transmission frequency as a function of time for (radio) tracking stations. This attribute is typically set automatically when loading tracking data files (e.g. ODF, IFMS, TNF, etc.)

Type:

TransmittingFrequencyCalculator

class GroundStationState

Object that performs computations of the current (body-fixed) position and frame conversions of the ground station.

Object that performs computations of the current (body-fixed) position and frame conversions of the ground station. In the simplest situation, only a Cartesian position is provided, which is then assumed constant. If time variations (for instance due to tides or plate motion) are present, their impact on station position is computed in this object.

get_cartesian_position(self: tudatpy.kernel.numerical_simulation.environment.GroundStationState, current_time: float, target_frame_origin: str = '') numpy.ndarray[numpy.float64[3, 1]]

This function computes the position of the station as a function of time.

This function computes the position of the station as a function of time, in a frame with body-fixed orientation. Some time-variations of the station position depend on the origin of the frame in which the computation is to be used. For instance, relativistic correction to the Earth-fixed position is different in a geocentric or barycentric frame. However, the output of this function is always given in the body-fixed, body-centered frame.

Parameters:
  • current_time (float) – Time (in seconds since J2000 TDB) at which the position is to be computed.

  • target_frame_origin (str, default = "") – Identifier for the frame origin w.r.t. which the computed position is to be used.

Returns:

Cartesian position of the station at the current epoch, in a body-centered, body-fixed frame

Return type:

numpy.ndarray

get_cartesian_state(self: tudatpy.kernel.numerical_simulation.environment.GroundStationState, current_time: float, target_frame_origin: str = '') numpy.ndarray[numpy.float64[6, 1]]
property cartesian_positon_at_reference_epoch

read-only

Cartesian position of the ground station, at the reference epoch, in a body-fixed, body-centered frame.

Type:

numpy.ndarray[numpy.float64[3, 1]]

property geodetic_position_at_reference_epoch

read-only

Geodetic position of the ground station (altitude w.r.t. body shape model, geodetic latitude, longitude), at the reference epoch, in a body-fixed, body-centered frame.

Type:

numpy.ndarray[numpy.float64[3, 1]]

property rotation_matrix_body_fixed_to_topocentric

read-only

Rotation matrix from the body-fixed frame (of the station’s central body) to the topocentric frame of the ground station. The body-fixed frame is defined by the rotation model of the body object (rotation_model). The axes of the topocentric frame are defined such that the x-axis is in East direction, the z-direction is upwards, perpendicular to the body’s surface sphere (with properties defined by the central body’s shape model shape_model). The y-axis completes the frame, and is in northern direction. For time-varying ground station positions, this function uses the station position at reference epoch for the computation of the axes.

Type:

numpy.ndarray[numpy.float64[3, 3]]

property spherical_positon_at_reference_epoch

read-only

Spherical position of the ground station (distance w.r.t. body center, latitude, longitude), at the reference epoch, in a body-fixed, body-centered frame.

Type:

numpy.ndarray[numpy.float64[3, 1]]

class FlightConditions

Object that calculates various state-derived quantities typically relevant for flight dynamics.

Object that calculates various state-derived quantities typically relevant for flight dynamics, such as latitude, longitude, altitude, etc. It also contains an AerodynamicAngleCalculator that computes derived angles (flight path, heading angle, etc.). This object is limited to non-atmospheric flight. For flight through Body objects endowed with an atmosphere model, the derived class AtmosphericFlightConditions is used. This object is stored inside a Body object, and represents the flight conditions of a single body w.r.t. a single central body.

update_conditions(self: tudatpy.kernel.numerical_simulation.environment.FlightConditions, current_time: float) None
property aerodynamic_angle_calculator

read-only

The object that is responsible for computing various relevant flight dynamics angles and frame rotations.

Type:

AerodynamicAngleCalculator

property altitude

read-only

The current time, at which this object was last updated

Type:

float

property body_centered_body_fixed_state

read-only

Cartesian translational state, expressed in a frame centered at, and fixed to, the central body. Note that, due to the rotation of the central body, the norm of the body-fixed, body-centered, velocity differs from the norm of the inertial body-centered velocity.

Type:

numpy.ndarray

property geodetic_latitude

read-only

The body-fixed geographic latitude of the body w.r.t. its central body.

Type:

float

property latitude

read-only

The body-fixed geographic latitude of the body w.r.t. its central body.

Type:

float

property longitude

read-only

The body-fixed longitude of the body w.r.t. its central body.

Type:

float

property time

read-only

The current time, at which this object was last updated

Type:

float

class AtmosphericFlightConditions

Object that calculates various state-derived quantities typically relevant for flight dynamics, for flight in an atmosphere.

Object that calculates various state-derived quantities typically relevant for flight dynamics, for flight in an atmosphere, such as latitude, longitude, altitude, density, Mach number etc. It also contains an AerodynamicAngleCalculator that computes derived angles (flight path, heading angle, etc.). This object is derived from FlightConditions, which performs computations for non-atmospheric flight only. This object is stored inside a Body object, and represents the flight conditions of a single body w.r.t. a single central body.

property aero_coefficient_independent_variables

read-only

List of current values of independent variables of aerodynamic coefficients. This list is only defined if the body has an AerodynamicCoefficientInterface that has dependencies on environmental variables (e.g. Mach number, angle of attack, etc.).

Type:

numpy.ndarray

property aerodynamic_coefficient_interface

read-only

Object extracted from the same Body object as this AtmosphericFlightConditions object, which defines the aerodynamic coefficients.

Type:

AerodynamicCoefficientInterface

property airspeed

read-only

The airspeed of the body w.r.t. the atmosphere.

Type:

float

property airspeed_velocity

read-only

The velocity vector of the body w.r.t. the freestream atmosphere (e.g. vectorial counterpart of airspeed).

Type:

numpy.ndarray

property control_surface_aero_coefficient_independent_variables

read-only

List of lists current values of independent variables of aerodynamic coefficients for control surfaces. The outer list defines the control surface, the inner list the values of the independent variables. This list is only defined if the body has an AerodynamicCoefficientInterface with control surfaces that have dependencies on environmental variables (e.g. Mach number, angle of attack, etc.).

Type:

numpy.ndarray

property density

read-only

The freestream atmospheric density at the body’s current location.

Type:

float

property dynamic_pressure

read-only

The freestream atmospheric dynamic pressure at the body’s current location.

Type:

float

property mach_number

read-only

The freestream Mach number of the body.

Type:

float

property pressure

read-only

The freestream atmospheric static pressure at the body’s current location.

Type:

float

property speed_of_sound

read-only

The freestream atmospheric speed of sound at the body’s current location.

Type:

float

property temperature

read-only

The freestream atmospheric temperature at the body’s current location.

Type:

float

class AerodynamicAngleCalculator

Object to calculate (aerodynamic) orientation angles, and frame transformations, from current vehicle state.

Object to calculate (aerodynamic) orientation angles (list given by the AerodynamicsReferenceFrameAngles enum) and transformations between frames (list given by the AerodynamicsReferenceFrames enum) from current vehicle state.

get_angle(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicAngleCalculator, angle_type: tudatpy.kernel.numerical_simulation.environment.AerodynamicsReferenceFrameAngles) float

Function to get a single orientation angle

Function to get a single orientation angle. This function is meant to be used only during a numerical propagation, in particular for the definition of a custom (e.g. guidance) model.

Parameters:

original_frame (AerodynamicsReferenceFrameAngles) – The identifier for the angle that is to be returned

Returns:

Value of requested angle

Return type:

double

get_rotation_matrix_between_frames(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicAngleCalculator, original_frame: tudatpy.kernel.numerical_simulation.environment.AerodynamicsReferenceFrames, target_frame: tudatpy.kernel.numerical_simulation.environment.AerodynamicsReferenceFrames) numpy.ndarray[numpy.float64[3, 3]]

Function to get the rotation matrix between two frames.

Function to get the rotation matrix between two frames. This function is meant to be used only during a numerical propagation, in particular for the definition of a custom (e.g. guidance) model.

Parameters:
Returns:

Rotation matrix \(\mathbf{R}^{B/A}\) from frame \(A\) to frame B

Return type:

numpy.ndarray

set_body_orientation_angle_functions(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicAngleCalculator, angle_of_attack_function: Callable[[], float] = None, angle_of_sideslip_function: Callable[[], float] = None, bank_angle_function: Callable[[], float] = None, angle_update_function: Callable[[float], None] = None, silence_warnings: bool = False) None
set_body_orientation_angles(self: tudatpy.kernel.numerical_simulation.environment.AerodynamicAngleCalculator, angle_of_attack: float = nan, angle_of_sideslip: float = nan, bank_angle: float = nan, silence_warnings: bool = False) None
class VehicleSystems

Object used to store physical (hardware) properties of a vehicle.

get_control_surface_deflection(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems, control_surface_id: str) float

Function to retrieve the current deflection of an aerodynamic control surface.

Function to retrieve the current deflection of an aerodynamic control surface, identified by its name. To extract the control surface deflection, the control surface has to exist. A control surface is created whenever control surfaces are defined in a body’s aerodynamic coefficient interface.

Parameters:

control_surface_id (str) – The identified (name) of the given control surface

Returns:

Current deflection (in radians) that the control surface

Return type:

float

get_engine_model(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems, engine_name: str) tudat::system_models::EngineModel

Function to retrieve an engine model from the vehicle

Parameters:

engine_name (str) – The identifier for the engine model that is to be retrieved

Returns:

Model for the engine that is requested

Return type:

EngineModel

set_control_surface_deflection(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems, control_surface_id: str, deflection_angle: float) None

Function to set the current deflection of an aerodynamic control surface.

Function to set the current deflection of an aerodynamic control surface, identified by its name. To set the control surface deflection, the control surface has to exist. A control surface is created whenever control surfaces are defined in a body’s aerodynamic coefficient interface.

Parameters:
  • control_surface_id (str) – The identified (name) of the given control surface

  • deflection_angle (float) – The deflection (in radians) that the control surface is to be set to. This will typically influence the aerodynamic coefficients of the vehicle

set_default_transponder_turnaround_ratio_function(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems) None

No documentation found.

set_reference_point(*args, **kwargs)

Overloaded function.

  1. set_reference_point(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems, reference_point: str, location: numpy.ndarray[numpy.float64[3, 1]], frame_origin: str = ‘’, frame_orientation: str = ‘’) -> None

No documentation found.

  1. set_reference_point(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems, reference_point: str, ephemeris: tudat::ephemerides::Ephemeris) -> None

No documentation found.

set_timing_system(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems, timing_system: tudat::system_models::TimingSystem) None
set_transponder_turnaround_ratio(self: tudatpy.kernel.numerical_simulation.environment.VehicleSystems, transponder_ratio_per_uplink_and_downlink_frequency_band: dict[tuple[tudat::observation_models::FrequencyBands, tudat::observation_models::FrequencyBands], float]) None

No documentation found.

class EngineModel
class Body

Object that stores the environment properties and current state of a single body.

Object that stores the environment properties and current state of a single celestial body (natural or artificial). Each separate environment model (gravity field, ephemeris, etc.) is stored as a member object in this class. During each time step, the Body gets updated to the current time/propagated state, and the current properties, in as much as they are time-dependent, can be extracted from this object

get_ground_station(self: tudatpy.kernel.numerical_simulation.environment.Body, station_name: str) tudatpy.kernel.numerical_simulation.environment.GroundStation

This function extracts a ground station object from the body.

This function extracts a ground station object, for a station of a given name, from the body. If no station of this name exists, an exception is thrown

Parameters:

station_name (str) – Name of the ground station that is to be retrieved.

Returns:

Ground station object of the station of requested name

Return type:

GroundStation

state_in_base_frame_from_ephemeris(self: tudatpy.kernel.numerical_simulation.environment.Body, time: float) numpy.ndarray[numpy.float64[6, 1]]

This function returns the body’s state, as computed from its ephemeris model (extracted from ephemeris) at the current time, and (if needed) translates this state to the global frame origin. For the case where the origin of the body’s ephemeris (extracted from frame_origin) is equal to the global frame origin of the system of bodies it is in (extracted from SystemOfBodies.global_frame_origin), this function is equal to Body.ephemeris.cartesian_state( time ). Where the global frame origin and ephemeris origin is not equal, other bodies’ ephemerides are queried as needed to provide this body’s state w.r.t. the global frame origin

Parameters:

time (float) – Time (in TDB seconds since J2000) at which the state is to be computed

Returns:

Cartesian state (position and velocity) of the body w.r.t. the global frame origin at the requested time.

Return type:

numpy.ndarray

property aerodynamic_coefficient_interface

Object defining the aerodynamic coefficients of the body (force-only, or force and moment) as a function of any number of independent variables. Depending on the selected type of model, the type of this attribute is of type AerodynamicCoefficientInterface, or a derived class thereof.

Type:

AerodynamicCoefficientInterface

property atmosphere_model

Object defining the atmosphere model of this body, used to calculate density, temperature, etc. at a given state/time. Depending on the selected type of model, the type of this attribute is of type AtmosphereModel, or a derived class thereof.

Type:

AtmosphereModel

property body_fixed_angular_velocity

read-only

Angular velocity vector of the body, expressed in body-fixed frame (see inertial_to_body_fixed_frame).

Type:

numpy.ndarray

property body_fixed_to_inertial_frame

read-only

The rotation from this Body’s body-fixed frame to inertial frame (see inertial_to_body_fixed_frame).

Type:

numpy.ndarray

property body_fixed_to_inertial_frame_derivative

read-only

Time derivative of rotation matrix from this Body’s body-fixed frame to inertial frame (see inertial_to_body_fixed_frame).

Type:

numpy.ndarray

property ephemeris

Object defining the ephemeris model of this body, used to calculate its current state as a function of time. Depending on the selected type of model, the type of this attribute is of type Ephemeris, or a derived class thereof.

Type:

Ephemeris

property flight_conditions

Object used to calculated and store the current flight conditions of a vehicle (altitude, latitude, longitude, flight-path angle, etc.) w.r.t. a central body. In case the central body contains an atmosphere, this object also stores current local density, Mach number, etc. This object is typically used for aerodynamic accelerations, guidance models or other central-body-related custom models.

Type:

FlightConditions

property gravitational_parameter

read-only

Attribute of convenience, equivalent to .gravity_field_model.gravitational_parameter

Type:

float

property gravity_field_model

Object defining the a gravity field model of this body, used to define the exterior gravitational potential, and its gradient(s). Depending on the selected type of model, the type of this attribute is of type GravityFieldModel, or a derived class thereof.

Type:

GravityFieldModel

property ground_station_list

Dictionary of all ground stations that exist in the body, with dictionary key being the name of the station, and the ground station object the key of the dictionary.

Type:

dict[str,GroundStation]

property inertia_tensor

The current inertia tensor \(\mathbf{I}\) of the vehicle, as used in the calculation of (for instance) the response to torques. This attribute is a shorthand for accessing the inertia tensor as computed/stored in the rigid_body_properties attribute. For certain types of rigid-body properties, this attribute cannot be used to (re)set the current mass.

Unlike the attributes containing the state, orientation, angular velocity of the Body, this attribute may be used to retrieve the state during the propagation and to define the mass of a vehicle.

Type:

numpy.ndarray

property inertial_angular_velocity

read-only

Angular velocity vector of the body, expressed in inertial frame (see inertial_to_body_fixed_frame).

Type:

numpy.ndarray

property inertial_to_body_fixed_frame

read-only

The rotation from inertial frame (with global frame orientation) to this Body’s body-fixed frame. The rotation is always returned here as a rotation matrix. If the body’s rotational state is numerically propagated, this property gets extracted from the propagated state vector. If it is not propagated, the state is extracted from this body’s rotational ephemeris.

Note

This function is only valid during the numerical propagation if any aspects of the dynamics or dependent variables require the body’s rotational state.

Type:

numpy.ndarray

property inertial_to_body_fixed_frame_derivative

read-only

Time derivative of rotation matrix from inertial frame to this Body’s body-fixed frame (see inertial_to_body_fixed_frame).

Type:

numpy.ndarray

property mass

The current mass \(m\) of the vehicle, as used in the calculation of non-conservative acceleration. This attribute is a shorthand for accessing the mass as computed/stored in the Body.rigid_body_properties attribute. For certain types of rigid-body properties, this attribute cannot be used to (re)set the current mass. If the body has no Body.rigid_body_properties, and this function is used to set a mass, a new object is automatically created, with settings analogous to the the constant_rigid_body_properties() setting.

Unlike the attributes containing the state, orientation, angular velocity of the Body, this attribute may be used to retrieve the state during the propagation and to define the mass of a vehicle.

Type:

float

property position

read-only

The translational position of the Body, as set during the current step of the numerical propagation (see state).

Type:

numpy.ndarray

property rigid_body_properties

Object defining the mass, center of mass and inertia tensor of the body. This object is distinct from the gravity field of a body (defined by the Body.gravity_field object). A body endowed with this property does not automatically have a gravity field created for it. However, the whenever a body is endowed with a gravity field, a rigid body properties attribute is created to be consistent with this gravity field (e.g. for a spherical harmonic gravity field the mass, center of mass and inertia tensor are created from the gravitational parameter, degree-1 coefficients, and degree-2 coefficients plus mean moment of inertia, respectively).

Type:

RigidBodyProperties

property rotation_model

Object defining the orientation of the body, used to calculate the rotation to/from a body-fixed frame (and its derivate). Depending on the selected type of model, the type of this attribute is of type RotationalEphemeris, or a derived class thereof.

Type:

RotationalEphemeris

property shape_model

Object defining the a shape model of this body, used to define the exterior shape of the body, for instance for the calculation of vehicle’s altitude. Depending on the selected type of model, the type of this attribute is of type BodyShapeModel, or a derived class thereof.

Type:

BodyShapeModel

property state

read-only

The translational state of the Body, as set during the current step of the numerical propagation. The translational state stored here is always in Cartesian elements, w.r.t. the global frame origin, with axes along the global frame orientation. If the body’s translational state is numerically propagated, this property gets extracted from the propagated state vector. If it is not propagated, the state is extracted from this body’s ephemeris. In both cases, any required state transformations are automatically applied. Note that this function is only valid during the numerical propagation if any aspects of the dynamics or dependent variables require the body’s state.

Type:

numpy.ndarray

property system_models

Object used to store physical (hardware) properties of a vehicle, such as engines, control surfaces, etc. This object is typically created automatically whenever such a hardware model needs to be assigned to a vehicle.

Type:

VehicleSystems

property velocity

read-only

The translational velocity of the Body, as set during the current step of the numerical propagation (see state).

Type:

numpy.ndarray

class SystemOfBodies

Object that contains a set of Body objects and associated frame information.

Object that contains a set of Body objects and associated frame information. This object stored the entire environment for a typical Tudat numerical simulation, and is fundamental for the overall Tudat architecture.

add_body(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies, body_to_add: tudatpy.kernel.numerical_simulation.environment.Body, body_name: str, process_body: bool = 1) None

This function adds an existing body, which the user has separately created, to the SystemOfBodies.

Parameters:
  • body_to_add (Body) – Body object that is to be added.

  • body_name (numpy.ndarray) – Name of the Body that is to be added.

  • process_body (bool, default=True) –

    Variable that defines whether this new Body will have its global frame origin/orientation set to conform to rest of the environment.

    Warning

    Only in very rare cases should this variable be anything other than True. Users are recommended to keep this default value intact.

create_empty_body(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies, body_name: str, process_body: bool = 1) None

This function creates a new empty body.

This function creates a new empty body, and adds it to the SystemOfBodies. Since the body is empty, it will not have any environment models defined. These must all be added manually by a user.

Parameters:
  • body_name (string) – Name of the Body that is to be added

  • process_body (bool, default=True) –

    Variable that defines whether this new Body will have its global frame origin/orientation set to conform to rest of the environment.

    Warning

    Only in very rare cases should this variable be anything other than True. Users are recommended to keep this default value intact.

Examples

This function is often used early on in the environment creation segment of a simulation, following the creation of a SystemOfBodies from the default settings for celestial bodies.

# Define string names for bodies to be created from default.
bodies_to_create = ["Sun", "Earth", "Moon", "Mars", "Venus"]

# Use "Earth"/"J2000" as global frame origin and orientation.
global_frame_origin = "Earth"
global_frame_orientation = "J2000"

# Create default body settings, usually from `spice`.
body_settings = environment_setup.get_default_body_settings(
    bodies_to_create,
    global_frame_origin,
    global_frame_orientation)

# Create system of selected celestial bodies
bodies = environment_setup.create_system_of_bodies(body_settings)

# Create vehicle objects.
bodies.create_empty_body("Delfi-C3")
does_body_exist(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies, body_name: str) bool

Function to check if a body with a given name exists in the SystemOfBodies

Parameters:

body_name (string) – Name of the Body whose existence is to be checked

Returns:

True if the body exists in this object, false if not

Return type:

bool

get(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies, body_name: str) tudatpy.kernel.numerical_simulation.environment.Body

This function extracts a single Body object from the SystemOfBodies.

Parameters:

body_name (str) – Name of the Body that is to be retrieved.

Returns:

Body object of the requested name

Return type:

Body

get_body(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies, body_name: str) tudatpy.kernel.numerical_simulation.environment.Body

Deprecated version of get()

global_frame_orientation(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies) str

Common global frame orientation for all bodies in this SystemOfBodies, described in more detail here.

global_frame_origin(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies) str

Common global frame origin for all bodies in this SystemOfBodies, described in more detail here.

list_of_bodies(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies) list[str]

List of names of bodies that are stored in this SystemOfBodies

remove_body(self: tudatpy.kernel.numerical_simulation.environment.SystemOfBodies, body_name: str) None

This function removes an existing body from the SystemOfBodies.

Warning

This function does not necessarily delete the Body object, it only removes it from this object. If any existing models in the simulation refer to this Body, it will persist in memory.

Parameters:

body_name (numpy.ndarray) – Name of the Body that is to be removed.