exo_k.util.radiation

Created on Jun 20 2019

@author: jeremy leconte

Library of useful functions for radiative transfer calculations

Module Contents

exo_k.util.radiation.PLANCK_CST1MKS
exo_k.util.radiation.PLANCK_CST1
exo_k.util.radiation.PLANCK_CST2
exo_k.util.radiation.PLANCK_CST3
exo_k.util.radiation.PLANCK_CST1_lamb
exo_k.util.radiation.Bnu(nu, T)[source]

Computes the Planck law in wavenumber domain.

Parameters:
  • nu (float) – The wavenumber in cm^-1.

  • T (float) – Temperature (K).

Returns:

  • float

    Bnu in W/m^2/str/cm^-1 (not SI units).

  • sigma*T^4/PI=Int(Bnu dnu) where nu is expressed in cm^-1

exo_k.util.radiation.Bnu_integral(nu_edges, T)[source]

Computes the integral of the Planck function in wavenumber bins.

Uses scipy.integrate.quad. Quite slow.

Parameters:
  • nu_edges (array, np.ndarray) – Edges of the wavenumber bins in cm^-1.

  • T (float) – Temperature (K).

Returns:

  • array

    Bnu_integral in W/m^2/str.

  • sigma*T^4/PI=sum of Bnu_integral

exo_k.util.radiation.Bnu_integral_num(nu_edges, T, n=30)[source]

Computes the integral of the Planck function in wavenumber bins.

Uses a series expansion of the integral derived in “PLANCK FUNCTIONS AND INTEGRALS; METHODS OF COMPUTATION by Thomas E. Michels Goddard Space Flight Center Greenbelt, Md. MARCH 1968”

Parameters:
  • nu_edges (array, np.ndarray) – Edges of the wavenumber bins in cm^-1.

  • T (float) – Temperature (K).

  • n (int) – number of terms to take into account in the black body calculation.

Returns:

  • Array of shape (T_array.size,nu_edges.size-1)

    Bnu_integral_num: Integral of the source function at temperatures T_array inside the bins in W/m^2/str

  • sigma*T^4/PI=sum of Bnu_integral_num

exo_k.util.radiation.Bnu_integral_array(nu_edges, T_array, Nw, Nt, n=30)[source]

Computes the integral of the Planck function in wavenumber bins.

Uses a series expansion of the integral derived in “PLANCK FUNCTIONS AND INTEGRALS; METHODS OF COMPUTATION by Thomas E. Michels Goddard Space Flight Center Greenbelt, Md. MARCH 1968”

Parameters:
  • nu_edges (array, np.ndarray) – Edges of the wavenumber bins in cm^-1.

  • T (array, np.ndarray) – Array of temperature (K).

  • Nw (int) – Number of spectral bins.

  • Nt (int) – Number of temperatures.

  • n (int, optional) – number of terms to take into account in the black body calculation.

Returns:

  • Array of shape (T_array.size,nu_edges.size-1)

    Bnu_integral_num: Integral of the source function at temperatures T_array inside the bins in W/m^2/str

  • sigma*T^4/PI=sum of Bnu_integral_num

exo_k.util.radiation.dBnudT_array(nu, T_array, Nw, Nt)[source]

Computes the derivative of the Planck function with respect to temperature.

Parameters:
  • nu (array, np.ndarray) – Wavenumbers in cm^-1.

  • T (array, np.ndarray) – Array of temperature (K).

  • Nw (int) – Number of spectral bins.

  • Nt (int) – Number of temperatures.

Returns:

Array of shape (T_array.size,nu.size)

dBnudT: derivative of the Planck function at temperatures T_array.

exo_k.util.radiation.Bmicron(lamb, T)[source]

Computes the Planck law in wavelength domain.

Parameters:
  • lamb (float) – The wavelength in microns.

  • T (float) – Temperature (K).

Returns:

float

Bmicron in W/m^2/str/micron (not SI units).

sigma*T^4/PI=Int(Bmicron d lamb) where lamb is expressed in microns

exo_k.util.radiation.BBnu_Stellar_Spectrum(nu, T, Flux)[source]

Computes the outgoing spectral flux of a black body at temperature T so that its bolometric flux is equal to the input ‘Flux’

Parameters:
  • nu (float) – The wavenumber in cm^-1.

  • T (float) – Temperature (K).

  • Flux (float) – bolometric flux for rescaling.

Returns:

  • float

    BBnu_Stellar_Spectrum is in W/m^2/cm^-1 (not SI units).

  • sigma*T^4=Int(Bnu dnu) where nu is expressed in cm^-1

exo_k.util.radiation.BBmicron_Stellar_Spectrum(lamb, T, Flux)[source]

Computes the outgoing spectral flux of a black body at temperature T so that its bolometric flux is equal to the input ‘Flux’

Parameters:
  • lamb (float) – The wavelength in microns.

  • T (float) – Temperature (K).

  • Flux (float) – bolometric flux for rescaling.

Returns:

  • float

    BBmicron_Stellar_Spectrum is in W/m^2/micron (not SI units).

  • sigma*T^4=Int(Bmicron d lamb) where lamb is expressed in microns.

exo_k.util.radiation.wavelength_grid_R(lambda_min, lambda_max, R)[source]

Creates a wavelength grid starting from lambda_min and ending at lambda_max with a resolution of R (roughly).

Parameters:
  • lambda_min/lambda_max (int or float) – Min/max wavelength (in micron).

  • R (int or float) – Resolution.

Returns:

array

Grid of wavelength.

exo_k.util.radiation.wavenumber_grid_R(nu_min, nu_max, R)[source]

Creates a wavenumber grid starting from nu_min and ending at nu_max with a resolution of R (roughly).

Parameters:
  • nu_min/nu_max (int or float) – Min/max wavenumber (in cm^-1).

  • R (int or float) – Resolution.

Returns:

array

Grid of wavenumber.

exo_k.util.radiation.rad_prop_corrk(dcol_density, opacity_prof, mu0)[source]

Computes the optical depth of each of the radiative layers for the opacity given for every wavelength and g point.

Parameters:
  • dcol_density (array, np.ndarray) – Column density per radiative layer (molecules/m^2/layer).

  • opacity_prof (array, np.ndarray) – effective cross section per molecule (m^2/molecule).

  • mu0 (cosine of the equivalent zenith angle for the rays.) – e.g. Cos(60deg)=0.5 is chosen (See, for example, chap 4 of Pierrehumbert 2010)

Returns:

tau: Array

cumulative optical depth from the top (with zeros at the top of the first radiative layer)

dtau: Array

optical depth of each radiative layer for each wavenumber.

exo_k.util.radiation.rad_prop_xsec(dcol_density, opacity_prof, mu0)[source]

Computes the optical depth of each of the radiative layers for the opacity given for every wavelength.

Parameters:
  • dcol_density (array, np.ndarray) – Column density per radiative layer (molecules/m^2/layer).

  • opacity_prof (array, np.ndarray) – effective cross section per molecule (m^2/molecule).

  • mu0 (cosine of the equivalent zenith angle for the rays.) – e.g. Cos(60deg)=0.5 is chosen (See, for example, chap 4 of Pierrehumbert 2010)

Returns:

  • tau: Array

    cumulative optical depth from the top (with zeros at the top of the first layer).

    dtau: Array

    optical depth of each radiative layer for each wavenumber.

  • The 1/mu0 factor is taken account so that these are the depths along the ray.

exo_k.util.radiation.path_integral_corrk(Nlay, Nw, Ng, tangent_path, density_prof, opacity_prof, weights)[source]

Computes the transmittance for each layer of the atmosphere for k-coefficients.

Parameters:
  • Nlay (int) – Number of layers.

  • Nw (int) – Number of wavenumber bins.

  • Ng (int) – Number of gauss points.

  • tangent_path (List of arrays) – Triangular array of tangent path for each layer.

  • density_prof (array, np.ndarray) – Array with the number density of molecules in the atmosphere.

  • opacity_prof (array, np.ndarray) – Effective cross section of the atmsophere for each (layer,wavenumber,g-point).

  • weights (array, np.ndarray) – Weights for the quadrature in g-space.

Returns:

transmittance: array, np.ndarray

Transmittance for each layer and wavenumber bin (Exp(-tau_sigma)).

exo_k.util.radiation.path_integral_xsec(Nlay, Nw, tangent_path, density_prof, opacity_prof)[source]

Computes the transmittance for each layer of the atmosphere from cross sections.

Parameters:
  • Nlay (int) – Number of layers.

  • Nw (int) – Number of wavenumber bins.

  • tangent_path (List of arrays) – Triangular array of tangent path for each layer.

  • density_prof (array, np.ndarray) – Array with the number density of molecules in the atmosphere.

  • opacity_prof (array, np.ndarray) – Effective cross section of the atmsophere for each (layer,wavenumber).

Returns:

transmittance: array, np.ndarray

Transmittance for each layer and wavenumber bin (Exp(-tau_sigma)).