yt.data_objects.static_output module

Dataset and related data structures.

class yt.data_objects.static_output.Dataset(filename, dataset_type=None, file_style=None, units_override=None, unit_system='cgs')[source]

Bases: object

add_deposited_particle_field(deposit_field, method, kernel_name='cubic', weight_field='particle_mass')[source]

Add a new deposited particle field

Creates a new deposited field based on the particle deposit_field.

Parameters:
  • deposit_field (tuple) – The field name tuple of the particle field the deposited field will be created from. This must be a field name tuple so yt can appropriately infer the correct particle type.
  • method (string) – This is the “method name” which will be looked up in the particle_deposit namespace as methodname_deposit. Current methods include simple_smooth, sum, std, cic, weighted_mean, mesh_id, and nearest.
  • kernel_name (string, default 'cubic') – This is the name of the smoothing kernel to use. It is only used for the simple_smooth method and is otherwise ignored. Current supported kernel names include cubic, quartic, quintic, wendland2, wendland4, and wendland6.
  • weight_field (string, default 'particle_mass') – Weighting field name for deposition method weighted_mean.
Returns:

Return type:

The field name tuple for the newly created field.

add_field(name, function=None, sampling_type=None, **kwargs)[source]

Dataset-specific call to add_field

Add a new field, along with supplemental metadata, to the list of available fields. This respects a number of arguments, all of which are passed on to the constructor for DerivedField.

Parameters:
  • name (str) – is the name of the field.
  • function (callable) – A function handle that defines the field. Should accept arguments (field, data)
  • units (str) – A plain text string encoding the unit. Powers must be in python syntax (** instead of ^).
  • take_log (bool) – Describes whether the field should be logged
  • validators (list) – A list of FieldValidator objects
  • particle_type (bool) – Is this a particle (1D) field?
  • vector_field (bool) – Describes the dimensionality of the field. Currently unused.
  • display_name (str) – A name used in the plots
  • force_override (bool) – Whether to override an existing derived field. Does not work with on-disk fields.
add_gradient_fields(input_field)[source]

Add gradient fields.

Creates four new grid-based fields that represent the components of the gradient of an existing field, plus an extra field for the magnitude of the gradient. Currently only supported in Cartesian geometries. The gradient is computed using second-order centered differences.

Parameters:input_field (tuple) – The field name tuple of the particle field the deposited field will be created from. This must be a field name tuple so yt can appropriately infer the correct field type.
Returns:
Return type:A list of field name tuples for the newly created fields.

Examples

>>> grad_fields = ds.add_gradient_fields(("gas","temperature"))
>>> print(grad_fields)
[('gas', 'temperature_gradient_x'),
 ('gas', 'temperature_gradient_y'),
 ('gas', 'temperature_gradient_z'),
 ('gas', 'temperature_gradient_magnitude')]
add_particle_filter(filter)[source]
add_particle_union(union)[source]
add_smoothed_particle_field(smooth_field, method='volume_weighted', nneighbors=64, kernel_name='cubic')[source]

Add a new smoothed particle field

Creates a new smoothed field based on the particle smooth_field.

Parameters:
  • smooth_field (tuple) – The field name tuple of the particle field the smoothed field will be created from. This must be a field name tuple so yt can appropriately infer the correct particle type.
  • method (string, default 'volume_weighted') – The particle smoothing method to use. Can only be ‘volume_weighted’ for now.
  • nneighbors (int, default 64) – The number of neighbors to examine during the process.
  • kernel_name (string, default cubic) – This is the name of the smoothing kernel to use. Current supported kernel names include cubic, quartic, quintic, wendland2, wendland4, and wendland6.
Returns:

Return type:

The field name tuple for the newly created field.

all_data(find_max=False, **kwargs)[source]

all_data is a wrapper to the Region object for creating a region which covers the entire simulation domain.

arr

Converts an array into a yt.units.yt_array.YTArray

The returned YTArray will be dimensionless by default, but can be cast to arbitrary units using the input_units keyword argument.

Parameters:
  • input_array (Iterable) – A tuple, list, or array to attach units to
  • input_units (String unit specification, unit symbol or astropy object) – The units of the array. Powers must be specified using python syntax (cm**3, not cm^3).
  • dtype (string or NumPy dtype object) – The dtype of the returned array data

Examples

>>> import yt
>>> import numpy as np
>>> ds = yt.load('IsolatedGalaxy/galaxy0030/galaxy0030')
>>> a = ds.arr([1, 2, 3], 'cm')
>>> b = ds.arr([4, 5, 6], 'm')
>>> a + b
YTArray([ 401.,  502.,  603.]) cm
>>> b + a
YTArray([ 4.01,  5.02,  6.03]) m

Arrays returned by this function know about the dataset’s unit system

>>> a = ds.arr(np.ones(5), 'code_length')
>>> a.in_units('Mpccm/h')
YTArray([ 1.00010449,  1.00010449,  1.00010449,  1.00010449,
         1.00010449]) Mpc
box(left_edge, right_edge, **kwargs)[source]

box is a wrapper to the Region object for creating a region without having to specify a center value. It assumes the center is the midpoint between the left_edge and right_edge.

checksum

Computes md5 sum of a dataset.

Note: Currently this property is unable to determine a complete set of files that are a part of a given dataset. As a first approximation, the checksum of parameter_file is calculated. In case parameter_file is a directory, checksum of all files inside the directory is calculated.

close()[source]
coordinates = None
create_field_info()[source]
default_field = ('gas', 'density')
default_fluid_type = 'gas'
define_unit(symbol, value, tex_repr=None, offset=None, prefixable=False)[source]

Define a new unit and add it to the dataset’s unit registry.

Parameters:
  • symbol (string) – The symbol for the new unit.
  • value (tuple or YTQuantity) – The definition of the new unit in terms of some other units. For example, one would define a new “mph” unit with (1.0, “mile/hr”)
  • tex_repr (string, optional) – The LaTeX representation of the new unit. If one is not supplied, it will be generated automatically based on the symbol string.
  • offset (float, optional) – The default offset for the unit. If not set, an offset of 0 is assumed.
  • prefixable (bool, optional) – Whether or not the new unit can use SI prefixes. Default: False

Examples

>>> ds.define_unit("mph", (1.0, "mile/hr"))
>>> two_weeks = YTQuantity(14.0, "days")
>>> ds.define_unit("fortnight", two_weeks)
derived_field_list
domain_center = None
domain_dimensions = None
domain_left_edge = None
domain_right_edge = None
domain_width = None
field_list
field_units = None
fields
find_field_values_at_point(fields, coords)[source]

Returns the values [field1, field2,...] of the fields at the given coordinates. Returns a list of field values in the same order as the input fields.

find_field_values_at_points(fields, coords)[source]

Returns the values [field1, field2,...] of the fields at the given [(x1, y1, z2), (x2, y2, z2),...] points. Returns a list of field values in the same order as the input fields.

find_max(field)[source]

Returns (value, location) of the maximum of a given field.

find_min(field)[source]

Returns (value, location) for the minimum of a given field.

fluid_types = ('gas', 'deposit', 'index')
geometry = 'cartesian'
get_smallest_appropriate_unit(v, quantity='distance', return_quantity=False)[source]

Returns the largest whole unit smaller than the YTQuantity passed to it as a string.

The quantity keyword can be equal to distance or time. In the case of distance, the units are: ‘Mpc’, ‘kpc’, ‘pc’, ‘au’, ‘rsun’, ‘km’, etc. For time, the units are: ‘Myr’, ‘kyr’, ‘yr’, ‘day’, ‘hr’, ‘s’, ‘ms’, etc.

If return_quantity is set to True, it finds the largest YTQuantity object with a whole unit and a power of ten as the coefficient, and it returns this YTQuantity.

get_unit_from_registry(unit_str)[source]

Creates a unit object matching the string expression, using this dataset’s unit registry.

Parameters:unit_str (str) – string that we can parse for a sympy Expr.
h
has_key(key)[source]

Checks units, parameters, and conversion factors. Returns a boolean.

hierarchy
hub_upload()[source]
index
ires_factor
known_filters = None
max_level
particle_fields_by_type
particle_type_counts
particle_types = ('io',)
particle_types_raw = ('io',)
particle_unions = None
particles_exist
print_key_parameters()[source]
print_stats()[source]
quan

Converts an scalar into a yt.units.yt_array.YTQuantity

The returned YTQuantity will be dimensionless by default, but can be cast to arbitrary units using the input_units keyword argument.

Parameters:
  • input_scalar (an integer or floating point scalar) – The scalar to attach units to
  • input_units (String unit specification, unit symbol or astropy object) – The units of the quantity. Powers must be specified using python syntax (cm**3, not cm^3).
  • dtype (string or NumPy dtype object) – The dtype of the array data.

Examples

>>> import yt
>>> ds = yt.load('IsolatedGalaxy/galaxy0030/galaxy0030')
>>> a = ds.quan(1, 'cm')
>>> b = ds.quan(2, 'm')
>>> a + b
201.0 cm
>>> b + a
2.01 m

Quantities created this way automatically know about the unit system of the dataset.

>>> a = ds.quan(5, 'code_length')
>>> a.in_cgs()
1.543e+25 cm
relative_refinement(l0, l1)[source]
set_code_units()[source]
set_units()[source]

Creates the unit registry for this dataset.

setup_deprecated_fields()[source]
storage_filename = None
class yt.data_objects.static_output.FieldNameContainer(ds, field_type)[source]

Bases: object

class yt.data_objects.static_output.FieldTypeContainer(ds)[source]

Bases: object

field_types
class yt.data_objects.static_output.IndexProxy(ds)[source]

Bases: object

class yt.data_objects.static_output.MutableAttribute[source]

Bases: object

A descriptor for mutable data

class yt.data_objects.static_output.ParticleDataset(filename, dataset_type=None, file_style=None, units_override=None, unit_system='cgs', n_ref=64, over_refine_factor=1)[source]

Bases: yt.data_objects.static_output.Dataset

add_deposited_particle_field(deposit_field, method, kernel_name='cubic', weight_field='particle_mass')

Add a new deposited particle field

Creates a new deposited field based on the particle deposit_field.

Parameters:
  • deposit_field (tuple) – The field name tuple of the particle field the deposited field will be created from. This must be a field name tuple so yt can appropriately infer the correct particle type.
  • method (string) – This is the “method name” which will be looked up in the particle_deposit namespace as methodname_deposit. Current methods include simple_smooth, sum, std, cic, weighted_mean, mesh_id, and nearest.
  • kernel_name (string, default 'cubic') – This is the name of the smoothing kernel to use. It is only used for the simple_smooth method and is otherwise ignored. Current supported kernel names include cubic, quartic, quintic, wendland2, wendland4, and wendland6.
  • weight_field (string, default 'particle_mass') – Weighting field name for deposition method weighted_mean.
Returns:

Return type:

The field name tuple for the newly created field.

add_field(name, function=None, sampling_type=None, **kwargs)

Dataset-specific call to add_field

Add a new field, along with supplemental metadata, to the list of available fields. This respects a number of arguments, all of which are passed on to the constructor for DerivedField.

Parameters:
  • name (str) – is the name of the field.
  • function (callable) – A function handle that defines the field. Should accept arguments (field, data)
  • units (str) – A plain text string encoding the unit. Powers must be in python syntax (** instead of ^).
  • take_log (bool) – Describes whether the field should be logged
  • validators (list) – A list of FieldValidator objects
  • particle_type (bool) – Is this a particle (1D) field?
  • vector_field (bool) – Describes the dimensionality of the field. Currently unused.
  • display_name (str) – A name used in the plots
  • force_override (bool) – Whether to override an existing derived field. Does not work with on-disk fields.
add_gradient_fields(input_field)

Add gradient fields.

Creates four new grid-based fields that represent the components of the gradient of an existing field, plus an extra field for the magnitude of the gradient. Currently only supported in Cartesian geometries. The gradient is computed using second-order centered differences.

Parameters:input_field (tuple) – The field name tuple of the particle field the deposited field will be created from. This must be a field name tuple so yt can appropriately infer the correct field type.
Returns:
Return type:A list of field name tuples for the newly created fields.

Examples

>>> grad_fields = ds.add_gradient_fields(("gas","temperature"))
>>> print(grad_fields)
[('gas', 'temperature_gradient_x'),
 ('gas', 'temperature_gradient_y'),
 ('gas', 'temperature_gradient_z'),
 ('gas', 'temperature_gradient_magnitude')]
add_particle_filter(filter)
add_particle_union(union)
add_smoothed_particle_field(smooth_field, method='volume_weighted', nneighbors=64, kernel_name='cubic')

Add a new smoothed particle field

Creates a new smoothed field based on the particle smooth_field.

Parameters:
  • smooth_field (tuple) – The field name tuple of the particle field the smoothed field will be created from. This must be a field name tuple so yt can appropriately infer the correct particle type.
  • method (string, default 'volume_weighted') – The particle smoothing method to use. Can only be ‘volume_weighted’ for now.
  • nneighbors (int, default 64) – The number of neighbors to examine during the process.
  • kernel_name (string, default cubic) – This is the name of the smoothing kernel to use. Current supported kernel names include cubic, quartic, quintic, wendland2, wendland4, and wendland6.
Returns:

Return type:

The field name tuple for the newly created field.

all_data(find_max=False, **kwargs)

all_data is a wrapper to the Region object for creating a region which covers the entire simulation domain.

arr

Converts an array into a yt.units.yt_array.YTArray

The returned YTArray will be dimensionless by default, but can be cast to arbitrary units using the input_units keyword argument.

Parameters:
  • input_array (Iterable) – A tuple, list, or array to attach units to
  • input_units (String unit specification, unit symbol or astropy object) – The units of the array. Powers must be specified using python syntax (cm**3, not cm^3).
  • dtype (string or NumPy dtype object) – The dtype of the returned array data

Examples

>>> import yt
>>> import numpy as np
>>> ds = yt.load('IsolatedGalaxy/galaxy0030/galaxy0030')
>>> a = ds.arr([1, 2, 3], 'cm')
>>> b = ds.arr([4, 5, 6], 'm')
>>> a + b
YTArray([ 401.,  502.,  603.]) cm
>>> b + a
YTArray([ 4.01,  5.02,  6.03]) m

Arrays returned by this function know about the dataset’s unit system

>>> a = ds.arr(np.ones(5), 'code_length')
>>> a.in_units('Mpccm/h')
YTArray([ 1.00010449,  1.00010449,  1.00010449,  1.00010449,
         1.00010449]) Mpc
box(left_edge, right_edge, **kwargs)

box is a wrapper to the Region object for creating a region without having to specify a center value. It assumes the center is the midpoint between the left_edge and right_edge.

checksum

Computes md5 sum of a dataset.

Note: Currently this property is unable to determine a complete set of files that are a part of a given dataset. As a first approximation, the checksum of parameter_file is calculated. In case parameter_file is a directory, checksum of all files inside the directory is calculated.

close()
coordinates = None
create_field_info()
default_field = ('gas', 'density')
default_fluid_type = 'gas'
define_unit(symbol, value, tex_repr=None, offset=None, prefixable=False)

Define a new unit and add it to the dataset’s unit registry.

Parameters:
  • symbol (string) – The symbol for the new unit.
  • value (tuple or YTQuantity) – The definition of the new unit in terms of some other units. For example, one would define a new “mph” unit with (1.0, “mile/hr”)
  • tex_repr (string, optional) – The LaTeX representation of the new unit. If one is not supplied, it will be generated automatically based on the symbol string.
  • offset (float, optional) – The default offset for the unit. If not set, an offset of 0 is assumed.
  • prefixable (bool, optional) – Whether or not the new unit can use SI prefixes. Default: False

Examples

>>> ds.define_unit("mph", (1.0, "mile/hr"))
>>> two_weeks = YTQuantity(14.0, "days")
>>> ds.define_unit("fortnight", two_weeks)
derived_field_list
domain_center = None
domain_dimensions = None
domain_left_edge = None
domain_right_edge = None
domain_width = None
field_list
field_units = None
fields
filter_bbox = False
find_field_values_at_point(fields, coords)

Returns the values [field1, field2,...] of the fields at the given coordinates. Returns a list of field values in the same order as the input fields.

find_field_values_at_points(fields, coords)

Returns the values [field1, field2,...] of the fields at the given [(x1, y1, z2), (x2, y2, z2),...] points. Returns a list of field values in the same order as the input fields.

find_max(field)

Returns (value, location) of the maximum of a given field.

find_min(field)

Returns (value, location) for the minimum of a given field.

fluid_types = ('gas', 'deposit', 'index')
geometry = 'cartesian'
get_smallest_appropriate_unit(v, quantity='distance', return_quantity=False)

Returns the largest whole unit smaller than the YTQuantity passed to it as a string.

The quantity keyword can be equal to distance or time. In the case of distance, the units are: ‘Mpc’, ‘kpc’, ‘pc’, ‘au’, ‘rsun’, ‘km’, etc. For time, the units are: ‘Myr’, ‘kyr’, ‘yr’, ‘day’, ‘hr’, ‘s’, ‘ms’, etc.

If return_quantity is set to True, it finds the largest YTQuantity object with a whole unit and a power of ten as the coefficient, and it returns this YTQuantity.

get_unit_from_registry(unit_str)

Creates a unit object matching the string expression, using this dataset’s unit registry.

Parameters:unit_str (str) – string that we can parse for a sympy Expr.
h
has_key(key)

Checks units, parameters, and conversion factors. Returns a boolean.

hierarchy
hub_upload()
index
ires_factor
known_filters = None
max_level
particle_fields_by_type
particle_type_counts
particle_types = ('io',)
particle_types_raw = ('io',)
particle_unions = None
particles_exist
print_key_parameters()
print_stats()
quan

Converts an scalar into a yt.units.yt_array.YTQuantity

The returned YTQuantity will be dimensionless by default, but can be cast to arbitrary units using the input_units keyword argument.

Parameters:
  • input_scalar (an integer or floating point scalar) – The scalar to attach units to
  • input_units (String unit specification, unit symbol or astropy object) – The units of the quantity. Powers must be specified using python syntax (cm**3, not cm^3).
  • dtype (string or NumPy dtype object) – The dtype of the array data.

Examples

>>> import yt
>>> ds = yt.load('IsolatedGalaxy/galaxy0030/galaxy0030')
>>> a = ds.quan(1, 'cm')
>>> b = ds.quan(2, 'm')
>>> a + b
201.0 cm
>>> b + a
2.01 m

Quantities created this way automatically know about the unit system of the dataset.

>>> a = ds.quan(5, 'code_length')
>>> a.in_cgs()
1.543e+25 cm
relative_refinement(l0, l1)
set_code_units()
set_units()

Creates the unit registry for this dataset.

setup_deprecated_fields()
storage_filename = None
class yt.data_objects.static_output.ParticleFile(ds, io, filename, file_id)[source]

Bases: object

count(selector)[source]
select(selector)[source]
class yt.data_objects.static_output.RegisteredDataset(name, b, d)[source]

Bases: type

mro() → list

return a type’s method resolution order

yt.data_objects.static_output.requires_index(attr_name)[source]