gstools_core/
lib.rs

1//! GSTools-Core
2//!
3//! `gstools_core` is a Rust implementation of the core algorithms of [GSTools].
4//! At the moment, it is a drop in replacement for the Cython code included in GSTools.
5//!
6//! This crate includes
7//! - [randomization methods](field) for the random field generation
8//! - the [matrix operations](krige) of the kriging methods
9//! - the [variogram estimation](variogram)
10//!
11//! [GSTools]: https://github.com/GeoStat-Framework/GSTools
12
13#![warn(missing_docs)]
14
15use numpy::{IntoPyArray, PyArray1, PyArray2, PyReadonlyArray1, PyReadonlyArray2};
16use pyo3::prelude::{pymodule, Bound, PyModule, PyResult, Python};
17
18use field::{summator, summator_fourier, summator_incompr};
19use krige::{calculator_field_krige, calculator_field_krige_and_variance};
20use variogram::{
21    variogram_directional, variogram_ma_structured, variogram_structured, variogram_unstructured,
22};
23
24pub mod field;
25pub mod krige;
26mod short_vec;
27pub mod variogram;
28
29#[pymodule]
30fn gstools_core(_py: Python<'_>, m: &Bound<'_, PyModule>) -> PyResult<()> {
31    m.add("__version__", env!("CARGO_PKG_VERSION"))?;
32
33    #[pyfn(m)]
34    #[pyo3(name = "summate")]
35    fn summate_py<'py>(
36        py: Python<'py>,
37        cov_samples: PyReadonlyArray2<f64>,
38        z1: PyReadonlyArray1<f64>,
39        z2: PyReadonlyArray1<f64>,
40        pos: PyReadonlyArray2<f64>,
41        num_threads: Option<usize>,
42    ) -> Bound<'py, PyArray1<f64>> {
43        let cov_samples = cov_samples.as_array();
44        let z1 = z1.as_array();
45        let z2 = z2.as_array();
46        let pos = pos.as_array();
47        summator(cov_samples, z1, z2, pos, num_threads).into_pyarray_bound(py)
48    }
49
50    #[pyfn(m)]
51    #[pyo3(name = "summate_incompr")]
52    fn summate_incompr_py<'py>(
53        py: Python<'py>,
54        cov_samples: PyReadonlyArray2<f64>,
55        z1: PyReadonlyArray1<f64>,
56        z2: PyReadonlyArray1<f64>,
57        pos: PyReadonlyArray2<f64>,
58        num_threads: Option<usize>,
59    ) -> Bound<'py, PyArray2<f64>> {
60        let cov_samples = cov_samples.as_array();
61        let z1 = z1.as_array();
62        let z2 = z2.as_array();
63        let pos = pos.as_array();
64        summator_incompr(cov_samples, z1, z2, pos, num_threads).into_pyarray_bound(py)
65    }
66
67    #[pyfn(m)]
68    #[pyo3(name = "summate_fourier")]
69    fn summate_fourier_py<'py>(
70        py: Python<'py>,
71        spectrum_factor: PyReadonlyArray1<f64>,
72        modes: PyReadonlyArray2<f64>,
73        z1: PyReadonlyArray1<f64>,
74        z2: PyReadonlyArray1<f64>,
75        pos: PyReadonlyArray2<f64>,
76        num_threads: Option<usize>,
77    ) -> Bound<'py, PyArray1<f64>> {
78        let spectrum_factor = spectrum_factor.as_array();
79        let modes = modes.as_array();
80        let z1 = z1.as_array();
81        let z2 = z2.as_array();
82        let pos = pos.as_array();
83        summator_fourier(spectrum_factor, modes, z1, z2, pos, num_threads).into_pyarray_bound(py)
84    }
85
86    #[pyfn(m)]
87    #[pyo3(name = "calc_field_krige_and_variance")]
88    fn calc_field_krige_and_variance_py<'py>(
89        py: Python<'py>,
90        krige_mat: PyReadonlyArray2<f64>,
91        krig_vecs: PyReadonlyArray2<f64>,
92        cond: PyReadonlyArray1<f64>,
93        num_threads: Option<usize>,
94    ) -> (Bound<'py, PyArray1<f64>>, Bound<'py, PyArray1<f64>>) {
95        let krige_mat = krige_mat.as_array();
96        let krig_vecs = krig_vecs.as_array();
97        let cond = cond.as_array();
98        let (field, error) =
99            calculator_field_krige_and_variance(krige_mat, krig_vecs, cond, num_threads);
100        let field = field.into_pyarray_bound(py);
101        let error = error.into_pyarray_bound(py);
102        (field, error)
103    }
104
105    #[pyfn(m)]
106    #[pyo3(name = "calc_field_krige")]
107    fn calc_field_krige_py<'py>(
108        py: Python<'py>,
109        krige_mat: PyReadonlyArray2<f64>,
110        krig_vecs: PyReadonlyArray2<f64>,
111        cond: PyReadonlyArray1<f64>,
112        num_threads: Option<usize>,
113    ) -> Bound<'py, PyArray1<f64>> {
114        let krige_mat = krige_mat.as_array();
115        let krig_vecs = krig_vecs.as_array();
116        let cond = cond.as_array();
117        calculator_field_krige(krige_mat, krig_vecs, cond, num_threads).into_pyarray_bound(py)
118    }
119
120    #[pyfn(m)]
121    #[pyo3(name = "variogram_structured")]
122    fn variogram_structured_py<'py>(
123        py: Python<'py>,
124        f: PyReadonlyArray2<f64>,
125        estimator_type: Option<char>,
126        num_threads: Option<usize>,
127    ) -> Bound<'py, PyArray1<f64>> {
128        let f = f.as_array();
129        let estimator_type = estimator_type.unwrap_or('m');
130        variogram_structured(f, estimator_type, num_threads).into_pyarray_bound(py)
131    }
132
133    #[pyfn(m)]
134    #[pyo3(name = "variogram_ma_structured")]
135    fn variogram_ma_structured_py<'py>(
136        py: Python<'py>,
137        f: PyReadonlyArray2<f64>,
138        mask: PyReadonlyArray2<bool>,
139        estimator_type: Option<char>,
140        num_threads: Option<usize>,
141    ) -> Bound<'py, PyArray1<f64>> {
142        let f = f.as_array();
143        let mask = mask.as_array();
144        let estimator_type = estimator_type.unwrap_or('m');
145        variogram_ma_structured(f, mask, estimator_type, num_threads).into_pyarray_bound(py)
146    }
147
148    #[pyfn(m)]
149    #[pyo3(name = "variogram_directional")]
150    #[allow(clippy::too_many_arguments)]
151    fn variogram_directional_py<'py>(
152        py: Python<'py>,
153        f: PyReadonlyArray2<f64>,
154        bin_edges: PyReadonlyArray1<f64>,
155        pos: PyReadonlyArray2<f64>,
156        direction: PyReadonlyArray2<f64>, //should be normed
157        angles_tol: Option<f64>,
158        bandwidth: Option<f64>,
159        separate_dirs: Option<bool>,
160        estimator_type: Option<char>,
161        num_threads: Option<usize>,
162    ) -> (Bound<'py, PyArray2<f64>>, Bound<'py, PyArray2<u64>>) {
163        let f = f.as_array();
164        let bin_edges = bin_edges.as_array();
165        let pos = pos.as_array();
166        let direction = direction.as_array();
167        let angles_tol = angles_tol.unwrap_or(std::f64::consts::PI / 8.0);
168        let bandwidth = bandwidth.unwrap_or(-1.0);
169        let separate_dirs = separate_dirs.unwrap_or(false);
170        let estimator_type = estimator_type.unwrap_or('m');
171        let (variogram, counts) = variogram_directional(
172            f,
173            bin_edges,
174            pos,
175            direction,
176            angles_tol,
177            bandwidth,
178            separate_dirs,
179            estimator_type,
180            num_threads,
181        );
182        let variogram = variogram.into_pyarray_bound(py);
183        let counts = counts.into_pyarray_bound(py);
184
185        (variogram, counts)
186    }
187
188    #[pyfn(m)]
189    #[pyo3(name = "variogram_unstructured")]
190    fn variogram_unstructured_py<'py>(
191        py: Python<'py>,
192        f: PyReadonlyArray2<f64>,
193        bin_edges: PyReadonlyArray1<f64>,
194        pos: PyReadonlyArray2<f64>,
195        estimator_type: Option<char>,
196        distance_type: Option<char>,
197        num_threads: Option<usize>,
198    ) -> (Bound<'py, PyArray1<f64>>, Bound<'py, PyArray1<u64>>) {
199        let f = f.as_array();
200        let bin_edges = bin_edges.as_array();
201        let pos = pos.as_array();
202        let estimator_type = estimator_type.unwrap_or('m');
203        let distance_type = distance_type.unwrap_or('e');
204        let (variogram, counts) = variogram_unstructured(
205            f,
206            bin_edges,
207            pos,
208            estimator_type,
209            distance_type,
210            num_threads,
211        );
212        let variogram = variogram.into_pyarray_bound(py);
213        let counts = counts.into_pyarray_bound(py);
214
215        (variogram, counts)
216    }
217
218    Ok(())
219}