source: trunk/UTIL/corrk_exo_k/chemistry.py @ 3609

Last change on this file since 3609 was 3605, checked in by afalco, 7 days ago

exo_k tools to generate corrk.
AF

File size: 2.0 KB
Line 
1### script to handle isotopologues
2### copy from https://github.com/bmorris3/shone/tree/refs/heads/opacity-downloads
3
4import re
5import numpy as np
6from periodictable import elements
7
8
9__all__ = [
10    'isotopologue_to_species',
11    'species_name_to_common_isotopologue_name'
12]
13
14
15def isotopologue_to_species(isotopologue):
16    """
17    Convert isotopologue name to common species name.
18
19    Example: Take 1H2-16O and turn it to H2O, or take 48Ti-16O and turn it to TiO.
20
21    Parameters
22    ----------
23    isotopologue : str
24        Isotopologue name, like "1H2-16O".
25
26    Returns
27    -------
28    common_name : str
29        Common name, like "H2O".
30    """
31    species = ""
32    for element in isotopologue.split('-'):
33        for s in re.findall(r'\D+\d*', element):
34            species += ''.join(s)
35    return species if len(species) > 0 else isotopologue
36
37
38def species_name_to_common_isotopologue_name(species):
39    """
40    Convert generic species name, like "H2O", to isotopologue name like "1H2-16O".
41
42    Parameters
43    ----------
44    species : str
45        Generic name, like "H2O".
46
47    Returns
48    -------
49    isotopologue_name : str
50        Isotopologue name, like "1H2-16O".
51    """
52    atoms = np.array(list(filter(
53        lambda x: len(x) > 0, re.split(r"(?<=[a-z])|(?=[A-Z])|\d", species)
54    )))
55
56    multipliers = np.array([
57        int(x) if len(x) > 0 else 1 for x in re.split(r'\D', species)
58    ])
59    lens = [len(''.join(atom)) for atom in atoms]
60    multipliers_skipped = np.array([multipliers[cs] for cs in np.cumsum(lens)])
61
62    masses = np.array([
63        round(getattr(elements, atom).mass) for atom, mult in zip(atoms, multipliers_skipped)
64    ])
65
66    if len(atoms) > 1:
67        correct_notation = '-'.join([
68            str(mass) + a + (str(mult) if mult > 1 else '')
69            for a, mult, mass in zip(atoms, multipliers_skipped, masses)
70        ])
71
72    # If single atom, give only the name of the atom:
73    else:
74        correct_notation = atoms[0]
75
76    return correct_notation
Note: See TracBrowser for help on using the repository browser.