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)[source]¶

Bases: object

molecular database of ExoMol

MdbExomol is a class for ExoMol.

Attributes:

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

QT_interp(T)[source]¶

interpolated partition function

Args:

T: temperature

Returns:

Q(T) interpolated in jnp.array

download(molec, extension, numtag=None)[source]¶

Downloading Exomol files

Args:

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”

Note:

The download URL is written in exojax.utils.url.

masking(mask)[source]¶

applying mask and (re)generate jnp.arrays

Args:

mask: mask to be applied. self.mask is updated.

qr_interp(T)[source]¶

interpolated partition function ratio

Args:

T: temperature

Returns:

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

to_df()[source]¶

Export Line Database to a RADIS-friendly Pandas DataFrame

to_partition_function_tabulator()[source]¶
fetch_exomol(molecule, database=None, local_databases='~/.radisdb/exomol/', databank_name='EXOMOL-{molecule}', isotope=None, load_wavenum_min=None, load_wavenum_max=None, cache=True, verbose=True, clean_cache_files=True, return_local_path=False, return_partition_function=False)[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 molecule

  • database (str) – database name. Ex:: POKAZATEL or BT2 for H2O. See KNOWN_EXOMOL_DATABASE_NAMES. If None and there is only one database available, use it.

  • local_databases (str) – where to create the RADIS HDF5 files. Default "~/.radisdb/exomol"

  • databank_name (str) – name of the databank in RADIS Configuration file Default "EXOMOL-{molecule}"

  • isotope (str) – load only certain isotopes : '2', '1,2', etc. Default 1.

  • load_wavenum_min, load_wavenum_max (float (cm-1)) – load only specific wavenumbers.

Other Parameters
  • cache (bool, or 'regen') – if True, use existing HDF5 file. If False or 'regen', rebuild it. Default True.

  • 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 a PartFuncExoMol object.

Returns

  • df (pd.DataFrame) – Line list A HDF5 file is also created in local_databases and referenced in the RADIS config file with name databank_name

  • local_path (str) – path of local database file if return_local_path

Notes

if using load_only_wavenum_above/below or isotope, the whole database is anyway downloaded and uncompressed to local_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 uses hdf2df()

See also

hdf2df()

get_exomol_database_list(molecule, isotope_full_name)[source]¶

Parse ExoMol website and return list of available databases, and recommended database

Parameters

Examples

::

databases, recommended = get_exomol_database_list(“CH4”, “12C-1H4”) >>> [‘xsec-YT10to10’, ‘YT10to10’, ‘YT34to10’], ‘YT34to10’

get_exomol_full_isotope_name(molecule, isotope)[source]¶

Get full isotope name for molecule and isotope number

Parameters
  • molecule (str)

  • isotope (int) – terrestrial abundance

Examples

::

get_exomol_full_isotope_name(“CH4”, 1) >>> ‘12C-1H4’