exo_k.ktable

@author: jeremy leconte

Module Contents

class exo_k.ktable.Ktable(*filename_filters, filename=None, xtable=None, path=None, p_unit='unspecified', file_p_unit='unspecified', kdata_unit='unspecified', file_kdata_unit='unspecified', remove_zeros=False, search_path=None, mol=None, **kwargs)[source]

Bases: exo_k.ktable_io.Ktable_io

A class that handles 4D tables of k-coefficients.

Based on the Data_table class that handles basic operations common to Xtable.

Based on the Ktable_io class that incorporates all io routines (read/write_xxx for all the xxx supported formats).

Initializes k coeff table and supporting data from various sources (see below by order of precedence)

Parameters:
  • filename (str, optional) – Relative or absolute name of the file to be loaded.

  • filename_filters (sequence of string) – As many strings as necessary to uniquely define a file in the global search path defined in Settings. This path will be searched for a file with all the filename_filters in the name. The filename_filters can contain ‘*’.

  • xtable (Xtable object) – If no filename nor filename_filters are provided, this xtable object will be used to create a ktable. In this case, wavenumber bins must be given with the wnedges keyword.

  • path (str) – If none of the above is specifed, path can point to a directory with a LMDZ type k coeff table. In this case, see read_LMDZ for the keywords to specify.

If none of the parameters above is specified, just creates an empty object to be filled later.

Parameters:
  • p_unit (str, optional) – String identifying the pressure units to convert to (e.g. ‘bar’, ‘Pa’, ‘mbar’, or any pressure unit recognized by the astropy.units library). If =’unspecified’, no conversion is done.

  • file_p_unit (str, optional) – String to specify the current pressure unit if it is unspecified or if you have reasons to believe it is wrong (e.g. you just read a file where you know that the pressure grid and the pressure unit do not correspond)

  • kdata_unit (str, optional) – String to identify the unit to convert to. Accepts ‘cm^2’, ‘m^2’ or any surface unit recognized by the astropy.units library. If =’unspecified’, no conversion is done. In general, kdata should be kept in ‘per number’ or ‘per volume’ units (as opposed to ‘per mass’ units) as composition will always be assumed to be a number or volume mixing ratio. Opacities per unit mass are not supported yet. Note that you do not need to specify the ‘/molec’ or ‘/molecule’ in the unit.

  • file_kdata_unit (str, optional) – String to specify the current kdata unit if it is unspecified or if you have reasons to believe it is wrong (e.g. you just read a file where you know that the kdata grid and the kdata unit do not correspond)

  • remove_zeros (boolean, optional) – If True, the zeros in the kdata table are replaced by a value 10 orders of magnitude smaller than the smallest positive value

  • mol (str, optional) – The name of the gas or molecule described by the Ktable

  • search_path (str, optional) – If search_path is provided, it locally overrides the global _search_path in Settings and only files in search_path are returned.

xtable_to_ktable(xtable=None, wnedges=None, weights=None, ggrid=None, quad='legendre', order=20, g_split=0.9, mid_dw=True, write=0, remove_zeros=False)[source]

Fills the Ktable object with a k-coeff table computed from a Xtable object (inplace).

The p and kcorr units are inherited from the Xtable object.

Parameters:
  • xtable (Xtable) – input Xtable object instance

  • wnedges (Array) – edges of the wavenumber bins to be used to compute the corrk

  • weights (array, np.ndarray, optional) – If weights are provided, they are used instead of the legendre quadrature.

  • quad (string, optional) – Type of quadrature used. Default is ‘legendre’. Also available: ‘split-legendre’ which uses half quadrature points between 0. and g_split and half between g_split and 1.

  • order (Integer, optional) – Order of the Gauss legendre quadrature used. Default is 20.

  • g_split (float, optional) – Used only if quad=’split-legendre’. See above.

  • mid_dw (boolean, optional) –

    • If True, the Xsec values in the high resolution xtable data are assumed to cover a spectral interval that is centered around the corresponding wavenumber value. The first and last Xsec values are discarded.

    • If False, each interval runs from the wavenumber value to the next one. The last Xsec value is dicarded.

hires_to_ktable(path=None, filename_grid=None, logpgrid=None, tgrid=None, wnedges=None, quad='legendre', order=20, g_split=0.9, weights=None, ggrid=None, mid_dw=True, write=0, mol=None, grid_p_unit='Pa', p_unit='unspecified', kdata_unit='unspecified', file_kdata_unit='unspecified', helios=False, remove_zeros=False, **kwargs)[source]

Computes a k coeff table from Hires_spectrum objects (inplace).

Warning

By default, log pressures are specified in Pa in logpgrid!!! If you want to use another unit, do not forget to specify it with the grid_p_unit keyword.

Parameters:
  • path (String) – directory with the input files

  • filename_grid (array, np.ndarray of str with shape (logpgrid.size,tgrid.size)) – Names of the input high-res spectra. If None, the files are assumed to follow Kspectrum/LMDZ convention, i.e. be of the type ‘k001’, ‘k002’, etc. See create_fname_grid_Kspectrum_LMDZ() for possible additional keyword arguments.

  • logpgrid (array, np.ndarray) – Grid in log(pressure) of the input. Default unit is Pa, but can be changed with the grid_p_unit keyword.

  • tgrid (array, np.ndarray) – Grid in temperature of the input.

  • wnedges (array, np.ndarray) – Edges of the wavenumber bins to be used to compute the corr-k

  • weights (array, np.ndarray, optional) – If weights are provided, they are used instead of the legendre quadrature.

  • quad (string, optional) – Type of quadrature used. Default is ‘legendre’. Also available: ‘split-legendre’ which uses half quadrature points between 0. and g_split and half between g_split and 1.

  • order (Integer, optional) – Order of the Gauss legendre quadrature used. Default is 20.

  • g_split (float, optional) – Used only if quad=’split-legendre’. See above.

  • mid_dw (boolean, optional) –

    • If True, the Xsec values in the high resolution xsec data are assumed to cover a spectral interval that is centered around the corresponding wavenumber value. The first and last Xsec values are discarded.

    • If False, each interval runs from the wavenumber value to the next one. The last Xsec value is dicarded.

  • mol (string, optional) – Give a name to the molecule. Useful when used later in a Kdatabase to track molecules.

  • p_unit (str, optional) – Pressure unit to convert to.

  • grid_p_unit (str, optional) – Unit of the specified logpgrid.

  • kdata_unit (str, optional) – Kdata unit to convert to.

  • file_kdata_unit (str, optional) – Kdata unit in input files.

  • helios (bool, optional) – Read HELIOS output .bin files.

See create_fname_grid_Kspectrum_LMDZ() or Hires_spectrum for a list of additional arguments that can be provided to those funtion through **kwargs.

copy(cp_kdata=True, ktab5d=False)[source]

Creates a new instance of Ktable object and (deep) copies data into it.

Parameters:
  • cp_kdata (bool, optional) – If false, the kdata table is not copied and only the structure and metadata are.

  • ktab5d (bool, optional) – If true, creates a Ktable5d object with the same structure. Data are not copied.

Returns:

Ktable or Ktable5d

A new Ktable or Ktable5d instance with the same structure as self.

gindex(g)[source]

Finds the index corresponding to the given g

spectrum_to_plot(p=1e-05, t=200.0, x=1.0, g=None)[source]

provide the spectrum for a given point to be plotted

Parameters:
  • p (float) – Pressure (Ktable pressure unit)

  • t (float) – Temperature(K)

  • g (float) – Gauss point

  • x (float) – Mixing ratio of the species

plot_distrib(ax, p=1e-05, t=200.0, wl=1.0, x=1.0, xscale=None, yscale='log', **kwarg)[source]

Plot the distribution for a given point

Parameters:
  • p (float) – Pressure (Ktable pressure unit)

  • t (float) – Temperature(K)

  • wl (float) – Wavelength (micron)

RandOverlap(other, x_self, x_other, write=0, use_rebin=False)[source]

Method to randomly mix the opacities of 2 species (self and other).

Parameters:
  • other (Ktable) – A Ktable object to be mixed with. Dimensions should be the same as self.

  • x_self (float or array) – Volume mixing ratio of the first species

  • x_other (float or array) – Volume mixing ratio of the species to be mixed with.

If one of these is None, the kcoeffs of the species in question are considered to be already normalized with respect to the mixing ratio.

Returns:

array

array of k-coefficients for the mix.

bin_down(wnedges=None, weights=None, ggrid=None, remove_zeros=False, num=300, use_rebin=False, write=0)[source]

Method to bin down a kcoeff table to a new grid of wavenumbers (inplace).

Parameters:
  • wnedges (array, np.ndarray) – Edges of the new bins of wavenumbers (cm-1) onto which the kcoeff should be binned down. if you want Nwnew bin in the end, wnedges.size must be Nwnew+1 wnedges[0] should be greater than self.wnedges[0] (JL20 not sure anymore) wnedges[-1] should be lower than self.wnedges[-1]

  • weights (array, np.ndarray, optional) – Desired weights for the resulting Ktable.

  • ggrid (array, np.ndarray, optional) – Desired g-points for the resulting Ktable. Must be consistent with provided weights. If not given, they are taken at the midpoints of the array given by the cumulative sum of the weights

remap_g(ggrid=None, weights=None)[source]

Method to resample a kcoeff table to a new g grid (inplace).

Parameters:
  • ggrid (array, np.ndarray) – New grid of abcissas for quadrature

  • weights (array, np.ndarray) – New grid of weights