smrt.core package

The core package contains the SMRT machinery. It provides the infrastructure that defines basic objects and orchestrates the “science” modules in the other packages (such as smrt.emmodel or smrt.rtsolver).

Amongst all, we suggest looking at the documentation of the Result object.

For developers

We strongly warn against changing anything in this directory. In principle this is not needed because no “science” is present and most objects and functions are generic enough to be extendable from outside (without affecting the core definition). Ask advice from the authors if you really want to change something here.

smrt.core.result module

The results of RT Solver are hold by the Result class. This class provides several functions to access to the Stokes Vector and Muller matrix in a simple way. Most notable ones are Result.TbV() and Result.TbH() for the passive mode calculations and Result.sigmaHH() and Result.sigmaVV(). Result.to_dataframe() is also very convenient for the sensors with a channel map (all specific satellite sensors have such a map, only generic sensors as smrt.sensor_list.active() and smrt.sensor_list.passive() does not provide a map by default).

In addition, the RT Solver stores some information in Result.other_data. Currently this includes the effective_permittivity, ks and ka for each layer. The data are accessed directly with e.g. result.other_data[‘ks’].

To save results of calculations in a file, simply use the pickle module or other serialization schemes. We may provide a unified and inter-operable solution in the future.

Under the hood, Result uses xarray module which provides multi-dimensional array with explicit, named, dimensions. Here the common dimensions are frequency, polarization, polarization_inc, theta_inc, theta, and phi. They are created by the RT Solver. The interest of using named dimension is that slice of the xarray (i.e. results) can be selected based on the dimension name whereas with numpy the order of the dimensions matters. Because this is very convenient, users may be interested in adding other dimensions specific to their context such as time, longitude, latitude, points, … To do so, smrt.core.model.Model.run() accepts a list of snowpack and optionally the parameter snowpack_dimension is used to specify the name and values of the new dimension to build.

Example:

times = [datetime(2012, 1, 1), datetime(2012, 1, 5), , datetime(2012, 1, 10)]
snowpacks = [snowpack_1jan, snowpack_5jan, snowpack_10jan]

res = model.run(sensor, snowpacks, snowpack_dimension=('time', times))

The res variable is a Result instance, so that for all the methods of this class that can be called, they will return a timeseries. For instance result.TbV(theta=53) returns a time-series of brightness temperature at V polarization and 53° incidence angle and the following code plots this timeseries:

plot(times, result.TbV(theta=53))
open_result(filename)

Reads a result save to disk. See Result.save() method.

make_result(sensor, *args, **kwargs)

Creates an active or passive result object according to the mode

class Result(intensity, coords=None, channel_map=None, other_data={}, mother_df=None)

Bases: object

Contains the results of a/many computations and provides convenience functions to access these results

property coords

Returns the coordinates of the result (theta, frequency, …). Note that the coordinates are also result attribute, so result.frequency works (and so on for all the coordinates).

save(filename)

Saves a result to disk. Under the hood, this is a netCDF file produced by xarray (http://xarray.pydata.org/en/stable/io.html).

to_series(**kwargs)

Returns the result as a series with the channels defined in the sensor as index. This requires that the sensor has declared a channel list.

optical_depth()

Returns the optical depth of each layer tau = ke * thickness, where ke = ka + ks calculated for each layer. This is useful to assess the e-folding depth (aka penetration depth), neglecting the layer transmittance.

For instance the direct incoming radiation (in active mode) or the radiation emanating directly (in passive mode) can be estimated as:

intensity = np.exp(-result.optical_depth().cumsum('layer'))
single_scattering_albedo()

Returns the single scattering albedo of each layer: ssalb = ks / (ks + ka). This is useful to assess if multiple scattering is significant (e.g. if ssalb > 0.2).

class PassiveResult(intensity, coords=None, channel_map=None, other_data={}, mother_df=None)

Bases: Result

Tb(channel=None, **kwargs)

Return brightness temperature. Any parameter can be added to slice the results (e.g. frequency=37e9 or polarization=’V’). See xarray slicing with sel method (to document). It is also posisble to select by channel if the sensor has a channel_map.

Parameters:
  • channel – channel to select

  • **kwargs – any parameter to slice the results.

Tb_as_dataframe(channel_axis=None, **kwargs)

See PassiveResult().to_dataframe

to_dataframe(channel_axis='auto', **kwargs)

Returns brightness temperature as a pandas.DataFrame. Any parameter can be added to slice the results (e.g. frequency=37e9 or polarization=’V’). See xarray slicing with sel method (to document). In addition channel_axis controls the format of the output. If set to None, the DataFrame has a multi-index with all the dimensions (frequency, polarization, …). If channel_axis is set to “column”, and if the sensor has a channel map, the channels are in columns and the other dimensions are in index. If set to “index”, the channel are in index with all the other dimensions.

The most conviennent is probably channel_axis=”column” while channel_axis=None (default) contains all the data even those not corresponding to a channel and applies to any sensor even those without channel_map. If set to “auto”, the channel_axis is “column” if the channel map exit, otherwise is None.

Parameters:

channel_axis – controls whether to use the sensor channel or not and if yes, as a column or index.

TbV(**kwargs)

Returns V polarization. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

TbH(**kwargs)

Returns H polarization. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

polarization_ratio(ratio='H_V', **kwargs)

Returns polarization ratio. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

class ActiveResult(intensity, coords=None, channel_map=None, other_data={}, mother_df=None)

Bases: Result

sigma(channel=None, **kwargs)

Returns backscattering coefficient. Any parameter can be added to slice the results (e.g. frequency=37e9 or polarization=’V’). See xarray slicing with sel method (to document). It is also posisble to select by channel if the sensor has a channel_map.

Args: channel: channel to select **kwargs: any parameter to slice the results.

sigma_dB(channel=None, **kwargs)

Returns backscattering coefficient. Any parameter can be added to slice the results (e.g. frequency=37e9, polarization_inc=’V’, polarization=’V’). See xarray slicing with sel method (to document)

sigma_as_dataframe(channel_axis=None, **kwargs)

Returns backscattering coefficient as a pandas.DataFrame. Any parameter can be added to slice the results (e.g. frequency=37e9 or polarization=’V’). See xarray slicing with sel method (to document). In addition channel_axis controls the format of the output. If set to None, the DataFrame has a multi-index formed with all the dimensions (frequency, polarization, …). If channel_axis is set to “column”, and if the sensor has named channels (channel_map in SMRT wording), the channel are in columns and the other dimensions are in index. If set to “index”, the channel are in index with all the other dimensions.

The most conviennent is probably channel_axis=”column” while channel_axis=None (default) contains all the data even those not corresponding to a channel and applies to any sensor even those without channel_map.

Parameters:

channel_axis – controls whether to use the sensor channel or not and if yes, as a column or index.

sigma_dB_as_dataframe(channel_axis=None, **kwargs)

See ActiveResult().to_dataframe

to_dataframe(channel_axis=None, **kwargs)

Return backscattering coefficient in dB as a pandas.DataFrame. Any parameter can be added to slice the results (e.g. frequency=37e9 or polarization=’V’). See xarray slicing with sel method (to document). In addition channel_axis controls the format of the output. If set to None, the DataFrame has a multi-index with all the dimensions (frequency, polarization, …). If channel_axis is set to “column”, and if the sensor has named channels (channel_map in SMRT wording), the channel are in columns and the other dimensions are in index. If set to “index”, the channel are in index with all the other dimensions.

If channel_axis is set to “column”, and if the sensor has a channel map, the channels are in columns and the other dimensions are in index. If set to “index”, the channel are in index with all the other dimensions.

The most conviennent is probably channel_axis=”column” while channel_axis=None (default) contains all the data even those not corresponding to a channel and applies to any sensor even those without channel_map. If set to “auto”, the channel_axis is “column” if the channel map exit, otherwise is None.

Parameters:

channel_axis – controls whether to use the sensor channel or not and if yes, as a column or index.

to_series(**kwargs)

Returns backscattering coefficients in dB as a series with the channels defined in the sensor as index. This requires that the sensor has declared a channel list.

sigmaVV(**kwargs)

Returns VV backscattering coefficient. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

sigmaVV_dB(**kwargs)

Returns VV backscattering coefficient in dB. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

sigmaHH(**kwargs)

Returns HH backscattering coefficient. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

sigmaHH_dB(**kwargs)

Returns HH backscattering coefficient in dB. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

sigmaHV(**kwargs)

Returns HV backscattering coefficient. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

sigmaHV_dB(**kwargs)

Returns HV backscattering coefficient in dB. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

sigmaVH(**kwargs)

Return VH backscattering coefficient. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

sigmaVH_dB(**kwargs)

Returns VH backscattering coefficient in dB. Any parameter can be added to slice the results (e.g. frequency=37e9). See xarray slicing with sel method (to document)

class AltimetryResult(intensity, coords=None, channel_map=None, other_data={}, mother_df=None)

Bases: ActiveResult

delay_doppler_map(**kwargs)

Returns the delay Doppler map

waveform(**kwargs)

Returns the waveform.

For simulations with return_contributions, this function returns the total only by default. Use explicit contribution=”all”” to get all the contributions or contribution=’…’ to access each contribution.

Parameters:

**kwargs – any dimension to select. See xarray.DataArray.sel.

contributions()

Returns the list of the contribution dimension. Raise an exception if the contribution does not exist.

concat_results(result_list, coord)

Concatenates several results from smrt.core.model.Model.run() (of type Result) into a single result (of type Result). This extends the number of dimension in the xarray hold by the instance. The new dimension is specified with coord

Parameters:
  • result_list – list of results returned by smrt.core.model.Model.run() or other functions.

  • coord – a tuple (dimension_name, dimension_values) for the new dimension. Dimension_values must be a sequence or

  • result_list. (array with the same length as)

Returns:

Result instance

prepare_kskaeps_profile_information(snowpack, emmodels, effective_permittivity=None, mu=0)

Return a dict with the profiles of ka, ks and effective permittivity. Can be directly used by Solver to insert data in other_data.

Parameters:
  • snowpack – the snowpack used for the calculation

  • emmodels – the list of instantiated emmodel

  • effective_permittivity – list of permittivity. If None, it is obtained from the emmodels

  • mu – cosine angles where to compute ks.

smrt.core.globalconstants module

Global constants used throughout the model are defined here and imported as needed. The constants are:

Parameter

Description

Value

DENSITY_OF_ICE

Density of pure ice at 273.15K

916.7 kg m -3

FREEZING_POINT

Freezing point of pure water

273.15 K

C_SPEED

Speed of light in a vacuum

2.99792458 x 10 8 ms -1

PERMITTIVITY_OF_AIR

Relative permittivity of air

1

Usage example:

from smrt.core.globalconstants import DENSITY_OF_ICE

smrt.core.atmosphere module

class AtmosphereResult(tb_down: numpy.ndarray, tb_up: numpy.ndarray, transmittance: numpy.ndarray)

Bases: object

smrt.core.error module

Definition of the Exception specific to SMRT.

exception SMRTError

Bases: Exception

Error raised by the model.

exception SMRTWarning

Bases: Warning

Warning raised by the model.

smrt.core.fresnel module

Fresnel coefficients formulae used in the packages smrt.interface and smrt.substrate.

fresnel_coefficients_old(eps_1, eps_2, mu1)

Computes the reflection in two polarizations (H and V). The equations are only valid for lossless media. Applying these equations for (strongly) lossy media results in (large) errors. Don’t use it. It is here for reference only. The returned reflection coefficients apply to the electric field. Use abs2(rv), abs2(rh) to obtain the power reflection coefficient.

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

  • mu1 – cosine zenith angle in medium 1.

Returns:

rv, rh, mu2 the cosine of the angle in medium 2

fresnel_coefficients_maezawa09_classical(eps_1, eps_2, mu1, full_output=False)

Computes the reflection in two polarizations (H and V) for lossly media with the “classical Fresnel” based on Maezawa, H., & Miyauchi, H. (2009). Rigorous expressions for the Fresnel equations at interfaces between absorbing media. Journal of the Optical Society of America A, 26(2), 330. https://doi.org/10.1364/josaa.26.000330

The classical derivation does not respect energy conservation, especially the transmittivity. Don’t use it. It is here for reference only. The returned reflection coefficients apply to the electric field. Use abs2(rv), abs2(rh) to obtain the power reflection coefficient.

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

  • mu1 – cosine zenith angle in medium 1.

  • full_output – (Default value = False)

Returns:

rv, rh, mu2 the cosine of the angle in medium 2

fresnel_coefficients_maezawa09_rigorous(eps_1, eps_2, mu1)

Computes the reflection in two polarizations (H and V) for lossly media with the “rigorous Fresnel” based on Maezawa, H., & Miyauchi, H. (2009). Rigorous expressions for the Fresnel equations at interfaces between absorbing media. Journal of the Optical Society of America A, 26(2), 330. https://doi.org/10.1364/josaa.26.000330

The ‘rigorous’ derivation respects the energy conservation even for strongly loosly media. The returned reflection coefficients apply to the electric field. Use abs2(rv), abs2(rh) to obtain the power reflection coefficient.

This function only returns the FIELD reflection coefficients and the cosine angle in the medium 2

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

  • mu1 – cosine zenith angle in medium 1.

Returns:

rv, rh, mu2 the cosine of the angle in medium 2

fresnel_coefficients_maezawa09_rigorous_full_output(eps_1, eps_2, mu1)

Computes the reflection in two polarizations (H and V) for lossly media with the “rigorous Fresnel” based on Maezawa, H., & Miyauchi, H. (2009). Rigorous expressions for the Fresnel equations at interfaces between absorbing media. Journal of the Optical Society of America A, 26(2), 330. https://doi.org/10.1364/josaa.26.000330

The ‘rigorous’ derivation respects the energy conservation even for strongly loosly media. The returned reflection coefficients apply to the electric field. Use abs2(rv), abs2(rh) to obtain the power reflection coefficient.

This function returns the FIELF and INTENSITY reflection and transmission coefficients and the cosine angle in the medium 2

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

  • mu1 – cosine zenith angle in medium 1.

  • full_output – (Default value = False)

Returns:

rv, rh, mu2 the cosine of the angle in medium 2

fresnel_coefficients(eps_1, eps_2, mu1)

Computes the reflection in two polarizations (H and V) for lossly media with the “rigorous Fresnel” based on Maezawa, H., & Miyauchi, H. (2009). Rigorous expressions for the Fresnel equations at interfaces between absorbing media. Journal of the Optical Society of America A, 26(2), 330. https://doi.org/10.1364/josaa.26.000330

The ‘rigorous’ derivation respects the energy conservation even for strongly loosly media. The returned reflection coefficients apply to the electric field. Use abs2(rv), abs2(rh) to obtain the power reflection coefficient.

This function only returns the FIELD reflection coefficients and the cosine angle in the medium 2

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

  • mu1 – cosine zenith angle in medium 1.

Returns:

rv, rh, mu2 the cosine of the angle in medium 2

snell_angle(eps_1, eps_2, mu1)

Computes mu2 the cos(angle) in the second medium according to Snell’s law.

Parameters:
  • eps_1

  • eps_2

  • mu1

brewster_angle(eps_1, eps_2)

Computes the brewster angle.

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

Returns:

angle in radians

fresnel_reflection_matrix(eps_1, eps_2, mu1, npol)

Computes the fresnel power reflection matrix for/in medium 1 laying above medium 2.

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

  • mu1 – cosine zenith angle in medium 1.

  • npol – number of polarizations to return.

Returns:

a matrix or the diagional depending on return_as_diagonal

fresnel_transmission_matrix(eps_1, eps_2, mu1, npol)

Computes the fresnel power transmission matrix for/in medium 1 lying above medium 2.

Parameters:
  • eps_1 – permittivity of medium 1.

  • eps_2 – permittivity of medium 2.

  • mu1 – cosine zenith angle in medium 1.

  • npol – number of polarizations to return.

Returns:

a matrix or the diagional depending on return_as_diagonal

smrt.core.interface module

This package implements the base class for all the substrate models. To create a substrate, it is recommended to use help functions such as make_soil() rather than the class constructor.

make_interface(inst_class_or_modulename, broadcast=True, **kwargs)

Returns an instance corresponding to the interface model with the provided arguments.

This function imports the interface module if necessary and returns an instance of the interface class with the provided arguments in **kwargs.

Parameters:
  • inst_class_or_modulename – a class, and instance or the name of the python module in smrt/interface

  • broadcast – (Default value = True)

  • **kwargs

class Interface(**kwargs)

Bases: object

Abstract class for interface between layer and substrate at the bottom of the snowpack. It provides argument handling.

class SubstrateBase(temperature=None, permittivity_model=None)

Bases: object

Abstract class for substrate at the bottom of the snowpack. It provides calculation of the permittivity constant for soil case. Argument handline is delegated to the instance of the interface

permittivity(frequency)

Computes the permittivity for the given frequency using permittivity_model. This method returns None when no permittivity model is available. This must be handled by the calling code and interpreted suitably.

Parameters:

frequency

substrate_from_interface(interface_cls)

this decorator transform an interface class into a substrate class with automatic method

Parameters:

interface_cls

class Substrate(temperature=None, permittivity_model=None, **kwargs)

Bases: SubstrateBase, Interface

get_substrate_model(substrate_model)

Returns the class corresponding to the substrate model called name.

Parameters:

substrate_model

Returns:

This function imports the correct module if possible and returns the class

smrt.core.layer module

Layer instance contains all the properties for a single snow layer (e.g. temperature, frac_volume, etc). It also contains a microstructure attribute that holds the microstructural properties (e.g. radius, corr_length, etc). The class of this attribute defines the microstructure model to use (see smrt.microstructure_model package).

To create a single layer, it is recommended to use the function make_snow_layer() rather than the class constructor. However, it is usually more convenient to create a snowpack using make_snowpack().

For developers

The Layer class should not be modified at all even if you need new properties to define the layer (e.g. brine concentration, humidity, …). If the property you need to add is related to geometric aspects, it is probably better to use an existing microstructure model or to create a new one. If the new parameter is not related to geometrical aspect, write a function similar to make_snow_layer() (choose an explicit name for your purpose). In this function, create the layer by calling the Layer constructor as in make_snow_layer() and then add your properties with lay.myproperty=xxx, … See the example of liquid water in make_snow_layer(). This approach avoids specialization of the Layer class. The new function can be in any file (inc. out of smrt directories), and should be added in make_medium if it is of general interest and written in a generic way, that is, covers many use cases for many users with default arguments, etc.

class Layer(thickness, microstructure_model=None, temperature=273.15, permittivity_model=None, inclusion_shape=None, **kwargs)

Bases: object

Contains the properties for a single layer including the microstructure attribute which holds the microstructure properties.

To create a layer, it is recommended to use the functions make_snow_layer or similar.

property ssa

Returns the SSA, computing it if necessary.

permittivity(i, frequency)

Returns the permittivity of the i-th medium depending on the frequency and internal layer properties.

Parameters:
  • i (int) – Number of the medium. 0 is reserved for the background.

  • frequency (float) – Frequency of the wave in Hz.

Returns:

Permittivity of the i-th medium.

Return type:

complex

Raises:

SMRTError – If the permittivity model is not defined.

basic_checks()

Performs basic input checks on the layer information.

Checks:
  • Temperature is between 100 and the freezing point (Kelvin units check).

  • Density is between 1 and DENSITY_OF_ICE (SI units check).

  • Layer thickness is above zero.

Raises:

SMRTError – If any of the checks fail.

inverted_medium()

Returns the layer with inverted autocorrelation and inverted permittivities.

Returns:

A new layer object with inverted properties.

Return type:

Layer

Raises:

SMRTError – If the microstructure model does not support inversion.

update(**kwargs)

Updates the attributes. This method is to be used when recalculation of the state of the object is necessary. See for instance SnowLayer.

Parameters:

**kwargs

get_microstructure_model(modulename, classname=None)

Returns the class corresponding to the microstructure_model defined in modulename.

This function imports the correct module if possible and returns the class. It is used internally and should not be needed for normal usage.

Parameters:
  • modulename – name of the python module in smrt/microstructure_model

  • classname – (Default value = None)

make_microstructure_model(modelname_or_class, **kwargs)

Creates a microstructure instance.

Parameters:
  • modelname_or_class (str or type) – Name of the module or directly the class.

  • **kwargs – Arguments needed for the specific autocorrelation.

  • modelname_or_class

  • **kwargs

Returns:

instance of the autocorrelation modelname with the parameters given in **kwargs

Example:

To import the StickyHardSpheres class with spheres radius of 1mm, stickiness of 0.5 and fractional_volume of 0.3:

shs = make_autocorrelation(“StickyHardSpheres”, radius=0.001, stickiness=0.5, frac_volume=0.3)

layer_properties(*required_arguments, optional_arguments=None, **kwargs)

This decorator is used for the permittivity functions (or any other functions) to inject layer’s attributes as arguments. The decorator declares the layer properties needed to call the function and the optional ones. This allows permittivity functions to use any property of the layer, as long as it is defined.

Parameters:
  • *required_arguments

  • optional_arguments – (Default value = None)

  • **kwargs

smrt.core.lib module

class_specializer(domain: str, cls: str | Type, **options) Type

Returns a subclass of cls (imported from the domain if cls is a string) that use the provided “options” for __init__. This is equivalent to functools.partial but for a class.

This is the same idea as: https://stackoverflow.com/questions/38911146/python-equivalent-of-functools-partial-for-a-class-constructor

class smrt_diag(arr)

Bases: object

Scipy.sparse is very slow for diagonal matrix and numpy has no good support for linear algebra. This diag class implements simple diagonal object without numpy subclassing (but without much features). It seems that proper subclassing numpy and overloading matmul is a very difficult problem.

class smrt_matrix(mat, mtype=None)

Bases: object

SMRT uses two formats of matrix: one most suitable to implement emmodel where equations are different for each polarization and another one suitable for DORT computation where stream and polarization are collapsed in one dimension to allow matrix operation. In addition, the reflection and transmission matrix are often diagonal matrix, which needs to be handled because it saves space and allow much faster operations. This class implemented all these features.

compress(mode=None, auto_reduce_npol=False)

Compresses a matrix. This comprises several actions: 1) select one mode, if relevant (dense5, and diagonal5). 2) reduce the number of polarization from 3 to 2 if mode==0 and auto_reduce_npol=True. 3) convert the format of the matrix to compressed numpy, involving a change of the dimension order (pola and streams are merged).

is_equal_zero(m)

Returns true if the smrt matrix is null

generic_ft_even_matrix(phase_function, m_max, nsamples=None)

Calculation of the Fourier decomposed of the phase or reflection or transmission matrix provided by the function.

This method calculates the 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]
Parameters:
  • phase_function – must be a function taking dphi as input. It is assumed that phi is symmetrical (it is in cos(phi))

  • m_max – maximum Fourier decomposition mode needed

set_max_numerical_threads(nthreads)

Sets the maximum number of threads for a few known library. This is useful to disable parallel computing in SMRT when using parallel computing to call multiple // SMRT runs. This avoid over-committing the CPUs and results in much better performance. Inspire from joblib.

cached_roots_legendre(n)

Cache roots_legendre results to speed up calls of the fixed_quad function.

smrt.core.model module

A model in SMRT is composed of the electromagnetic scattering theory (smrt.emmodel) and the radiative transfer solver (smrt.rtsolver).

The smrt.emmodel computes the scattering and absorption coefficients and the phase function of a layer. It is applied to each layer, and it is even possible to choose different emmodels for each layer (e.g., for a complex medium made of different materials: snow, soil, water, atmosphere, …). The smrt.rtsolver propagates the incident or emitted energy through the layers, up to the surface, and eventually through the atmosphere.

To build a model, use the make_model() function with the type of emmodel and type of rtsolver as arguments. Then call the Model.run() method of the model instance by specifying the sensor (smrt.core.sensor.Sensor), snowpack (smrt.core.snowpack.Snowpack) and optionally atmosphere (see smrt.atmosphere). The results are returned as a Result which can then been interrogated to retrieve brightness temperature, backscattering coefficient, and other information.

Example:

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

result = m.run(sensor, snowpack)  # sensor and snowpack are created before

print(result.TbV())

The model can be run on a list of snowpacks or even more conveniently on a pandas.Series or pandas.DataFrame including snowpacks. The first advantage is that by setting parallel_computation=True, the Model.run() method performs the simulation in parallel on all the available cores of your machine and even possibly remotely on a high performance cluster using dask. The second advantage is that the returned Result object contains all the simulations and provide an easier way to plot the results or compute statistics.

If a list of snowpacks is provided, it is recommended to also set the snowpack_dimension argument. It takes the form of a tuple (list of snowpack_dimension values, dimension name). The name and values are used to define the coordinates in the Result object. This is useful with timeseries or sensitivity analysis for instance.

Example:

snowpacks = []
times = []
for file in filenames:
    #  create a snowpack for each time series
    sp = ...
    snowpacks.append(sp)
    times.append(sp)

# now run the model

res = m.run(sensor, snowpacks, snowpack_dimension=('time', times))

The res variable has now a coordinate time and res.TbV() returns a timeseries.

Using pandas.Series offers an even more elegant way to run SMRT and assemble the results of all the simulations.

Example:

thickness_list = np.arange(0, 10, 1)
snowpacks = pd.Series([make_snowpack(thickness=t, ........) for t in thickness_list], index=thickness_list)
# snowpacks is a pandas Series of snowpack objects with the thickness as index

# now run the model

res = m.run(sensor, snowpacks, parallel_computation=True)

# convert the result into a datframe
res = res.to_dataframe()

The res variable is a dataframe with the thickness as index and the channels of the sensor as column.

Using pandas.DataFrame is similar. One column must contain Snowpack objects (see snowpack_column argument). The results of the simulations are automatically joined with this dataframe and returned by to_dataframe() or to_dataframe().

Example:

# df is a DataFrame with several parameters in each row.

# add a snowpack object for each row
df['snowpack'] = [make_snowpack(thickness=row['thickness'], ........) for i, row in df.iterrows()]]

# now run the model
res = m.run(sensor, snowpacks, parallel_computation=True)

# convert the result into a datframe
res = res.to_dataframe()

The res variable is a pandas.DataFrame equal to df + the results at all sensor’s channel added.

Most rtsolvers and some emmodels take arguments (usually optional but still useful) that can be specified in two ways in make_model. Either using the rtsolver_options and emmodel_options arguments of that function or using the functions make_rtsolver() and make_emmodel() to build a new class where the prescribed options are applied by default.

Examples of usage:

make_model("iba", "dort", rtsolver_options=dict(n_max_stream=128))   # original approach to specify options

make_model("iba", make_rtsolver("dort", n_max_stream=128))                # newer approach that is more readible

Both are equivalent. There is no plan to depreciate the original approach that has some nice use-cases.

make_model(emmodel, rtsolver=None, emmodel_options=None, rtsolver_options=None, emmodel_kwargs=None, rtsolver_kwargs=None)

Creates a new model with a given EM model and RT solver. The model is then ready to be run using the Model.run() method. This function is the privileged way to create models compared to class instantiation. It supports automatic import of the emmodel and rtsolver modules.

Parameters:

emmodel (string or class; or list of strings or classes; or dict of strings or classes.) – type of emmodel to use. Can be given by the name of a file/module in the emmodel directory (as a string) or a class.

List (and dict, respectively) can be provided when a different emmodel is needed for every layer (or every kind of layer medium). If a list of emmodels is given, the size must be the same as the number of layers in the snowpack. If a dict is given, the keys are the kinds of medium and the values are the associated emmodels to each sort of medium. The layer attribute ‘medium’ is used to determine the emmodel to use for each layer.

rtsolver (string or class., optional): type of RT solver to use. Can be given by the name of a file/module in the rtsolver directeory (as a string)

or a class. This argument is optional when only the computation of the layer electromagnetic properties is needed. (Default value = None)

emmodel_options (dict or a list of dict. In the latter case, the size of the list must be the same as, optional): arguments used to create the emmodel instance of each layer. Valid arguments depend on the

selected emmodel (refer to the documentation of the selected emmodel). The function emmodel() provides an alternative to setting emmodel_options. the number of layers in the snowpack. (Default value = None)

rtsolver_options (dict, optional): arguments used to create the rtsolver instance (refer to the documentation of the rtsolvers). The function

rtsolver() provides an alternative to setting rtsolver_options. (Default value = None)

emmodel_kwargs: (Default value = None) rtsolver_kwargs: (Default value = None)

Returns:

a model instance

make_rtsolver(rtsolver_class: str | Type, **options) Type

Returns a rtsolver subclass of cls (either given as a string or a class) where the provided options are applied to __init__.

Parameters:
  • rtsolver_class (Union[str, Type])

  • **options

Returns:

This function provides an alternative to setting rtsolver_options in make_model()).

Return type:

Type

Example:

make_model(..., make_rtsolver("dort", n_max_stream=128))
make_emmodel(emmodel_class: str | Type, **options) Type

Returns a emmodel subclass of cls (either given as a string or a class) where the provided options are applied to __init__.

Parameters:
  • emmodel_class (Union[str, Type])

  • **options

Returns:

This function provides an alternative to setting emmodel_options in make_model()).

Return type:

Type

Example:

make_model(make_emmodel("iba", dense_snow_correction=True), ...)
get_emmodel(emmodel)

Returns an emmodel class from the file name ‘emmodel’.

Parameters:

emmodel

make_emmodel_instance(emmodel, sensor, layer, **emmodel_options)

Creates a new emmodel instance based on the emmodel class or string. This function used to be called make_emmodel but has been renamed from SMRT v1.4 and will soon be depreciated. It is recommended to use instead:

em = make_emmodel(emmodel)(sensor, layer, **emmodel_options)

or:

emmodel_class = make_emmodel(emmodel)
em = emodel_class(sensor, layer, **emmodel_options)
Parameters:
  • emmodel – type of emmodel to use. Can be given by the name of a file/module in the emmodel directory (as a string) or a class.

  • sensor – sensor to use for the calculation.

  • layer – layer to use for the calculation

  • **emmodel_options

class Model(emmodel, rtsolver, emmodel_options=None, rtsolver_options=None)

Bases: object

Drives the whole calculation.

set_rtsolver_options(options=None, **kwargs)

Sets the option for the rtsolver.

Parameters:
  • options – (Default value = None)

  • **kwargs

set_emmodel_options(options=None, **kwargs)

Sets the options for the emmodel.

Parameters:
  • options – (Default value = None)

  • **kwargs

run(sensor, snowpack, atmosphere=None, snowpack_dimension=None, snowpack_column='snowpack', progressbar=False, parallel_computation=False, runner=None)

Runs the model for the given sensor configuration and returns the results.

Parameters:

sensor – sensor to use for the calculation. Can be a list of the same size as the snowpack list.

In this case, the computation is performed for each pair (sensor, snowpack).

snowpack: snowpack to use for the calculation. Can be a single snowpack, a list of snowpack, a dict of snowpack or

a SensitivityStudy object.

atmosphere: (Default value = None) snowpack_dimension: name and values (as a tuple) of the dimension to create for the results when a list of snowpack

is provided. E.g. time, point, longitude, latitude. By default the dimension is called ‘snowpack’ and the values are rom 1 to the number of snowpacks.

snowpack_column: when snowpack is a DataFrame this argument is used to specify which column contians the Snowpack objects (Default value = ‘snowpack’) progressbar: if True, display a progress bar during multi-snowpacks computation (Default value = False) parallel_computation: if True, use the joblib library to run the simulations of many snowpacks in parallel.

Otherwise, the simulations are run sequentially, one after one. See ‘runner’ for a more advanced control on parallel computations. Note for users seeking performances: numpy and scipy usually also perform low- level parallel computations that may (inefficiently) interact with the high-level parallelism activated by parallel_computation. For this reason joblib and other parallel runners try to desactivate numpy and scipy low-level parallelism (see set_max_numerical_threads()) to maximize performances. Conversely it means that when parallel_computation is False, the simulations are run sequentially, but numpy and scipy parallelism is NOT disabled. If you really want to use a single core for the simulations, you must first call set_max_numerical_threads() with 1 as argument and then call Model.run with parallel_computation=False. (Default value = False)

runner: a ‘runner’ is a function (or more likely a class with a __call__ method) that takes a function and a

list/generator of simulations, executes the function on each simulation and returns a list of results. ‘parallel_computation’ allows to select between two default (basic) runners (sequential and joblib). Use ‘runner’ for more advanced parallel distributed computations. To develop a costum runner, see the implementation of JoblibParallelRunner for instance.

Returns:

result of the calculation(s) as a Results instance

class SequentialRunner(progressbar)

Bases: object

Runs the simulations sequentially on a single (local) core. This is the most simple way to run smrt simulations, but the efficiency is poor.

class JoblibParallelRunner(progressbar, backend='loky', n_jobs=-1, max_numerical_threads=1)

Bases: object

Runs the simulations on the local machine on all the cores, using the joblib library for parallelism.

smrt.core.plugin module

import_class(scope: str, modulename: str, classname: str | None = None) Type

Imports the modulename and return either the class named “classname” or the first class defined in the module if classname is None.

Args: scope: scope where to search for the module. modulename: name of the module to load. classname: name of the class to read from the module.

smrt.core.sensitivity_study module

SensitivityStudy is used to easily conduct sensitivity studies. This class may be depreciated in the future. A more modern alternative is to use a pandas.DataFrame of snowpacks.

Example:

times = [datetime(2012, 1, 1), datetime(2012, 1, 5), , datetime(2012, 1, 10)]
snowpacks = SensitivityStudy("time", times, [snowpack_1jan, snowpack_5jan, snowpack_10jan])

res = model.run(sensor, snowpacks)

The res variable is a Result instance, so that for all the methods of this class that can be called, they will return a timeseries. For instance result.TbV(theta=53) returns a time-series of brightness temperature at V polarization and 53° incidence angle and the following code plots this timeseries:

plot(times, result.TbV(theta=53))
sensitivity_study(name, values, snowpacks)

Creates a sensitivity study

Parameters:
  • name – name of the variable to investigate

  • values – values taken by the variable

  • snowpacks – list of snowpacks. Can be a sequence or a function that takes one argument and return a snowpack.

In the latter case, the function is called for each values to build the list of snowpacks

smrt.core.sensor module

The sensor configuration includes all the information describing the sensor viewing geometry (incidence, …) and operating parameters (frequency, polarization, …). The easiest and recommended way to create a Sensor instance is to use one of the convenience functions such as passive(), active(), amsre(), etc. Adding a function for a new or unlisted sensor can be done in sensor_list if the sensor is common and of general interest. Otherwise, we recommend to add these functions in your own files (outside of smrt directories).

passive(frequency, theta, polarization=None, channel_map=None, name=None)

Generic configuration for passive microwave sensor.

Return a Sensor for a microwave radiometer with given frequency, incidence angle and polarization

Parameters:
  • frequency – frequency in Hz

  • theta – viewing angle or list of viewing angles in degrees from vertical. Note that some RT solvers compute all

viewing angles whatever this configuration because it is internally needed part of the multiple scattering calculation. It it therefore often more efficient to call the model once with many viewing angles instead of calling it many times with a single angle.

polarization (list of characters, optional): H and/or V polarizations. Both polarizations is the default. Note that most RT solvers compute all

the polarizations whatever this configuration because the polarizations are coupled in the RT equation.

channel_map (dict, optional): map channel names (keys) to configuration (values). A configuration is a dict with frequency, polarization and other

such parameters to be used by Result to select the results. (Default value = None)

name (string, optional): name of the sensor (Default value = None)

Returns:

py:class:Sensor instance

Usage example:

from smrt import sensor_list radiometer = sensor_list.passive(18e9, 50) radiometer = sensor_list.passive(18e9, 50, “V”) radiometer = sensor_list.passive([18e9,36.5e9], [50,55], [“V”,”H”])

channel_map_for_radar(frequency=None, polarization='HV', order='fp')
Parameters:
  • frequency – (Default value = None)

  • polarization – (Default value = “HV”)

  • order – (Default value = “fp”)

Returns:

in GHz with leading 0 if necessary. The polarization is after the frequency if order is ‘fp’ and before if order is ‘pf’.

active(frequency, theta_inc, theta=None, phi=None, polarization_inc=None, polarization=None, channel_map=None, name=None)

Configuration for active microwave sensor.

Return a Sensor for a radar with given frequency, incidence and viewing angles and polarization

If polarizations are not specified, quad-pol is the default (VV, VH, HV and HH). If the angle of incident radiation is not specified, backscatter will be simulated

Parameters:
  • frequency – frequency in Hz

  • theta_inc – incident angle in degrees from the vertical

  • theta – viewing zenith angle in degrees from the vertical. By default, it is equal to theta_inc which corresponds

to the backscatter direction

phi: viewing azimuth angle in degrees from the incident direction. By default, it is pi which corresponds

to the backscatter direction

polarization_inc (list of 1-character strings, optional): list of polarizations of the incidence wave (‘H’ or ‘V’ or both.) (Default value = None) polarization (list of 1-character strings, optional): list of viewing polarizations (‘H’ or ‘V’ or both) (Default value = None) channel_map (dict, optional): map channel names (keys) to configuration (values). A configuration is a dict with frequency, polarization and other

such parameters to be used by Result to select the results. (Default value = None)

name (string, optional): name of the sensor (Default value = None)

Returns:

py:class:Sensor instance

Usage example:

from smrt import sensor_list scatterometer = sensor_list.active(frequency=18e9, theta_inc=50) scatterometer = sensor_list.active(18e9, 50, 50, 0, “V”, “V”) scatterometer = sensor_list.active([18e9,36.5e9], theta=50, theta_inc=50, polarization_inc=[“V”, “H”], polarization=[“V”, “H”])

class Sensor(frequency=None, theta_inc_deg=None, theta_deg=None, phi_deg=None, polarization_inc=None, polarization=None, channel_map=None, name=None, wavelength=None)

Bases: SensorBase

Configuration for sensor. Use of the functions passive(), active(), or the sensor specific functions e.g. amsre() are recommended to access this class.

property mode

“A” for active or “P” for passive.

Type:

Returns the mode of observation

iterate(axis)

Iterates over the configuration for the given axis.

Parameters:

axis – one of the attribute of the sensor (frequency, …) to iterate along

class Altimeter(frequency, altitude, beamwidth, pulse_bandwidth, sigma_p=None, antenna_gain=1, off_nadir_angle=0, beam_asymmetry=0, ngate=1024, nominal_gate=40, theta_inc_deg=0.0, polarization_inc=None, polarization=None, channel=None)

Bases: Sensor

Configuration for altimeter. Use of the functions altimeter(), or the sensor specific functions e.g. envisat_ra2() are recommended to access this class.

smrt.core.snowpack module

Snowpack instance contains the description of the snowpack, including a list of layers and interfaces between the layers, and the substrate (soil, ice, …).

To create a snowpack, it is recommended to use the make_snowpack() function which avoids the complexity of creating each layer and then the snowpack from the layers. For more complex media (like lake ice or sea ice), it may be necessary to directly call the functions to create the different layers (such as make_snow_layer()).

Example:

# create a 10-m thick snowpack with a single layer,
# density is 350 kg/m3. The exponential autocorrelation function is
# used to describe the snow and the "size" parameter is therefore
# the correlation length which is given as an optional
# argument of this function (but is required in practice)

sp = make_snowpack([10], "exponential", [350], corr_length=[3e-3])
class Snowpack(layers=None, interfaces=None, substrate=None, atmosphere=None, terrain_info=None)

Bases: object

Holds the description of the snowpack, including the layers, interfaces, and the substrate.

property nlayer

Returns the number of layers.

property layer_thicknesses

Returns the thickness of each layer.

property layer_depths

Returns the depth of the bottom of each layer.

property bottom_layer_depths

Returns the depth of the bottom of each layer.

property top_layer_depths

Returns the depth of the bottom of each layer.

property mid_layer_depths

Returns the depth of the bottom of each layer.

property z

Returns the depth of each interface, that is, 0 and the depths of the bottom of each layer.

property layer_densities

Returns the density of each layer.

profile(property_name, where='all', raise_attributeerror=False)

Returns the vertical profile of property_name. The property is searched either in the layer, microstructure or interface.

Parameters:
  • property_name (str) – Name of the property.

  • where (str) – Where to search the property. Can be ‘all’, ‘layer’, ‘microstructure’, or ‘interface’.

  • raise_attributeerror (bool) – Raise an attribute error if the attribute is not found.

Returns:

Array of the property values.

Return type:

np.ndarray

append(layer, interface=None)

Appends a new layer at the bottom of the stack of layers. The interface is that at the top of the appended layer.

Parameters:
  • layer (Layer) – Instance of Layer.

  • interface – Type of interface. By default, flat surface (Flat) is considered meaning the coefficients are calculated with Fresnel coefficient and using the effective permittivity of the surrounding layers.

delete(ilayer)

Deletes a layer and the upper interface.

Parameters:

ilayer (int) – Index of the layer.

delete_bottom(ilayer)

Deletes the bottom of the snowpack from layer n. Deletes also the substrate.

Parameters:

ilayer (int) – Index of the first layer to delete.

copy(cut_bottom=None)

Makes a shallow copy of a snowpack by copying the list of layers and interfaces but not the layers and interfaces themselves which are still shared with the original snowpack. This method allows the user to create a new snowpack and remove, append or replace some layers or interfaces afterward. It does not allow to alter the layers or interfaces without changing the original snowpack. See deepcopy.

Parameters:

cut_bottom (int, optional) – If cut_bottom is a number, all layers below the layer indexed by ‘cut_bottom’ are removed as well as the substrate.

Returns:

The shallow copy of the snowpack.

Return type:

Snowpack

deepcopy()

Makes a deep copy of a snowpack.

Returns:

The deep copy of the snowpack.

Return type:

Snowpack