krigekit.kriging#

Python wrapper for the Fortran kriging module via ISO C Binding.

Build the shared library first:

gfortran -O2 -fPIC -fdefault-real-8 -fopenmp -shared: common.f90 utils.F90 rotation.f90 variogram.f90 kriging.F90 kriging_capi.f90 -o libkriging.so

Then use this module:

from krigekit import Kriging import numpy as np

# Spatial ordinary kriging k = Kriging(ndim=2, nvar=1) k.set_obs(ivar=1, coord=coord, value=value, nmax=20) k.set_grid(coord=grid_coord) k.set_vgm(ivar=1, jvar=1, vtype=”sph”, nugget=0, sill=1.0, a_major=1000, a_minor1=500, a_minor2=500) k.set_search(ivar=1) k.solve() est, var = k.get_results() del k # release memory

# Space-time kriging — same entry point, st=True k = Kriging(st=True, nvar=1) k.set_st_model(model=’sum_metric’, transform=’bounded’, at=5.0) …

Kriging is a factory function that returns either a _SpatialKriging or SpaceTimeKriging instance depending on the st keyword.

Classes#

Kriging

Python interface to the Fortran t_kriging spatial kriging/simulation engine.

Functions#

`omp_info`(→ dict)	Return a dict with OpenMP thread counts as seen by the Fortran runtime.
`get_omp_info`()
`ordinary_kriging`(→ tuple[numpy.ndarray, numpy.ndarray])	One-shot ordinary kriging with a single isotropic (or anisotropic) variogram.
`cokriging`(→ tuple[numpy.ndarray, numpy.ndarray])	One-shot ordinary co-kriging with multiple variables.
`sequential_gaussian_simulation`(→ numpy.ndarray)	Sequential Gaussian Simulation.

Module Contents#

krigekit.kriging.omp_info() → dict#

Return a dict with OpenMP thread counts as seen by the Fortran runtime.

Keys#

max_threadsint: Value of omp_get_max_threads() — the number of threads that will be used in the next parallel region (respects OMP_NUM_THREADS and any omp_set_num_threads() calls). Returns 1 when OpenMP is not compiled in.
num_threadsint: Value of omp_get_num_threads() measured inside an actual parallel region — the number of threads that are actually running. Returns 1 when OpenMP is not compiled in.
openmpbool: True when the library was compiled with OpenMP support.

Example

>>> import os; os.environ["OMP_NUM_THREADS"] = "4"
>>> from _kriging import omp_info
>>> omp_info()
{'max_threads': 4, 'num_threads': 4, 'openmp': True}

krigekit.kriging.get_omp_info()#

class krigekit.kriging.Kriging(ndim: int = 2, nvar: int = 1, ndrift: int = 0, unbias: int = 1, nsim: int = 0, anisotropic_search: bool = False, weight_correction: bool = False, use_old_weight: bool = False, store_weight: bool = False, cross_validation: bool = False, write_mat: bool = False, neglect_error: bool = True, varying_vgm: bool = False, std_ck: bool = False, verbose: bool = False, pf_cache: bool = False, weight_file: str = '', bounds: tuple | None = None, seed: int | None = None)#

Python interface to the Fortran t_kriging spatial kriging/simulation engine.

Array convention#

All coordinate arrays use (nobs, ndim) shape — rows are points, columns are spatial dimensions. This matches NumPy, pandas, and scikit-learn conventions. The wrapper transparently transposes to Fortran’s (ndim, nobs) before calling the library.

Typical workflow#

>>> k = Kriging(ndim=2, nvar=1)
>>> k.set_obs(ivar=1, coord=obs_coord, value=obs_value, nmax=20)
>>> k.set_grid(coord=grid_coord)
>>> k.set_vgm(ivar=1, jvar=1, vtype="sph", nugget=0, sill=1.0, a_major=1000, a_minor1=500, a_minor2=50)
>>> k.set_search(ivar=1)
>>> k.solve()
>>> estimate, variance = k.get_results()
>>> del k    # release memory

For sequential Gaussian simulation add nsim=N to the constructor and call set_sim() after set_grid().

set_obs(ivar: int, coord: numpy.ndarray, value: numpy.ndarray, variance: numpy.ndarray | None = None, nmax: int | None = None, maxdist: float | None = None, sk_mean: float = 0.0)#

Set observations for variable ivar.

Drift values are set separately via set_obs_drift() after this call, when ndrift > 0.

Parameters:

ivar (int) – Variable index, 1-based.
coord (ndarray, shape (nobs, ndim)) – Duplicate coordinate tuples within the same variable are rejected. Observation coordinates. Rows are points, columns are spatial dimensions — standard Python/NumPy convention. The wrapper transposes to Fortran’s (ndim, nobs) internally.
value (ndarray, shape (nobs,)) – Observed values.
variance (ndarray, shape (nobs,), optional) – Per-observation measurement error variance added to the diagonal of the covariance matrix. Defaults to zeros (no measurement error).
nmax (int, optional) – Maximum number of neighbours. Default: use all observations.
maxdist (float, optional) – Maximum search distance. Default: unlimited.
sk_mean (float) – Global mean for simple kriging (unbias=0). Default 0.0.

set_obs_drift(ivar: int, drift: numpy.ndarray)#

Set external drift values at observation locations for variable ivar.

Call after set_obs() for the same ivar. Only needed when ndrift > 0 was passed to the constructor.

Parameters:

ivar (int) – Variable index, 1-based.
drift (ndarray, shape (nobs, ndrift)) – Drift values. Rows are observations, columns are drift functions. Transposed to (ndrift, nobs) internally before calling Fortran.

update_obs_value(ivar: int, value: numpy.ndarray)#

Replace observation values for variable ivar in-place.

Coordinates and the kd-tree are unchanged. The primary use case is weight reuse: after solving once (with store_weight=True or use_old_weight=True), call this method with new observed values at the same locations and call solve() again to get updated estimates without recomputing search neighbourhoods or the LHS factorization.

Parameters:

ivar (int) – Variable index, 1-based.
value (ndarray, shape (nobs,)) – New observed values. Length must match the nobs passed to the previous set_obs() call for this variable.

set_vgm(ivar: int, jvar: int, vtype: str, nugget: float = 0.0, sill: float = 1.0, a_major: float = 1.0, a_minor1: float | None = None, a_minor2: float | None = None, azimuth: float = 0.0, dip: float = 0.0, plunge: float = 0.0, append: bool = True, product: bool = False)#

Add one nested variogram structure for the (ivar, jvar) pair. Call multiple times with append=True to build a nested (multi-structure) model. Pass append=False to clear any previously set structures for the pair before adding this one (useful when reusing a Kriging object with a different variogram).

Parameters:

ivar (int) – Variable indices (1-based). Use ivar=jvar for auto-variograms, ivar≠jvar for cross-variograms. The LMC constraint b12² ≤ b11 × b22 must be satisfied for each nested structure.
jvar (int) – Variable indices (1-based). Use ivar=jvar for auto-variograms, ivar≠jvar for cross-variograms. The LMC constraint b12² ≤ b11 × b22 must be satisfied for each nested structure.
vtype (str) – Variogram type: one of sph, exp, gau, pow, lin, hol, bsq, cir, nug, cyc, dco. cyc is the GP periodic (exp-sine-squared) kernel where a_major is the period; dco is a damped cosine.
nugget (float) – Nugget contribution of this structure (default 0).
sill (float) – Partial sill of this structure (default 1).
a_major (float) – Range along the major axis (default 1).
a_minor1 (float, optional) – Range along the first minor axis. Defaults to a_major (isotropic in the horizontal plane).
a_minor2 (float, optional) – Range along the second minor axis. Defaults to a_minor1.
azimuth (float) – Rotation angles in degrees (default 0).
dip (float) – Rotation angles in degrees (default 0).
plunge (float) – Rotation angles in degrees (default 0).
append (bool) – True (default) appends this structure to any existing ones. False clears the current model before adding.
product (bool) –
False (default) adds this structure to the previous ones (standard additive nesting). True multiplies this structure with the immediately preceding one, forming a product group. The Schur product of two positive-definite covariances is also positive-definite, so products are always valid.

Example — damped cosine with independent decay and period:
```
k.set_vgm(1, 1, vtype="exp", sill=1.0, a_major=5.0)   # decay over 5 yr
k.set_vgm(1, 1, vtype="hol", sill=1.0, a_major=0.5,   # 1-yr period
          product=True)
```

Example

>>> k.set_vgm(1, 1, vtype="sph", nugget=0.0, sill=1.0, a_major=500.0)
>>> k.set_vgm(1, 1, vtype="nug", nugget=0.1, sill=0.0, a_major=1.0)
>>> k.set_vgm(1, 1, vtype="sph", nugget=0.0, sill=0.9, a_major=500.0)

set_vgm_block(ib: int, ivar: int, jvar: int, vtype: str, nugget: float = 0.0, sill: float = 1.0, a_major: float = 1.0, a_minor1: float | None = None, a_minor2: float | None = None, azimuth: float = 0.0, dip: float = 0.0, plunge: float = 0.0)#

Add one nested variogram structure for a specific block ib.

Requires varying_vgm=True in the constructor and set_grid() to have been called first (because the number of blocks must be known before the per-block variogram array can be allocated in Fortran).

Call multiple times for the same ib to build a nested model.

Parameters:

ib (int) – Block index (1-based).
ivar (int) – Variable indices (1-based).
jvar (int) – Variable indices (1-based).
vtype (str) – Variogram type: sph, exp, gau, pow, lin, hol, bsq, cir, nug, cyc, or dco.
nugget (float) – Nugget contribution (default 0).
sill (float) – Partial sill (default 1).
a_major (float) – Range along the major axis (default 1).
a_minor1 (float, optional) – First minor-axis range (defaults to a_major).
a_minor2 (float, optional) – Second minor-axis range (defaults to a_minor1).
azimuth (float) – Rotation angles in degrees (default 0).
dip (float) – Rotation angles in degrees (default 0).
plunge (float) – Rotation angles in degrees (default 0).

set_grid(coord: numpy.ndarray | None = None, rangescale: numpy.ndarray | None = None, localnugget: numpy.ndarray | None = None)#

Set the estimation grid for point kriging (one node per block).

For block kriging use set_grid_block(). For cross-validation use set_grid_cv(). Drift is set separately via set_grid_drift() when ndrift > 0.

Parameters:

coord (ndarray, shape (ngrid, ndim)) – Grid coordinates. Rows are grid nodes, columns are spatial dimensions.
rangescale (ndarray, shape (ngrid,), optional) – Per-block variogram range scaling factor. Values > 1 increase the effective range, useful to account for data sparsity. Default: 1.0 for all blocks.
localnugget (ndarray, shape (ngrid,), optional) – Additional nugget added per block to model local uncertainty. Default: 0.0 for all blocks.

set_grid_block(coord: numpy.ndarray, block_type: int, nblockpnt: numpy.ndarray, pointweight: numpy.ndarray | None = None, blocksize: numpy.ndarray | None = None, rangescale: numpy.ndarray | None = None, localnugget: numpy.ndarray | None = None)#

Set the estimation grid for block kriging.

Drift is set separately via set_grid_drift() when ndrift > 0.

Parameters:

coord (ndarray, shape (ngrid, ndim)) – Sub-node coordinates across all blocks (total ngrid = sum(nblockpnt)).
block_type (int) – -4 = Gaussian quadrature nodes (auto-generated); >0 = user-supplied sub-nodes (coord contains sub-node positions).
nblockpnt (ndarray of int, shape (nblock,)) – Number of sub-nodes per block.
pointweight (ndarray, shape (sum(nblockpnt),), optional) – Weight of each sub-node. Uniform weights (1/nblockpnt) used if omitted.
blocksize (ndarray, shape (nblock,ndim), optional) – Block size in each dimension when block_type == -4.
rangescale (ndarray, shape (nblock,), optional) – Per-block variogram range scaling. Default: 1.0.
localnugget (ndarray, shape (nblock,), optional) – Per-block additional nugget. Default: 0.0.

set_grid_cv()#

Set up the grid for cross-validation mode.

No coordinate argument is needed — Fortran derives the grid from the observation coordinates automatically. Call instead of set_grid() when cross_validation=True was passed to the constructor.

set_grid_drift(drift: numpy.ndarray, ivar: int | None = None)#

Set external drift values at grid/block locations.

Call after set_grid(), set_grid_block(), or set_grid_cv(). Only needed when ndrift > 0.

Parameters:

drift (ndarray, shape (nblocks, ndrift)) – Drift values. Rows are blocks, columns are drift functions. Note: use nblocks (number of blocks), not ngrid (number of sub-nodes), even for block kriging. Transposed to (ndrift, nblocks) internally before calling Fortran.
ivar (int, optional) – Target-variable index (1-based) whose RHS receives this drift. None (default) broadcasts the same drift to all target variables — the usual case when external drift is independent of which variable is being estimated.

Note

ivar here refers to the target variable (which variable’s estimate uses this drift in its RHS), not the source variable. This is the opposite end from set_obs_drift(), whose ivar identifies the source variable (whose observations form the F-matrix column).

set_sim(randpath: numpy.ndarray | None = None, sample: numpy.ndarray | None = None)#

Set up Sequential Gaussian Simulation parameters.

Call after set_grid() and before set_search(). Only needed when nsim > 0.

Parameters:

randpath (ndarray of int, shape (nblocks,), optional) – 1-based random visiting order for the block loop. When None, Fortran generates a random permutation.
sample (ndarray, shape (nblocks, nvar, nsim), optional) – Pre-drawn samples. When None, Fortran generates them via the object’s set_sim override: N(0,1) for Kriging, U(0,1) for IndicatorKriging.

set_search(ivar: int = 1, anis1: float = 1.0, anis2: float = 1.0, azimuth: float = 0.0, dip: float = 0.0, plunge: float = 0.0, sector_search: bool = False)#

Build the KD-tree and configure the search ellipse for variable ivar. Call once per variable after set_obs() (and set_sim() for SGSIM).

Parameters:

ivar (int) – Variable index (1-based).
anis1 (float) – Horizontal anisotropy ratio (minor / major range). 1.0 = isotropic.
anis2 (float) – Vertical anisotropy ratio (vertical / major range). 1.0 = isotropic.
azimuth (float) – Azimuth of the major axis (degrees, clockwise from North).
dip (float) – Dip angle (degrees, positive downward).
plunge (float) – Plunge angle (degrees).
sector_search (bool) – Enable sector (quadrant in 2D, octant in 3D) search limiting candidates per sector. If True, the search space is divided into quadrants/octants centered on the prediction location. Candidates are sorted by distance, and at most nmax (from set_obs()) are selected from each sector. This ensures a balanced spatial distribution of neighbours and prevents clustering artifacts. The maximum total neighbours selected is 4 * nmax in 2D or 8 * nmax in 3D. If search anisotropy is enabled, coordinates are rotated and scaled according to the anisotropy parameters before sector assignment.

set_grad(coord1: numpy.ndarray, coord2: numpy.ndarray, grad_value: numpy.ndarray, ivar: int = 1, variance: numpy.ndarray | None = None, drift_ext: numpy.ndarray | None = None)#

Set gradient observation pairs (Delhomme 1979: “Kriging in hydrology”).

Each pair (coord1[i], coord2[i]) approximates the directional gradient at a boundary as a finite difference. The constraint Z(xs1) - Z(xs2) = grad_value[i] is enforced as a hard equality. For a no-flow (zero normal gradient) boundary use grad_value = 0.

Call set_grad() after set_search() and before solve().

Parameters:

coord1 (ndarray, shape (ngrad, ndim)) – Positive-side virtual node coordinates.
coord2 (ndarray, shape (ngrad, ndim)) – Negative-side virtual node coordinates.
grad_value (ndarray, shape (ngrad,)) – Known gradient values. Use 0 for no-flow boundaries.
ivar (int, default 1) – Variable index (1-based) the gradient pairs constrain. For cokriging, specifies which variable’s gradient is observed.
variance (ndarray, shape (ngrad,), optional) – Gradient observation variance (default 0 = exact constraint). A non-zero value relaxes the constraint, analogous to obs nugget.
drift_ext (ndarray, shape (ngrad, ndrift), optional) – External drift differences f_ext(xs1) - f_ext(xs2) for each pair. Required when ndrift > 0; omit for ordinary kriging.

solve(nthread: int = 0, ncache: int | None = None)#

Run the kriging or SGSIM loop over all blocks. Calls prepare(), then the parallel block loop internally.

Parameters:

nthread (int, optional) – Number of OpenMP threads to use for this call. 0 (default) leaves the OMP runtime setting unchanged. 1 forces single-threaded execution (useful for reproducible results or when calling solve() from inside another parallel region).
ncache (int, optional) – Number of per-thread multi-slot hcache entries to use for this solve call. None keeps the compiled/object default, 0 disables the hcache, and 1 gives a one-slot hcache for cache-overhead comparisons. The single-entry ctx%cache and optional persistent factor cache are unaffected.

property solver_stats: dict#

Solver statistics from the most recent solve() call.

Returns a dict with five integer counts, reset to zero at the start of every solve():

fail: Blocks where both Cholesky and SSYTRF failed; solution set to NaN (only possible when neglect_error is True).
chol_fact: Fresh Cholesky factorizations performed — O(n³) each, one per unique neighbourhood on a cache miss.
chol_reuse: Blocks solved via a cached Cholesky factorization — O(n²) each. A large value relative to chol_fact means the neighbourhood cache is working effectively.
ssytrf_fact: SSYTRF (Bunch-Kaufman LDL^T) factorizations performed — O(n³) each, once per unique neighbourhood. Non-zero means Cholesky failed for at least one neighbourhood (e.g. a non-SPD system).
ssytrf_reuse: Blocks solved by a cached SSYTRF via SSYTRS — O(n²) each.

get_results(copy: bool = False, squeeze: bool = True)#

Retrieve the kriging estimates and variances after solve().

Fortran fills estimate(nsim, nblocks) directly into a Fortran-contiguous Python-owned buffer.

Parameters:

copy (bool, default False) – If True, return C-contiguous copies for downstream NumPy/Pandas use. If False, return views / Fortran-order arrays when possible.
squeeze (bool, default True) – If True, return a 1-D estimate when nsim == 1.

Returns:

estimate (ndarray) – (ngrid, nvar, nsim). Shape (ngrid,) when (nsim == 1 or ==1) and squeeze; otherwise shape (ngrid, nvar, nsim).
variance (ndarray, shape (nblocks,)) – Kriging variance at each block.

Example

>>> est, var = k.get_results()
>>> kriging_estimate = est[0]          # shape (nblocks,)
>>> sim_realisation1 = est[0]          # same for nsim=1

get_result_array() → numpy.ndarray#

Return all results as a NumPy structured (record) array — one row per block.

The array contains block centroid coordinates alongside every estimate, simulation realization, and variance produced by the last solve().

Fields#

Coordinates (always present):: x, y [, z] — block centroid coordinates.
Estimates / simulations:: Kriging (nsim == 0), nvar == 1: estimate Kriging (nsim == 0), nvar > 1: est_v1, est_v2, … SGSIM (nsim > 0), nvar == 1: sim_1, sim_2, …, sim_{nsim} SGSIM (nsim > 0), nvar > 1: v1_s1, v1_s2, …, v{nvar}_s{nsim}
Variances (always present, diagonal of conditional covariance):: nvar == 1: variance nvar > 1: var_v1, var_v2, …

returns:: np.ndarray with named fields (structured array / record array) – Shape (nblocks,). Access a column with arr['estimate'], convert to a plain 2-D array with np.column_stack([arr[f] for f in arr.dtype.names]).

Example

>>> k.solve()
>>> ra = k.get_result_array()
>>> ra.dtype.names
('x', 'y', 'estimate', 'variance')
>>> ra['estimate']          # 1-D array, shape (nblocks,)
>>> import pandas as pd
>>> df = pd.DataFrame(ra)   # direct conversion to DataFrame

get_result_df() → pandas.DataFrame#

Return all results as a pandas.DataFrame.

Wraps get_result_array(); column names and field descriptions are identical to those documented there.

Returns:: pandas.DataFrame – One row per block; columns match the fields of the structured array returned by get_result_array().

Example

>>> k.solve()
>>> df = k.get_result_df()
>>> df.columns.tolist()
['x', 'y', 'estimate', 'variance']

get_factor() → dict#

Return the persistent LHS factorization cached after the last solve().

The Cholesky factorization of the kriging covariance matrix K and the related Schur-complement matrices are computed once per solve() call (or reused across blocks with the same neighbour set). Starting from the second solve() call on unchanged observations and variogram, the cached factors allow the Fortran engine to skip kriging_setup entirely.

This method exposes those matrices and the assembled linear system so that users can inspect or verify the factorization inputs.

Returns:

dict with keys
``valid`` (bool) – True if a persistent factor exists (i.e. at least one solve() has been completed and observations/variogram have not changed since).
``npp`` (int) – Number of neighbours in the LHS matrix (size of K).
``p`` (int) – Number of drift + unbiasedness columns (size of the Schur complement).
``L`` (ndarray, shape (npp, npp)) – Lower-triangular Cholesky factor of K (stored column-major by Fortran, returned as a C-contiguous array).
``kinv_drift`` (ndarray, shape (npp, max(1, p))) – K⁻¹ F (K inverse applied to the drift matrix F).
``schur`` (ndarray, shape (max(1, p), max(1, p))) – Cholesky factor of the Schur complement F’ K⁻¹ F.
``matA`` (ndarray, shape (npp + p, npp + p)) – Assembled linear-system LHS before factorization.
``rhsB`` (ndarray, shape (nvar, npp + p)) – Assembled linear-system RHS before solving.

Example

>>> k = Kriging(ndim=2, nvar=1, ndrift=1, unbias=0)
>>> # ... set_obs, set_vgm, set_grid, set_obs_drift, set_grid_drift ...
>>> k.solve()
>>> f = k.get_factor()
>>> if f['valid']:
...     L = f['L']            # Cholesky factor of covariance matrix
...     kinv = f['kinv_drift']  # K^{-1} F
...     schur = f['schur']    # Cholesky of Schur complement

Notes

The factor is invalidated (valid = False) whenever set_obs() or set_vgm() is called. Call solve() again to repopulate.

get_estimate_all(copy: bool = False)#

Return multivariable estimates / simulations for all variables.

Populated when nvar > 1. For co-kriging without simulation, the leading dimension is 1.

Parameters:: copy (bool, default False) – If True, return a C-contiguous copy. If False, return the Fortran-contiguous output buffer filled by the Fortran core.
Returns:: np.ndarray, shape (nblock, nvar, max(nsim, 1)) – Values of all variables. out[ib, kvar, isim] is the value at block ib+1 for variable kvar+1 in realization isim+1.

get_variance_all(copy: bool = False)#

Return the conditional covariance matrix for all variables.

Returns:: np.ndarray, shape (nblock, nvar, nvar) – Conditional covariance matrix at each block. The diagonal contains each variable’s kriging variance, and out[:, 0, 0] matches the variance returned by get_results().

free_weight_store()#: Release the in-memory weight store, freeing its memory.

set_weights(weights: dict) → None#

Load kriging weights into the in-memory store so solve() reuses them.

Activated when use_old_weight=True and no weight_file is given. solve() then applies the supplied neighbour indices and kriging weights directly — skipping the kriging-system solve — and restores the stored variance. This is the in-memory equivalent of the use_old_weight=True + factor-file workflow.

Typical workflow#

>>> # First run: solve and capture weights + variance
>>> k1 = Kriging(ndim=2, nvar=1, store_weight=True)
>>> k1.set_obs(...); k1.set_grid(...); k1.set_vgm(...); k1.set_search(...)
>>> k1.solve()
>>> w = k1.get_weights()           # {'nnear', 'inear', 'weight', 'variance'}
>>>
>>> # Second run: same grid/vgm, new obs values, reuse weights
>>> k2 = Kriging(ndim=2, nvar=1, use_old_weight=True)   # no weight_file
>>> k2.set_obs(...new_values...)
>>> k2.set_grid(...); k2.set_vgm(...); k2.set_search(...)
>>> k2.set_weights(w)              # populate the in-memory store
>>> k2.solve()                     # fast: skips kriging system solve
>>> est, var = k2.get_results()

param weights:

Dict as returned by get_weights(), with keys:

nnearndarray (nblock, ngroups), int32: Number of active neighbours per block and group.
inearndarray (nblock, ngroups, nmax), int32: 1-based neighbour indices.
weightndarray (nblock, nvar, ngroups, nmax), float64: Kriging weights. For nvar==1 the array may be 3-D (nblock, ngroups, nmax) (the shape returned by get_weights()).
orderndarray (nblock,), optional: Random visiting order for SGSIM. Defaults to None.
variancendarray (nblock,) or (nblock, nvar, nvar), optional: Per-block conditional variance. Included automatically when the dict was produced by get_weights(). Defaults to zeros if absent.

type weights:

dict

Notes

Call after set_obs(), set_grid(), set_vgm(), and set_search(). The Fortran-side use_old_weight flag is set automatically; you may also pass use_old_weight=True to the constructor to declare intent explicitly.

get_weights() → dict#

Return the stored kriging weights and neighbour indices.

alloc_weight_store() must have been called before solve().

Returns:

dict with keys
``nnear`` (ndarray, shape (nblock, ngroups), dtype int32) – Number of active neighbours for each block and group. ngroups = ngroups_base when set_grad has not been called, or ngroups_base + nvar when gradient data is present. ngroups_base = nvar (kriging) or 2*nvar (SGSIM). Group layout:

indices 0..nvar-1 — obs groups (variable 1..nvar) indices nvar..2*nvar-1 — sim groups (SGSIM only) indices ngroups_base..ngroups-1 — grad groups (present only when set_grad called)
``inear`` (ndarray, shape (nblock, ngroups, nmax), dtype int32) – 1-based neighbour indices. Entries beyond nnear[ib, ig] are zero.
``weight`` (ndarray, shape (nblock, nvar, ngroups, nmax), dtype float64) – Kriging weights. Entries beyond nnear[ib, ig] are zero. Shape is (nblock, ngroups, nmax) when nvar == 1.
``variance`` (ndarray, shape (nblock,) for nvar==1, else (nblock, nvar, nvar)) – Per-block conditional kriging variance stored alongside the weights. Present only when the compiled library supports the variance store (i.e. built with the current source). Pass this dict directly to set_weights() to get a full round-trip.

get_info()#

krigekit.kriging.ordinary_kriging(obs_coord: numpy.ndarray, obs_value: numpy.ndarray, grid_coord: numpy.ndarray, vgm_spec: dict | list[dict], nmax: int | None = None, maxdist: float | None = None, search_anis1: float = 1.0, search_anis2: float = 1.0, search_azimuth: float = 0.0, rangescale: float | None = None, localnugget: float | None = None, nthread=0, ncache: int | None = None) → tuple[numpy.ndarray, numpy.ndarray]#

One-shot ordinary kriging with a single isotropic (or anisotropic) variogram.

Parameters:

obs_coord (ndarray, shape (nobs, ndim)) – Observation coordinates. Rows are points, columns are spatial dimensions.
obs_value (ndarray, shape (nobs,)) – Observation values.
grid_coord (ndarray, shape (ngrid, ndim)) – Grid coordinates to estimate.
vgm_spec (dict or list of dict) – One variogram structure dict, or a list of dicts for nested models. Each dict is passed as keyword arguments to Kriging.set_vgm() (keys: vtype, nugget, sill, a_major, and optionally a_minor1, a_minor2, azimuth, dip, plunge).
nmax (int) – Maximum number of neighbours.
maxdist (float, optional) – Maximum search distance.
search_anis1 (float) – Anisotropy ratios for search ellipse (1.0 = isotropic).
search_anis2 (float) – Anisotropy ratios for search ellipse (1.0 = isotropic).
search_azimuth (float) – Azimuth of search ellipse major axis (degrees from North).
nthread (int) – max OMP threads for this call (0 or absent = OMP default)
ncache (int, optional) – Per-thread hcache slots for this solve. None uses the default.

Returns:

estimate (ndarray, shape (ngrid,))
variance (ndarray, shape (ngrid,))

Example

>>> est, var = ordinary_kriging(
...     obs_coord, obs_value, grid_coord,
...     vgm_spec=dict(vtype="sph", nugget=100, sill=900, a_major=1000, a_minor1=500),
...     nmax=20)

krigekit.kriging.cokriging(obs_coords: list[numpy.ndarray], obs_values: list[numpy.ndarray], grid_coord: numpy.ndarray, vgm_spec: dict, nmax: int | None = None, rangescale: float | None = None, localnugget: float | None = None, nthread: int = 0, ncache: int | None = None, std_ck: bool = False) → tuple[numpy.ndarray, numpy.ndarray]#

One-shot ordinary co-kriging with multiple variables.

Parameters:

obs_coords (list of ndarray, each shape (nobs_i, ndim)) – Observation coordinates per variable. Rows are points.
obs_values (list of ndarray, each shape (nobs_i,)) – Observation values per variable.
grid_coord (ndarray, shape (ngrid, ndim)) – Grid coordinates.
vgm_spec (dict) – Mapping (ivar, jvar) to a variogram dict or list of dicts. Each dict is passed as keyword arguments to Kriging.set_vgm(). Both (i,j) and (j,i) can be provided; if only (i,j) is given, (j,i) will mirror it automatically (handled inside Fortran set_vgm).
nmax (int) – Maximum neighbours per variable.
nthread (int) – max OMP threads for this call (0 or absent = OMP default)
ncache (int, optional) – Per-thread hcache slots for this solve. None uses the default.
std_ck (bool) – Use standard Ordinary Kriging.

Returns:

estimate (ndarray, shape (ngrid,))
variance (ndarray, shape (ngrid,))

Example

>>> est, var = cokriging(
...     obs_coords=[coord1, coord2],
...     obs_values=[val1, val2],
...     grid_coord=grid,
...     vgm_spec={
...         (1,1): dict(vtype="sph", nugget=100, sill=900, a_major=1000, a_minor1=500),
...         (2,2): dict(vtype="sph", nugget=50,  sill=450, a_major=1000, a_minor1=500),
...         (1,2): dict(vtype="sph", nugget=0,   sill=600, a_major=1000, a_minor1=500),
...     })

krigekit.kriging.sequential_gaussian_simulation(obs_coord: numpy.ndarray, obs_value: numpy.ndarray, grid_coord: numpy.ndarray, vgm_spec: str, nsim: int, nmax: int | None = None, randpath: numpy.ndarray | None = None, sample: numpy.ndarray | None = None, seed: int | None = None, rangescale: float | None = None, localnugget: float | None = None, nthread: int = 0, ncache: int | None = None) → numpy.ndarray#

Sequential Gaussian Simulation.

Parameters:

obs_coord (ndarray, shape (nobs, ndim)) – Observation coordinates. Rows are points, columns are spatial dimensions.
obs_value (ndarray, shape (nobs,)) – Observation values.
grid_coord (ndarray, shape (ngrid, ndim)) – Grid coordinates.
vgm_spec (dict or list of dict) – One or more nested variogram structure dicts, each passed as keyword arguments to Kriging.set_vgm().
nsim (int) – Number of realisations.
nmax (int) – Maximum neighbours (includes previously simulated nodes).
seed (int, optional) – Random seed for reproducibility.
nthread (int) – max OMP threads for this call (0 or absent = OMP default)
ncache (int, optional) – Per-thread hcache slots for this solve. None uses the default.

Returns:

simulations (ndarray, shape (nsim, ngrid)) – Each row is one realisation in the original (non-randomised) block order.