radis.io.exomol module¶
Molecular database (MDB) class
MdbExomol is the MDB for ExoMol
Initial code borrowed from the Exojax code (which you should also have a look at !), by @HajimeKawahara, under MIT License.
- class MdbExomol(path, nurange=[- inf, inf], margin=1.0, crit=- inf, bkgdatm='H2', broadf=True, engine='vaex', verbose=True)[source]¶
Bases:
object
molecular database of ExoMol
MdbExomol is a class for ExoMol.
- Parameters
nu_lines (nd array) (line center (cm-1))
Sij0 (nd array) (line strength at T=Tref (cm))
dev_nu_lines (np array) (line center in device (cm-1))
logsij0 (np array) (log line strength at T=Tref)
A (np array) (Einstein A coefficient)
elower (np array) (the lower state energy (cm-1))
gpp (np array) (statistical weight)
jlower (np array) (J_lower)
jupper (np array) (J_upper)
n_Tref (np array) (temperature exponent)
alpha_ref (np array) (alpha_ref (gamma0))
n_Tref_def (default temperature exponent in .def file, used for jlower not given in .broad)
alpha_ref_def (default alpha_ref (gamma0) in .def file, used for jlower not given in .broad)
.. minigallery:: radis.fetch_exomol
- QT_interp(T)[source]¶
interpolated partition function
- Parameters
T (temperature)
- Return type
Q(T) interpolated in jnp.array
- download(molec, extension, numtag=None)[source]¶
Downloading Exomol files
- Parameters
molec (like “12C-16O__Li2015”)
extension (extension list e.g. [“.pf”,”.def”,”.trans.bz2”,”.states.bz2”,”.broad”])
numtag (number tag of transition file if exists. e.g. “11100-11200”)
Notes
The download URL is written in exojax.utils.url.
- masking(mask, mask_needed=True)[source]¶
applying mask and (re)generate jnp.arrays
- Args:
mask: mask to be applied. self.mask is updated. mask_needed: whether mask needs to be applied or not
- qr_interp(T)[source]¶
interpolated partition function ratio
- Parameters
T (temperature)
- Return type
qr(T)=Q(T)/Q(Tref) interpolated in jnp.array
- set_broadening(alpha_ref_def=None, n_Texp_def=None)[source]¶
setting broadening parameters
- Args:
alpha_ref: set default alpha_ref and apply it. None=use self.alpha_ref_def n_Texp_def: set default n_Texp and apply it. None=use self.n_Texp_def
- fetch_exomol(molecule, database=None, local_databases=None, databank_name='EXOMOL-{molecule}', isotope='1', load_wavenum_min=None, load_wavenum_max=None, columns=None, cache=True, verbose=True, clean_cache_files=True, return_local_path=False, return_partition_function=False, engine='default')[source]¶
Stream ExoMol file from EXOMOL website. Unzip and build a HDF5 file directly.
Returns a Pandas DataFrame containing all lines.
- Parameters
molecule (
str
) – ExoMol moleculedatabase (
str
) – database name. Ex::POKAZATEL
orBT2
forH2O
. SeeKNOWN_EXOMOL_DATABASE_NAMES
. IfNone
and there is only one database available, use it.local_databases (
str
) – where to create the RADIS HDF5 files. Default"~/.radisdb/exomol"
. Can be changed inradis.config["DEFAULT_DOWNLOAD_PATH"]
or in ~/radis.json config filedatabank_name (
str
) – name of the databank in RADIS Configuration file Default"EXOMOL-{molecule}"
isotope (
str
orint
) – load only certain isotopes, sorted by terrestrial abundances :'1'
,'2'
, etc. Default1
.Note
In RADIS, isotope abundance is included in the line intensity calculation. However, the terrestrial abundances used may not be relevant to non-terrestrial applications. By default, the abundance is given reading HITRAN data. If the molecule does not exist in the HITRAN database, the abundance is read from the
radis/radis_default.json
configuration file, which can be modified by editingradis.config
after import or directly by editing the user~/radis.json
user configuration file (overwritesradis_default.json
). In theradis/radis_default.json
file, values were calculated with a simple model based on the terrestrial isotopic abundance of each element.load_wavenum_min, load_wavenum_max (float (cm-1)) – load only specific wavenumbers.
columns (list of str) – list of columns to load. If
None
, returns all columns in the file.
- Other Parameters
cache (bool, or
'regen'
) – ifTrue
, use existing HDF5 file. IfFalse
or'regen'
, rebuild it. DefaultTrue
.verbose (bool)
clean_cache_files (bool) – if
True
clean downloaded cache files after HDF5 are created.return_local_path (bool) – if
True
, also returns the path of the local database file.return_partition_function (bool) – if
True
, also returns aPartFuncExoMol
object.engine (‘vaex’, ‘feather’) – which memory-mapping library to use. If ‘default’ use the value from ~/radis.json
- Returns
df (pd.DataFrame) – Line list A HDF5 file is also created in
local_databases
and referenced in the RADIS config file with namedatabank_name
local_path (str) – path of local database file if
return_local_path
Examples
Notes
if using
load_only_wavenum_above/below
orisotope
, the whole database is anyway downloaded and uncompressed tolocal_databases
fast access .HDF5 files (which will take a long time on first call). Only the expected wavenumber range & isotopes are returned. The .HFD5 parsing useshdf2df()
See also
- get_exomol_database_list(molecule, isotope_full_name)[source]¶
Parse ExoMol website and return list of available databases, and recommended database
- Parameters
molecule (str)
isotope_full_name (str) – isotope full name (ex.
12C-1H4
for CH4,1). Get it fromradis.io.exomol.get_exomol_full_isotope_name()
Examples
Get CH4 from ExoMol :
databases, recommended = get_exomol_database_list("CH4", "12C-1H4") >>> ['xsec-YT10to10', 'YT10to10', 'YT34to10'], 'YT34to10'
Or combine with
get_exomol_full_isotope_name()
to get the isopologue (sorted by terrestrial abundance)from radis.io.exomol import get_exomol_database_list, get_exomol_full_isotope_name databases, recommended = get_exomol_database_list("CH4", get_exomol_full_isotope_name("CH4", 1)) >>> ['xsec-YT10to10', 'YT10to10', 'YT34to10'], 'YT34to10'
See also