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(). Other methods could be developed for specific needs.

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)

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

class Result(intensity, coords=None)

Bases: object

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

coords

Return 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).

Tb(**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)

Tb_as_dataframe(**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)

TbV(**kwargs)

Return 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)

Return 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)

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

sigma(**kwargs)

Return 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_dB(**kwargs)

Return 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(**kwargs)

Return 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_dB_as_dataframe(**kwargs)

Return backscattering coefficient in dB. 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)

sigmaVV(**kwargs)

Return 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)

Return 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)

Return 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)

Return 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)

Return 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)

Return 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)

Return 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)

save(filename)

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

concat_results(result_list, coord)

Concatenate 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 array with the same length as result_list.
Returns:

Result instance