Contributions (taurex.contributions)

Classes related to the computation of the optical depth

Base Contribution

Base contribution classes and functions for computing optical depth

class Contribution(name)[source]

Bases: taurex.data.fittable.Fittable, taurex.log.logger.Logger, taurex.output.writeable.Writeable, taurex.data.citation.Citable

Abstract class

The base class for modelling contributions to the optical depth. By default this handles contributions from cross-sections. If the type of contribution being implemented is a sigma-type like the form given in contribute_tau() then To function in Taurex3, it only requires the concrete implementation of:

Different forms may require reimplementing contribute() as well as prepare()

Parameters

name (str) – Identifier of the contribution.

build(model)[source]

Called during forward model build phase Does nothing by default

Parameters

model (ForwardModel) – Forward model

contribute(model, start_layer, end_layer, density_offset, layer, density, tau, path_length=None)[source]

Computes an integral for a single layer for the optical depth.

Parameters
  • model (ForwardModel) – A forward model

  • start_layer (int) – Lowest layer limit for integration

  • end_layer (int) – Upper layer limit of integration

  • density_offset (int) – offset in density layer

  • layer (int) – atmospheric layer being computed

  • density (array) – density profile of atmosphere

  • tau (array) – optical depth to store result

  • path_length (array) – integration length

finalize(model, tau)[source]

Called in the last phase of the calculation, after the optical depth has be completely computed.

classmethod input_keywords()[source]
property name

Name of the contribution. Identifier for plots

property order

Computational order. Lower numbers are given higher priority and are computed first.

Returns

Order of computation

Return type

int

prepare(model, wngrid)[source]

Used to prepare the contribution for the calculation. Called before the forward model performs the main optical depth calculation. Default behaviour is to loop through prepare_each() and sum all results into a single cross-section.

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

prepare_each(model, wngrid)[source]

Requires implementation

Used to prepare each component of the contribution. For context when the main taurex program is run with the option each spectra is the component for the contribution. For cross-section based contributions, the components are each molecule Should yield the name of the component and the component itself

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

Yields

component (tuple of type (str, array)) – Name of component and component itself

property sigma
write(output)[source]

Writes contribution class and arguments to file

Parameters

output (Output) –

contribute_tau[source]

Generic cross-section integration function for tau, numba-fied for performance.

This has the form:

\[\tau_{\lambda}(z) = \int_{z_{0}}^{z_{1}} \sigma(z') \rho(z') dz',\]

where \(z\) is the layer, \(z_0\) and \(z_1\) are startK and endK respectively. \(\sigma\) is the weighted cross-section sigma. \(rho\) is the density and \(dz'\) is the integration path length path

Parameters
  • startK (int) – starting layer in integration

  • endK (int) – last layer in integration

  • density_offset (int) – Which part of the density profile to start from

  • sigma (array) – cross-section

  • density (array_like) – density profile of atmosphere

  • path (array_like) – path-length or altitude gradient

  • nlayers (int) – Total number of layers (unused)

  • ngrid (int) – total number of grid points

  • layer (int) – Which layer we currently on

Returns

tau – optical depth (well almost you still need to do exp(-tau) yourself)

Return type

array_like

Absorption

class AbsorptionContribution[source]

Bases: taurex.contributions.contribution.Contribution

build(model)[source]

Called during forward model build phase Does nothing by default

Parameters

model (ForwardModel) – Forward model

contribute(model, start_horz_layer, end_horz_layer, density_offset, layer, density, tau, path_length=None)[source]

Computes an integral for a single layer for the optical depth.

Parameters
  • model (ForwardModel) – A forward model

  • start_layer (int) – Lowest layer limit for integration

  • end_layer (int) – Upper layer limit of integration

  • density_offset (int) – offset in density layer

  • layer (int) – atmospheric layer being computed

  • density (array) – density profile of atmosphere

  • tau (array) – optical depth to store result

  • path_length (array) – integration length

finalize(model)[source]

Called in the last phase of the calculation, after the optical depth has be completely computed.

classmethod input_keywords()[source]
prepare(model, wngrid)[source]

Used to prepare the contribution for the calculation. Called before the forward model performs the main optical depth calculation. Default behaviour is to loop through prepare_each() and sum all results into a single cross-section.

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

prepare_each(model, wngrid)[source]

Requires implementation

Used to prepare each component of the contribution. For context when the main taurex program is run with the option each spectra is the component for the contribution. For cross-section based contributions, the components are each molecule Should yield the name of the component and the component itself

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

Yields

component (tuple of type (str, array)) – Name of component and component itself

property sigma

CIA

class CIAContribution(cia_pairs=None)[source]

Bases: taurex.contributions.contribution.Contribution

Computes the contribution to the optical depth occuring from collisionally induced absorption.

Parameters

cia_pairs (list of str) – list of molecule pairs of the form mol1-mol2 e.g. H2-He

property ciaPairs

Returns list of molecular pairs involved

Returns

Return type

list of str

contribute(model, start_layer, end_layer, density_offset, layer, density, tau, path_length=None)[source]

Computes an integral for a single layer for the optical depth.

Parameters
  • model (ForwardModel) – A forward model

  • start_layer (int) – Lowest layer limit for integration

  • end_layer (int) – Upper layer limit of integration

  • density_offset (int) – offset in density layer

  • layer (int) – atmospheric layer being computed

  • density (array) – density profile of atmosphere

  • tau (array) – optical depth to store result

  • path_length (array) – integration length

classmethod input_keywords()[source]
prepare_each(model, wngrid)[source]

Computes and weighs cross-section for a single pair of molecules

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

Yields

component (tuple of type (str, array)) – Molecular pair and the weighted cia opacity.

write(output)[source]

Writes contribution class and arguments to file

Parameters

output (Output) –

contribute_cia[source]

Collisionally induced absorption integration function

This has the form:

\[\tau_{\lambda}(z) = \int_{z_{0}}^{z_{1}} \sigma(z') \rho(z')^{2} dz',\]

where \(z\) is the layer, \(z_0\) and \(z_1\) are startK and endK respectively. \(\sigma\) is the weighted cross-section sigma. \(rho\) is the density and \(dz'\) is the integration path length path

Parameters
  • startK (int) – starting layer in integration

  • endK (int) – last layer in integration

  • density_offset (int) – Which part of the density profile to start from

  • sigma (array) – cross-section

  • density (array_like) – density profile of atmosphere

  • path (array_like) – path-length or altitude gradient

  • nlayers (int) – Total number of layers (unused)

  • ngrid (int) – total number of grid points

  • layer (int) – Which layer we currently on

Returns

tau – optical depth (well almost you still need to do exp(-tau) yourself)

Return type

array_like

Rayleigh

class RayleighContribution[source]

Bases: taurex.contributions.contribution.Contribution

Computes contribution from Rayleigh scattering

BIBTEX_ENTRIES = ['\n @book{cox_allen_rayleigh,\n title={Allen’s astrophysical quantities},\n author={Cox, Arthur N},\n year={2015},\n publisher={Springer}\n }\n ']
classmethod input_keywords()[source]
prepare_each(model, wngrid)[source]

Computes the weighted opacity due to rayleigh scattering for any possible molecules within atmosphere.

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

Yields

component (tuple of type (str, array)) – Name of scattering molecule and the weighted rayeligh opacity.

SimpleClouds

class SimpleCloudsContribution(clouds_pressure=1000.0)[source]

Bases: taurex.contributions.contribution.Contribution

Optically thick cloud deck up to a certain height

These have the form:

\[\begin{split}\tau(\lambda,z) = \begin{cases} \infty & \quad \text{if } P(z) >= P_{0}\\ 0 & \quad \text{if } P(z) < P_{0} \end{cases}\end{split}\]

Where \(P_{0}\) is the pressure at the top of the cloud-deck

Parameters

clouds_pressure (float) – Pressure at top of cloud deck

property cloudsPressure

Cloud top pressure in Pascal

contribute(model, start_layer, end_layer, density_offset, layer, density, tau, path_length=None)[source]

Computes an integral for a single layer for the optical depth.

Parameters
  • model (ForwardModel) – A forward model

  • start_layer (int) – Lowest layer limit for integration

  • end_layer (int) – Upper layer limit of integration

  • density_offset (int) – offset in density layer

  • layer (int) – atmospheric layer being computed

  • density (array) – density profile of atmosphere

  • tau (array) – optical depth to store result

  • path_length (array) – integration length

classmethod input_keywords()[source]
property order

Computational order. Lower numbers are given higher priority and are computed first.

Returns

Order of computation

Return type

int

prepare_each(model, wngrid)[source]

Returns an absorbing cross-section that is infinitely absorping up to a certain height

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

Yields

component (tuple of type (str, array)) – Clouds and opacity array.

write(output)[source]

Writes contribution class and arguments to file

Parameters

output (Output) –

Mie Scattering (BH)

Warning

This is no longer available in the base TauREx 3 since version 3.1. To use this you must install the taurex_bhmie plugin.

Mie Scattering (Lee)

class LeeMieContribution(lee_mie_radius=0.01, lee_mie_q=40, lee_mie_mix_ratio=1e-10, lee_mie_bottomP=-1, lee_mie_topP=-1)[source]

Bases: taurex.contributions.contribution.Contribution

Computes Mie scattering contribution to optica depth Formalism taken from: Lee et al. 2013, ApJ, 778, 97

Parameters
  • lee_mie_radius (float) – Particle radius in um

  • lee_mie_q (float) – Extinction coefficient

  • lee_mie_mix_ratio (float) – Mixing ratio in atmosphere

  • lee_mie_bottomP (float) – Bottom of cloud deck in Pa

  • lee_mie_topP (float) – Top of cloud deck in Pa

BIBTEX_ENTRIES = ['\n @article{Lee_2013,\n doi = {10.1088/0004-637x/778/2/97},\n url = {https://doi.org/10.1088%2F0004-637x%2F778%2F2%2F97},\n year = 2013,\n month = {nov},\n publisher = {{IOP} Publishing},\n volume = {778},\n number = {2},\n pages = {97},\n author = {Jae-Min Lee and Kevin Heng and Patrick G. J. Irwin},\n title = {{ATMOSPHERIC} {RETRIEVAL} {ANALYSIS} {OF} {THE} {DIRECTLY} {IMAGED} {EXOPLANET} {HR} 8799b},\n journal = {The Astrophysical Journal},\n abstract = {Directly imaged exoplanets are unexplored laboratories for the application of the spectral and temperature retrieval method, where the chemistry and composition of their atmospheres are inferred from inverse modeling of the available data. As a pilot study, we focus on the extrasolar gas giant HR\xa08799b, for which more than 50 data points are available. We upgrade our non-linear optimal estimation retrieval method to include a phenomenological model of clouds that requires the cloud optical depth and monodisperse particle size to be specified. Previous studies have focused on forward models with assumed values of the exoplanetary properties; there is no consensus on the best-fit values of the radius, mass, surface gravity, and effective temperature of HR\xa08799b. We show that cloud-free models produce reasonable fits to the data if the atmosphere is of super-solar metallicity and non-solar elemental abundances. Intermediate cloudy models with moderate values of the cloud optical depth and micron-sized particles provide an equally reasonable fit to the data and require a lower mean molecular weight. We report our best-fit values for the radius, mass, surface gravity, and effective temperature of HR\xa08799b. The mean molecular weight is about 3.8, while the carbon-to-oxygen ratio is about unity due to the prevalence of carbon monoxide. Our study emphasizes the need for robust claims about the nature of an exoplanetary atmosphere to be based on analyses involving both photometry and spectroscopy and inferred from beyond a few photometric data points, such as are typically reported for hot Jupiters.}\n }\n ']
classmethod input_keywords()[source]
property mieBottomPressure

Pressure at bottom of cloud deck in Pa

property mieMixing

Mixing ratio in atmosphere

property mieQ

Extinction coefficient

property mieRadius

Particle radius in um

property mieTopPressure

Pressure at top of cloud deck in Pa

prepare_each(model, wngrid)[source]

Computes and weights the mie opacity for the pressure regions given

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

Yields

component (tuple of type (str, array)) – Lee and the weighted mie opacity.

write(output)[source]

Writes contribution class and arguments to file

Parameters

output (Output) –

Mie Scattering (Flat)

class FlatMieContribution(flat_mix_ratio=1e-10, flat_bottomP=-1, flat_topP=-1)[source]

Bases: taurex.contributions.contribution.Contribution

Computes a flat absorption contribution across all wavelengths to the optical depth

Parameters
  • flat_mix_ratio (float) – Opacity value

  • flat_bottomP (float) – Bottom of absorbing region in Pa

  • flat_topP (float) – Top of absorbing region in Pa

classmethod input_keywords()[source]
property mieBottomPressure

Pressure at bottom of absorbing region in Pa

property mieMixing

Opacity of absorbing region in m2

property mieTopPressure

Pressure at top of absorbing region in Pa

prepare_each(model, wngrid)[source]

Computes and flat absorbing opacity for the pressure regions given

Parameters
  • model (ForwardModel) – Forward model

  • wngrid (array) – Wavenumber grid

Yields

component (tuple of type (str, array)) – Flat and the weighted mie opacity.

write(output)[source]

Writes contribution class and arguments to file

Parameters

output (Output) –