retrieval module¶
Module: skretrieval.retrieval
-
class
ForwardModel
¶ A ForwardModel is an object which is capable of calculating a radiance. This serves as the primary interface to the retrieval, along with the RetrievalTarget.
-
abstract
calculate_radiance
()¶
-
abstract
-
class
Minimizer
¶ A class which performs minimization between some aspect of measurement level1 data and the forward model simulations.
-
abstract
retrieve
(measurement_l1: skretrieval.core.radianceformat.RadianceBase, forward_model: skretrieval.retrieval.ForwardModel, retrieval_target: skretrieval.retrieval.RetrievalTarget)¶ - Parameters
measurement_l1 (RadianceBase) – The data we are trying to match, either from a real instrument or simulations.
forward_model (ForwardModel) – A model for the data in measurement_l1
retrieval_target (RetrievalTarget) – What we are trying to retrieve
- Returns
Various parameters specific to the minimizer
- Return type
dict
-
abstract
-
class
RetrievalTarget
¶ The retrieval target defines the parameter that is to be retrieved, and also what measurements are going to be used to retrieve it. Notation is similar to that of Rodgers.
-
adjust_parameters
(forward_model, y_dict, chi_sq, chi_sq_linear, iter_idx, predicted_delta_y)¶
-
abstract
apriori_state
() → numpy.array¶ - Returns
Apriori state vector, x_a. If no apriori is used return None
- Return type
np.array
-
abstract
inverse_apriori_covariance
()¶ - Returns
Inverse of the apriori covariance matrix. If no apriori is used return None.
- Return type
np.array
-
abstract
measurement_vector
(l1_data: skretrieval.core.radianceformat.RadianceBase)¶ - Parameters
l1_data (RadianceBase) – Radiance data. Usually this is an instrument specific instance of RadianceBase, and the RetrievalTarget only works with specific formats.
- Returns
Keys ‘y’ for the measurement vector, ‘jacobian’ for the jacobian of the measurement vector (if weighting functions are in l1_data, ‘y_error’ the covariance of ‘y’ (if error information is provided in l1_data)
- Return type
dict
-
static
measurement_vector_allowed_to_change
()¶ - Returns
True if the measurement_vector may change shape between iterations, False otherwise.
- Return type
bool
-
abstract
state_vector
()¶ - Returns
The state vector, x
- Return type
np.array
-
static
state_vector_allowed_to_change
()¶ - Returns
True if the state vector/apriori may change shape between iterations, False otherwise.
- Return type
bool
-
abstract
update_state
(x: numpy.ndarray)¶ Updates the state for the new state vector. Note that this change has to propagate backwards to the ForwardModel somehow. Typically this is done by passing a climatology into the RetrievalTarget at initiliazation which is used in the ForwardModel.
- Parameters
x (np.array) – New state vector
-
ForwardModel¶
Minimizer¶
-
class
Minimizer
¶ A class which performs minimization between some aspect of measurement level1 data and the forward model simulations.
-
abstract
retrieve
(measurement_l1: skretrieval.core.radianceformat.RadianceBase, forward_model: skretrieval.retrieval.ForwardModel, retrieval_target: skretrieval.retrieval.RetrievalTarget)¶ - Parameters
measurement_l1 (RadianceBase) – The data we are trying to match, either from a real instrument or simulations.
forward_model (ForwardModel) – A model for the data in measurement_l1
retrieval_target (RetrievalTarget) – What we are trying to retrieve
- Returns
Various parameters specific to the minimizer
- Return type
dict
-
abstract
AirDensityRetrieval¶
-
class
AirDensityRetrieval
(air_species: sasktran.species.Species, ret_wavel=350, high_alt_normalize=False, tikh_factor=1)¶ Limb retrieval for air number density. Radiances should be supplied in the RadianceGridded format.
-
apriori_state
()¶ - Returns
Apriori state vector, x_a. If no apriori is used return None
- Return type
np.array
-
inverse_apriori_covariance
()¶ - Returns
Inverse of the apriori covariance matrix. If no apriori is used return None.
- Return type
np.array
-
measurement_vector
(l1_data: skretrieval.core.radianceformat.RadianceBase)¶ The measurement vector is the logarithm of radiance at ret_wavel
-
state_vector
()¶ State vector is the logarithm of air number density
-
temperature
(hires_spacing=100, T0=200, earth_radius=6372000)¶ Integrates the air number density using hydrostatic balance and the ideal gas law.
- Parameters
hires_spacing (float, optional) – Internal high resolution spacing to use for the integral. Default 100 m
T0 (float, optional) – Upper altitude pin temperature, this is at self._retrieval_altitudes[-1]. Default 200 K
earth_radius (float, optional) – Earth radius in m, necessary to approximate gravity. Default 6372000
- Returns
Temperature in K on self._retrieval_altitudes
- Return type
np.array
-
update_state
(x: numpy.ndarray)¶ Updates the state for the new state vector. Note that this change has to propagate backwards to the ForwardModel somehow. Typically this is done by passing a climatology into the RetrievalTarget at initiliazation which is used in the ForwardModel.
- Parameters
x (np.array) – New state vector
-
OzoneRetrieval¶
-
class
OzoneRetrieval
(ozone_species: sasktran.species.Species)¶ -
apriori_state
()¶ - Returns
Apriori state vector, x_a. If no apriori is used return None
- Return type
np.array
-
inverse_apriori_covariance
()¶ - Returns
Inverse of the apriori covariance matrix. If no apriori is used return None.
- Return type
np.array
-
measurement_vector
(l1_data: skretrieval.core.radianceformat.RadianceBase)¶ - Parameters
l1_data (RadianceBase) – Radiance data. Usually this is an instrument specific instance of RadianceBase, and the RetrievalTarget only works with specific formats.
- Returns
Keys ‘y’ for the measurement vector, ‘jacobian’ for the jacobian of the measurement vector (if weighting functions are in l1_data, ‘y_error’ the covariance of ‘y’ (if error information is provided in l1_data)
- Return type
dict
-
state_vector
()¶ - Returns
The state vector, x
- Return type
np.array
-
update_state
(x: numpy.ndarray)¶ Updates the state for the new state vector. Note that this change has to propagate backwards to the ForwardModel somehow. Typically this is done by passing a climatology into the RetrievalTarget at initiliazation which is used in the ForwardModel.
- Parameters
x (np.array) – New state vector
-
Rodgers¶
-
class
Rodgers
(max_iter=10, lm_damping=0, iterative_update_lm=False, retreat_lm=False, lm_change_factor=1.5, convergence_factor=1)¶ -
retrieve
(measurement_l1, forward_model: skretrieval.retrieval.ForwardModel, retrieval_target: skretrieval.retrieval.RetrievalTarget)¶ - Parameters
measurement_l1 (RadianceBase) – The data we are trying to match, either from a real instrument or simulations.
forward_model (ForwardModel) – A model for the data in measurement_l1
retrieval_target (RetrievalTarget) – What we are trying to retrieve
- Returns
Various parameters specific to the minimizer
- Return type
dict
-
Tikhonov Functions¶
-
two_dim_vertical_second_deriv
(numangle, numalt, factor=1, sparse=False)¶ Calculates the second derivatvie Tikhonov regularization matrix for a two dimensional uniform grid. The matrix is calculated assuming that the measurement vector is constructed with altitude being the leading dimension
- Parameters
numangle (scalar) – The number of angular grid points
numalt (scalar) – The number of altitude grid points
factor (scalar or length numalt, optional) – If scalar, the resulting matrix is multiplied by this value. If a vector of length numalt, then each altitude level is multiplied by its corresponding factor. Default is 1