SMRT

smrt.emmodel package

Submodules

smrt.emmodel.common module

rayleigh_scattering_matrix_and_angle_tsang00(mu_s, mu_i, dphi, npol=2)

compute the Rayleigh matrix and half scattering angle. Based on Tsang theory and application p271 Eq 7.2.16

phase_matrix_from_scattering_amplitude(fvv, fvh, fhv, fhh, npol=2)

compute the phase function according to the scattering amplitude. This follows Tsang’s convention.

extinction_matrix(sigma_V, sigma_H=None, npol=2, mu=None)

compute the extinction matrix from the extinction in V and in H-pol. If sigma_V or sigma_H are a scalar, they are expanded in a diagonal matrix provided mu is given. If sigma_H is None, sigma_V is used.

rayleigh_scattering_matrix_and_angle_maetzler06(mu_s, mu_i, dphi, npol=2)

compute the Rayleigh matrix and half scattering angle. Based on Mätzler 2006 book p111. This version is relatively slow because it uses phase matrix rotations which is unnecessarily complex for the Rayleigh phase matrix but would be of interest for other phase matrices.

Lmatrix(cos_phi, sin_phi_sign, npol)
rayleigh_scattering_matrix_and_angle(mu_s, mu_i, dphi, npol=2)

compute the Rayleigh matrix and half scattering angle. Based on Tsang theory and application p271 Eq 7.2.16

class AdjustableEffectivePermittivityMixins

Bases: object

Mixin that allow an EM model to have the effective permittivity model defined by the user instead of by the theory of the EM Model.

The EM model must declare a default effective permittivity model.

effective_permittivity()

Calculation of complex effective permittivity of the medium.

Returns effective_permittivity:
 complex effective permittivity of the medium
derived_EMModel(base_class, effective_permittivity_model)

return a new IBA model with variant from the default IBA.

Parameters:effective_permittivity_model – permittivity mixing formula.

:returns a new class inheriting from IBA but with patched methods

smrt.emmodel.commontest module

test_energy_conservation(em, tolerance_pc, npol=None, subset=16)

test energy conservation

Parameters:
  • em – the electromagnetic model that has been set up
  • tolerance_pc – relative tolerance

smrt.emmodel.dmrt_qca_shortrange module

Compute scattering with DMRT QCA Short range. Short range means that it is accurate only for small and weakly sticky spheres (high stickiness value). It diverges (increasing scattering coefficient) if these conditions are not met. Numerically the size conditions can be evaluated with the ratio radius/wavelength as for Rayleigh scatterers. For the stickiness, it is more difficult as this depends on the size of the scatterers and the fractional volume. In any case, it is dangerous to use too small a stickiness value, especially if the grains are big.

This model is only compatible with the SHS microstructure model.

Examples:

from smrt import make_snowpack, make_sensor

density = [345.0] temperature = [260.0] thickness = [70] radius = [750e-6] stickiness = [0.1]

snowpack = make_snowpack(thickness, “sticky_hard_spheres”,
density=density, temperature=temperature, radius=radius, stickiness=stickiness)

# create the EM Model - Equivalent DMRTML m = make_model(“dmrt_shortrange”, “dort”)

# create the sensor theta = np.arange(5.0, 80.0, 5.0) radiometer = sensor.amsre()

class DMRT_QCA_ShortRange(sensor, layer, dense_snow_correction='auto')

Bases: smrt.emmodel.rayleigh.Rayleigh

DMRT electromagnetic model in the short range limit (grains AND aggregates are small) as implemented in DMRTML

Parameters:
  • sensor – sensor instance
  • layer – layer instance
Dense_snow_correction:
 

set how snow denser than half the ice density (ie. fractional volume larger than 0.5 is handled).

“auto” means that snow is modeled as air bubble in ice instead of ice spheres in air. “bridging” should be developed in the future.

basic_check()

smrt.emmodel.dmrt_qcacp_shortrange module

Compute scattering with DMRT QCACP Short range like in DMRT-ML. Short range means that it is accurate only for small and weakly sticky spheres (high stickiness value). It diverges (increasing scattering coefficient) if these conditions are not met. Numerically the size conditions can be evaluated with the ratio radius/wavelength as for Rayleigh scatterers. For the stickiness, it is more difficult as this depends on the size of the scatterers and the fractional volume. In any case, it is dangerous to use too small a stickiness value, especially if the grains are big.

This model is only compatible with the SHS microstructure model.

Examples:

from smrt import make_snowpack, make_sensor

density = [345.0] temperature = [260.0] thickness = [70] radius = [750e-6] stickiness = [0.1]

snowpack = make_snowpack(thickness, “sticky_hard_spheres”,
density=density, temperature=temperature, radius=radius, stickiness=stickiness)

# create the EM Model - Equivalent DMRTML m = make_model(“dmrt_shortrange”, “dort”)

# create the sensor theta = np.arange(5.0, 80.0, 5.0) radiometer = sensor.amsre()

class DMRT_QCACP_ShortRange(sensor, layer, dense_snow_correction='auto')

Bases: smrt.emmodel.rayleigh.Rayleigh

DMRT electromagnetic model in the short range limit (grains AND aggregates are small) as implemented in DMRTML

Parameters:
  • sensor – sensor instance
  • layer – layer instance
Dense_snow_correction:
 

set how snow denser than half the ice density (ie. fractional volume larger than 0.5 is handled).

“auto” means that snow is modeled as air bubble in ice instead of ice spheres in air. “bridging” should be developed in the future.

basic_check()

smrt.emmodel.iba module

Compute scattering from Improved Born Approximation theory as described in Mätzler 1998 and Mätzler and Wiesman 1999, except the absorption coefficient which is computed with Polden von Staten formulation instead of the Eq 24 in Mätzler 1998. See iba_original.py for a fully conforming IBA version.

This model allows for different microstructural models provided that the Fourier transform of the correlation function

may be performed. All properties relate to a single layer.

derived_IBA(effective_permittivity_model=<function polder_van_santen>)

return a new IBA model with variant from the default IBA.

Parameters:effective_permittivity_model – permittivity mixing formula.

:returns a new class inheriting from IBA but with patched methods

class IBA(sensor, layer)

Bases: smrt.emmodel.common.AdjustableEffectivePermittivityMixins

Improved Born Approximation electromagnetic model class.

As with all electromagnetic modules, this class is used to create an electromagnetic object that holds information about the effective permittivity, extinction coefficient and phase function for a particular snow layer. Due to the frequency dependence, information about the sensor is required. Passive and active sensors also have different requirements on the size of the phase matrix as redundant information is not calculated for the passive case.

Parameters:
  • sensor – object containing sensor characteristics
  • layer – object containing snow layer characteristics (single layer)

Usage Example:

This class is not normally accessed directly by the user, but forms part of the smrt model, together with the radiative solver (in this example, dort) i.e.:

from smrt import make_model
model = make_model("iba", "dort")

iba does not need to be imported by the user due to autoimport of electromagnetic model modules

static effective_permittivity_model(frac_volume, e0=1, eps=3.185, depol_xyz=None, length_ratio=None, inclusion_shape=None, mixing_ratio=1)

Calculates effective permittivity of snow by solution of quadratic Polder Van Santen equation for spherical inclusion.

Parameters:
  • frac_volume – Fractional volume of inclusions
  • e0 – Permittivity of background (default is 1)
  • eps – Permittivity of scattering material (default is 3.185 to compare with MEMLS)
  • depol_xyz – [Optional] Depolarization factors, spherical isotropy is default. It is not taken into account here.
  • length_ratio – Length_ratio. Used to estimate depolarization factors when they are not given.
  • inclusion_shape – Assumption for shape(s) of brine inclusions. Can be a string for single shape, or a list/tuple/dict of strings for mixture of shapes. So far, we have the following shapes: “spheres” and “random_needles” (i.e. randomly-oriented elongated ellipsoidal inclusions). If the argument is a dict, the keys are the shapes and the values are the mixing ratio. If it is a list, the mixing_ratio argument is required.
  • mixing_ratio – The mixing ratio of the shapes. This is only relevant when inclusion_shape is a list/tuple. Mixing ratio must be a sequence with length len(inclusion_shape)-1. The mixing ratio of the last shapes is deduced as the sum of the ratios must equal to 1.
Returns:

Effective permittivity

Usage example:

from smrt.permittivity.generic_mixing_formula import polder_van_santen
effective_permittivity = polder_van_santen(frac_volume, e0, eps)

# for a mixture of 30% spheres and 70% needles
effective_permittivity = polder_van_santen(frac_volume, e0, eps, inclusion_shape={"spheres": 0.3, "random_needles": 0.7})
# or
effective_permittivity = polder_van_santen(frac_volume, e0, eps, inclusion_shape=("spheres", "random_needles"), mixing_ratio=0.3)

Todo

Extend Polder Van Santen model to account for ellipsoidal inclusions

compute_iba_coeff()

Calculate angular independent IBA coefficient: used in both scattering coefficient and phase function calculations

Note

Requires mean squared field ratio (uses mean_sq_field_ratio method)

mean_sq_field_ratio(e0, eps)

Mean squared field ratio calculation

Uses layer effective permittivity

Parameters:
  • e0 – background relative permittivity
  • eps – scattering constituent relative permittivity
basic_check()
compute_ks()

Calculate scattering coefficient: integrate p11+p12 over mu

ks_integrand(mu)

This is the scattering function for the IBA model.

It uses the phase matrix in the 1-2 frame. With incident angle chosen to be 0, the scattering angle becomes the scattering zenith angle:

\Theta = \theta

Scattering coefficient is determined by integration over the scattering angle (0 to pi)

Parameters:mu – cosine of the scattering angle (single angle)

ks\_int = p11 + p22

The integration is performed outside this method.

phase(mu_s, mu_i, dphi, npol=2)

IBA Phase function (not decomposed).

ft_even_phase(mu_s, mu_i, m_max, npol=None)

Calculation of the Fourier decomposed IBA phase function.

This method calculates the Improved Born Approximation phase matrix for all Fourier decomposition modes and return the output.

Coefficients within the phase function are

Passive case (m = 0 only) and active (m = 0)

M  = [Pvvp  Pvhp]
     [Phvp  Phhp]

Active case (m > 0):

M =  [Pvvp Pvhp Pvup]
     [Phvp Phhp Phup]
     [Puvp Puhp Puup]

The IBA phase function is given in Mätzler, C. (1998). Improved Born approximation for scattering of radiation in a granular medium. Journal of Applied Physics, 83(11), 6111-6117. Here, calculation of the phase matrix is based on the phase matrix in the 1-2 frame, which is then rotated according to the incident and scattering angles, as described in e.g. Thermal Microwave Radiation: Applications for Remote Sensing, Mätzler (2006). Fourier decomposition is then performed to separate the azimuthal dependency from the incidence angle dependency.

Parameters:
  • mu_s – 1-D array of cosine of viewing radiation stream angles (set by solver)
  • mu_i – 1-D array of cosine of incident radiation stream angles (set by solver)
  • m_max – maximum Fourier decomposition mode needed
  • npol – number of polarizations considered (set from sensor characteristics)
compute_ka()

IBA absorption coefficient calculated from the low-loss assumption of a general lossy medium.

Calculates ka from wavenumber in free space (determined from sensor), and effective permittivity of the medium (snow layer property)

Return ka:absorption coefficient [m -1]

Note

This may not be suitable for high density material

ke(mu, npol=2)

IBA extinction coefficient matrix

The extinction coefficient is defined as the sum of scattering and absorption coefficients. However, the radiative transfer solver requires this in matrix form, so this method is called by the solver.

param mu:1-D array of cosines of radiation stream incidence angles
param npol:number of polarization
returns ke:extinction coefficient matrix [m -1]

Note

Spherical isotropy assumed (all elements in matrix are identical).

Size of extinction coefficient matrix depends on number of radiation streams, which is set by the radiative transfer solver.

class IBA_MM(sensor, layer)

Bases: smrt.emmodel.iba.IBA

smrt.emmodel.iba_original module

Compute scattering from Improved Born Approximation theory. This model allows for different microstructural models provided that the Fourier transform of the correlation function may be performed. All properties relate to a single layer. The absorption is calculated with the original formula in Mätzler 1998

class IBA_original(sensor, layer)

Bases: smrt.emmodel.iba.IBA

Original Improved Born Approximation electromagnetic model class.

As with all electromagnetic modules, this class is used to create an electromagnetic object that holds information about the effective permittivity, extinction coefficient and phase function for a particular snow layer. Due to the frequency dependence, information about the sensor is required. Passive and active sensors also have different requirements on the size of the phase matrix as redundant information is not calculated for the passive case.

Parameters:
  • sensor – object containing sensor characteristics
  • layer – object containing snow layer characteristics (single layer)
compute_ka()

IBA absorption coefficient calculated from the low-loss assumption of a general lossy medium.

Calculates ka from wavenumber in free space (determined from sensor), and effective permittivity of the medium (snow layer property)

Return ka:absorption coefficient [m -1]

Note

This may not be suitable for high density material

smrt.emmodel.nonscattering module

Non-scattering medium can be applied to medium without heteoreneity (like water or pure ice layer).

class NonScattering(sensor, layer)

Bases: object

basic_check()
ft_even_phase(mu_s, mu_i, m_max, npol=None)

Non-scattering phase matrix.

Returns : null phase matrix

phase(mu_s, mu_i, dphi, npol=2)

Non-scattering phase matrix.

Returns : null phase matrix

ke(mu, npol=2)
effective_permittivity()

smrt.emmodel.prescribed_kskaeps module

Use prescribed scattering ks and absorption ka coefficients and effective permittivity in the layer.

The phase matrix has the Rayleigh form with prescribed scattering coefficient

This model is compatible with any microstructure but requires that ks, ka, and optionally effective permittivity to be set in the layer

Example:

m = make_model("prescribed_kskaeps", "dort")
snowpack.layers[0].ks = ks
snowpack.layers[0].ka = ka
snowpack.layers[0].effective_permittivity = eff_eps
class Prescribed_KsKaEps(sensor, layer)

Bases: smrt.emmodel.rayleigh.Rayleigh

smrt.emmodel.rayleigh module

Compute Rayleigh scattering. This theory requires the scatterers to be smaller than the wavelength and the medium to be sparsely populated (eq. very low density in the case of snow).

This model is only compatible with the Independent Sphere microstructure model

class Rayleigh(sensor, layer)

Bases: object

basic_check()
ft_even_phase_baseonUlaby(mu_s, mu_i, m_max, npol=None)

# # Equations are from pg 1188-1189 Ulaby, Moore, Fung. Microwave Remote Sensing Vol III. # See also pg 157 of Tsang, Kong and Shin: Theory of Microwave Remote Sensing (1985) - can be used to derive # the Ulaby equations.

ft_even_phase_basedonJin(mu_s, mu_i, m_max, npol=None)

Rayleigh phase matrix.

These are the Fourier decomposed phase matrices for modes m = 0, 1, 2. It is based on Y.Q. Jin

Coefficients within the phase function are:

M = [Pvvp Pvhp]
[Phvp Phhp]

Inputs are: :param m: mode for decomposed phase matrix (0, 1, 2) :param mu: vector of cosines of incidence angle

Returns P: phase matrix

ft_even_phase_tsang(mu_s, mu_i, m_max, npol=None)

Rayleigh phase matrix.

These are the Fourier decomposed phase matrices for modes m = 0, 1, 2. Equations are from p128 Tsang Application and Theory 200 and sympy calculations

Coefficients within the phase function are:

M = [P[v, v] P[v, h] -P[v, u]]
[P[h, v] P[h, h] -P[h, u]] [P[u, v] P[u, h] P[u, u]]

Inputs are: :param m: mode for decomposed phase matrix (0, 1, 2) :param mu: vector of cosines of incidence angle

Returns P: phase matrix

ft_even_phase(mu_s, mu_i, m_max, npol=None)

# # Equations are from pg 1188-1189 Ulaby, Moore, Fung. Microwave Remote Sensing Vol III. # See also pg 157 of Tsang, Kong and Shin: Theory of Microwave Remote Sensing (1985) - can be used to derive # the Ulaby equations.

phase(mu_s, mu_i, dphi, npol=2)
ke(mu, npol=2)

return the extinction coefficient matrix

The extinction coefficient is defined as the sum of scattering and absorption coefficients. However, the radiative transfer solver requires this in matrix form, so this method is called by the solver.

param mu:1-D array of cosines of radiation stream incidence angles
param npol:number of polarizations
returns ke:extinction coefficient matrix [m -1]

Note

Spherical isotropy assumed (all elements in matrix are identical).

Size of extinction coefficient matrix depends on number of radiation streams, which is set by the radiative transfer solver.

effective_permittivity()

smrt.emmodel.sft_rayleigh module

Compute Strong Fluctuation Theory scattering. This theory requires the scatterers to be smaller than the wavelength

This model is only compatible with the Exponential autocorrelation function only

class SFT_Rayleigh(sensor, layer)

Bases: smrt.emmodel.rayleigh.Rayleigh

smrt.emmodel.test_iba module

setup_func_sp()
setup_func_indep(radius=0.0005)
setup_func_shs()
setup_func_pc(pc)
setup_func_em(testpack=None)
setup_func_active(testpack=None)
setup_func_rayleigh()
setup_mu(stepsize, bypass_exception=None)
test_ks_pc_is_0p3_mm()
test_ks_pc_is_0p25_mm()
test_ks_pc_is_0p15_mm()
test_ks_pc_is_0p1_mm()
test_ks_pc_is_0p2_mm()
test_energy_conservation_exp()
test_energy_conservation_indep()
test_energy_conservation_shs()
test_npol_passive_is_2()
test_npol_active_is_3()
test_energy_conservation_exp_active()
test_energy_conservation_indep_active()
test_energy_conservation_shs_active()
test_iba_vs_rayleigh_passive_m0()
test_iba_vs_rayleigh_active_m0()
test_iba_vs_rayleigh_active_m1()
test_iba_vs_rayleigh_active_m2()
test_permittivity_model()
test_iba_raise_exception_mu_is_1()

smrt.emmodel.test_iba_original module

setup_func_sp()
setup_func_indep(radius=0.0005)
setup_func_shs()
setup_func_pc(pc)
setup_func_em(testpack=None)
setup_func_active(testpack=None)
setup_func_rayleigh()
setup_mu(stepsize, bypass_exception=None)
test_ks_pc_is_0p3_mm()
test_ks_pc_is_0p25_mm()
test_ks_pc_is_0p15_mm()
test_ks_pc_is_0p1_mm()
test_ks_pc_is_0p2_mm()
test_energy_conservation_exp()
test_energy_conservation_indep()
test_energy_conservation_shs()
test_npol_passive_is_2()
test_npol_active_is_3()
test_energy_conservation_exp_active()
test_energy_conservation_indep_active()
test_energy_conservation_shs_active()
test_iba_vs_rayleigh_passive_m0()
test_iba_vs_rayleigh_active_m0()
test_iba_vs_rayleigh_active_m1()
test_iba_vs_rayleigh_active_m2()
test_iba_raise_exception_mu_is_1()

smrt.emmodel.test_prescribed_kskaeps module

setup_func_sp()
setup_func_em(testpack=None)
test_energy_conservation()

smrt.emmodel.test_rayleigh module

setup_func_sp()
setup_func_em(testpack=None)
test_energy_conservation()
test_energy_conservation_tsang()
test_energy_conservation_jin()

smrt.emmodel.test_sft_rayleigh module

setup_func_sp()
setup_func_em(testpack=None)
test_energy_conservation()

Module contents

This directory contains the different electromagnetic (EM) models that compute the scattering and absorption coefficients and the phase function in a _given_ _layer_. The computation of the inter-layer propagation is done by the rtsolver package.

The EM models differ in many aspects, one of which is the constraint on the microstructure model they can be used with. The smrt.emmodel.iba model can use any microstructure model that defines autocorrelation functions (or its FT). In contrast others such as smrt.emmodel.dmrt_shortrange is bound to the smrt.microstructuremodel.sticky_hard_spheres microstructure for theoretical reasons.

The selection of the EM model is done with the smrt.core.model.make_model() function

For developers

To implement a new scattering formulation / phase function, we recommend to start from an existing module, probably rayleigh.py is the simplest. Copy this file to myscatteringtheory.py or any meaningful name. It can be directly used with make_model() function as follows:

m = make_model("myscatteringtheory", "dort")

Note that if the file is not in the emmodels directory, you must explicitly import the module and pass it to make_model as a module object (instead of a string).

An emmodel model must define:
  • ks and ka attributes/properties
  • ke() and effective_permittivity() methods
  • at least one of the phase and ft_even_phase methods (both is better).

For the details it is recommended to contact the authors as the calling arguments and required methods may change time to time.