diff --git a/tests/test_line_database.py b/tests/test_line_database.py index ff742e874..c09ad4ded 100644 --- a/tests/test_line_database.py +++ b/tests/test_line_database.py @@ -62,6 +62,14 @@ def test_line_database_from_input(): assert ld.lines_all[0].identifier == HI.identifier print(ld) +def test_line_database_add_line_linetools(): + ld = LineDatabase() + HI = Line('H', 'I', 1216, 626500000.0, 2.3, identifier='Ly a') + ld.add_line('H', 'I', 1216, use_linetools=True, identifier='Ly a') + np.testing.assert_allclose( ld.lines_all[0].gamma, HI.gamma ) + assert ld.lines_all[0].identifier == HI.identifier + print(ld) + def test_select_lines_from_line_database(): ld = LineDatabase('lines.txt') assert len(ld.select_lines('Ne')) == 8 # 8 listed Ne lines diff --git a/trident/constants.py b/trident/constants.py new file mode 100644 index 000000000..6941057e5 --- /dev/null +++ b/trident/constants.py @@ -0,0 +1,49 @@ +""" +Constants and Tables + +""" + +#----------------------------------------------------------------------------- +# Copyright (c) 2016, Trident Development Team. +# +# Distributed under the terms of the Modified BSD License. +# +# The full license is in the file LICENSE, distributed with this software. +#----------------------------------------------------------------------------- + +# Taken from Cloudy documentation. +solar_abundance = { + 'H' : 1.00e+00, 'He': 1.00e-01, 'Li': 2.04e-09, + 'Be': 2.63e-11, 'B' : 6.17e-10, 'C' : 2.45e-04, + 'N' : 8.51e-05, 'O' : 4.90e-04, 'F' : 3.02e-08, + 'Ne': 1.00e-04, 'Na': 2.14e-06, 'Mg': 3.47e-05, + 'Al': 2.95e-06, 'Si': 3.47e-05, 'P' : 3.20e-07, + 'S' : 1.84e-05, 'Cl': 1.91e-07, 'Ar': 2.51e-06, + 'K' : 1.32e-07, 'Ca': 2.29e-06, 'Sc': 1.48e-09, + 'Ti': 1.05e-07, 'V' : 1.00e-08, 'Cr': 4.68e-07, + 'Mn': 2.88e-07, 'Fe': 2.82e-05, 'Co': 8.32e-08, + 'Ni': 1.78e-06, 'Cu': 1.62e-08, 'Zn': 3.98e-08} + +atomic_mass = { + 'H' : 1.00794, 'He': 4.002602, 'Li': 6.941, + 'Be': 9.012182, 'B' : 10.811, 'C' : 12.0107, + 'N' : 14.0067, 'O' : 15.9994, 'F' : 18.9984032, + 'Ne': 20.1797, 'Na': 22.989770, 'Mg': 24.3050, + 'Al': 26.981538, 'Si': 28.0855, 'P' : 30.973761, + 'S' : 32.065, 'Cl': 35.453, 'Ar': 39.948, + 'K' : 39.0983, 'Ca': 40.078, 'Sc': 44.955910, + 'Ti': 47.867, 'V' : 50.9415, 'Cr': 51.9961, + 'Mn': 54.938049, 'Fe': 55.845, 'Co': 58.933200, + 'Ni': 58.6934, 'Cu': 63.546, 'Zn': 65.409} + +atomic_number = { + 'H' : 1, 'He': 2, 'Li': 3, + 'Be': 4, 'B' : 5, 'C' : 6, + 'N' : 7, 'O' : 8, 'F' : 9, + 'Ne': 10, 'Na': 11, 'Mg': 12, + 'Al': 13, 'Si': 14, 'P' : 15, + 'S' : 16, 'Cl': 17, 'Ar': 18, + 'K' : 19, 'Ca': 20, 'Sc': 21, + 'Ti': 22, 'V' : 23, 'Cr': 24, + 'Mn': 25, 'Fe': 26, 'Co': 27, + 'Ni': 28, 'Cu': 29, 'Zn': 30} diff --git a/trident/ion_balance.py b/trident/ion_balance.py index 39a9c062d..8e2adfa73 100644 --- a/trident/ion_balance.py +++ b/trident/ion_balance.py @@ -29,6 +29,8 @@ uniquify from trident.roman import \ from_roman +from trident.constants import \ + solar_abundance, atomic_mass, atomic_number H_mass_fraction = 0.76 to_nH = H_mass_fraction / mh @@ -245,7 +247,8 @@ def add_ion_fields(ds, ions, ftype='gas', # If line_database is set, then use the underlying file as the line list # to select ions from. if line_database is not None: - line_database = LineDatabase(line_database) + if not isinstance( line_database, LineDatabase ): + line_database = LineDatabase(line_database) ion_list = line_database.parse_subset_to_ions(ions) # Otherwise, any ion can be selected (not just ones in the line list). @@ -861,41 +864,3 @@ def _alias_field(ds, alias_name, name): ds.field_info.alias(alias_name, name) ds.derived_field_list.append(alias_name) return - - -# Taken from Cloudy documentation. -solar_abundance = { - 'H' : 1.00e+00, 'He': 1.00e-01, 'Li': 2.04e-09, - 'Be': 2.63e-11, 'B' : 6.17e-10, 'C' : 2.45e-04, - 'N' : 8.51e-05, 'O' : 4.90e-04, 'F' : 3.02e-08, - 'Ne': 1.00e-04, 'Na': 2.14e-06, 'Mg': 3.47e-05, - 'Al': 2.95e-06, 'Si': 3.47e-05, 'P' : 3.20e-07, - 'S' : 1.84e-05, 'Cl': 1.91e-07, 'Ar': 2.51e-06, - 'K' : 1.32e-07, 'Ca': 2.29e-06, 'Sc': 1.48e-09, - 'Ti': 1.05e-07, 'V' : 1.00e-08, 'Cr': 4.68e-07, - 'Mn': 2.88e-07, 'Fe': 2.82e-05, 'Co': 8.32e-08, - 'Ni': 1.78e-06, 'Cu': 1.62e-08, 'Zn': 3.98e-08} - -atomic_mass = { - 'H' : 1.00794, 'He': 4.002602, 'Li': 6.941, - 'Be': 9.012182, 'B' : 10.811, 'C' : 12.0107, - 'N' : 14.0067, 'O' : 15.9994, 'F' : 18.9984032, - 'Ne': 20.1797, 'Na': 22.989770, 'Mg': 24.3050, - 'Al': 26.981538, 'Si': 28.0855, 'P' : 30.973761, - 'S' : 32.065, 'Cl': 35.453, 'Ar': 39.948, - 'K' : 39.0983, 'Ca': 40.078, 'Sc': 44.955910, - 'Ti': 47.867, 'V' : 50.9415, 'Cr': 51.9961, - 'Mn': 54.938049, 'Fe': 55.845, 'Co': 58.933200, - 'Ni': 58.6934, 'Cu': 63.546, 'Zn': 65.409} - -atomic_number = { - 'H' : 1, 'He': 2, 'Li': 3, - 'Be': 4, 'B' : 5, 'C' : 6, - 'N' : 7, 'O' : 8, 'F' : 9, - 'Ne': 10, 'Na': 11, 'Mg': 12, - 'Al': 13, 'Si': 14, 'P' : 15, - 'S' : 16, 'Cl': 17, 'Ar': 18, - 'K' : 19, 'Ca': 20, 'Sc': 21, - 'Ti': 22, 'V' : 23, 'Cr': 24, - 'Mn': 25, 'Fe': 26, 'Co': 27, - 'Ni': 28, 'Cu': 29, 'Zn': 30} diff --git a/trident/line_database.py b/trident/line_database.py index 1f3f45e04..8a76c6106 100644 --- a/trident/line_database.py +++ b/trident/line_database.py @@ -12,12 +12,16 @@ #----------------------------------------------------------------------------- import os +import warnings from yt.funcs import \ mylog from trident.config import \ trident_path from trident.roman import \ - from_roman + from_roman, to_roman +from trident.constants import \ + atomic_number +import astropy.units as u def uniquify(list): @@ -145,10 +149,10 @@ def __init__(self, input_file=None): if input_file is not None: self.load_line_list_from_file(input_file) else: - self.input_file = 'Manually Entered' + self.load_line_list_from_linetools() - def add_line(self, element, ion_state, wavelength, gamma, - f_value, field=None, identifier=None): + def add_line(self, element, ion_state, wavelength, gamma=None, + f_value=None, field=None, identifier=None, use_linetools=False): """ Manually add a line to the :class:`~trident.LineDatabase`. @@ -188,6 +192,10 @@ def add_line(self, element, ion_state, wavelength, gamma, An optional identifier for the transition Example: 'Ly a' for Lyman alpha + :use_linetools: bool + + If True retrieve the values for wavelength, gamma, and f_value + from the linetools python package (linetools.readthedocs.io). **Example** @@ -199,9 +207,67 @@ def add_line(self, element, ion_state, wavelength, gamma, >>> ldb.add_line('H', 'I', 1215.67, 469860000, 0.41641, 'Ly a') >>> print(ldb.lines_all) """ + + if use_linetools: + # Retrieve + def linetools_id( wrest ): + return '%s%s %d' % (element, ion_state, round(float(wrest), 0)) + line = self.linetools_linelist[ linetools_id(wavelength) ] + + # Check for rounding differences + if line is None: + warnings.warn( 'Line ' + linetools_id(wavelength) + \ + ' not found in linetools database. Checking adjacent wavelengths.') + line = self.linetools_linelist[ linetools_id(wavelength-1) ] + if line is None: + line = self.linetools_linelist[ linetools_id(wavelength+1) ] + if line is None: + raise KeyError( 'No line found in linetools database.' ) + + # Format + wavelength = line['wrest'].value + gamma = line['gamma'].value + f_value = line['f'] + + else: + assert gamma is not None + assert f_value is not None + self.lines_all.append(Line(element, ion_state, wavelength, gamma, f_value, field, identifier)) + def add_line_from_linetools(self, name): + """ + name is the linetools name for a feature: "ion wrest" + e.g., HI 1215 + """ + line = self.linetools_linelist[name] + reversed_atomic_number = {value : key for (key, value) in atomic_number.items()} + if line['Z'] not in reversed_atomic_number.keys(): + # Not a line < Zn so not a line we care about. + return + element = reversed_atomic_number[line['Z']] + ion_state = to_roman(line['ion']) + wavelength = line['wrest'].value + gamma = line['gamma'].value + f_value = line['f'] + identifier = line['name'] + self.add_line(element, ion_state, wavelength, gamma, f_value, + identifier=identifier) + + def load_line_list_from_linetools(self): + """ + Create a line list that incorporates all of the "ISM" lines between 100 and + 1 micron angstroms from the external linetools code. + """ + bounds = [100, 10000] * u.AA + name_list = self.linetools_linelist.available_transitions(bounds)['name'] + for name in name_list: + # Don't include fine structure lines or deuterium lines + if "*" not in name and \ + "DI" not in name: + self.add_line_from_linetools(name) + def load_line_list_from_file(self, filename): """ Load a line list from a file into the LineDatabase. @@ -451,6 +517,21 @@ def parse_subset_to_ions(self, subsets=None): ions = uniquify(ions) return ions + @property + def linetools_linelist(self): + '''Handler for a linetools LineList class. + ''' + + if not hasattr( self, '_linetools_linelist' ): + + # Only doing the import here so rest of the class doesn't + # break if the feature isn't enabled. + from linetools.lists.linelist import LineList + + self._linetools_linelist = LineList('ISM') + + return self._linetools_linelist + def __repr__(self): disp = "" for line in self.lines_all: diff --git a/trident/ray_generator.py b/trident/ray_generator.py index 8442b7fd4..44ec82400 100644 --- a/trident/ray_generator.py +++ b/trident/ray_generator.py @@ -24,7 +24,7 @@ from_roman from yt.data_objects.static_output import \ Dataset -from trident.ion_balance import \ +from trident.constants import \ atomic_number def make_simple_ray(dataset_file, start_position, end_position, diff --git a/trident/spectrum_generator.py b/trident/spectrum_generator.py index a7cb5a6eb..d2782a6ca 100644 --- a/trident/spectrum_generator.py +++ b/trident/spectrum_generator.py @@ -35,7 +35,8 @@ from trident.instrument import \ Instrument from trident.ion_balance import \ - add_ion_number_density_field, \ + add_ion_number_density_field +from trident.constants import \ atomic_mass from trident.line_database import \ LineDatabase