:py:mod:`exo_k.ktable5d` ======================== .. py:module:: exo_k.ktable5d .. autoapi-nested-parse:: @author: jeremy leconte Module Contents --------------- .. py:class:: Ktable5d(*filename_filters, filename=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) Bases: :py:obj:`exo_k.data_table.Data_table` A class that handles tables of k-coefficients with a variable gas. Based on the Data_table class that handles basic operations common to Ktable and Xtable. This class is specifically designed to deal with LMDZ type ktable where there is a variable gas. Initializes k coeff table with variable gas and supporting data from various sources (see below by order of precedence) :param filename: Relative or absolute name of the file to be loaded. :type filename: :class:`str (optional)` :param path: 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. :type path: :class:`str` If there is no input, just creates an empty object to be filled later See :class:`~exo_k.ktable.Ktable` __init__ mehthod for documentation on `p_unit`, `file_p_unit`, `kdata_unit`, `file_kdata_unit`, `remove_zeros`, `search_path`, and `mol` arguments. .. py:method:: read_hdf5(filename=None, mol=None, wn_range=None, wl_range=None) Initializes k coeff table and supporting data from an Exomol hdf5 file :param filename: Name of the input hdf5 file :type filename: :class:`str` .. py:method:: write_hdf5(filename, compression='gzip', compression_level=9, kdata_unit=None, p_unit=None) Saves data in a hdf5 format :param filename: Name of the file to be created and saved :type filename: :class:`str` .. py:method:: read_LMDZ(path=None, res=None, band=None, mol=None) Initializes k coeff table and supporting data from a .dat file in a gcm friendly format. Units are assumed to be cm^2 for kdata and mbar for pressure. :param path: Name of the directory with the various input files :type path: :class:`str` :param res: "IRxVI" where IR and VI are the numbers of bands in the infrared and visible of the k table to load. :type res: :class:`str` :param band: "IR" or "VI" to specify which band to load. :type band: :class:`str` .. py:method:: write_LMDZ(path, band='IR', fmt='%22.15e', write_only_metadata=False) Saves data in a LMDZ friendly format. The gcm requires p in mbar and kdata in cm^2/molec. The conversion is done automatically. :param path: Name of the directory to be created and saved, the one that will contain all the necessary files :type path: :class:`str` :param band: The band you are computing: 'IR' or 'VI' :type band: :class:`str` :param fmt: Fortran format for the corrk file. :type fmt: :class:`str` :param write_only_metadata: If `True`, only supporting files are written (T.dat, p.dat, etc.) :type write_only_metadata: :class:`bool`, *optional* .. py:method:: hires_to_ktable(path=None, filename_grid=None, logpgrid=None, tgrid=None, xgrid=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', remove_zeros=False, **kwargs) Computes a k coeff table from :class:`~exo_k.util.hires_spectrum.Hires_spectrum` objects. see :func:`exo_k.ktable.Ktable.hires_to_ktable` method for details on the arguments and options. :param xgrid: Input grid in vmr of the variable gas. Needed for a Ktable5d. :type xgrid: :class:`array`, :class:`np.ndarray` .. 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. .. py:method:: setup_interpolation(log_interp=None) Creates interpolating functions to be called later on. and loads it as attribute (inplace). .. py:method:: set_kdata(new_kdata) Changes kdata (inplace). this is preferred to directly accessing kdata because one could forget to run setup_interpolation(). :param new_kdata: New array of kdata. :type new_kdata: :class:`array`, :class:`np.ndarray` .. py:method:: interpolate_kdata(logp_array=None, t_array=None, x_array=None, log_interp=None, logp_interp=True, wngrid_limit=None, **kwargs) interpolate_kdata interpolates the kdata at on a given temperature and log pressure profile. :param logp_array: log 10 pressure array to interpolate to :type logp_array: :class:`Array` :param t_array: Temperature array to interpolate to :type t_array: :class:`Array`, :class:`same size a logp_array` :param x_array: vmr of variable gas array to interpolate to :type x_array: :class:`Array`, :class:`same size a logp_array` :param If floats are given: :param they are interpreted as arrays of size 1.: :param wngrid_limit: if an array is given, interpolates only within this array :type wngrid_limit: :class:`list` or :class:`array`, *optional* :param log_interp: Dummy variable to be consistent with interpolate_kdata in data_table. Whether the interpolation is linear in kdata or in log(kdata) is actually controlled by self._settings._log_interp but only when the ktable is loaded. If you change that after the loading, you should rerun setup_interpolation(). :type log_interp: :class:`bool`, :class:`dummy` :returns: array of shape (logp_array.size, self.Nw , self.Ng) The interpolated kdata. .. py:method:: remap_logPT(logp_array=None, t_array=None, x_array=None) remap_logPT re-interpolates the kdata on a new temprature and log pressure grid (inplace). If an `x_array` is specified, the interpolation occurs along the vmr axis as well. :param logp_array: log 10 pressure array to interpolate to :type logp_array: :class:`Array` :param t_array: temperature array to interpolate to :type t_array: :class:`Array` :param x_array: vmr of variable gas array to interpolate to :type x_array: :class:`Array` .. warning:: Whether the interpolation is linear in kdata or in log(kdata) is controlled by self._settings._log_interp but only when the ktable is loaded. If you change that after the loading, you should rerun setup_interpolation(). .. py:method:: copy(cp_kdata=True) Creates a new instance of :class:`Ktable5d` object and (deep) copies data into it :param cp_kdata: If false, the kdata table is not copied and only the structure and metadata are. :type cp_kdata: :class:`bool`, *optional* :returns: :class:`Ktable` A new :class:`Ktable5d` instance with the same structure as self. .. py:method:: gindex(g) Finds the index corresponding to the given g .. py:method:: xindex(x) Finds the index corresponding to the given x .. py:method:: spectrum_to_plot(p=1e-05, t=200.0, x=1.0, g=None) provide the spectrum for a given point to be plotted :param p: Pressure (Ktable pressure unit) :type p: :class:`float` :param t: Temperature(K) :type t: :class:`float` :param g: Gauss point :type g: :class:`float` :param x: Mixing ratio of the species :type x: :class:`float` .. py:method:: plot_distrib(ax, p=1e-05, t=200.0, wl=1.0, x=1.0, xscale=None, yscale='log', **kwarg) Plot the distribution for a given point :param p: Pressure (Ktable pressure unit) :type p: :class:`float` :param t: Temperature(K) :type t: :class:`float` :param wl: Wavelength (micron) :type wl: :class:`float` .. py:method:: combine_with(other, x_self=None, x_other=None, **kwargs) Method to create a new :class:`Ktable5d` where the kdata of 'self' are randomly mixed with 'other' (that must be a :class:`Ktable`). The main purpose is to add the opacity of a trace species to the background gas of the :class:`Ktable5d` instance. .. warning:: Because: * the opacity from the background and variable gases cannot be isolated, * The values of the array for the vmr of the variable gas (self.xgrid) are not modified (diluted), the treatment here is valid only if `x_other` << 1. For this reason, `x_self` should be either left to None, or 1-`x_other` depending exactly what you want to do. But if you put `x_self`<<1, you are on your own. :param other: A :class:`Ktable` object to be mixed with. Dimensions should be the same as self (except for xgrid). :type other: :class:`~exo_k.ktable.Ktable` :param x_self: Volume mixing ratio of self. :type x_self: :class:`float only`, *optional* :param x_other: Volume mixing ratio of the species to be mixed with (other). :type x_other: :class:`float` or :class:`array`, *optional* If either x_self or x_other are set to `None` (default), the cross section of the species in question are considered to be already normalized with respect to the mixing ratio. :returns: :class:`Ktable5d` A new Ktable5d with the opacity of the new species added. .. py:method:: bin_down(wnedges=None, weights=None, ggrid=None, remove_zeros=False, num=300, use_rebin=False, write=0) Method to bin down a kcoeff table to a new grid of wavenumbers (inplace). :param wnedges: 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] :type wnedges: :class:`array`, :class:`np.ndarray` :param weights: Desired weights for the resulting Ktable. :type weights: :class:`array`, :class:`np.ndarray`, *optional* :param ggrid: 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 :type ggrid: :class:`array`, :class:`np.ndarray`, *optional* .. py:method:: clip_spectral_range(wn_range=None, wl_range=None) Limits the data to the provided spectral range (inplace): * Wavenumber in cm^-1 if using wn_range argument * Wavelength in micron if using wl_range .. py:method:: extend_spectral_range(**kwargs) Extends the spectral range of an existing table (inplace). The new bins are filled with zeros (except if remove_zeros=True). See :func:`exo_k.data_table.Data_table.extend_spectral_range` method for details on the arguments and options. .. py:method:: remap_g(ggrid=None, weights=None) Method to resample a kcoeff table to a new g grid (inplace). :param ggrid: New grid of abcissas for quadrature :type ggrid: :class:`array`, :class:`np.ndarray` :param weights: New grid of weights :type weights: :class:`array`, :class:`np.ndarray` .. py:function:: read_Qdat(filename) Reads Q.dat files LMDZ style and extract the vmr grid. :param filename: Path to file to read. :type filename: :class:`str` :returns: background_mol_names: list list of names of molecules in background gas var_mol: str Name of variable molecule Nx: int Size of xgrid xgrid: array, np.ndarray grid of vmr for variable gas