gstools.krige.Krige
- class gstools.krige.Krige(model, cond_pos, cond_val, drift_functions=None, ext_drift=None, mean=None, normalizer=None, trend=None, unbiased=True, exact=False, cond_err='nugget', pseudo_inv=True, pseudo_inv_type='pinv', fit_normalizer=False, fit_variogram=False)[source]
Bases:
Field
A Swiss Army knife for kriging.
A Kriging class enabling the basic kriging routines: Simple-, Ordinary-, Universal-, External Drift- and detrended/regression-Kriging as well as Kriging the Mean [Wackernagel2003].
- Parameters:
model (
CovModel
) – Covariance Model used for kriging.cond_pos (
list
) – tuple, containing the given condition positions (x, [y, z])cond_val (
numpy.ndarray
) – the values of the conditions (nan values will be ignored)drift_functions (
list
ofcallable
,str
orint
) – Either a list of callable functions, an integer representing the polynomial order of the drift or one of the following strings:“linear” : regional linear drift (equals order=1)
“quadratic” : regional quadratic drift (equals order=2)
ext_drift (
numpy.ndarray
orNone
, optional) – the external drift values at the given cond. positions.mean (
float
, optional) – mean value used to shift normalized conditioning data. Could also be a callable. The default is None.normalizer (
None
orNormalizer
, optional) – Normalizer to be applied to the input data to gain normality. The default is None.trend (
None
orfloat
orcallable
, optional) – A callable trend function. Should have the signature: f(x, [y, z, …]) This is used for detrended kriging, where the trended is subtracted from the conditions before kriging is applied. This can be used for regression kriging, where the trend function is determined by an external regression algorithm. If no normalizer is applied, this behaves equal to ‘mean’. The default is None.unbiased (
bool
, optional) – Whether the kriging weights should sum up to 1, so the estimator is unbiased. If unbiased is False and no drifts are given, this results in simple kriging. Default: Trueexact (
bool
, optional) – Whether the interpolator should reproduce the exact input values. If False, cond_err is interpreted as measurement error at the conditioning points and the result will be more smooth. Default: Falsecond_err (
str
, :classfloat
orlist
, optional) – The measurement error at the conditioning points. Either “nugget” to apply the model-nugget, a single value applied to all points or an array with individual values for each point. The “exact=True” variant only works with “cond_err=’nugget’”. Default: “nugget”pseudo_inv (
bool
, optional) – Whether the kriging system is solved with the pseudo inverted kriging matrix. If True, this leads to more numerical stability and redundant points are averaged. But it can take more time. Default: Truepseudo_inv_type (
str
orcallable
, optional) –Here you can select the algorithm to compute the pseudo-inverse matrix:
“pinv”: use pinv from scipy which uses SVD
“pinvh”: use pinvh from scipy which uses eigen-values
If you want to use another routine to invert the kriging matrix, you can pass a callable which takes a matrix and returns the inverse. Default: “pinv”
fit_normalizer (
bool
, optional) – Whether to fit the data-normalizer to the given conditioning data. Default: Falsefit_variogram (
bool
, optional) – Whether to fit the given variogram model to the data. Directional variogram fitting is triggered by setting any anisotropy factor of the model to anything unequal 1 but the main axes of correlation are taken from the model rotation angles. If the model is a spatio-temporal latlon model, this will raise an error. This assumes the sill to be the data variance and with standard bins provided by thestandard_bins
routine. Default: False
Notes
If you have changed any properties in the class, you can update the kriging setup by calling
Krige.set_condition
without any arguments.References
[Wackernagel2003]Wackernagel, H., “Multivariate geostatistics”, Springer, Berlin, Heidelberg (2003)
- Attributes:
all_fields
list
: All fields as stacked list.cond_err
list
: The measurement errors at the condition points.cond_ext_drift
numpy.ndarray
: The ext. drift at the conditions.cond_mean
numpy.ndarray
: Trend at the conditions.cond_no
int
: The number of the conditions.cond_pos
list
: The position tuple of the conditions.cond_trend
numpy.ndarray
: Trend at the conditions.cond_val
list
: The values of the conditions.dim
int
: Dimension of the field.drift_functions
drift_no
int
: Number of drift values per point.exact
bool
: Whether the interpolator is exact.ext_drift_no
int
: Number of external drift values per point.field_names
list
: Names of present fields.field_shape
tuple
: The shape of the field.has_const_mean
bool
: Whether the field has a constant mean or not.int_drift_no
int
: Number of internal drift values per point.krige_size
int
: Size of the kriging system.latlon
bool
: Whether the field depends on geographical coords.mean
mesh_type
str
: The mesh type of the field.model
CovModel
: The covariance model of the field.name
str
: The name of the class.normalizer
Normalizer
: Normalizer of the field.pos
tuple
: The position tuple of the field.pseudo_inv
bool
: Whether pseudo inverse matrix is used.pseudo_inv_type
str
: Method selector for pseudo inverse calculation.temporal
bool
: Whether the field depends on time.trend
unbiased
bool
: Whether the kriging is unbiased or not.value_type
str
: Type of the field values (scalar, vector).
Methods
__call__
([pos, mesh_type, ext_drift, ...])Generate the kriging field.
delete_fields
([select])Delete selected fields.
get_mean
([post_process])Calculate the estimated mean of the detrended field.
get_store_config
(store[, default, fld_cnt])Get storage configuration from given selection.
mesh
(mesh[, points, direction, name])Generate a field on a given meshio, ogs5py or PyVista mesh.
plot
([field, fig, ax])Plot the spatial random field.
post_field
(field[, name, process, save])Postprocessing field values.
pre_pos
([pos, mesh_type, info])Preprocessing positions and mesh_type.
set_condition
([cond_pos, cond_val, ...])Set the conditions for kriging.
set_drift_functions
([drift_functions])Set the drift functions for universal kriging.
set_pos
(pos[, mesh_type, info])Set positions and mesh_type.
structured
(*args, **kwargs)Generate a field on a structured mesh.
to_pyvista
([field_select, fieldname])Create a VTK/PyVista grid of the stored field.
transform
(method[, field, store, process])Apply field transformation.
unstructured
(*args, **kwargs)Generate a field on an unstructured mesh.
vtk_export
(filename[, field_select, fieldname])Export the stored field to vtk.
- __call__(pos=None, mesh_type='unstructured', ext_drift=None, chunk_size=None, only_mean=False, return_var=True, post_process=True, store=True)[source]
Generate the kriging field.
The field is saved as self.field and is also returned. The error variance is saved as self.krige_var and is also returned.
- Parameters:
pos (
list
, optional) – the position tuple, containing main direction and transversal directions (x, [y, z])mesh_type (
str
, optional) – ‘structured’ / ‘unstructured’ext_drift (
numpy.ndarray
orNone
, optional) – the external drift values at the given positions (only for EDK)chunk_size (
int
, optional) – Chunk size to cut down the size of the kriging system to prevent memory errors. Default: Noneonly_mean (
bool
, optional) – Whether to only calculate the mean of the kriging field. Default: Falsereturn_var (
bool
, optional) – Whether to return the variance along with the field. Default: Truepost_process (
bool
, optional) – Whether to apply mean, normalizer and trend to the field. Default: Truestore (
str
orbool
orlist
, optional) – Whether to store kriging fields (True/False) with default name or with specified names. The default isTrue
for default names [“field”, “krige_var”] or “mean_field” if only_mean=True.
- Returns:
field (
numpy.ndarray
) – the kriged field or mean_fieldkrige_var (
numpy.ndarray
, optional) – the kriging error variance (if return_var is True and only_mean is False)
- delete_fields(select=None)
Delete selected fields.
- get_mean(post_process=True)[source]
Calculate the estimated mean of the detrended field.
- Parameters:
post_process (
bool
, optional) – Whether to apply field-mean and normalizer. Default: True- Returns:
mean – Mean of the Kriging System.
- Return type:
Notes
Only not
None
if the Kriging System has a constant mean. This means, no drift is given and the given field-mean is constant. The result is neglecting a potential given trend.
- get_store_config(store, default=None, fld_cnt=None)
Get storage configuration from given selection.
- Parameters:
store (
str
orbool
orlist
, optional) – Whether to store fields (True/False) with default names or with specified names. The default isTrue
for default names.default (
str
orlist
, optional) – Default field names. The default is “field”.fld_cnt (
None
orint
, optional) – Number of fields when using lists. The default is None.
- Returns:
- mesh(mesh, points='centroids', direction='all', name='field', **kwargs)
Generate a field on a given meshio, ogs5py or PyVista mesh.
- Parameters:
mesh (meshio.Mesh or ogs5py.MSH or PyVista mesh) – The given mesh
points (
str
, optional) – The points to evaluate the field at. Either the “centroids” of the mesh cells (calculated as mean of the cell vertices) or the “points” of the given mesh. Default: “centroids”direction (
str
orlist
, optional) – Here you can state which direction should be chosen for lower dimension. For example, if you got a 2D mesh in xz direction, you have to pass “xz”. By default, all directions are used. One can also pass a list of indices. Default: “all”name (
str
orlist
ofstr
, optional) – Name(s) to store the field(s) in the given mesh as point_data or cell_data. If to few names are given, digits will be appended. Default: “field”**kwargs – Keyword arguments forwarded to
__call__
.
Notes
This will store the field in the given mesh under the given name, if a meshio or PyVista mesh was given.
- See:
- plot(field='field', fig=None, ax=None, **kwargs)
Plot the spatial random field.
- Parameters:
field (
str
, optional) – Field that should be plotted. Default: “field”fig (
Figure
orNone
) – Figure to plot the axes on. If None, a new one will be created. Default: Noneax (
Axes
orNone
) – Axes to plot on. If None, a new one will be added to the figure. Default: None**kwargs – Forwarded to the plotting routine.
- post_field(field, name='field', process=True, save=True)
Postprocessing field values.
- Parameters:
field (
numpy.ndarray
) – Field values.name (
str
, optional) – Name. to store the field. The default is “field”.process (
bool
, optional) – Whether to process field to apply mean, normalizer and trend. The default is True.save (
bool
, optional) – Whether to store the field under the given name. The default is True.
- Returns:
field – Processed field values.
- Return type:
- pre_pos(pos=None, mesh_type='unstructured', info=False)
Preprocessing positions and mesh_type.
- Parameters:
- Returns:
iso_pos ((d, n),
numpy.ndarray
) – Isometrized position tuple.shape (
tuple
) – Shape of the resulting field.info (
dict
, optional) – Information about settings.
Warning
When setting a new position tuple that differs from the present one, all stored fields will be deleted.
- set_condition(cond_pos=None, cond_val=None, ext_drift=None, cond_err=None, fit_normalizer=False, fit_variogram=False)[source]
Set the conditions for kriging.
This method could also be used to update the kriging setup, when properties were changed. Then you can call it without arguments.
- Parameters:
cond_pos (
list
, optional) – the position tuple of the conditions (x, [y, z]). Default: current.cond_val (
numpy.ndarray
, optional) – the values of the conditions (nan values will be ignored). Default: current.ext_drift (
numpy.ndarray
orNone
, optional) – the external drift values at the given conditions (only for EDK) For multiple external drifts, the first dimension should be the index of the drift term. When passing None, the extisting external drift will be used.cond_err (
str
, :classfloat
,list
, optional) – The measurement error at the conditioning points. Either “nugget” to apply the model-nugget, a single value applied to all points or an array with individual values for each point. The measurement error has to be <= nugget. The “exact=True” variant only works with “cond_err=’nugget’”. Default: “nugget”fit_normalizer (
bool
, optional) – Whether to fit the data-normalizer to the given conditioning data. Default: Falsefit_variogram (
bool
, optional) – Whether to fit the given variogram model to the data. Directional variogram fitting is triggered by setting any anisotropy factor of the model to anything unequal 1 but the main axes of correlation are taken from the model rotation angles. If the model is a spatio-temporal latlon model, this will raise an error. This assumes the sill to be the data variance and with standard bins provided by thestandard_bins
routine. Default: False
- set_drift_functions(drift_functions=None)[source]
Set the drift functions for universal kriging.
- Parameters:
drift_functions (
list
ofcallable
,str
orint
) – Either a list of callable functions, an integer representing the polynomial order of the drift or one of the following strings:“linear” : regional linear drift (equals order=1)
“quadratic” : regional quadratic drift (equals order=2)
- Raises:
ValueError – If the given drift functions are not callable.
- set_pos(pos, mesh_type='unstructured', info=False)
Set positions and mesh_type.
- Parameters:
- Returns:
info – Information about settings.
- Return type:
dict
, optional
Warning
When setting a new position tuple that differs from the present one, all stored fields will be deleted.
- to_pyvista(field_select='field', fieldname='field')
Create a VTK/PyVista grid of the stored field.
- transform(method, field='field', store=True, process=False, **kwargs)
Apply field transformation.
- Parameters:
method (
str
) – Method to use. Seegstools.transform
for available transformations.field (
str
, optional) – Name of field to be transformed. The default is “field”.store (
str
orbool
, optional) – Whether to store field inplace (True/False) or under a given name. The default is True.process (
bool
, optional) – Whether to process in/out fields with trend, normalizer and mean of given Field instance. The default is False.**kwargs – Keyword arguments forwarded to selected method.
- Raises:
ValueError – When method is unknown.
- Returns:
Transformed field.
- Return type:
- vtk_export(filename, field_select='field', fieldname='field')
Export the stored field to vtk.
- Parameters:
filename (
str
) – Filename of the file to be saved, including the path. Note that an ending (.vtr or .vtu) will be added to the name.field_select (
str
, optional) – Field that should be stored. Can be: “field”, “raw_field”, “krige_field”, “err_field” or “krige_var”. Default: “field”fieldname (
str
, optional) – Name of the field in the VTK file. Default: “field”
- property cond_ext_drift
The ext. drift at the conditions.
- Type:
- property cond_mean
Trend at the conditions.
- Type:
- property cond_trend
Trend at the conditions.
- Type:
- property normalizer
Normalizer of the field.
- Type: