trident.ion_balance.add_ion_fields
- trident.ion_balance.add_ion_fields(ds, ions, ftype='gas', ionization_table=None, field_suffix=False, line_database=None, sampling_type='local', particle_type=None)[source]
Preferred method for adding ion fields to a yt dataset.
Select ions based on the selection indexing set up in
parse_subset_to_ions
function, that is, by specifying a list of strings where each string represents an ion or line. Strings are of one of three forms:<element>
<element> <ion state>
<element> <ion state> <line_wavelength>
If a line_database is selected, then the ions chosen will be a subset of the ions present in the equivalent
LineDatabase
, nominally located intrident.__path__/data/line_lists
.For each ion species selected, four fields will be added (example for Mg II):
Ion fraction field. e.g. (“gas”, ‘Mg_p1_ion_fraction’)
Number density field. e.g. (“gas”, ‘Mg_p1_number_density’)
Density field. e.g. (“gas”, ‘Mg_p1_density’)
Mass field. e.g. (“gas”, ‘Mg_p1_mass’)
This function is the preferred method for adding ion fields to one’s dataset, but for more fine-grained control, one can also employ the
add_ion_fraction_field
,add_ion_number_density_field
,add_ion_density_field
,add_ion_mass_field
functions individually.Fields are added assuming collisional ionization equilibrium and photoionization in the optically thin limit from a redshift-dependent metagalactic ionizing background using the ionization_table specified.
Parameters
- Ds
yt dataset object
This is the dataset to which the ion fraction field will be added.
- Ions
list of strings
List of strings matching possible lines. Strings can be of the form: * Atom - Examples: “H”, “C”, “Mg” * Ion - Examples: “H I”, “H II”, “C IV”, “Mg II” * Line - Examples: “H I 1216”, “C II 1336”, “Mg II 1240”
If set to ‘all’, creates all ions for the first 30 elements: (ie hydrogen to zinc). If set to ‘all’ with
line_database
keyword set, then creates all ions associated with the lines specified in the equivalentLineDatabase
.- Ionization_table
string, optional
Path to an appropriately formatted HDF5 table that can be used to compute the ion fraction as a function of density, temperature, metallicity, and redshift. When set to None, it uses the table specified in ~/.trident/config Default: None
- Field_suffix
boolean, optional
Determines whether or not to append a suffix to the field name that indicates what ionization table was used. Useful when using generating ion_fields that already exist in a dataset.
- Line_database
string, optional
Ions are selected out of the set of ions present in the line_database constructed from the line list filename specified here. See
LineDatabase
for more information.- Ftype
string, optional
This is deprecated and no longer necessary since all relevant fields are aliased to the ‘gas’ ftype. Default: ‘gas’
- Sampling_type
string, optional
This is deprecated and no longer necessary. Default: ‘local’
- Particle_type
boolean, optional
This is deprecated and no longer necessary. Default: ‘auto’
Example
To add ionized hydrogen, doubly-ionized Carbon, and all of the Magnesium species fields to a dataset, you would run:
>>> import yt >>> import trident >>> ds = yt.load('path/to/file') >>> trident.add_ion_fields(ds, ions=['H II', 'C III', 'Mg'])