Optical Properties

Module: sasktran.opticalproperty

Classes that describe the optical properties of some atmospheric Species. These objects are used to create Species objects.

See our Useful Shorthands for optical properties that we use a lot.

Base Classes

class sasktran.OpticalProperty(name: str)

Bases: object

Generic optical property container that is applicable to a wide variety of optical properties.

Parameters:

name (str) – Name of the optical property to create

Examples

>>> from sasktran import OpticalProperty, Climatology
>>> opt_prop = OpticalProperty('O3_OSIRISRES')
>>> print(opt_prop)
SasktranIF Optical Property: O3_OSIRISRES
>>> atmospheric_state = Climatology('msis90')
>>> opt_prop.calculate_cross_sections(atmospheric_state, latitude=0, longitude=0, altitude=20000, mjd=54372,                                          wavelengths=[350, 600])
CrossSections(wavelengths=array([350, 600]), absorption=array([  7.09857544e-23,   5.27719024e-21]), scattering=array([ 0.,  0.]), total=array([  7.09857544e-23,   5.27719024e-21]))
calculate_cross_sections(atmospheric_state: Climatology, latitude: float, longitude: float, altitude: float, mjd: float, wavelengths: array)

Calculates absorption and scattering cross sections for the given optical property with a specific atmospheric state and location.

Parameters:
  • atmospheric_state (sasktran.Climatology) – Atmospheric state (sometime called background climatology, corresponds to sasktran.Atmosphere.atmospheric_state). Typically the background climatology must support temperature and pressure, e.g. msis90. Some optical properties may not need a background climatology, for example Mie aerosols, however one must still be passed in for the function to operate.

  • latitude (float) – Latitude in degrees.

  • longitude (float) – Longitude in degrees.

  • altitude (float) – Altitude in meters.

  • mjd (float) – Modified Julian Date.

  • wavelengths (np.array) – Wavelengths in nm.

Returns:

NamedTuple with fields wavelengths, absorption, scattering, total. wavelengths field matches the input wavelengths, absorption and scattering are the absorption and scattering cross sections respectively in \(\mathrm{cm^2}\) and total is the sum of the two.

Return type:

CrossSections

Examples

>>> from sasktran import OpticalProperty, Climatology
>>> opt_prop = OpticalProperty('O3_OSIRISRES')
>>> print(opt_prop)
SasktranIF Optical Property: O3_OSIRISRES
>>> atmospheric_state = Climatology('msis90')
>>> opt_prop.calculate_cross_sections(atmospheric_state, latitude=0, longitude=0, altitude=20000, mjd=54372,                                              wavelengths=[350, 600])
CrossSections(wavelengths=array([350, 600]), absorption=array([  7.09857544e-23,   5.27719024e-21]), scattering=array([ 0.,  0.]), total=array([  7.09857544e-23,   5.27719024e-21]))
calculate_phase_matrix(atmospheric_state: Climatology, latitude: float, longitude: float, altitude: float, mjd: float, wavelengths: array, cosine_scatter_angles: array)

Calculates the scattering phase matrix for the given optical property with a specific atmospheric state and location.

Parameters:
  • atmospheric_state (sasktran.Climatology) – Atmospheric state (sometime called background climatology, corresponds to sasktran.Atmosphere.atmospheric_state). Typically the background climatology must support temperature and pressure, e.g. msis90. Some optical properties may not need a background climatology, for example Mie aerosols, however one must still be passed in for the function to operate.

  • latitude (float) – Latitude in degrees.

  • longitude (float) – Longitude in degrees.

  • altitude (float) – Altitude in meters.

  • mjd (float) – Modified Julian Date.

  • wavelengths (np.array) – Wavelengths in nm. Shape (N_wavel,)

  • cosine_scatter_angles (np.array) – Array of cosine of the scattering angle. Shape (N_angle,)

Returns:

phase_matrix – Array of shape (N_wavel, N_angle, 4, 4) of phase matrices.

Return type:

np.ndarray

Examples

>>> import sasktran as sk
>>> opt_prop = sk.Rayleigh()
>>> atmospheric_state = sk.MSIS90()
>>> opt_prop.calculate_phase_matrix(atmospheric_state, latitude=0, longitude=0, altitude=20000, mjd=54372,                                              wavelengths=[350], cosine_scatter_angles=[1])
array([[[[ 1.47818063, -0.        ,  0.        ,  0.        ],
         [-0.        ,  1.4345419 ,  0.        ,  0.        ],
         [ 0.        ,  0.        ,  1.4345419 ,  0.        ],
         [ 0.        ,  0.        ,  0.        ,  1.39018959]]]])
class sasktran.OpticalPropertyConvolved(optical_property: OpticalProperty, psf_wavelength: ndarray, psf: ndarray, wavel_spacing: float = 0.01, num_stdev: float = 4, output_spacing: ndarray | None = None)

Bases: OpticalProperty

Generic optical property container that is applicable to a wide variety of optical properties.

Parameters:
  • name (str) – Name of the optical property to create

  • psf_wavelength (np.ndarray) – array of wavelengths in nanometers where the point spread function is defined.

  • psf (np.ndarray) – array of point spread function values in nanometers. Same length as psf_wavelength

  • wavel_spacing (scalar, optional) – Spacing to calculate the high resolution wavelengths on in nanometers, default is 0.01nm

  • num_stdev (scalar, optional) – Number of standard deviations to use in the gaussian for convolution, default is 4

  • output_spacing (scalar, string optional) – Spacing to calculate the convolved cross section on. Either a scalar to use a fixed spacing, or None to use the wavelengths of the point spread function (wavel parameter). Default is none.

Examples

>>> from sasktran import OpticalPropertyConvolved, MSIS90, O3DBM
>>> import numpy as np
>>> convolved_dbm = OpticalPropertyConvolved(O3DBM(), psf_wavelength=np.linspace(250, 1000, 50),                                                 psf=np.linspace(0.5, 10, 50))
>>> print(convolved_dbm)
SasktranIF Optical Property: O3_DBM_Convolved
>>> convolved_dbm.calculate_cross_sections(MSIS90(), latitude=0, longitude=0, altitude=20000, mjd=54372,                                               wavelengths=[350, 600])
CrossSections(wavelengths=array([350, 600]), absorption=array([  2.83084701e-22,   5.01481379e-21]), scattering=array([ 0.,  0.]), total=array([  2.83084701e-22,   5.01481379e-21]))
static convolved_optical_property(hires_optprop: OpticalProperty, wavel: ndarray, psf: ndarray, wavel_spacing: float = 0.01, num_stdev: float = 4, output_spacing: ndarray = None)

Convolves down the hi-resolution cross sections and creates a user defined optical property. Convolution assumes the hi-resolution version has infinite resolution

Parameters:
  • hires_optprop (sasktran.OpticalProperty) – sasktran OpticalProperty that will be convolved down to desired resolution.

  • wavel (numpy array) – Wavelengths to calculate the new cross section for in nm

  • psf (numpy array) – Standard deviations of the gaussian to convolve by for each wavelength Can be the same size as wavel, or size 1 in which case the psf is assumed to be the same for all wavelengths

  • wavel_spacing (scalar, optional) – Spacing to calculate the high resolution wavelengths on in nanometers, default is 0.01nm

  • num_stdev (scalar, optional) – Number of standard deviations to use in the gaussian for convolution, default is 4

  • output_spacing (scalar,string optional) – Spacing to calculate the convolved cross section on. Either a scalar to use a fixed spacing, or None to use the wavelengths of the point spread function (wavel parameter). Default is none.

Returns:

optprop – Convolved optical property

Return type:

ISKOpticalProperty

Useful Shorthands

class sasktran.Rayleigh

Bases: OpticalProperty

An optical property object that provides Rayleigh molecular scattering in dry-air. The code closely follows the algorithm published by Bates 1984 and exactly replicates his cross-section calculations to the 4 significant digits in his Table 1. The cross-section is weighted to account for the different gas ratios in standard atmospheric composition. No attempt is made to track changes in CO2 composition. Note that water vapour effects are implicitly ignored as it only considers dry air. The only difference with Bates is that this object takes into account the tiny fraction of gas that is not N2, O2, Argon or CO2. Bates ignores this component while this object assumes it the residual gas with properties similar to Argon.

class sasktran.O3DBM

Bases: OpticalProperty

Tabulated high resolution cross-sections of O3 measured by Daumont, Brion and Malicet in the early 1990’s [1]. The wavelength range slightly varies with temperature but covers the entire UV to NIR region, from 194.50 nm to 830.00 nm at 0.01 to 0.02 nm resolution. The cross-section data were collected at 0.01-0.02 nm resolution and each wavelength/cross-section table varies in size from 22,052 to 63,501 entries. The data consists of 5 tables of wavelength versus cross-section for 5 temperatures.

Notes

Temerature Range

Measurements are provided at 5 temperatures covering typical stratospheric and tropospheric conditions:

218 K
228 K
243 K
273 K
295 K
Wavelength Range

The wavelength range of each temperature table is slightly different and is given below. Note that most of the temperature variation occurs in the Huggins band between 315 and 360 nm:

218K -> 194.50nm to 650.01nm
228K -> 194.50nm to 520.01nm
243K -> 194.50nm to 519.01nm
273K -> 299.50nm to 520.01nm
295K -> 195.00nm to 830.00nm

We looked into temperature interpolation and while DBM suggest that a quadratic interpolation scheme [3] they do not indicate an explicit technique. We tested several quadratic fitting routines and found that a truncated linear fit in temperature was visually more appealing than any of the quadratic fits and had none of the undesirable artifacts (excessive curvature etc.) that naturally arise with quadratic curve fitting. Consequently this object uses a truncated linear fit in temperature.

Data Source

These data are an exact replication of the data files:

O3_CRS_BDM_218K.dat
O3_CRS_BDM_228K.dat
O3_CRS_BDM_243K.dat
O3_CRS_BDM_273K.dat
O3_CRS_BDM_295K.dat

Data is from the IGACO site, http://igaco-o3.fmi.fi/ACSO/files/cross_sections. The files were copied on July 16-July 25 2012.

References

class sasktran.O3OSIRISRes

Bases: OpticalProperty

” These are the cross-sections used in the OSIRIS level 2 analysis for Saskmart V5.07. This table is based upon the cross-sections of Bogumil Orphal and Burrows and has been reduced to the nominal resolution of the OSIRIS instrument (approx. 1.0 nm). The cross-sections range from 270nm to 820 nm in 0.1 nm steps.

class sasktran.NO2Vandaele1998

Bases: OpticalProperty

Calculates the absorption cross section of NO2 molecules from 230 nm to 1000 nm at 220 K to 294 K following [1]

class sasktran.NO2OSIRISRes

Bases: OpticalProperty

Calculates the absorption cross section of NO2 molecules from 230 nm to 795 nm and 221K to 293K. The cross-sections have been reduced to the resolution of OSIRIS and these cross-sections have been used in the OSIRIS level 2 MART retrievals.

class sasktran.HITRANChemical(chemical_name: str, line_tolerance=None, max_line_strength=None, isotope_filter=None, use_cache=True)

Bases: OpticalProperty

Calculates the optical absorption and extinction of various atmospheric molecules using the Voigt line-shape and the HITRAN spectral line database. The object supports all of the HITRAN species specified in the HITRAN database file molparam.txt. This optical property requires additional configuration, see Configuration for more information.

Parameters:
  • chemical_name (str) – Chemical abbreviation of the molecule of interest.

  • isotope_filter (int, optional) – Allows the HITRAN object to load in just one isotope of the requested molecule. The value set must match one of the isotope labels used for the given molecule in the HITRAN database file, molparam.txt. Note that the code does not adjust the line strength but uses the line strength value as written in the HITRAN database. This means you may have to account for and/or remove the abundance automatically built into the HITRAN database line strength values.

  • line_tolerance (float, optional) – Allows the user to set the tolerance used to reject weak lines from the current micro-window as part of a speed optimization strategy. The default value is 0.0 which disables the optimization. Larger values speed up calculation of spectra but may result in choppy spectra at the smaller intensities, especially in extinction/absorption spectra which typically follow the log of the cross-section. A smaller value will reduce choppiness but increase computational speed. Only values greater than or equal to zero are acceptable. A value of 1.0E-09 is a good starting value for users wishing to investigate the effectiveness of the speed optimization.

  • max_line_strength (float, optional) – Allows the user to manually set the maximum line strength within a micro-window. By default the object will take this value from the strongest line in the micro-window. The value is used with the line tolerance to reject weak lines from spectral calculations. A value of zero will disable the manual setting and reinstate usage of the default.

  • use_cache (bool, optional) – If true then cross sections will be cached when used internally inside a radiative transfer calculation. The default is True, this should be left as True unless trying to calculate millions of wavelengths where memory starts to become an issue.

skif_object(**kwargs)
Returns:

Underlying ISKOpticalProperty object

Return type:

skif.ISKOpticalProperty

class sasktran.MieAerosol(particlesize_climatology: Climatology, species: str)

Bases: OpticalProperty

Specialized OpticalProperty which supports Mie Aerosol calculations.

Parameters:
  • particlesize_climatology (sasktran.Climatology) –

  • species (str) – Molecule to use, one of [‘H2SO4’, ‘ICE’, ‘WATER’]

class sasktran.SimpleRayleigh

Bases: OpticalProperty

An optical property object that provides Rayleigh molecular scattering without any corrections. This is very similar to the sasktran.Rayleigh optical property except that the Bates correction factor is not included. For most calculations it is preferred to use the sasktran.Rayleigh optical property instead.

class sasktran.BaumIceCrystal(effective_size_microns, use_delta_eddington=True)

Bases: OpticalProperty

Scattering non-spherical ice crystals based upon the database from Baum.

Parameters:
  • effective_size_microns (float) – Size of the particles, typically from 10 microns to 50 microns

  • use_delta_eddington (bool, optional) – True if the delta eddington approximation is to be used. Should be True for use in either the HR or DO engines. Default: True.

class sasktran.SO2Vandaele2009

Bases: OpticalProperty

Calculates the absorption cross section of SO2 molecules from 227 nm to 420 nm at 298 K, 318 K, 338 K, and 358 K following [1] and [2].

class sasktran.O2O2Fally2000

Bases: OpticalProperty

The collision induced absorption cross-sections measured by Fally et al. 2000. The cross-sections were measured at room temperature with a Fourier Transform spectrometer in the 15000 to 29800 cm−1 region (335-667 nm) at a maximal optical path difference of 0.45 cm (resolution 2 cm-1). Note that an appropriate climatology is the square of the O2 number density. This is provided by climatologies that support SKCLIMATOLOGY_O2_O2_CM6.

class sasktran.O2O2HITRAN2016

Bases: OpticalProperty

The O2-O2 collision induced absorption cross-sections distributed in Hitran 2016 as described by Karman et al. 2019. It is composed of 8 spectral regions measured by several researchers that extend from 335nm to over 8 microns. Some of the regions are temperature dependent and others are not

class sasktran.O2O2Thalman2013

Bases: OpticalProperty

The collision induced absorption cross-sections measured by Thalman et al. 2013. The cross-sections were measured at room temperature with a Fourier Transform spectrometer in the 15000 to 29800 cm−1 region (335-667 nm) at a maximal optical path difference of 0.45 cm (resolution 2 cm-1). Note that an appropriate climatology is the square of the O2 number density. This is provided by climatologies that support SKCLIMATOLOGY_O2_O2_CM6. For example, see MSIS90 and ECMWF:

class sasktran.UserDefinedAbsorptionPressure(broadener_vmr: array, pressures: array, temperatures: ndarray, wavel_nm: array, xs: ndarray, broadener_clim: Climatology | None = None, broadener_handle: str | None = None)

Bases: OpticalProperty

A user defined absorption optical property that is a function of temperature, pressure, wavelength, and VMR of a secondary species that broadens the profile.

Parameters:
  • broadener_vmr (np.array) – VMR of the broadening species. Must be in ascending order. Shape (Nv)

  • pressures (np.array) – Pressure in Pa, must be in ascending order. Shape (Np)

  • temperatures (np.array) – Temperatures in K at each pressure level. Must be in ascending order. Shape (Np, Nt)

  • wavel_nm (np.array) – Wavelengths in nm. Must be in ascending order. Shape (Nw)

  • xs (np.ndarray) – Cross sections in cm^2. Shape (Nv, Np, Nt, Nw)

  • broadener_clim (Optional[Climatology]) – A climatology for the broadening species, can be None if there is no broadener. Default None.

  • broadener_handle (Optional[str]) – The climatology handle for the broadener. Can be None if there is no broadener. Default None.

class sasktran.UserDefinedScatterConstantHeight(wavelengths: array, xs_scat: array, xs_abs: array, legendre_moments: ndarray | None = None, lm_a1: ndarray | None = None, lm_a2: ndarray | None = None, lm_a3: ndarray | None = None, lm_a4: ndarray | None = None, lm_b1: ndarray | None = None, lm_b2: ndarray | None = None, delta_m_order: int | None = None)

Bases: OpticalProperty

A user defined optical property for a scattering species that has constant optical properties in location space.

Currently this optical property only supports specifying the phase function in terms of the Legendre expansion, this means that it is only useful for the SASKTRAN-DO engine.

The Legendre expansion is defined in terms of “greek” coefficients, where the phase matrix has the standard form

\[\begin{split}P(\Theta) = \begin{pmatrix} a1(\Theta) & b1(\Theta) & 0 & 0 \\ b1(\Theta) & a2(\Theta) & 0 & 0 \\ 0 & 0 & a3(\Theta) & b2(\Theta) \\ 0 & 0 & -b2(\Theta) & a4(\Theta) \end{pmatrix}\end{split}\]

The phase function is then expanded in terms of Wigner functions (generalized spherical functions related to Legendre polynomials/associated Legendre functions)

\[ \begin{align}\begin{aligned}a_1(\Theta) = \sum_l \alpha_{1,l} \, d^{l}_{00}(\Theta)\\a_2(\Theta) + a_3(\Theta) = \sum_l (\alpha_{2,l} + \alpha_{3,l}) \, d^{l}_{2,2}(\Theta)\\a_2(\Theta) - a_3(\Theta) = \sum_l (\alpha_{2,l} - \alpha_{3,l}) \, d^{l}_{2,-2}(\Theta)\\a_4(\Theta) = \sum_l \alpha_{4,l} \, d^{l}_{00}(\Theta)\\b_1(\Theta) = \sum_l \beta_{1,l} \, d^{l}_{02}(\Theta)\\b_2(\Theta) = -\sum_l \beta_{2,l} \, d^{l}_{02}(\Theta)\end{aligned}\end{align} \]
Parameters:
  • wavelengths (np.array) – Array of wavelengths in nm. Shape (n,)

  • xs_scat (np.array) – Scattering cross section in cm2 at wavelengths. Shape (n,)

  • xs_abs (np.array) – Absorption cross section in cm2 at wavelengths. Shape (n,)

  • lm_a1 (np.ndarray) – Legendre coefficients \(\alpha_{1,l}\). Shape (n, m).

  • lm_a2 (np.ndarray) – Legendre coefficients \(\alpha_{2,l}\). Shape (n, m).

  • lm_a3 (np.ndarray) – Legendre coefficients \(\alpha_{3,l}\). Shape (n, m).

  • lm_a4 (np.ndarray) – Legendre coefficients \(\alpha_{4,l}\). Shape (n, m).

  • lm_b1 (np.ndarray) – Legendre coefficients \(\beta_{1,l}\). Shape (n, m).

  • lm_b2 (np.ndarray) – Legendre coefficients \(\beta_{2,l}\). Shape (n, m).

  • delta_m_order (int) – Delta-m order to truncaate the optical property by. None implies no scaling

See Also

Species

Representation of a atmospheric constituent.

Climatology

A profile quantifying something in/about the atmosphere. For example climatologies could be used to define profiles for temperature, pressure, gas number densities, etc.

Atmosphere

An object that describes the atmospheric constituents, and surface type to the engine. Atmospheric constituents are defined by a list of Species, and surfaces are defined by a given BRDF.