From 8d23e30294a37002744f366788765f416064c808 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Wed, 15 Jan 2020 18:55:56 +0100
Subject: [PATCH 01/49] Add new data encapsulation and refactor

2 new classes are added in order to separate data and methods more
cleanly. This is very much WIP, but I'd like to get a discussion
started before I put in too much work.
The `FieldData` class contains information about the mesh and it
holds a dict containing the individual field values, which are
defined on that mesh, e.g. "krige" and "cond".
---
 gstools/field/base.py      |  23 ++-----
 gstools/field/condition.py |  16 ++---
 gstools/field/data.py      | 119 +++++++++++++++++++++++++++++++++++++
 gstools/field/srf.py       |  37 ++++++------
 gstools/krige/ordinary.py  |   6 +-
 gstools/krige/simple.py    |   9 +--
 gstools/tools/export.py    |   4 +-
 gstools/transform/field.py |  80 +++++++++++++++----------
 tests/test_krige.py        |  40 ++++++-------
 tests/test_srf.py          |  16 +++--
 10 files changed, 237 insertions(+), 113 deletions(-)
 create mode 100644 gstools/field/data.py

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 9e2ff36e..ec7be97c 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -15,6 +15,7 @@
 
 import numpy as np
 
+from gstools.field.data import FieldData
 from gstools.covmodel.base import CovModel
 from gstools.tools.export import to_vtk, vtk_export
 from gstools.field.tools import _get_select
@@ -22,24 +23,19 @@
 __all__ = ["Field"]
 
 
-class Field:
+class Field(FieldData):
     """A field base class for random and kriging fields ect.
 
     Parameters
     ----------
     model : :any:`CovModel`
         Covariance Model related to the field.
-    mean : :class:`float`, optional
-        Mean value of the field.
     """
 
-    def __init__(self, model, mean=0.0):
+    def __init__(self, model):
         # initialize attributes
-        self.pos = None
-        self.mesh_type = None
-        self.field = None
+        super().__init__(self)
         # initialize private attributes
-        self._mean = mean
         self._model = None
         self.model = model
         self._value_type = None
@@ -299,15 +295,6 @@ def plot(self, field="field", fig=None, ax=None):  # pragma: no cover
 
         return r
 
-    @property
-    def mean(self):
-        """:class:`float`: The mean of the field."""
-        return self._mean
-
-    @mean.setter
-    def mean(self, mean):
-        self._mean = float(mean)
-
     @property
     def model(self):
         """:any:`CovModel`: The covariance model of the field."""
@@ -333,7 +320,7 @@ def __str__(self):
 
     def __repr__(self):
         """Return String representation."""
-        return "Field(model={0}, mean={1})".format(self.model, self.mean)
+        return "Field(model={0})".format(self.model)
 
 
 if __name__ == "__main__":  # pragma: no cover
diff --git a/gstools/field/condition.py b/gstools/field/condition.py
index 176be302..c2e3b09b 100644
--- a/gstools/field/condition.py
+++ b/gstools/field/condition.py
@@ -38,7 +38,7 @@ def ordinary(srf):
     krige_ok = Ordinary(
         model=srf.model, cond_pos=srf.cond_pos, cond_val=srf.cond_val
     )
-    krige_field, krige_var = krige_ok(srf.pos, srf.mesh_type)
+    krige_field = krige_ok(srf.pos, srf.mesh_type)
 
     # evaluate the field at the conditional points
     x, y, z = pos2xyz(srf.cond_pos, max_dim=srf.model.dim)
@@ -50,10 +50,10 @@ def ordinary(srf):
     err_ok = Ordinary(
         model=srf.model, cond_pos=srf.cond_pos, cond_val=err_data
     )
-    err_field, __ = err_ok(srf.pos, srf.mesh_type)
-    cond_field = srf.raw_field + krige_field - err_field
+    err_field = err_ok(srf.pos, srf.mesh_type)
+    cond_field = srf["raw"] + krige_field.values - err_field.values
     info = {"mean": krige_ok.mean}
-    return cond_field, krige_field, err_field, krige_var, info
+    return cond_field, krige_field, err_field, krige_field.krige_var, info
 
 
 def simple(srf):
@@ -81,7 +81,7 @@ def simple(srf):
         cond_pos=srf.cond_pos,
         cond_val=srf.cond_val,
     )
-    krige_field, krige_var = krige_sk(srf.pos, srf.mesh_type)
+    krige_field = krige_sk(srf.pos, srf.mesh_type)
 
     # evaluate the field at the conditional points
     x, y, z = pos2xyz(srf.cond_pos, max_dim=srf.model.dim)
@@ -96,7 +96,7 @@ def simple(srf):
         cond_pos=srf.cond_pos,
         cond_val=err_data,
     )
-    err_field, __ = err_ok(srf.pos, srf.mesh_type)
-    cond_field = srf.raw_field + krige_field - err_field + srf.mean
+    err_field = err_ok(srf.pos, srf.mesh_type)
+    cond_field = srf["raw"] + krige_field.values - err_field.values + srf.mean
     info = {}
-    return cond_field, krige_field, err_field, krige_var, info
+    return cond_field, krige_field, err_field, krige_field.krige_var, info
diff --git a/gstools/field/data.py b/gstools/field/data.py
new file mode 100644
index 00000000..4bb9b6c0
--- /dev/null
+++ b/gstools/field/data.py
@@ -0,0 +1,119 @@
+# -*- coding: utf-8 -*-
+"""
+GStools subpackage providing a data class for spatial fields.
+
+.. currentmodule:: gstools.field.data
+
+The following classes are provided
+
+.. autosummary::
+   FieldData
+"""
+# pylint: disable=C0103
+
+__all__ = ["FieldData"]
+
+
+from typing import List, Dict
+import numpy as np
+
+
+class Data:
+    """A data class mainly storing the specific values of an individual field.
+    """
+
+    def __init__(self, values: np.ndarray, mean: float = 0.0, value_type: str = "scalar"):
+        self.values = values
+        self.mean = mean
+        self.value_type = value_type
+
+
+class FieldData:
+    """A base class encapsulating field data.
+
+    It holds a position array, which define the spatial locations of the
+    field values.
+    It can hold multiple fields in the :any:`self.values` list. This assumes
+    that each field is defined on the same positions.
+    The mesh type must also be specified.
+
+    Parameters
+    ----------
+    pos : :class:`numpy.ndarray`, optional
+        positions of the field values
+    values : :any:`list`, optional
+        a list of the values of the fields
+
+    Examples
+    --------
+    >>> import numpy as np
+    >>> pos = np.random.random((100, 100))
+    >>> z = np.random.random((100, 100))
+    >>> z2 = np.random.random((100, 100))
+    >>> field_data = FieldData(pos)
+    >>> field_data.add_field("test_field1", z)
+    >>> field_data.add_field("test_field2", z2)
+    >>> field_data.set_default_field("test_field2")
+    >>> print(field.values)
+
+    """
+
+    def __init__(
+        self, pos: np.ndarray = None, mesh_type: str = "unstructured",
+    ):
+        # initialize attributes
+        self.pos = pos
+        self.fields: Dict[str, np.ndarray] = {}
+        self.mesh_type = mesh_type
+
+    def add_field(
+        self,
+        name: str,
+        values: np.ndarray,
+        mean: float = 0.0,
+        *,
+        value_type: str = "scalar",
+        default_field: bool = False,
+    ):
+        self.fields[name] = Data(values, mean, value_type)
+        # set the default_field to the first field added
+        if len(self.fields) == 1 or default_field:
+            self.default_field = name
+
+    def get_data(self, key):
+        """:class:`Data`: The field data class."""
+        return self.fields[key]
+
+    def set_default_field(self, default_field):
+        self.default_field = default_field
+
+    def __setitem__(self, key, value):
+        self.fields[key].values = value
+
+    def __getitem__(self, key):
+        """:any:`numpy.ndarray`: The values of the field."""
+        return self.fields[key].values
+
+    @property
+    def field(self):
+        return self.fields[self.default_field]
+
+    @property
+    def values(self):
+        return self.fields[self.default_field].values
+
+    @values.setter
+    def values(self, values):
+        self.fields[self.default_field].values = values
+
+    @property
+    def value_type(self):
+        return self.fields[self.default_field].value_type
+
+    @property
+    def mean(self):
+        return self.fields[self.default_field].mean
+
+    @mean.setter
+    def mean(self, value):
+        self.fields[self.default_field].mean = mean
diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index bbd7e320..c1851b5e 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -89,8 +89,9 @@ def __init__(
         generator="RandMeth",
         **generator_kwargs
     ):
-        super().__init__(model, mean)
+        super().__init__(model)
         # initialize private attributes
+        self._mean = mean
         self._generator = None
         self._upscaling = None
         self._upscaling_func = None
@@ -99,10 +100,6 @@ def __init__(
         self._cond_val = None
         self._krige_type = None
         # initialize attributes
-        self.raw_field = None
-        self.krige_field = None
-        self.err_field = None
-        self.krige_var = None
         self.set_generator(generator, **generator_kwargs)
         self.upscaling = upscaling
         if self._value_type is None:
@@ -160,13 +157,13 @@ def __call__(
         y, z = make_isotropic(self.model.dim, self.model.anis, y, z)
 
         # generate the field
-        self.raw_field = self.generator.__call__(x, y, z, mesh_type)
+        self.add_field("raw", self.generator.__call__(x, y, z, mesh_type), self._mean)
 
         # reshape field if we got an unstructured mesh
         if mesh_type_changed:
             mesh_type = mesh_type_old
-            self.raw_field = reshape_field_from_unstruct_to_struct(
-                self.model.dim, self.raw_field, axis_lens
+            self["raw"] = reshape_field_from_unstruct_to_struct(
+                self.model.dim, self["raw"], axis_lens
             )
 
         # apply given conditions to the field
@@ -179,23 +176,23 @@ def __call__(
                 info,
             ) = self.cond_func(self)
             # store everything in the class
-            self.field = cond_field
-            self.krige_field = krige_field
-            self.err_field = err_field
-            self.krige_var = krigevar
-            if "mean" in info:  # ordinary krging estimates mean
-                self.mean = info["mean"]
+            self.add_field("cond", cond_field, default_field=True)
+            self.add_field("krige", krige_field)
+            self.add_field("error", err_field)
+            self.get_data("krige").krige_var = krigevar
+            if "mean" in info:  # ordinary kriging estimates mean
+                self.get_data("cond").mean = info["mean"]
         else:
-            self.field = self.raw_field + self.mean
+            self.add_field("srf", self["raw"] + self.mean)
 
         # upscaled variance
         if not np.isscalar(point_volumes) or not np.isclose(point_volumes, 0):
             scaled_var = self.upscaling_func(self.model, point_volumes)
-            self.field -= self.mean
-            self.field *= np.sqrt(scaled_var / self.model.sill)
-            self.field += self.mean
+            self["srf"] -= self.mean
+            self["srf"] *= np.sqrt(scaled_var / self.model.sill)
+            self["srf"] += self.mean
 
-        return self.field
+        return self.values
 
     def set_condition(
         self, cond_pos=None, cond_val=None, krige_type="ordinary"
@@ -309,7 +306,7 @@ def upscaling(self, upscaling):
     def __repr__(self):
         """Return String representation."""
         return "SRF(model={0}, mean={1}, generator={2}".format(
-            self.model, self.mean, self.generator
+            self.model, self._mean, self.generator
         )
 
 
diff --git a/gstools/krige/ordinary.py b/gstools/krige/ordinary.py
index 5e3ff3d2..bc8d8fa9 100644
--- a/gstools/krige/ordinary.py
+++ b/gstools/krige/ordinary.py
@@ -121,9 +121,9 @@ def __call__(self, pos, mesh_type="unstructured"):
                 self.model.dim, krige_var, axis_lens
             )
         # save the field
-        self.krige_var = krige_var
-        self.field = field
-        return self.field, self.krige_var
+        self.add_field("krige", field, default_field=True)
+        self.field.krige_var = krige_var
+        return self.field
 
     def _get_krig_mat(self, pos1, pos2):
         size = pos1[0].size
diff --git a/gstools/krige/simple.py b/gstools/krige/simple.py
index 2a62eedf..18bf95e1 100644
--- a/gstools/krige/simple.py
+++ b/gstools/krige/simple.py
@@ -117,11 +117,12 @@ def __call__(self, pos, mesh_type="unstructured"):
             krige_var = reshape_field_from_unstruct_to_struct(
                 self.model.dim, krige_var, axis_lens
             )
-        # calculate the kriging error
-        self.krige_var = self.model.sill - krige_var
         # add the given mean
-        self.field = field + self.mean
-        return self.field, self.krige_var
+        self.add_field("krige", field + self.mean, default_field=True)
+
+        # calculate the kriging error
+        self.field.krige_var = self.model.sill - krige_var
+        return self.field
 
     def _get_cov_mat(self, pos1, pos2):
         return self.model.cov_nugget(
diff --git a/gstools/tools/export.py b/gstools/tools/export.py
index 94210b51..77c42913 100644
--- a/gstools/tools/export.py
+++ b/gstools/tools/export.py
@@ -51,7 +51,7 @@ def _vtk_structured_helper(pos, fields):
         z = np.array([0])
     # need fortran order in VTK
     for field in fields:
-        fields[field] = fields[field].reshape(-1, order="F")
+        fields[field] = fields[field].values.reshape(-1, order="F")
         if len(fields[field]) != len(x) * len(y) * len(z):
             raise ValueError(
                 "gstools.vtk_export_structured: "
@@ -119,7 +119,7 @@ def _vtk_unstructured_helper(pos, fields):
     if z is None:
         z = np.zeros_like(x)
     for field in fields:
-        fields[field] = fields[field].reshape(-1)
+        fields[field] = fields[field].values.reshape(-1)
         if (
             len(fields[field]) != len(x)
             or len(fields[field]) != len(y)
diff --git a/gstools/transform/field.py b/gstools/transform/field.py
index c3a90e43..488e1afb 100644
--- a/gstools/transform/field.py
+++ b/gstools/transform/field.py
@@ -49,7 +49,7 @@ def binary(fld, divide=None, upper=None, lower=None):
         Field will be transformed inplace.
     divide : :class:`float`, optional
         The dividing value.
-        Default: ``fld.mean``
+        Default: ``fld.field.mean``
     upper : :class:`float`, optional
         The resulting upper value of the field.
         Default: ``mean + sqrt(fld.model.sill)``
@@ -57,14 +57,14 @@ def binary(fld, divide=None, upper=None, lower=None):
         The resulting lower value of the field.
         Default: ``mean - sqrt(fld.model.sill)``
     """
-    if fld.field is None:
+    if fld.values is None:
         print("binary: no field stored in SRF class.")
     else:
-        divide = fld.mean if divide is None else divide
-        upper = fld.mean + np.sqrt(fld.model.sill) if upper is None else upper
-        lower = fld.mean - np.sqrt(fld.model.sill) if lower is None else lower
-        fld.field[fld.field > divide] = upper
-        fld.field[fld.field <= divide] = lower
+        divide = fld.field.mean if divide is None else divide
+        upper = fld.field.mean + np.sqrt(fld.model.sill) if upper is None else upper
+        lower = fld.field.mean - np.sqrt(fld.model.sill) if lower is None else lower
+        fld.values[fld.values > divide] = upper
+        fld.values[fld.values <= divide] = lower
 
 
 def boxcox(fld, lmbda=1, shift=0):
@@ -90,16 +90,16 @@ def boxcox(fld, lmbda=1, shift=0):
         The field will be shifted by that value before transformation.
         Default: ``0``
     """
-    if fld.field is None:
+    if fld.values is None:
         print("Box-Cox: no field stored in SRF class.")
     else:
-        fld.mean += shift
-        fld.field += shift
+        fld.field.mean += shift
+        fld.field.values += shift
         if np.isclose(lmbda, 0):
             normal_to_lognormal(fld)
-        if np.min(fld.field) < -1 / lmbda:
+        if np.min(fld.values) < -1 / lmbda:
             warn("Box-Cox: Some values will be cut off!")
-        fld.field = (np.maximum(lmbda * fld.field + 1, 0)) ** (1 / lmbda)
+        fld.values = (np.maximum(lmbda * fld.values + 1, 0)) ** (1 / lmbda)
 
 
 def zinnharvey(fld, conn="high"):
@@ -117,10 +117,12 @@ def zinnharvey(fld, conn="high"):
         Desired connectivity. Either "low" or "high".
         Default: "high"
     """
-    if fld.field is None:
+    if fld.values is None:
         print("zinnharvey: no field stored in SRF class.")
     else:
-        fld.field = _zinnharvey(fld.field, conn, fld.mean, fld.model.sill)
+        fld.values = _zinnharvey(
+            fld.values, conn, fld.field.mean, fld.model.sill
+        )
 
 
 def normal_force_moments(fld):
@@ -135,10 +137,12 @@ def normal_force_moments(fld):
         Spatial Random Field class containing a generated field.
         Field will be transformed inplace.
     """
-    if fld.field is None:
+    if fld.values is None:
         print("normal_force_moments: no field stored in SRF class.")
     else:
-        fld.field = _normal_force_moments(fld.field, fld.mean, fld.model.sill)
+        fld.values = _normal_force_moments(
+            fld.values, fld.field.mean, fld.model.sill
+        )
 
 
 def normal_to_lognormal(fld):
@@ -153,10 +157,10 @@ def normal_to_lognormal(fld):
         Spatial Random Field class containing a generated field.
         Field will be transformed inplace.
     """
-    if fld.field is None:
+    if fld.values is None:
         print("normal_to_lognormal: no field stored in SRF class.")
     else:
-        fld.field = _normal_to_lognormal(fld.field)
+        fld.values = _normal_to_lognormal(fld.values)
 
 
 def normal_to_uniform(fld):
@@ -171,10 +175,12 @@ def normal_to_uniform(fld):
         Spatial Random Field class containing a generated field.
         Field will be transformed inplace.
     """
-    if fld.field is None:
+    if fld.values is None:
         print("normal_to_uniform: no field stored in SRF class.")
     else:
-        fld.field = _normal_to_uniform(fld.field, fld.mean, fld.model.sill)
+        fld.values = _normal_to_uniform(
+            fld.values, fld.field.mean, fld.model.sill
+        )
 
 
 def normal_to_arcsin(fld, a=None, b=None):
@@ -197,15 +203,15 @@ def normal_to_arcsin(fld, a=None, b=None):
         Parameter b of the arcsin distribution (upper bound).
         Default: keep mean and variance
     """
-    if fld.field is None:
+    if fld.values is None:
         print("normal_to_arcsin: no field stored in SRF class.")
     else:
-        a = fld.mean - np.sqrt(2.0 * fld.model.sill) if a is None else a
-        b = fld.mean + np.sqrt(2.0 * fld.model.sill) if b is None else b
-        fld.field = _normal_to_arcsin(
-            fld.field, fld.mean, fld.model.sill, a, b
+        a = fld.field.mean - np.sqrt(2.0 * fld.model.sill) if a is None else a
+        b = fld.field.mean + np.sqrt(2.0 * fld.model.sill) if b is None else b
+        fld.values = _normal_to_arcsin(
+            fld.values, fld.field.mean, fld.model.sill, a, b
         )
-        fld.mean = (a + b) / 2.0
+        fld.field.mean = (a + b) / 2.0
 
 
 def normal_to_uquad(fld, a=None, b=None):
@@ -228,13 +234,23 @@ def normal_to_uquad(fld, a=None, b=None):
         Parameter b of the U-quadratic distribution (upper bound).
         Default: keep mean and variance
     """
-    if fld.field is None:
+    if fld.values is None:
         print("normal_to_uquad: no field stored in SRF class.")
     else:
-        a = fld.mean - np.sqrt(5.0 / 3.0 * fld.model.sill) if a is None else a
-        b = fld.mean + np.sqrt(5.0 / 3.0 * fld.model.sill) if b is None else b
-        fld.field = _normal_to_uquad(fld.field, fld.mean, fld.model.sill, a, b)
-        fld.mean = (a + b) / 2.0
+        a = (
+            fld.field.mean - np.sqrt(5.0 / 3.0 * fld.model.sill)
+            if a is None
+            else a
+        )
+        b = (
+            fld.field.mean + np.sqrt(5.0 / 3.0 * fld.model.sill)
+            if b is None
+            else b
+        )
+        fld.values = _normal_to_uquad(
+            fld.values, fld.field.mean, fld.model.sill, a, b
+        )
+        fld.field.mean = (a + b) / 2.0
 
 
 # low level functions
@@ -442,5 +458,5 @@ def _uniform_to_uquad(field, a=0, b=1):
     y_raw = 3 * field / al + ga
     out = np.zeros_like(y_raw)
     out[y_raw > 0] = y_raw[y_raw > 0] ** (1 / 3)
-    out[y_raw < 0] = -(-y_raw[y_raw < 0]) ** (1 / 3)
+    out[y_raw < 0] = -((-y_raw[y_raw < 0]) ** (1 / 3))
     return out + be
diff --git a/tests/test_krige.py b/tests/test_krige.py
index 335084e4..8e4edfbb 100644
--- a/tests/test_krige.py
+++ b/tests/test_krige.py
@@ -37,13 +37,13 @@ def test_simple(self):
             simple = krige.Simple(
                 model, self.mean, self.cond_pos[0], self.cond_val
             )
-            field_1, __ = simple.unstructured(self.pos[0])
-            field_2, __ = simple.structured(self.pos[0])
+            field_1 = simple.unstructured(self.pos[0])
+            field_2 = simple.structured(self.pos[0])
             for i, val in enumerate(self.cond_val):
-                self.assertAlmostEqual(val, field_1[i], places=2)
-                self.assertAlmostEqual(val, field_2[(i,)], places=2)
-            self.assertAlmostEqual(self.mean, field_1[-1], places=2)
-            self.assertAlmostEqual(self.mean, field_2[(-1,)], places=2)
+                self.assertAlmostEqual(val, field_1.values[i], places=2)
+                self.assertAlmostEqual(val, field_2.values[(i,)], places=2)
+            self.assertAlmostEqual(self.mean, field_1.values[-1], places=2)
+            self.assertAlmostEqual(self.mean, field_2.values[(-1,)], places=2)
 
             for dim in self.dims[1:]:
                 model = Model(
@@ -56,14 +56,14 @@ def test_simple(self):
                 simple = krige.Simple(
                     model, self.mean, self.cond_pos[:dim], self.cond_val
                 )
-                field_1, __ = simple.unstructured(self.pos[:dim])
-                field_2, __ = simple.structured(self.pos[:dim])
+                field_1 = simple.unstructured(self.pos[:dim])
+                field_2 = simple.structured(self.pos[:dim])
                 for i, val in enumerate(self.cond_val):
-                    self.assertAlmostEqual(val, field_1[i], places=2)
-                    self.assertAlmostEqual(val, field_2[dim * (i,)], places=2)
-                self.assertAlmostEqual(self.mean, field_1[-1], places=2)
+                    self.assertAlmostEqual(val, field_1.values[i], places=2)
+                    self.assertAlmostEqual(val, field_2.values[dim * (i,)], places=2)
+                self.assertAlmostEqual(self.mean, field_1.values[-1], places=2)
                 self.assertAlmostEqual(
-                    self.mean, field_2[dim * (-1,)], places=2
+                    self.mean, field_2.values[dim * (-1,)], places=2
                 )
 
     def test_ordinary(self):
@@ -72,11 +72,11 @@ def test_ordinary(self):
                 dim=1, var=0.5, len_scale=2, anis=[0.1, 1], angles=[0.5, 0, 0]
             )
             ordinary = krige.Ordinary(model, self.cond_pos[0], self.cond_val)
-            field_1, __ = ordinary.unstructured(self.pos[0])
-            field_2, __ = ordinary.structured(self.pos[0])
+            field_1 = ordinary.unstructured(self.pos[0])
+            field_2 = ordinary.structured(self.pos[0])
             for i, val in enumerate(self.cond_val):
-                self.assertAlmostEqual(val, field_1[i], places=2)
-                self.assertAlmostEqual(val, field_2[(i,)], places=2)
+                self.assertAlmostEqual(val, field_1.values[i], places=2)
+                self.assertAlmostEqual(val, field_2.values[(i,)], places=2)
 
             for dim in self.dims[1:]:
                 model = Model(
@@ -89,11 +89,11 @@ def test_ordinary(self):
                 ordinary = krige.Ordinary(
                     model, self.cond_pos[:dim], self.cond_val
                 )
-                field_1, __ = ordinary.unstructured(self.pos[:dim])
-                field_2, __ = ordinary.structured(self.pos[:dim])
+                field_1 = ordinary.unstructured(self.pos[:dim])
+                field_2 = ordinary.structured(self.pos[:dim])
                 for i, val in enumerate(self.cond_val):
-                    self.assertAlmostEqual(val, field_1[i], places=2)
-                    self.assertAlmostEqual(val, field_2[dim * (i,)], places=2)
+                    self.assertAlmostEqual(val, field_1.values[i], places=2)
+                    self.assertAlmostEqual(val, field_2.values[dim * (i,)], places=2)
 
 
 if __name__ == "__main__":
diff --git a/tests/test_srf.py b/tests/test_srf.py
index a30e30bd..88d12a29 100644
--- a/tests/test_srf.py
+++ b/tests/test_srf.py
@@ -12,9 +12,12 @@
 
 class TestSRF(unittest.TestCase):
     def setUp(self):
-        self.cov_model = Gaussian(dim=2, var=1.5, len_scale=4.0, mode_no=100)
         self.mean = 0.3
+        self.var = 1.5
         self.mode_no = 100
+        self.cov_model = Gaussian(
+            dim=2, var=self.var, len_scale=4.0, mode_no=self.mode_no
+        )
 
         self.seed = 825718662
         self.x_grid = np.linspace(0.0, 12.0, 48)
@@ -219,7 +222,7 @@ def test_calls(self):
         srf = SRF(self.cov_model, mean=self.mean, mode_no=self.mode_no)
         field = srf((self.x_tuple, self.y_tuple), seed=self.seed)
         field2 = srf.unstructured((self.x_tuple, self.y_tuple), seed=self.seed)
-        self.assertAlmostEqual(field[0], srf.field[0])
+        self.assertAlmostEqual(field[0], field[0])
         self.assertAlmostEqual(field[0], field2[0])
         field = srf(
             (self.x_tuple, self.y_tuple),
@@ -227,16 +230,17 @@ def test_calls(self):
             mesh_type="structured",
         )
         field2 = srf.structured((self.x_tuple, self.y_tuple), seed=self.seed)
-        self.assertAlmostEqual(field[0, 0], srf.field[0, 0])
+        self.assertAlmostEqual(field[0, 0], field[0, 0])
         self.assertAlmostEqual(field[0, 0], field2[0, 0])
 
     def test_transform(self):
         self.cov_model.dim = 2
-        srf = SRF(self.cov_model, mean=self.mean, mode_no=self.mode_no)
+        srf = SRF(self.cov_model, mode_no=self.mode_no)
         srf((self.x_grid, self.y_grid), seed=self.seed, mesh_type="structured")
+        srf.field.mean = self.mean
         tf.normal_force_moments(srf)  # force ergodicity of the given field
-        self.assertAlmostEqual(srf.field.mean(), srf.mean)
-        self.assertAlmostEqual(srf.field.var(), srf.model.var)
+        self.assertAlmostEqual(srf.field.mean, self.mean)
+        self.assertAlmostEqual(srf.model.var, srf.model.var)
         tf.zinnharvey(srf)  # make high values mostly connected
         tf.normal_force_moments(srf)  # force ergodicity of the given field
         tf.normal_to_lognormal(srf)  # log-normal

From 0bff5deed1b27cc34099b19735d99cbbbb4b8740 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Mon, 20 Jan 2020 17:00:57 +0100
Subject: [PATCH 02/49] Blackened

---
 gstools/field/data.py | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/gstools/field/data.py b/gstools/field/data.py
index 4bb9b6c0..71f4ec59 100644
--- a/gstools/field/data.py
+++ b/gstools/field/data.py
@@ -22,7 +22,9 @@ class Data:
     """A data class mainly storing the specific values of an individual field.
     """
 
-    def __init__(self, values: np.ndarray, mean: float = 0.0, value_type: str = "scalar"):
+    def __init__(
+        self, values: np.ndarray, mean: float = 0.0, value_type: str = "scalar"
+    ):
         self.values = values
         self.mean = mean
         self.value_type = value_type

From 89e76b4b9028b9be4d0b534a34c3c9ab8ee74af1 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 23 Jan 2020 10:49:25 +0100
Subject: [PATCH 03/49] Fix bug in setter method

---
 gstools/field/data.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/gstools/field/data.py b/gstools/field/data.py
index 71f4ec59..e75689ee 100644
--- a/gstools/field/data.py
+++ b/gstools/field/data.py
@@ -118,4 +118,4 @@ def mean(self):
 
     @mean.setter
     def mean(self, value):
-        self.fields[self.default_field].mean = mean
+        self.fields[self.default_field].mean = value

From bb6741175480820356b1b33dede8447f92f6a390 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 23 Jan 2020 10:50:14 +0100
Subject: [PATCH 04/49] [WIP] Add methods for data checking

---
 gstools/field/data.py | 13 +++++++++++++
 1 file changed, 13 insertions(+)

diff --git a/gstools/field/data.py b/gstools/field/data.py
index e75689ee..2ca6478c 100644
--- a/gstools/field/data.py
+++ b/gstools/field/data.py
@@ -66,6 +66,8 @@ def __init__(
         # initialize attributes
         self.pos = pos
         self.fields: Dict[str, np.ndarray] = {}
+        if mesh_type != "unstructured" and mesh_type != "structured":
+            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
         self.mesh_type = mesh_type
 
     def add_field(
@@ -77,6 +79,8 @@ def add_field(
         value_type: str = "scalar",
         default_field: bool = False,
     ):
+        values = np.array(values)
+        self._check_field(values)
         self.fields[name] = Data(values, mean, value_type)
         # set the default_field to the first field added
         if len(self.fields) == 1 or default_field:
@@ -119,3 +123,12 @@ def mean(self):
     @mean.setter
     def mean(self, value):
         self.fields[self.default_field].mean = value
+
+    def _check_field(self, values: np.ndarray):
+        # TODO
+        if self.mesh_type == "unstructured":
+            pass
+        elif self.mesh_type == "structured":
+            pass
+        else:
+            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type)

From 9e53dc9f237eb6c9574cd04d00372047c057a4ab Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 23 Jan 2020 10:50:57 +0100
Subject: [PATCH 05/49] Move mean calc. after krige field is added

This way, we do not need a temporary field to store the mean value.
---
 gstools/krige/ordinary.py | 12 +++++++-----
 1 file changed, 7 insertions(+), 5 deletions(-)

diff --git a/gstools/krige/ordinary.py b/gstools/krige/ordinary.py
index bc8d8fa9..b60c21e8 100644
--- a/gstools/krige/ordinary.py
+++ b/gstools/krige/ordinary.py
@@ -105,11 +105,6 @@ def __call__(self, pos, mesh_type="unstructured"):
         krig_vecs = self._get_vario_mat((c_x, c_y, c_z), (x, y, z), add=True)
         # generate the kriged field
         field, krige_var = krigesum(krig_mat, krig_vecs, cond)
-        # calculate the estimated mean (kriging field at infinity)
-        mean_est = np.concatenate(
-            (np.full_like(self.cond_val, self.model.sill), [1])
-        )
-        self.mean = np.einsum("i,ij,j", cond, krig_mat, mean_est)
 
         # reshape field if we got an unstructured mesh
         if mesh_type_changed:
@@ -122,7 +117,14 @@ def __call__(self, pos, mesh_type="unstructured"):
             )
         # save the field
         self.add_field("krige", field, default_field=True)
+        # calculate the estimated mean (kriging field at infinity)
+        mean_est = np.concatenate(
+            (np.full_like(self.cond_val, self.model.sill), [1])
+        )
+        self.mean = np.einsum("i,ij,j", cond, krig_mat, mean_est)
+
         self.field.krige_var = krige_var
+
         return self.field
 
     def _get_krig_mat(self, pos1, pos2):

From 32a64fa7fd6e169037298d725cf761541c0ee6f4 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 23 Jan 2020 10:52:07 +0100
Subject: [PATCH 06/49] Add todo note

---
 tests/test_condition.py | 1 +
 1 file changed, 1 insertion(+)

diff --git a/tests/test_condition.py b/tests/test_condition.py
index e7c58879..3658df89 100644
--- a/tests/test_condition.py
+++ b/tests/test_condition.py
@@ -46,6 +46,7 @@ def test_simple(self):
             )
             srf = SRF(model, self.mean, seed=19970221)
             srf.set_condition(self.cond_pos[0], self.cond_val, "simple")
+            # TODO should this be possible?!?
             field_1 = srf.unstructured(self.pos[0])
             field_2 = srf.structured(self.pos[0])
             for i, val in enumerate(self.cond_val):

From e18138b8a17082bdc9e13ba0cfa76e9c20a5b48b Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 23 Jan 2020 16:51:16 +0100
Subject: [PATCH 07/49] [WIP] Move field base classes and add docstrings

---
 gstools/field/base.py | 166 +++++++++++++++++++++++++++++++++++++++++-
 gstools/field/data.py | 134 ----------------------------------
 2 files changed, 164 insertions(+), 136 deletions(-)
 delete mode 100644 gstools/field/data.py

diff --git a/gstools/field/base.py b/gstools/field/base.py
index ec7be97c..b0b5d75b 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -8,19 +8,181 @@
 
 .. autosummary::
    Field
+   FieldData
 """
 # pylint: disable=C0103
 
 from functools import partial
+from typing import List, Dict
 
 import numpy as np
 
-from gstools.field.data import FieldData
 from gstools.covmodel.base import CovModel
 from gstools.tools.export import to_vtk, vtk_export
 from gstools.field.tools import _get_select
 
-__all__ = ["Field"]
+__all__ = ["Field", "FieldData"]
+
+
+class Data:
+    """A data class mainly storing the specific values of an individual field.
+
+    Instances of this class are stored in a dictionary owned by FieldBase.
+    The positions on which these field values are defined are also saved in
+    FieldBase.
+    """
+
+    def __init__(
+        self,
+        values: np.ndarray = None,
+        mean: float = 0.0,
+        value_type: str = "scalar",
+    ):
+        self.values = values
+        self.mean = mean
+        self.value_type = value_type
+
+
+class FieldData:
+    """A base class encapsulating field data.
+
+    It holds a position array, which define the spatial locations of the
+    field values.
+    It can hold multiple fields in the :any:`self.fields` dict. This assumes
+    that each field is defined on the same positions.
+    The mesh type must also be specified.
+
+    Parameters
+    ----------
+    pos : :class:`numpy.ndarray`, optional
+        positions of the field values
+    name : :any:`str`, optional
+        key of the field values
+    values : :any:`list`, optional
+        a list of the values of the fields
+    mean : :any:`float`, optional
+        mean of the field
+    value_type : :any:`str`, optional
+        the value type of the default field, can be
+        * scalar
+        * vector
+    mesh_type : :any:`str`, optional
+        the type of mesh on which the field is defined on, can be
+        * unstructured
+        * structured
+
+    Examples
+    --------
+    >>> import numpy as np
+    >>> pos = np.random.random((100, 100))
+    >>> z = np.random.random((100, 100))
+    >>> z2 = np.random.random((100, 100))
+    >>> field_data = FieldData(pos)
+    >>> field_data.add_field("test_field1", z)
+    >>> field_data.add_field("test_field2", z2)
+    >>> field_data.set_default_field("test_field2")
+    >>> print(field.values)
+
+    """
+
+    def __init__(
+        self,
+        pos: np.ndarray = None,
+        name: str = "default_field",
+        values: np.ndarray = None,
+        *,
+        mean: float = 0.0,
+        value_type: str = "scalar",
+        mesh_type: str = "unstructured",
+    ):
+        # initialize attributes
+        self.pos = pos
+        self.fields: Dict[str, np.ndarray] = {}
+        self._default_field = "default_field"
+        if mesh_type != "unstructured" and mesh_type != "structured":
+            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
+        self.mesh_type = mesh_type
+        self.add_field(name, values, mean=mean, value_type=value_type)
+
+    def add_field(
+        self,
+        name: str,
+        values: np.ndarray,
+        *,
+        mean: float = 0.0,
+        value_type: str = "scalar",
+        default_field: bool = False,
+    ):
+        values = np.array(values)
+        self._check_field(values)
+        self.fields[name] = Data(values, mean, value_type)
+        # set the default_field to the first field added
+        if len(self.fields) == 1 or default_field:
+            self._default_field = name
+
+    def get_data(self, key):
+        """:class:`Data`: The field data class."""
+        return self.fields[key]
+
+    def __getitem__(self, key):
+        """:any:`numpy.ndarray`: The values of the field."""
+        return self.fields[key].values
+
+    def __setitem__(self, key, value):
+        self.fields[key].values = value
+
+    @property
+    def default_field(self):
+        """:any:`str`: The key of the default field."""
+        return self._default_field
+
+    @default_field.setter
+    def default_field(self, value):
+        self._default_field = value
+
+    @property
+    def field(self):
+        """:any:`Data`: The Data instance of the default field."""
+        return self.fields[self.default_field]
+
+    @property
+    def values(self):
+        """:any:`numpy.ndarray`: The values of the default field."""
+        return self.fields[self.default_field].values
+
+    @values.setter
+    def values(self, values):
+        self.fields[self.default_field].values = values
+
+    @property
+    def value_type(self):
+        """:any:`str`: The value type of the default field."""
+        return self.fields[self.default_field].value_type
+
+    @property
+    def mean(self):
+        """:any:`float`: The mean of the default field."""
+        return self.fields[self.default_field].mean
+
+    @mean.setter
+    def mean(self, value):
+        self.fields[self.default_field].mean = value
+
+    def _check_field(self, values: np.ndarray):
+        """Compare field shape to pos shape.
+
+        Parameters
+        ----------
+        values : :class:`numpy.ndarray`
+            the values of the field to be checked
+        """
+        # TODO
+        if self.mesh_type == "unstructured":
+            pass
+        elif self.mesh_type == "structured":
+            pass
+        else:
+            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
 
 
 class Field(FieldData):
diff --git a/gstools/field/data.py b/gstools/field/data.py
deleted file mode 100644
index 2ca6478c..00000000
--- a/gstools/field/data.py
+++ /dev/null
@@ -1,134 +0,0 @@
-# -*- coding: utf-8 -*-
-"""
-GStools subpackage providing a data class for spatial fields.
-
-.. currentmodule:: gstools.field.data
-
-The following classes are provided
-
-.. autosummary::
-   FieldData
-"""
-# pylint: disable=C0103
-
-__all__ = ["FieldData"]
-
-
-from typing import List, Dict
-import numpy as np
-
-
-class Data:
-    """A data class mainly storing the specific values of an individual field.
-    """
-
-    def __init__(
-        self, values: np.ndarray, mean: float = 0.0, value_type: str = "scalar"
-    ):
-        self.values = values
-        self.mean = mean
-        self.value_type = value_type
-
-
-class FieldData:
-    """A base class encapsulating field data.
-
-    It holds a position array, which define the spatial locations of the
-    field values.
-    It can hold multiple fields in the :any:`self.values` list. This assumes
-    that each field is defined on the same positions.
-    The mesh type must also be specified.
-
-    Parameters
-    ----------
-    pos : :class:`numpy.ndarray`, optional
-        positions of the field values
-    values : :any:`list`, optional
-        a list of the values of the fields
-
-    Examples
-    --------
-    >>> import numpy as np
-    >>> pos = np.random.random((100, 100))
-    >>> z = np.random.random((100, 100))
-    >>> z2 = np.random.random((100, 100))
-    >>> field_data = FieldData(pos)
-    >>> field_data.add_field("test_field1", z)
-    >>> field_data.add_field("test_field2", z2)
-    >>> field_data.set_default_field("test_field2")
-    >>> print(field.values)
-
-    """
-
-    def __init__(
-        self, pos: np.ndarray = None, mesh_type: str = "unstructured",
-    ):
-        # initialize attributes
-        self.pos = pos
-        self.fields: Dict[str, np.ndarray] = {}
-        if mesh_type != "unstructured" and mesh_type != "structured":
-            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
-        self.mesh_type = mesh_type
-
-    def add_field(
-        self,
-        name: str,
-        values: np.ndarray,
-        mean: float = 0.0,
-        *,
-        value_type: str = "scalar",
-        default_field: bool = False,
-    ):
-        values = np.array(values)
-        self._check_field(values)
-        self.fields[name] = Data(values, mean, value_type)
-        # set the default_field to the first field added
-        if len(self.fields) == 1 or default_field:
-            self.default_field = name
-
-    def get_data(self, key):
-        """:class:`Data`: The field data class."""
-        return self.fields[key]
-
-    def set_default_field(self, default_field):
-        self.default_field = default_field
-
-    def __setitem__(self, key, value):
-        self.fields[key].values = value
-
-    def __getitem__(self, key):
-        """:any:`numpy.ndarray`: The values of the field."""
-        return self.fields[key].values
-
-    @property
-    def field(self):
-        return self.fields[self.default_field]
-
-    @property
-    def values(self):
-        return self.fields[self.default_field].values
-
-    @values.setter
-    def values(self, values):
-        self.fields[self.default_field].values = values
-
-    @property
-    def value_type(self):
-        return self.fields[self.default_field].value_type
-
-    @property
-    def mean(self):
-        return self.fields[self.default_field].mean
-
-    @mean.setter
-    def mean(self, value):
-        self.fields[self.default_field].mean = value
-
-    def _check_field(self, values: np.ndarray):
-        # TODO
-        if self.mesh_type == "unstructured":
-            pass
-        elif self.mesh_type == "structured":
-            pass
-        else:
-            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type)

From d5053ada8096b01f8ab896401b31ff83d86b6283 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 23 Jan 2020 17:18:20 +0100
Subject: [PATCH 08/49] [WIP] 'mean' back to c'tor and refactor add_field

---
 gstools/field/base.py | 4 ++--
 gstools/field/srf.py  | 6 +++++-
 2 files changed, 7 insertions(+), 3 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index b0b5d75b..036099d7 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -194,9 +194,9 @@ class Field(FieldData):
         Covariance Model related to the field.
     """
 
-    def __init__(self, model):
+    def __init__(self, model, mean=0.0):
         # initialize attributes
-        super().__init__(self)
+        super().__init__(self, mean=mean)
         # initialize private attributes
         self._model = None
         self.model = model
diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index c1851b5e..b5ac5e46 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -157,7 +157,11 @@ def __call__(
         y, z = make_isotropic(self.model.dim, self.model.anis, y, z)
 
         # generate the field
-        self.add_field("raw", self.generator.__call__(x, y, z, mesh_type), self._mean)
+        self.add_field(
+                "raw",
+                self.generator.__call__(x, y, z, mesh_type),
+                mean=self._mean
+            )
 
         # reshape field if we got an unstructured mesh
         if mesh_type_changed:

From 0dcd4bf6e526a7ba3bc308bfb201d33cf6943117 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 10:28:06 +0100
Subject: [PATCH 09/49] SRF __call__ is now using correct default field

---
 gstools/field/srf.py | 10 +++++-----
 1 file changed, 5 insertions(+), 5 deletions(-)

diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index b5ac5e46..bfb8575a 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -180,21 +180,21 @@ def __call__(
                 info,
             ) = self.cond_func(self)
             # store everything in the class
-            self.add_field("cond", cond_field, default_field=True)
+            self.add_field("default_field", cond_field, default_field=True)
             self.add_field("krige", krige_field)
             self.add_field("error", err_field)
             self.get_data("krige").krige_var = krigevar
             if "mean" in info:  # ordinary kriging estimates mean
                 self.get_data("cond").mean = info["mean"]
         else:
-            self.add_field("srf", self["raw"] + self.mean)
+            self.add_field("default_field", self["raw"] + self.mean)
 
         # upscaled variance
         if not np.isscalar(point_volumes) or not np.isclose(point_volumes, 0):
             scaled_var = self.upscaling_func(self.model, point_volumes)
-            self["srf"] -= self.mean
-            self["srf"] *= np.sqrt(scaled_var / self.model.sill)
-            self["srf"] += self.mean
+            self.values -= self.mean
+            self.values *= np.sqrt(scaled_var / self.model.sill)
+            self.values += self.mean
 
         return self.values
 

From 60491a3817d85f665cec2af070eab75d123dcb0b Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 11:03:31 +0100
Subject: [PATCH 10/49] Remove a getter which is already defined in parent

---
 gstools/field/base.py | 5 -----
 1 file changed, 5 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 036099d7..44c0664a 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -471,11 +471,6 @@ def model(self, model):
                 "Field: 'model' is not an instance of 'gstools.CovModel'"
             )
 
-    @property
-    def value_type(self):
-        """:class:`str`: Type of the field values (scalar, vector)."""
-        return self._value_type
-
     def __str__(self):
         """Return String representation."""
         return self.__repr__()

From 12357f44fbfd6b7d27ae1eba952b1686fb0ca64a Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 11:17:37 +0100
Subject: [PATCH 11/49] Remove (now) wrong doc line

---
 gstools/field/srf.py | 2 --
 1 file changed, 2 deletions(-)

diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index bfb8575a..7fb98b2b 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -113,8 +113,6 @@ def __call__(
     ):
         """Generate the spatial random field.
 
-        The field is saved as `self.field` and is also returned.
-
         Parameters
         ----------
         pos : :class:`list`

From dfc8a38283b2a267867c90520cf626fd56d961d0 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 11:22:23 +0100
Subject: [PATCH 12/49] Add pos setter and getter

---
 gstools/field/base.py | 11 ++++++++++-
 1 file changed, 10 insertions(+), 1 deletion(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 44c0664a..478e8a5a 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -96,7 +96,7 @@ def __init__(
         mesh_type: str = "unstructured",
     ):
         # initialize attributes
-        self.pos = pos
+        self._pos = pos
         self.fields: Dict[str, np.ndarray] = {}
         self._default_field = "default_field"
         if mesh_type != "unstructured" and mesh_type != "structured":
@@ -131,6 +131,15 @@ def __getitem__(self, key):
     def __setitem__(self, key, value):
         self.fields[key].values = value
 
+    @property
+    def pos(self):
+        """:any:`numpy.ndarray`: The pos. on which the field is defined."""
+        return self._pos
+
+    @pos.setter
+    def pos(self, value):
+        self._pos = value
+
     @property
     def default_field(self):
         """:any:`str`: The key of the default field."""

From f0f341fea396eb656c283273544f9d667291ab11 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 11:33:29 +0100
Subject: [PATCH 13/49] Amend to commit 60491a3

---
 gstools/field/base.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 478e8a5a..8ad87a07 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -209,7 +209,6 @@ def __init__(self, model, mean=0.0):
         # initialize private attributes
         self._model = None
         self.model = model
-        self._value_type = None
 
     def __call__(*args, **kwargs):
         """Generate the field."""

From 4eee22b243040126824459c90e79e44ccb9406a5 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 12:19:46 +0100
Subject: [PATCH 14/49] Add more getters, setters

---
 gstools/field/base.py | 25 +++++++++++++++++++++----
 1 file changed, 21 insertions(+), 4 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 8ad87a07..630b46c6 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -88,7 +88,7 @@ class FieldData:
     def __init__(
         self,
         pos: np.ndarray = None,
-        name: str = "default_field",
+        name: str = "field",
         values: np.ndarray = None,
         *,
         mean: float = 0.0,
@@ -98,10 +98,10 @@ def __init__(
         # initialize attributes
         self._pos = pos
         self.fields: Dict[str, np.ndarray] = {}
-        self._default_field = "default_field"
+        self._default_field = "field"
         if mesh_type != "unstructured" and mesh_type != "structured":
             raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
-        self.mesh_type = mesh_type
+        self._mesh_type = mesh_type
         self.add_field(name, values, mean=mean, value_type=value_type)
 
     def add_field(
@@ -138,7 +138,11 @@ def pos(self):
 
     @pos.setter
     def pos(self, value):
+        """
+        Warning: setting new positions deletes all previously stored fields.
+        """
         self._pos = value
+        self.fields = {}
 
     @property
     def default_field(self):
@@ -177,6 +181,19 @@ def mean(self):
     def mean(self, value):
         self.fields[self.default_field].mean = value
 
+    @property
+    def mesh_type(self):
+        """:any:`str`: The mesh type of the fields."""
+        return self._mesh_type
+
+    @mesh_type.setter
+    def mesh_type(self, value):
+        """
+        Warning: setting a new mesh type deletes all previously stored fields.
+        """
+        self._mesh_type = value
+        self.fields = {}
+
     def _check_field(self, values: np.ndarray):
         """Compare field shape to pos shape.
 
@@ -210,7 +227,7 @@ def __init__(self, model, mean=0.0):
         self._model = None
         self.model = model
 
-    def __call__(*args, **kwargs):
+    def __call__(self, *args, **kwargs):
         """Generate the field."""
         pass
 

From c75a45c1edf5d5735e0e43cd1df63e586f2c63ec Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 13:09:30 +0100
Subject: [PATCH 15/49] Add reset fct. and checks to FieldData class

---
 gstools/field/base.py | 14 ++++++++++++--
 1 file changed, 12 insertions(+), 2 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 630b46c6..c58a98bf 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -102,6 +102,11 @@ def __init__(
         if mesh_type != "unstructured" and mesh_type != "structured":
             raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
         self._mesh_type = mesh_type
+        if value_type != "scalar" and value_type != "vector":
+            raise ValueError(
+                "Unknown field value type, "
+                + "specify 'scalar' or 'vector'."
+            )
         self.add_field(name, values, mean=mean, value_type=value_type)
 
     def add_field(
@@ -142,7 +147,7 @@ def pos(self, value):
         Warning: setting new positions deletes all previously stored fields.
         """
         self._pos = value
-        self.fields = {}
+        self._reset()
 
     @property
     def default_field(self):
@@ -192,7 +197,12 @@ def mesh_type(self, value):
         Warning: setting a new mesh type deletes all previously stored fields.
         """
         self._mesh_type = value
-        self.fields = {}
+        self._reset()
+
+    def _reset(self):
+        self._default_field = "field"
+        self._mesh_type = "unstructured"
+        self.add_field(self._default_field, None, mean=0.0, value_type="scalar")
 
     def _check_field(self, values: np.ndarray):
         """Compare field shape to pos shape.

From b6bbb9472225a359bb13515edff161d9d4a011d5 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Jan 2020 13:10:36 +0100
Subject: [PATCH 16/49] [WIP] updating SRF class

---
 gstools/field/srf.py | 30 ++++++++++++++++++------------
 1 file changed, 18 insertions(+), 12 deletions(-)

diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index 7fb98b2b..9a90445b 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -89,9 +89,8 @@ def __init__(
         generator="RandMeth",
         **generator_kwargs
     ):
-        super().__init__(model)
+        super().__init__(model, mean)
         # initialize private attributes
-        self._mean = mean
         self._generator = None
         self._upscaling = None
         self._upscaling_func = None
@@ -102,14 +101,14 @@ def __init__(
         # initialize attributes
         self.set_generator(generator, **generator_kwargs)
         self.upscaling = upscaling
-        if self._value_type is None:
-            raise ValueError(
-                "Unknown field value type, "
-                + "specify 'scalar' or 'vector' before calling SRF."
-            )
 
     def __call__(
-        self, pos, seed=np.nan, point_volumes=0.0, mesh_type="unstructured"
+        self,
+        pos,
+        name="field",
+        seed=np.nan,
+        point_volumes=0.0,
+        mesh_type="unstructured",
     ):
         """Generate the spatial random field.
 
@@ -118,6 +117,11 @@ def __call__(
         pos : :class:`list`
             the position tuple, containing main direction and transversal
             directions
+        name : :class:`str`, optional
+            the name of the dictionary key, in which the values will be saved.
+            If this method is called multiple times, the name determines, if
+            the field is overwritten or if it is added to the values dict.
+            Default: "field"
         seed : :class:`int`, optional
             seed for RNG for reseting. Default: keep seed from generator
         point_volumes : :class:`float` or :class:`numpy.ndarray`
@@ -136,8 +140,10 @@ def __call__(
         """
         # internal conversation
         x, y, z = pos2xyz(pos, max_dim=self.model.dim)
+        mean_temp = self.mean
         self.pos = xyz2pos(x, y, z)
         self.mesh_type = mesh_type
+        self.mean = mean_temp
         # update the model/seed in the generator if any changes were made
         self.generator.update(self.model, seed)
         # format the positional arguments of the mesh
@@ -158,7 +164,7 @@ def __call__(
         self.add_field(
                 "raw",
                 self.generator.__call__(x, y, z, mesh_type),
-                mean=self._mean
+                mean=self.mean
             )
 
         # reshape field if we got an unstructured mesh
@@ -178,14 +184,14 @@ def __call__(
                 info,
             ) = self.cond_func(self)
             # store everything in the class
-            self.add_field("default_field", cond_field, default_field=True)
+            self.add_field(name, cond_field, default_field=True)
             self.add_field("krige", krige_field)
             self.add_field("error", err_field)
             self.get_data("krige").krige_var = krigevar
             if "mean" in info:  # ordinary kriging estimates mean
                 self.get_data("cond").mean = info["mean"]
         else:
-            self.add_field("default_field", self["raw"] + self.mean)
+            self.add_field(name, self["raw"] + self.mean)
 
         # upscaled variance
         if not np.isscalar(point_volumes) or not np.isclose(point_volumes, 0):
@@ -308,7 +314,7 @@ def upscaling(self, upscaling):
     def __repr__(self):
         """Return String representation."""
         return "SRF(model={0}, mean={1}, generator={2}".format(
-            self.model, self._mean, self.generator
+            self.model, self.mean, self.generator
         )
 
 

From fe257c8ccfa5de6ea2c20e911c2db0f4eea574e4 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Mon, 10 Feb 2020 17:39:04 +0100
Subject: [PATCH 17/49] [WIP] Add new Mesh class for discussion

---
 tests/test_mesh.py | 374 +++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 374 insertions(+)
 create mode 100644 tests/test_mesh.py

diff --git a/tests/test_mesh.py b/tests/test_mesh.py
new file mode 100644
index 00000000..291a1470
--- /dev/null
+++ b/tests/test_mesh.py
@@ -0,0 +1,374 @@
+import unittest
+from typing import Dict, Tuple
+import numpy as np
+
+
+def value_type(mesh_type, shape):
+    """Determine the value type ("scalar" or "vector")"""
+    r = "scalar"
+    if mesh_type == "unstructured":
+        # TODO is this the right place for doing these checks?!
+        if len(shape) == 2 and 2 <= shape[0] <= 3:
+            r = "vector"
+    else:
+        # for very small (2-3 points) meshes, this could break
+        # a 100% solution would require dim. information.
+        if len(shape) == shape[0] + 1:
+            r = "vector"
+    return r
+
+
+class Mesh:
+    def __init__(
+        self,
+        pos=None,
+        name: str = "field",
+        values: np.ndarray = None,
+        *,
+        mesh_type: str = "unstructured",
+    ) -> None:
+        # the pos/ points of the mesh
+        self._pos = pos
+
+        # data stored at each pos/ point, the "fields"
+        self.point_data: Dict[str, np.ndarray] = {name: values}
+
+        # data valid for the global field
+        self.field_data = {}
+
+        self.set_field_data("default_field", name)
+
+        # mesh_type needs a special setter, therefore, `set_field_data` is not
+        # used here
+        self.mesh_type = mesh_type
+        self.field_data["mesh_type"] = mesh_type
+
+    def set_field_data(self, name: str, value) -> None:
+        """Add an attribute to this instance and add it the `field_data`
+
+        This helper method is used to create attributes for easy access
+        while using this class, but it also adds an entry to the dictionary
+        `field_data`, which is used for exporting the data.
+        """
+        setattr(self, name, value)
+        self.field_data[name] = value
+
+    def add_field(
+        self,
+        values: np.ndarray,
+        name: str = "field",
+        *,
+        default_field: bool = False,
+    ) -> None:
+        """Add a field (point data) to the mesh
+
+        .. note::
+            If no field has existed before calling this method,
+            the given field will be set to the default field.
+
+        .. warning::
+            If point data with the same `name` already exists, it will be
+            overwritten.
+
+        Parameters
+        ----------
+        values : :class:`numpy.ndarray`
+            the point data, has to be the same shape as the mesh
+        name : :class:`str`, optional
+            the name of the point data
+        default_field : :class:`bool`, optional
+            is this the default field of the mesh?
+
+        """
+
+        values = np.array(values)
+        self._check_point_data(values)
+        self.point_data[name] = values
+        # set the default_field to the first field added
+        if len(self.point_data) == 1 or default_field:
+            self._default_field = name
+
+    def __getitem__(self, key: str) -> np.ndarray:
+        """:any:`numpy.ndarray`: The values of the field."""
+        return self.point_data[key]
+
+    def __setitem__(self, key: str, value):
+        self.point_data[key] = value
+
+    @property
+    def pos(self) -> Tuple[np.ndarray]:
+        """:any:`numpy.ndarray`: The pos. on which the field is defined."""
+        return self._pos
+
+    @pos.setter
+    def pos(self, value: Tuple[np.ndarray]):
+        """
+        Warning: setting new positions deletes all previously stored fields.
+        """
+        # TODO what is the best way to handle different fields?
+        # * delete all except for default field
+        # * delete them all
+        # * only delete the values, not the keys
+        self.point_data = {self.default_field: None}
+        self._pos = value
+
+    @property
+    def field(self) -> np.ndarray:
+        """:class:`numpy.ndarray`: The point data of the default field."""
+        return self.point_data[self.default_field]
+
+    @field.setter
+    def field(self, values: np.ndarray):
+        self._check_point_data(values)
+        self.point_data[self.default_field] = values
+
+    @property
+    def value_type(self, field="field") -> str:
+        """:any:`str`: The value type of the default field."""
+        if self.point_data[field] is None:
+            r = None
+        else:
+            r = value_type(self.mesh_type, self.point_data[field].shape)
+        return r
+
+    @property
+    def mesh_type(self) -> str:
+        """:any:`str`: The mesh type of the fields."""
+        return self._mesh_type
+
+    @mesh_type.setter
+    def mesh_type(self, value: str):
+        """
+        Warning: setting a new mesh type deletes all previously stored fields.
+        """
+        self._check_mesh_type(value)
+        self.point_data = {self.default_field: None}
+        self._mesh_type = value
+
+    def _check_mesh_type(self, mesh_type: str) -> None:
+        if mesh_type != "unstructured" and mesh_type != "structured":
+            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
+
+    def _check_point_data(self, values: np.ndarray):
+        """Compare field shape to pos shape.
+
+        Parameters
+        ----------
+        values : :class:`numpy.ndarray`
+            the values of the field to be checked
+        """
+        err = True
+        if self.mesh_type == "unstructured":
+            # scalar
+            if len(values.shape) == 1:
+                if values.shape[0] == len(self.pos[0]):
+                    err = False
+            # vector
+            elif len(values.shape) == 2:
+                if (
+                    values.shape[1] == len(self.pos[0])
+                    and 2 <= values.shape[0] <= 3
+                ):
+                    err = False
+            if err:
+                raise ValueError(
+                    "Wrong field shape: {0} does not match mesh shape ({1},)".format(
+                        values.shape, len(self.pos[0])
+                    )
+                )
+        else:
+            # scalar
+            if len(values.shape) == len(self.pos):
+                if all(
+                    [
+                        values.shape[i] == len(self.pos[i])
+                        for i in range(len(self.pos))
+                    ]
+                ):
+                    err = False
+            # vector
+            elif len(values.shape) == len(self.pos) + 1:
+                if all(
+                    [
+                        values.shape[i + 1] == len(self.pos[i])
+                        for i in range(len(self.pos))
+                    ]
+                ) and values.shape[0] == len(self.pos):
+                    err = False
+            if err:
+                raise ValueError(
+                    "Wrong field shape: {0} does not match mesh shape [0/2/3]{1}".format(
+                        list(values.shape),
+                        [len(self.pos[i]) for i in range(len(self.pos))],
+                    )
+                )
+
+
+class TestMesh(unittest.TestCase):
+    def setUp(self):
+        self.x_grid = np.linspace(0.0, 12.0, 48)
+        self.y_grid = np.linspace(0.0, 10.0, 46)
+        self.z_grid = np.linspace(0.0, 10.0, 40)
+
+        self.f1_grid = self.x_grid
+        self.f2_grid = self.x_grid.reshape(
+            (len(self.x_grid), 1)
+        ) * self.y_grid.reshape((1, len(self.y_grid)))
+        self.f3_grid = (
+            self.x_grid.reshape((len(self.x_grid), 1, 1))
+            * self.y_grid.reshape((1, len(self.y_grid), 1))
+            * self.z_grid.reshape((1, 1, len(self.z_grid)))
+        )
+
+        self.rng = np.random.RandomState(123018)
+        self.x_tuple = self.rng.uniform(0.0, 10, 100)
+        self.y_tuple = self.rng.uniform(0.0, 10, 100)
+        self.z_tuple = self.rng.uniform(0.0, 10, 100)
+
+        self.f1_tuple = self.x_tuple
+        self.f2_tuple = self.x_tuple * self.y_tuple
+        self.f3_tuple = self.x_tuple * self.y_tuple * self.z_tuple
+
+        self.m1_grid = Mesh((self.x_grid,), mesh_type="structured")
+        self.m2_grid = Mesh((self.x_grid, self.y_grid), mesh_type="structured")
+        self.m3_grid = Mesh(
+            (self.x_grid, self.y_grid, self.z_grid), mesh_type="structured"
+        )
+        self.m1_tuple = Mesh((self.x_tuple,))
+        self.m2_tuple = Mesh((self.x_tuple, self.y_tuple))
+        self.m3_tuple = Mesh((self.x_tuple, self.y_tuple, self.z_tuple))
+
+    def test_item_getter_setter(self):
+        self.m3_grid.add_field(256.0 * self.f3_grid, name="2nd")
+        self.m3_grid.add_field(512.0 * self.f3_grid, name="3rd")
+        self.assertEqual(
+            self.m3_grid["2nd"].all(), (256.0 * self.f3_grid).all()
+        )
+        self.assertEqual(
+            self.m3_grid["3rd"].all(), (512.0 * self.f3_grid).all()
+        )
+
+        self.m3_tuple["tmp_field"] = 2.0 * self.f3_tuple
+        self.assertEqual(
+            self.m3_tuple["tmp_field"].all(), (2.0 * self.f3_tuple).all()
+        )
+
+    def test_pos_getter(self):
+        self.assertEqual(self.m1_grid.pos, (self.x_grid,))
+        self.assertEqual(self.m1_tuple.pos, (self.x_tuple,))
+        self.assertEqual(self.m2_grid.pos, (self.x_grid, self.y_grid))
+        self.assertEqual(
+            self.m3_grid.pos, (self.x_grid, self.y_grid, self.z_grid)
+        )
+        self.assertEqual(
+            self.m3_tuple.pos, (self.x_tuple, self.y_tuple, self.z_tuple)
+        )
+
+    def test_default_value_type(self):
+        # value_type getter with no field set
+        self.assertEqual(self.m1_tuple.value_type, None)
+        self.assertEqual(self.m2_tuple.value_type, None)
+        self.assertEqual(self.m3_tuple.value_type, None)
+        self.assertEqual(self.m1_grid.value_type, None)
+        self.assertEqual(self.m2_grid.value_type, None)
+        self.assertEqual(self.m3_grid.value_type, None)
+
+    def test_field_data_setter(self):
+        # attribute creation by adding field_data
+        self.m2_tuple.set_field_data("mean", 3.14)
+        self.assertEqual(self.m2_tuple.field_data["mean"], 3.14)
+        self.assertEqual(self.m2_tuple.mean, 3.14)
+
+    def test_new_pos(self):
+        # set new pos. (which causes reset)
+        x_tuple2 = self.rng.uniform(0.0, 10, 100)
+        y_tuple2 = self.rng.uniform(0.0, 10, 100)
+        self.m2_tuple.add_field(self.f2_tuple)
+
+        self.m2_tuple.pos = (x_tuple2, y_tuple2)
+
+        self.assertEqual(self.m2_tuple.pos, (x_tuple2, y_tuple2))
+
+        # previous field has to be deleted
+        self.assertEqual(self.m2_tuple.field, None)
+
+    def test_add_field(self):
+        # structured
+        self.m1_grid.add_field(self.f1_grid)
+        self.assertEqual(self.m1_grid.field.all(), self.f1_grid.all())
+        self.m2_grid.add_field(self.f2_grid)
+        self.assertEqual(self.m2_grid.field.all(), self.f2_grid.all())
+        self.m3_grid.add_field(self.f3_grid)
+        self.assertEqual(self.m3_grid.field.all(), self.f3_grid.all())
+
+        # unstructured
+        self.m1_tuple.add_field(self.f1_tuple)
+        self.assertEqual(self.m1_tuple.field.all(), self.f1_tuple.all())
+        self.m2_tuple.add_field(self.f2_tuple)
+        self.assertEqual(self.m2_tuple.field.all(), self.f2_tuple.all())
+        self.m3_tuple.add_field(self.f3_tuple)
+        self.assertEqual(self.m3_tuple.field.all(), self.f3_tuple.all())
+
+        # multiple fields
+        new_field = 10.0 * self.f1_grid
+        self.m1_grid.add_field(new_field, name="2nd")
+        # default field
+        self.assertEqual(self.m1_grid.field.all(), self.f1_grid.all())
+        self.assertEqual(self.m1_grid["2nd"].all(), new_field.all())
+        # overwrite default field
+        newer_field = 100.0 * self.f1_grid
+        self.m1_grid.add_field(newer_field, name="3rd", default_field=True)
+        self.assertEqual(self.m1_grid.field.all(), newer_field.all())
+
+    def test_point_data_check(self):
+        self.assertRaises(ValueError, self.m1_tuple.add_field, self.f1_grid)
+        self.assertRaises(ValueError, self.m1_tuple.add_field, self.f2_grid)
+        self.assertRaises(ValueError, self.m1_tuple.add_field, self.f3_grid)
+        self.assertRaises(ValueError, self.m2_tuple.add_field, self.f1_grid)
+        self.assertRaises(ValueError, self.m2_tuple.add_field, self.f2_grid)
+        self.assertRaises(ValueError, self.m2_tuple.add_field, self.f3_grid)
+        self.assertRaises(ValueError, self.m3_tuple.add_field, self.f1_grid)
+        self.assertRaises(ValueError, self.m3_tuple.add_field, self.f2_grid)
+        self.assertRaises(ValueError, self.m3_tuple.add_field, self.f3_grid)
+
+        self.assertRaises(ValueError, self.m1_grid.add_field, self.f2_grid)
+        self.assertRaises(ValueError, self.m1_grid.add_field, self.f3_grid)
+        self.assertRaises(ValueError, self.m2_grid.add_field, self.f1_grid)
+        self.assertRaises(ValueError, self.m2_grid.add_field, self.f3_grid)
+        self.assertRaises(ValueError, self.m3_grid.add_field, self.f1_grid)
+        self.assertRaises(ValueError, self.m3_grid.add_field, self.f2_grid)
+        self.assertRaises(ValueError, self.m1_grid.add_field, self.f1_tuple)
+        self.assertRaises(ValueError, self.m1_grid.add_field, self.f2_tuple)
+        self.assertRaises(ValueError, self.m1_grid.add_field, self.f3_tuple)
+        self.assertRaises(ValueError, self.m2_grid.add_field, self.f1_tuple)
+        self.assertRaises(ValueError, self.m2_grid.add_field, self.f2_tuple)
+        self.assertRaises(ValueError, self.m2_grid.add_field, self.f3_tuple)
+        self.assertRaises(ValueError, self.m3_grid.add_field, self.f1_tuple)
+        self.assertRaises(ValueError, self.m3_grid.add_field, self.f2_tuple)
+        self.assertRaises(ValueError, self.m3_grid.add_field, self.f3_tuple)
+
+        x_tuple2 = self.rng.uniform(0.0, 10, 100)
+        y_tuple2 = self.rng.uniform(0.0, 10, 100)
+        f_tuple2 = np.vstack((x_tuple2, y_tuple2))
+        x_tuple3 = self.rng.uniform(0.0, 10, (3, 100))
+        y_tuple3 = self.rng.uniform(0.0, 10, (3, 100))
+        z_tuple3 = self.rng.uniform(0.0, 10, (3, 100))
+        f_tuple3 = np.vstack((x_tuple2, y_tuple2, z_tuple3))
+
+        m2_tuple = Mesh((x_tuple2, y_tuple2))
+        m3_tuple = Mesh((x_tuple3, y_tuple3, z_tuple3))
+
+        self.assertRaises(ValueError, m2_tuple.add_field, f_tuple3)
+        self.assertRaises(ValueError, m3_tuple.add_field, f_tuple2)
+
+        f_grid2 = np.zeros((2, len(self.x_grid), len(self.y_grid)))
+        f_grid3 = np.zeros(
+            (3, len(self.x_grid), len(self.y_grid), len(self.z_grid))
+        )
+
+        self.assertRaises(ValueError, self.m2_grid.add_field, f_grid3)
+        self.assertRaises(ValueError, self.m3_grid.add_field, f_grid2)
+
+
+if __name__ == "__main__":
+    unittest.main()

From 5e54b66b73b1f35e73e4abb5d04078becf6893a2 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 31 Mar 2020 18:33:16 +0200
Subject: [PATCH 18/49] Rename argument name for better distinction

---
 tests/test_mesh.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index 291a1470..d681ada3 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -58,7 +58,7 @@ def add_field(
         values: np.ndarray,
         name: str = "field",
         *,
-        default_field: bool = False,
+        is_default_field: bool = False,
     ) -> None:
         """Add a field (point data) to the mesh
 
@@ -76,7 +76,7 @@ def add_field(
             the point data, has to be the same shape as the mesh
         name : :class:`str`, optional
             the name of the point data
-        default_field : :class:`bool`, optional
+        is_default_field : :class:`bool`, optional
             is this the default field of the mesh?
 
         """
@@ -85,7 +85,7 @@ def add_field(
         self._check_point_data(values)
         self.point_data[name] = values
         # set the default_field to the first field added
-        if len(self.point_data) == 1 or default_field:
+        if len(self.point_data) == 1 or is_default_field:
             self._default_field = name
 
     def __getitem__(self, key: str) -> np.ndarray:
@@ -317,7 +317,7 @@ def test_add_field(self):
         self.assertEqual(self.m1_grid["2nd"].all(), new_field.all())
         # overwrite default field
         newer_field = 100.0 * self.f1_grid
-        self.m1_grid.add_field(newer_field, name="3rd", default_field=True)
+        self.m1_grid.add_field(newer_field, name="3rd", is_default_field=True)
         self.assertEqual(self.m1_grid.field.all(), newer_field.all())
 
     def test_point_data_check(self):

From 2dab1a04d7440563afbae08c80e4267c33bcea92 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 31 Mar 2020 18:34:09 +0200
Subject: [PATCH 19/49] Add getter setter

---
 tests/test_mesh.py | 9 +++++++++
 1 file changed, 9 insertions(+)

diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index d681ada3..97a9d55a 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -95,6 +95,15 @@ def __getitem__(self, key: str) -> np.ndarray:
     def __setitem__(self, key: str, value):
         self.point_data[key] = value
 
+    @property
+    def default_field(self) -> str:
+        """:class:`str`: the name of the default field."""
+        return self._default_field
+
+    @default_field.setter
+    def default_field(self, value: str):
+        self._default_field = value
+
     @property
     def pos(self) -> Tuple[np.ndarray]:
         """:any:`numpy.ndarray`: The pos. on which the field is defined."""

From ef3945fcf3ca9451ff0175dde3e977ae0c961872 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 31 Mar 2020 18:34:24 +0200
Subject: [PATCH 20/49] Delete comments

---
 tests/test_mesh.py | 4 ----
 1 file changed, 4 deletions(-)

diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index 97a9d55a..14f37df6 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -114,10 +114,6 @@ def pos(self, value: Tuple[np.ndarray]):
         """
         Warning: setting new positions deletes all previously stored fields.
         """
-        # TODO what is the best way to handle different fields?
-        # * delete all except for default field
-        # * delete them all
-        # * only delete the values, not the keys
         self.point_data = {self.default_field: None}
         self._pos = value
 

From b2e848df2a3f324aac5a13232efdcb2b0157d65a Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 31 Mar 2020 19:37:54 +0200
Subject: [PATCH 21/49] [WIP] Move Mesh class from tests to code base

Unittest tests/test_mesh.py working
---
 gstools/__init__.py       |   2 +
 gstools/field/__init__.py |   4 +-
 gstools/field/base.py     | 293 ++++++++++++++++++++++----------------
 tests/test_mesh.py        | 221 ++--------------------------
 4 files changed, 193 insertions(+), 327 deletions(-)

diff --git a/gstools/__init__.py b/gstools/__init__.py
index 5cdcf262..115e7ee6 100644
--- a/gstools/__init__.py
+++ b/gstools/__init__.py
@@ -33,6 +33,7 @@
 
 .. autosummary::
    SRF
+   Mesh
 
 Covariance Base-Class
 ^^^^^^^^^^^^^^^^^^^^^
@@ -103,6 +104,7 @@
 
 from gstools._version import __version__
 from gstools import field, variogram, random, covmodel, tools, krige, transform
+from gstools.field.base import Mesh
 from gstools.field import SRF
 from gstools.tools.export import (
     vtk_export_structured,
diff --git a/gstools/field/__init__.py b/gstools/field/__init__.py
index 6884e377..bba7f848 100644
--- a/gstools/field/__init__.py
+++ b/gstools/field/__init__.py
@@ -17,10 +17,12 @@
 
 .. autosummary::
    SRF
+   Mesh
 
 ----
 """
 
+from gstools.field.base import Mesh
 from gstools.field.srf import SRF
 
-__all__ = ["SRF"]
+__all__ = ["SRF", "Mesh"]
diff --git a/gstools/field/base.py b/gstools/field/base.py
index c58a98bf..6c74364c 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -8,12 +8,12 @@
 
 .. autosummary::
    Field
-   FieldData
+   Mesh
 """
 # pylint: disable=C0103
 
 from functools import partial
-from typing import List, Dict
+from typing import List, Dict, Tuple
 
 import numpy as np
 
@@ -21,34 +21,30 @@
 from gstools.tools.export import to_vtk, vtk_export
 from gstools.field.tools import _get_select
 
-__all__ = ["Field", "FieldData"]
+__all__ = ["Field", "Mesh"]
 
 
-class Data:
-    """A data class mainly storing the specific values of an individual field.
+def value_type(mesh_type, shape):
+    """Determine the value type ("scalar" or "vector")"""
+    r = "scalar"
+    if mesh_type == "unstructured":
+        # TODO is this the right place for doing these checks?!
+        if len(shape) == 2 and 2 <= shape[0] <= 3:
+            r = "vector"
+    else:
+        # for very small (2-3 points) meshes, this could break
+        # a 100% solution would require dim. information.
+        if len(shape) == shape[0] + 1:
+            r = "vector"
+    return r
 
-    Instances of this class are stored in a dictionary owned by FieldBase.
-    The positions on which these field values are defined are also saved in
-    FieldBase.
-    """
-
-    def __init__(
-        self,
-        values: np.ndarray = None,
-        mean: float = 0.0,
-        value_type: str = "scalar",
-    ):
-        self.values = values
-        self.mean = mean
-        self.value_type = value_type
 
-
-class FieldData:
+class Mesh:
     """A base class encapsulating field data.
 
     It holds a position array, which define the spatial locations of the
     field values.
-    It can hold multiple fields in the :any:`self.fields` dict. This assumes
+    It can hold multiple fields in the :any:`self.point_data` dict. This assumes
     that each field is defined on the same positions.
     The mesh type must also be specified.
 
@@ -60,12 +56,6 @@ class FieldData:
         key of the field values
     values : :any:`list`, optional
         a list of the values of the fields
-    mean : :any:`float`, optional
-        mean of the field
-    value_type : :any:`str`, optional
-        the value type of the default field, can be
-        * scalar
-        * vector
     mesh_type : :any:`str`, optional
         the type of mesh on which the field is defined on, can be
         * unstructured
@@ -74,137 +64,148 @@ class FieldData:
     Examples
     --------
     >>> import numpy as np
-    >>> pos = np.random.random((100, 100))
-    >>> z = np.random.random((100, 100))
-    >>> z2 = np.random.random((100, 100))
-    >>> field_data = FieldData(pos)
-    >>> field_data.add_field("test_field1", z)
-    >>> field_data.add_field("test_field2", z2)
-    >>> field_data.set_default_field("test_field2")
-    >>> print(field.values)
+    >>> from gstools import Mesh
+    >>> pos = np.linspace(0., 100., 40)
+    >>> z = np.random.random(40)
+    >>> z2 = np.random.random(40)
+    >>> mesh = Mesh((pos,))
+    >>> mesh.add_field(z, "test_field1")
+    >>> mesh.add_field(z2, "test_field2")
+    >>> mesh.set_default_field("test_field2")
+    >>> print(mesh.values)
 
     """
 
     def __init__(
         self,
-        pos: np.ndarray = None,
+        pos=None,
         name: str = "field",
         values: np.ndarray = None,
         *,
-        mean: float = 0.0,
-        value_type: str = "scalar",
         mesh_type: str = "unstructured",
-    ):
-        # initialize attributes
+    ) -> None:
+        # mesh_type needs a special setter, therefore, `set_field_data` is not
+        # used here
+        self.mesh_type = mesh_type
+
+        # the pos/ points of the mesh
         self._pos = pos
-        self.fields: Dict[str, np.ndarray] = {}
-        self._default_field = "field"
-        if mesh_type != "unstructured" and mesh_type != "structured":
-            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
-        self._mesh_type = mesh_type
-        if value_type != "scalar" and value_type != "vector":
-            raise ValueError(
-                "Unknown field value type, "
-                + "specify 'scalar' or 'vector'."
-            )
-        self.add_field(name, values, mean=mean, value_type=value_type)
+
+        # data stored at each pos/ point, the "fields"
+        if values is not None:
+            self.point_data: Dict[str, np.ndarray] = {name: values}
+        else:
+            self.point_data: Dict[str, np.ndarray] = {}
+
+        # data valid for the global field
+        self.field_data = {}
+
+        self.set_field_data("default_field", name)
+
+        self.field_data["mesh_type"] = mesh_type
+
+    def set_field_data(self, name: str, value) -> None:
+        """Add an attribute to this instance and add it the `field_data`
+
+        This helper method is used to create attributes for easy access
+        while using this class, but it also adds an entry to the dictionary
+        `field_data`, which is used for exporting the data.
+        """
+        setattr(self, name, value)
+        self.field_data[name] = value
 
     def add_field(
         self,
-        name: str,
         values: np.ndarray,
+        name: str = "field",
         *,
-        mean: float = 0.0,
-        value_type: str = "scalar",
-        default_field: bool = False,
-    ):
+        is_default_field: bool = False,
+    ) -> None:
+        """Add a field (point data) to the mesh
+
+        .. note::
+            If no field has existed before calling this method,
+            the given field will be set to the default field.
+
+        .. warning::
+            If point data with the same `name` already exists, it will be
+            overwritten.
+
+        Parameters
+        ----------
+        values : :class:`numpy.ndarray`
+            the point data, has to be the same shape as the mesh
+        name : :class:`str`, optional
+            the name of the point data
+        is_default_field : :class:`bool`, optional
+            is this the default field of the mesh?
+
+        """
         values = np.array(values)
-        self._check_field(values)
-        self.fields[name] = Data(values, mean, value_type)
+        self._check_point_data(values)
+        self.point_data[name] = values
         # set the default_field to the first field added
-        if len(self.fields) == 1 or default_field:
-            self._default_field = name
+        if len(self.point_data) == 1 or is_default_field:
+            self.set_field_data("default_field", name)
 
-    def get_data(self, key):
-        """:class:`Data`: The field data class."""
-        return self.fields[key]
-
-    def __getitem__(self, key):
+    def __getitem__(self, key: str) -> np.ndarray:
         """:any:`numpy.ndarray`: The values of the field."""
-        return self.fields[key].values
+        return self.point_data[key]
 
-    def __setitem__(self, key, value):
-        self.fields[key].values = value
+    def __setitem__(self, key: str, value):
+        self.point_data[key] = value
 
     @property
-    def pos(self):
+    def pos(self) -> Tuple[np.ndarray]:
         """:any:`numpy.ndarray`: The pos. on which the field is defined."""
         return self._pos
 
     @pos.setter
-    def pos(self, value):
+    def pos(self, value: Tuple[np.ndarray]):
         """
         Warning: setting new positions deletes all previously stored fields.
         """
+        self.point_data = {self.default_field: None}
         self._pos = value
-        self._reset()
 
     @property
-    def default_field(self):
-        """:any:`str`: The key of the default field."""
-        return self._default_field
-
-    @default_field.setter
-    def default_field(self, value):
-        self._default_field = value
+    def field(self) -> np.ndarray:
+        """:class:`numpy.ndarray`: The point data of the default field."""
+        return self.point_data[self.default_field]
 
-    @property
-    def field(self):
-        """:any:`Data`: The Data instance of the default field."""
-        return self.fields[self.default_field]
+    @field.setter
+    def field(self, values: np.ndarray):
+        self._check_point_data(values)
+        self.point_data[self.default_field] = values
 
     @property
-    def values(self):
-        """:any:`numpy.ndarray`: The values of the default field."""
-        return self.fields[self.default_field].values
-
-    @values.setter
-    def values(self, values):
-        self.fields[self.default_field].values = values
-
-    @property
-    def value_type(self):
+    def value_type(self, field="field") -> str:
         """:any:`str`: The value type of the default field."""
-        return self.fields[self.default_field].value_type
-
-    @property
-    def mean(self):
-        """:any:`float`: The mean of the default field."""
-        return self.fields[self.default_field].mean
-
-    @mean.setter
-    def mean(self, value):
-        self.fields[self.default_field].mean = value
+        if self.point_data[field] is None:
+            r = None
+        else:
+            r = value_type(self.mesh_type, self.point_data[field].shape)
+        return r
 
     @property
-    def mesh_type(self):
+    def mesh_type(self) -> str:
         """:any:`str`: The mesh type of the fields."""
         return self._mesh_type
 
     @mesh_type.setter
-    def mesh_type(self, value):
+    def mesh_type(self, value: str):
         """
         Warning: setting a new mesh type deletes all previously stored fields.
         """
+        self._check_mesh_type(value)
+        self.point_data = {}
         self._mesh_type = value
-        self._reset()
 
-    def _reset(self):
-        self._default_field = "field"
-        self._mesh_type = "unstructured"
-        self.add_field(self._default_field, None, mean=0.0, value_type="scalar")
+    def _check_mesh_type(self, mesh_type: str) -> None:
+        if mesh_type != "unstructured" and mesh_type != "structured":
+            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
 
-    def _check_field(self, values: np.ndarray):
+    def _check_point_data(self, values: np.ndarray):
         """Compare field shape to pos shape.
 
         Parameters
@@ -212,17 +213,55 @@ def _check_field(self, values: np.ndarray):
         values : :class:`numpy.ndarray`
             the values of the field to be checked
         """
-        # TODO
+        err = True
         if self.mesh_type == "unstructured":
-            pass
-        elif self.mesh_type == "structured":
-            pass
+            # scalar
+            if len(values.shape) == 1:
+                if values.shape[0] == len(self.pos[0]):
+                    err = False
+            # vector
+            elif len(values.shape) == 2:
+                if (
+                    values.shape[1] == len(self.pos[0])
+                    and 2 <= values.shape[0] <= 3
+                ):
+                    err = False
+            if err:
+                raise ValueError(
+                    "Wrong field shape: {0} does not match mesh shape ({1},)".format(
+                        values.shape, len(self.pos[0])
+                    )
+                )
         else:
-            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
+            # scalar
+            if len(values.shape) == len(self.pos):
+                if all(
+                    [
+                        values.shape[i] == len(self.pos[i])
+                        for i in range(len(self.pos))
+                    ]
+                ):
+                    err = False
+            # vector
+            elif len(values.shape) == len(self.pos) + 1:
+                if all(
+                    [
+                        values.shape[i + 1] == len(self.pos[i])
+                        for i in range(len(self.pos))
+                    ]
+                ) and values.shape[0] == len(self.pos):
+                    err = False
+            if err:
+                raise ValueError(
+                    "Wrong field shape: {0} does not match mesh shape [0/2/3]{1}".format(
+                        list(values.shape),
+                        [len(self.pos[i]) for i in range(len(self.pos))],
+                    )
+                )
 
 
-class Field(FieldData):
-    """A field base class for random and kriging fields ect.
+class Field(Mesh):
+    """A field base class for random and kriging fields, etc.
 
     Parameters
     ----------
@@ -230,9 +269,23 @@ class Field(FieldData):
         Covariance Model related to the field.
     """
 
-    def __init__(self, model, mean=0.0):
+    def __init__(
+        self,
+        model,
+        *,
+        pos=None,
+        name: str = "field",
+        values: np.ndarray = None,
+        mesh_type: str = "unstructured",
+    ) -> None:
         # initialize attributes
-        super().__init__(self, mean=mean)
+        super().__init__(
+            self,
+            pos=pos,
+            name=name,
+            values=values,
+            mesh_type=mesh_type
+        )
         # initialize private attributes
         self._model = None
         self.model = model
diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index 14f37df6..5beeb680 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -1,212 +1,12 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+"""
+This is the unittest of Mesh class.
+"""
 import unittest
-from typing import Dict, Tuple
 import numpy as np
 
-
-def value_type(mesh_type, shape):
-    """Determine the value type ("scalar" or "vector")"""
-    r = "scalar"
-    if mesh_type == "unstructured":
-        # TODO is this the right place for doing these checks?!
-        if len(shape) == 2 and 2 <= shape[0] <= 3:
-            r = "vector"
-    else:
-        # for very small (2-3 points) meshes, this could break
-        # a 100% solution would require dim. information.
-        if len(shape) == shape[0] + 1:
-            r = "vector"
-    return r
-
-
-class Mesh:
-    def __init__(
-        self,
-        pos=None,
-        name: str = "field",
-        values: np.ndarray = None,
-        *,
-        mesh_type: str = "unstructured",
-    ) -> None:
-        # the pos/ points of the mesh
-        self._pos = pos
-
-        # data stored at each pos/ point, the "fields"
-        self.point_data: Dict[str, np.ndarray] = {name: values}
-
-        # data valid for the global field
-        self.field_data = {}
-
-        self.set_field_data("default_field", name)
-
-        # mesh_type needs a special setter, therefore, `set_field_data` is not
-        # used here
-        self.mesh_type = mesh_type
-        self.field_data["mesh_type"] = mesh_type
-
-    def set_field_data(self, name: str, value) -> None:
-        """Add an attribute to this instance and add it the `field_data`
-
-        This helper method is used to create attributes for easy access
-        while using this class, but it also adds an entry to the dictionary
-        `field_data`, which is used for exporting the data.
-        """
-        setattr(self, name, value)
-        self.field_data[name] = value
-
-    def add_field(
-        self,
-        values: np.ndarray,
-        name: str = "field",
-        *,
-        is_default_field: bool = False,
-    ) -> None:
-        """Add a field (point data) to the mesh
-
-        .. note::
-            If no field has existed before calling this method,
-            the given field will be set to the default field.
-
-        .. warning::
-            If point data with the same `name` already exists, it will be
-            overwritten.
-
-        Parameters
-        ----------
-        values : :class:`numpy.ndarray`
-            the point data, has to be the same shape as the mesh
-        name : :class:`str`, optional
-            the name of the point data
-        is_default_field : :class:`bool`, optional
-            is this the default field of the mesh?
-
-        """
-
-        values = np.array(values)
-        self._check_point_data(values)
-        self.point_data[name] = values
-        # set the default_field to the first field added
-        if len(self.point_data) == 1 or is_default_field:
-            self._default_field = name
-
-    def __getitem__(self, key: str) -> np.ndarray:
-        """:any:`numpy.ndarray`: The values of the field."""
-        return self.point_data[key]
-
-    def __setitem__(self, key: str, value):
-        self.point_data[key] = value
-
-    @property
-    def default_field(self) -> str:
-        """:class:`str`: the name of the default field."""
-        return self._default_field
-
-    @default_field.setter
-    def default_field(self, value: str):
-        self._default_field = value
-
-    @property
-    def pos(self) -> Tuple[np.ndarray]:
-        """:any:`numpy.ndarray`: The pos. on which the field is defined."""
-        return self._pos
-
-    @pos.setter
-    def pos(self, value: Tuple[np.ndarray]):
-        """
-        Warning: setting new positions deletes all previously stored fields.
-        """
-        self.point_data = {self.default_field: None}
-        self._pos = value
-
-    @property
-    def field(self) -> np.ndarray:
-        """:class:`numpy.ndarray`: The point data of the default field."""
-        return self.point_data[self.default_field]
-
-    @field.setter
-    def field(self, values: np.ndarray):
-        self._check_point_data(values)
-        self.point_data[self.default_field] = values
-
-    @property
-    def value_type(self, field="field") -> str:
-        """:any:`str`: The value type of the default field."""
-        if self.point_data[field] is None:
-            r = None
-        else:
-            r = value_type(self.mesh_type, self.point_data[field].shape)
-        return r
-
-    @property
-    def mesh_type(self) -> str:
-        """:any:`str`: The mesh type of the fields."""
-        return self._mesh_type
-
-    @mesh_type.setter
-    def mesh_type(self, value: str):
-        """
-        Warning: setting a new mesh type deletes all previously stored fields.
-        """
-        self._check_mesh_type(value)
-        self.point_data = {self.default_field: None}
-        self._mesh_type = value
-
-    def _check_mesh_type(self, mesh_type: str) -> None:
-        if mesh_type != "unstructured" and mesh_type != "structured":
-            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
-
-    def _check_point_data(self, values: np.ndarray):
-        """Compare field shape to pos shape.
-
-        Parameters
-        ----------
-        values : :class:`numpy.ndarray`
-            the values of the field to be checked
-        """
-        err = True
-        if self.mesh_type == "unstructured":
-            # scalar
-            if len(values.shape) == 1:
-                if values.shape[0] == len(self.pos[0]):
-                    err = False
-            # vector
-            elif len(values.shape) == 2:
-                if (
-                    values.shape[1] == len(self.pos[0])
-                    and 2 <= values.shape[0] <= 3
-                ):
-                    err = False
-            if err:
-                raise ValueError(
-                    "Wrong field shape: {0} does not match mesh shape ({1},)".format(
-                        values.shape, len(self.pos[0])
-                    )
-                )
-        else:
-            # scalar
-            if len(values.shape) == len(self.pos):
-                if all(
-                    [
-                        values.shape[i] == len(self.pos[i])
-                        for i in range(len(self.pos))
-                    ]
-                ):
-                    err = False
-            # vector
-            elif len(values.shape) == len(self.pos) + 1:
-                if all(
-                    [
-                        values.shape[i + 1] == len(self.pos[i])
-                        for i in range(len(self.pos))
-                    ]
-                ) and values.shape[0] == len(self.pos):
-                    err = False
-            if err:
-                raise ValueError(
-                    "Wrong field shape: {0} does not match mesh shape [0/2/3]{1}".format(
-                        list(values.shape),
-                        [len(self.pos[i]) for i in range(len(self.pos))],
-                    )
-                )
+from gstools import Mesh
 
 
 class TestMesh(unittest.TestCase):
@@ -325,6 +125,15 @@ def test_add_field(self):
         self.m1_grid.add_field(newer_field, name="3rd", is_default_field=True)
         self.assertEqual(self.m1_grid.field.all(), newer_field.all())
 
+    def test_default_field(self):
+        # first field added, should be default_field
+        f2_grid = 5.0 * self.f1_grid
+        self.m1_grid.add_field(self.f1_grid, name="test_field1")
+        self.m1_grid.add_field(f2_grid, name="test_field2")
+        self.assertEqual(self.m1_grid.default_field, "test_field1")
+        self.m1_grid.default_field = "test_field2"
+        self.assertEqual(self.m1_grid.default_field, "test_field2")
+
     def test_point_data_check(self):
         self.assertRaises(ValueError, self.m1_tuple.add_field, self.f1_grid)
         self.assertRaises(ValueError, self.m1_tuple.add_field, self.f2_grid)

From 72cbeea4d211aae11b75b8eab273bc05985561b4 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 31 Mar 2020 19:45:07 +0200
Subject: [PATCH 22/49] [WIP] Fix docstring

---
 gstools/field/base.py | 2 +-
 tests/test_mesh.py    | 1 +
 2 files changed, 2 insertions(+), 1 deletion(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 6c74364c..c3e805f1 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -72,7 +72,7 @@ class Mesh:
     >>> mesh.add_field(z, "test_field1")
     >>> mesh.add_field(z2, "test_field2")
     >>> mesh.set_default_field("test_field2")
-    >>> print(mesh.values)
+    >>> print(mesh.field)
 
     """
 
diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index 5beeb680..be160ba7 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -133,6 +133,7 @@ def test_default_field(self):
         self.assertEqual(self.m1_grid.default_field, "test_field1")
         self.m1_grid.default_field = "test_field2"
         self.assertEqual(self.m1_grid.default_field, "test_field2")
+        self.assertEqual(self.m1_grid.field[5], self.m1_grid["test_field2"][5])
 
     def test_point_data_check(self):
         self.assertRaises(ValueError, self.m1_tuple.add_field, self.f1_grid)

From 453c9e099f26b7ec0cca6919744c72003c75970f Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 31 Mar 2020 19:45:43 +0200
Subject: [PATCH 23/49] [WIP] Fix bug

---
 gstools/field/base.py | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index c3e805f1..111b2d7b 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -181,10 +181,10 @@ def field(self, values: np.ndarray):
     @property
     def value_type(self, field="field") -> str:
         """:any:`str`: The value type of the default field."""
-        if self.point_data[field] is None:
-            r = None
-        else:
+        if field in self.point_data:
             r = value_type(self.mesh_type, self.point_data[field].shape)
+        else:
+            r = None
         return r
 
     @property

From 49d8bb3cb14789ad8c4fc256795005db84267e6a Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 3 Apr 2020 19:26:37 +0200
Subject: [PATCH 24/49] [WIP] Fix bug

---
 gstools/field/base.py | 1 -
 gstools/field/srf.py  | 2 +-
 2 files changed, 1 insertion(+), 2 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index b2235ca2..538d5e5f 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -288,7 +288,6 @@ def __init__(
     ) -> None:
         # initialize attributes
         super().__init__(
-            self,
             pos=pos,
             name=name,
             values=values,
diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index b00b41ee..33a89439 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -82,7 +82,7 @@ def __init__(
         generator="RandMeth",
         **generator_kwargs
     ):
-        super().__init__(model, mean)
+        super().__init__(model, mean=mean)
         # initialize private attributes
         self._generator = None
         self._upscaling = None

From ccfab457c4da471351c3be2f0981534dc5d89b52 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 7 Apr 2020 12:40:10 +0200
Subject: [PATCH 25/49] [WIP] Update Krige class

---
 gstools/krige/base.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/gstools/krige/base.py b/gstools/krige/base.py
index 366d02cd..99c38a64 100755
--- a/gstools/krige/base.py
+++ b/gstools/krige/base.py
@@ -69,7 +69,7 @@ def __init__(
         drift_functions=None,
         trend_function=None,
     ):
-        super().__init__(model, mean)
+        super().__init__(model, mean=mean)
         self.krige_var = None
         # initialize private attributes
         self._unbiased = True

From 074012d74fd1c9985a9fa98926a21f46f1de5e86 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 7 Apr 2020 12:44:38 +0200
Subject: [PATCH 26/49] [WIP] Update SRF class

---
 gstools/field/srf.py | 81 +++++++++++++++-----------------------------
 1 file changed, 27 insertions(+), 54 deletions(-)

diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index 33a89439..8c5a9a17 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -98,27 +98,24 @@ def __init__(
         self.krige_var = None
         self.set_generator(generator, **generator_kwargs)
         self.upscaling = upscaling
+        if self._value_type is None:
+            raise ValueError(
+                "Unknown field value type, "
+                + "specify 'scalar' or 'vector' before calling SRF."
+            )
 
     def __call__(
-        self,
-        pos,
-        name="field",
-        seed=np.nan,
-        point_volumes=0.0,
-        mesh_type="unstructured",
+        self, pos, seed=np.nan, point_volumes=0.0, mesh_type="unstructured"
     ):
         """Generate the spatial random field.
 
+        The field is saved as `self.field` and is also returned.
+
         Parameters
         ----------
         pos : :class:`list`
             the position tuple, containing main direction and transversal
             directions
-        name : :class:`str`, optional
-            the name of the dictionary key, in which the values will be saved.
-            If this method is called multiple times, the name determines, if
-            the field is overwritten or if it is added to the values dict.
-            Default: "field"
         seed : :class:`int`, optional
             seed for RNG for reseting. Default: keep seed from generator
         point_volumes : :class:`float` or :class:`numpy.ndarray`
@@ -135,42 +132,20 @@ def __call__(
         field : :class:`numpy.ndarray`
             the SRF
         """
-        # internal conversation
-        x, y, z = pos2xyz(pos, max_dim=self.model.dim)
-        mean_temp = self.mean
-        self.pos = xyz2pos(x, y, z)
         self.mesh_type = mesh_type
-        self.mean = mean_temp
         # update the model/seed in the generator if any changes were made
         self.generator.update(self.model, seed)
-        # format the positional arguments of the mesh
-        check_mesh(self.model.dim, x, y, z, mesh_type)
-        mesh_type_changed = False
-        if self.model.do_rotation:
-            if mesh_type == "structured":
-                mesh_type_changed = True
-                mesh_type_old = mesh_type
-                mesh_type = "unstructured"
-                x, y, z, axis_lens = reshape_axis_from_struct_to_unstruct(
-                    self.model.dim, x, y, z
-                )
-            x, y, z = unrotate_mesh(self.model.dim, self.model.angles, x, y, z)
-        y, z = make_isotropic(self.model.dim, self.model.anis, y, z)
-
+        # internal conversation
+        x, y, z, self.pos, mt_gen, mt_changed, axis_lens = self._pre_pos(
+            pos, mesh_type
+        )
         # generate the field
-        self.add_field(
-                "raw",
-                self.generator.__call__(x, y, z, mesh_type),
-                mean=self.mean
-            )
-
+        self.raw_field = self.generator.__call__(x, y, z, mt_gen)
         # reshape field if we got an unstructured mesh
-        if mesh_type_changed:
-            mesh_type = mesh_type_old
-            self["raw"] = reshape_field_from_unstruct_to_struct(
-                self.model.dim, self["raw"], axis_lens
+        if mt_changed:
+            self.raw_field = reshape_field_from_unstruct_to_struct(
+                self.model.dim, self.raw_field, axis_lens
             )
-
         # apply given conditions to the field
         if self.condition:
             (
@@ -181,23 +156,21 @@ def __call__(
                 info,
             ) = self.cond_func(self)
             # store everything in the class
-            self.add_field(name, cond_field, default_field=True)
-            self.add_field("krige", krige_field)
-            self.add_field("error", err_field)
-            self.get_data("krige").krige_var = krigevar
-            if "mean" in info:  # ordinary kriging estimates mean
-                self.get_data("cond").mean = info["mean"]
+            self.field = cond_field
+            self.krige_field = krige_field
+            self.err_field = err_field
+            self.krige_var = krigevar
+            if "mean" in info:  # ordinary krging estimates mean
+                self.mean = info["mean"]
         else:
-            self.add_field(name, self["raw"] + self.mean)
-
+            self.field = self.raw_field + self.mean
         # upscaled variance
         if not np.isscalar(point_volumes) or not np.isclose(point_volumes, 0):
             scaled_var = self.upscaling_func(self.model, point_volumes)
-            self.values -= self.mean
-            self.values *= np.sqrt(scaled_var / self.model.sill)
-            self.values += self.mean
-
-        return self.values
+            self.field -= self.mean
+            self.field *= np.sqrt(scaled_var / self.model.sill)
+            self.field += self.mean
+        return self.field
 
     def set_condition(
         self, cond_pos=None, cond_val=None, krige_type="ordinary"

From 4570a42dd264d38d25aa9f21e203f40722e8e6be Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 7 Apr 2020 12:49:34 +0200
Subject: [PATCH 27/49] [WIP] Update export functions

---
 gstools/tools/export.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/gstools/tools/export.py b/gstools/tools/export.py
index f9b79e19..9174830b 100644
--- a/gstools/tools/export.py
+++ b/gstools/tools/export.py
@@ -51,7 +51,7 @@ def _vtk_structured_helper(pos, fields):
         z = np.array([0])
     # need fortran order in VTK
     for field in fields:
-        fields[field] = fields[field].values.reshape(-1, order="F")
+        fields[field] = fields[field].reshape(-1, order="F")
         if len(fields[field]) != len(x) * len(y) * len(z):
             raise ValueError(
                 "gstools.vtk_export_structured: "
@@ -119,7 +119,7 @@ def _vtk_unstructured_helper(pos, fields):
     if z is None:
         z = np.zeros_like(x)
     for field in fields:
-        fields[field] = fields[field].values.reshape(-1)
+        fields[field] = fields[field].reshape(-1)
         if (
             len(fields[field]) != len(x)
             or len(fields[field]) != len(y)

From fafa5d837d88f19a1700e7e081862521235b9ea1 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 7 Apr 2020 14:47:55 +0200
Subject: [PATCH 28/49] Remove type hints for Py3.5 support

---
 gstools/field/base.py | 61 +++++++++++++++++--------------------------
 1 file changed, 24 insertions(+), 37 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 538d5e5f..d26eadfa 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -13,7 +13,6 @@
 # pylint: disable=C0103
 
 from functools import partial
-from typing import List, Dict, Tuple
 
 import numpy as np
 
@@ -84,13 +83,8 @@ class Mesh:
     """
 
     def __init__(
-        self,
-        pos=None,
-        name: str = "field",
-        values: np.ndarray = None,
-        *,
-        mesh_type: str = "unstructured",
-    ) -> None:
+        self, pos=None, name="field", values=None, *, mesh_type="unstructured",
+    ):
         # mesh_type needs a special setter, therefore, `set_field_data` is not
         # used here
         self.mesh_type = mesh_type
@@ -100,9 +94,9 @@ def __init__(
 
         # data stored at each pos/ point, the "fields"
         if values is not None:
-            self.point_data: Dict[str, np.ndarray] = {name: values}
+            self.point_data = {name: values}
         else:
-            self.point_data: Dict[str, np.ndarray] = {}
+            self.point_data = {}
 
         # data valid for the global field
         self.field_data = {}
@@ -111,7 +105,7 @@ def __init__(
 
         self.field_data["mesh_type"] = mesh_type
 
-    def set_field_data(self, name: str, value) -> None:
+    def set_field_data(self, name, value):
         """Add an attribute to this instance and add it the `field_data`
 
         This helper method is used to create attributes for easy access
@@ -122,12 +116,8 @@ def set_field_data(self, name: str, value) -> None:
         self.field_data[name] = value
 
     def add_field(
-        self,
-        values: np.ndarray,
-        name: str = "field",
-        *,
-        is_default_field: bool = False,
-    ) -> None:
+        self, values, name="field", *, is_default_field=False,
+    ):
         """Add a field (point data) to the mesh
 
         .. note::
@@ -155,20 +145,20 @@ def add_field(
         if len(self.point_data) == 1 or is_default_field:
             self.set_field_data("default_field", name)
 
-    def __getitem__(self, key: str) -> np.ndarray:
+    def __getitem__(self, key):
         """:any:`numpy.ndarray`: The values of the field."""
         return self.point_data[key]
 
-    def __setitem__(self, key: str, value):
+    def __setitem__(self, key, value):
         self.point_data[key] = value
 
     @property
-    def pos(self) -> Tuple[np.ndarray]:
+    def pos(self):
         """:any:`numpy.ndarray`: The pos. on which the field is defined."""
         return self._pos
 
     @pos.setter
-    def pos(self, value: Tuple[np.ndarray]):
+    def pos(self, value):
         """
         Warning: setting new positions deletes all previously stored fields.
         """
@@ -176,17 +166,17 @@ def pos(self, value: Tuple[np.ndarray]):
         self._pos = value
 
     @property
-    def field(self) -> np.ndarray:
+    def field(self):
         """:class:`numpy.ndarray`: The point data of the default field."""
         return self.point_data[self.default_field]
 
     @field.setter
-    def field(self, values: np.ndarray):
+    def field(self, values):
         self._check_point_data(values)
         self.point_data[self.default_field] = values
 
     @property
-    def value_type(self, field="field") -> str:
+    def value_type(self, field="field"):
         """:any:`str`: The value type of the default field."""
         if field in self.point_data:
             r = value_type(self.mesh_type, self.point_data[field].shape)
@@ -195,12 +185,12 @@ def value_type(self, field="field") -> str:
         return r
 
     @property
-    def mesh_type(self) -> str:
+    def mesh_type(self):
         """:any:`str`: The mesh type of the fields."""
         return self._mesh_type
 
     @mesh_type.setter
-    def mesh_type(self, value: str):
+    def mesh_type(self, value):
         """
         Warning: setting a new mesh type deletes all previously stored fields.
         """
@@ -208,11 +198,11 @@ def mesh_type(self, value: str):
         self.point_data = {}
         self._mesh_type = value
 
-    def _check_mesh_type(self, mesh_type: str) -> None:
+    def _check_mesh_type(self, mesh_type):
         if mesh_type != "unstructured" and mesh_type != "structured":
             raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
 
-    def _check_point_data(self, values: np.ndarray):
+    def _check_point_data(self, values):
         """Compare field shape to pos shape.
 
         Parameters
@@ -281,17 +271,14 @@ def __init__(
         model,
         *,
         pos=None,
-        name: str = "field",
-        values: np.ndarray = None,
-        mesh_type: str = "unstructured",
-        mean: float = 0.0,
-    ) -> None:
+        name="field",
+        values=None,
+        mesh_type="unstructured",
+        mean=0.0,
+    ):
         # initialize attributes
         super().__init__(
-            pos=pos,
-            name=name,
-            values=values,
-            mesh_type=mesh_type
+            pos=pos, name=name, values=values, mesh_type=mesh_type
         )
         # initialize private attributes
         self._model = None

From 557fe8ac6690626d86702bbadb8d9c1e89992a21 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 7 Apr 2020 14:59:21 +0200
Subject: [PATCH 29/49] Remove positional only arguments due to Py3.5

---
 gstools/field/base.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index d26eadfa..0ba59d56 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -83,7 +83,7 @@ class Mesh:
     """
 
     def __init__(
-        self, pos=None, name="field", values=None, *, mesh_type="unstructured",
+        self, pos=None, name="field", values=None, mesh_type="unstructured",
     ):
         # mesh_type needs a special setter, therefore, `set_field_data` is not
         # used here
@@ -116,7 +116,7 @@ def set_field_data(self, name, value):
         self.field_data[name] = value
 
     def add_field(
-        self, values, name="field", *, is_default_field=False,
+        self, values, name="field", is_default_field=False,
     ):
         """Add a field (point data) to the mesh
 

From 3ac633e5a7f320922c357570e48b5ecc9fada109 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Tue, 7 Apr 2020 16:01:24 +0200
Subject: [PATCH 30/49] Amend

---
 gstools/field/base.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 0ba59d56..dab0b71d 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -269,7 +269,6 @@ class Field(Mesh):
     def __init__(
         self,
         model,
-        *,
         pos=None,
         name="field",
         values=None,

From 27603ec52b144cb34829e5e45556533db337dcaa Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 17 Apr 2020 18:29:16 +0200
Subject: [PATCH 31/49] Rename helper fcts and add deprecation warnings

---
 gstools/tools/export.py | 29 +++++++++++++++++++++++------
 1 file changed, 23 insertions(+), 6 deletions(-)

diff --git a/gstools/tools/export.py b/gstools/tools/export.py
index 9174830b..3d69e5d3 100644
--- a/gstools/tools/export.py
+++ b/gstools/tools/export.py
@@ -16,6 +16,7 @@
 """
 # pylint: disable=C0103, E1101
 
+import warnings
 import numpy as np
 from pyevtk.hl import gridToVTK, pointsToVTK
 
@@ -39,7 +40,7 @@
 # export routines #############################################################
 
 
-def _vtk_structured_helper(pos, fields):
+def _vtk_structured_reshape(pos, fields):
     """An internal helper to extract what is needed for the vtk rectilinear grid
     """
     if not isinstance(fields, dict):
@@ -79,7 +80,11 @@ def to_vtk_structured(pos, fields):  # pragma: no cover
         A PyVista rectilinear grid of the structured field data. Data arrays
         live on the point data of this PyVista dataset.
     """
-    x, y, z, fields = _vtk_structured_helper(pos=pos, fields=fields)
+    warnings.warn(
+        "Export routines are deprecated, use Mesh export methods instead",
+        DeprecationWarning
+    )
+    x, y, z, fields = _vtk_structured_reshape(pos=pos, fields=fields)
     try:
         import pyvista as pv
 
@@ -106,11 +111,15 @@ def vtk_export_structured(filename, pos, fields):  # pragma: no cover
         Either a single numpy array as returned by SRF,
         or a dictionary of fields with theirs names as keys.
     """
-    x, y, z, fields = _vtk_structured_helper(pos=pos, fields=fields)
+    warnings.warn(
+        "Export routines are deprecated, use Mesh export methods instead",
+        DeprecationWarning
+    )
+    x, y, z, fields = _vtk_structured_reshape(pos=pos, fields=fields)
     return gridToVTK(filename, x, y, z, pointData=fields)
 
 
-def _vtk_unstructured_helper(pos, fields):
+def _vtk_unstructured_reshape(pos, fields):
     if not isinstance(fields, dict):
         fields = {"field": fields}
     x, y, z = pos2xyz(pos)
@@ -152,7 +161,11 @@ def to_vtk_unstructured(pos, fields):  # pragma: no cover
         live on the point data of this PyVista dataset. This is essentially
         a point cloud with no topology.
     """
-    x, y, z, fields = _vtk_unstructured_helper(pos=pos, fields=fields)
+    warnings.warn(
+        "Export routines are deprecated, use Mesh export methods instead",
+        DeprecationWarning
+    )
+    x, y, z, fields = _vtk_unstructured_reshape(pos=pos, fields=fields)
     try:
         import pyvista as pv
 
@@ -179,7 +192,11 @@ def vtk_export_unstructured(filename, pos, fields):  # pragma: no cover
         Either a single numpy array as returned by SRF,
         or a dictionary of fields with theirs names as keys.
     """
-    x, y, z, fields = _vtk_unstructured_helper(pos=pos, fields=fields)
+    warnings.warn(
+        "Export routines are deprecated, use Mesh export methods instead",
+        DeprecationWarning
+    )
+    x, y, z, fields = _vtk_unstructured_reshape(pos=pos, fields=fields)
     return pointsToVTK(filename, x, y, z, data=fields)
 
 

From d272b0fd579aa4725008b3921d9f8523c16061be Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 17 Apr 2020 18:30:02 +0200
Subject: [PATCH 32/49] Move export methods from Field to Mesh class

---
 gstools/field/base.py | 220 +++++++++++++++++++++---------------------
 1 file changed, 110 insertions(+), 110 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index dab0b71d..a4f0c6f3 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -198,6 +198,116 @@ def mesh_type(self, value):
         self.point_data = {}
         self._mesh_type = value
 
+    def _to_vtk_helper(
+        self, filename=None, field_select="field", fieldname="field"
+    ):  # pragma: no cover
+        """Create a VTK/PyVista grid of the field or save it as a VTK file.
+
+        This is an internal helper that will handle saving or creating objects
+
+        Parameters
+        ----------
+        filename : :class:`str`
+            Filename of the file to be saved, including the path. Note that an
+            ending (.vtr or .vtu) will be added to the name. If ``None`` is
+            passed, a PyVista dataset of the appropriate type will be returned.
+        field_select : :class:`str`, optional
+            Field that should be stored. Can be:
+            "field", "raw_field", "krige_field", "err_field" or "krige_var".
+            Default: "field"
+        fieldname : :class:`str`, optional
+            Name of the field in the VTK file. Default: "field"
+        """
+        if self.value_type is None:
+            raise ValueError(
+                "Field value type not set! "
+                + "Specify 'scalar' or 'vector' before plotting."
+            )
+        elif self.value_type == "vector":
+            if hasattr(self, field_select):
+                field = getattr(self, field_select)
+            else:
+                field = None
+            if not (
+                self.pos is None or field is None or self.mesh_type is None
+            ):
+                suf = ["_X", "_Y", "_Z"]
+                fields = {}
+                for i in range(self.model.dim):
+                    fields[fieldname + suf[i]] = field[i]
+                if filename is None:
+                    return to_vtk(self.pos, fields, self.mesh_type)
+                else:
+                    return vtk_export(
+                        filename, self.pos, fields, self.mesh_type
+                    )
+        elif self.value_type == "scalar":
+            if hasattr(self, field_select):
+                field = getattr(self, field_select)
+            else:
+                field = None
+            if not (
+                self.pos is None or field is None or self.mesh_type is None
+            ):
+                if filename is None:
+                    return to_vtk(self.pos, {fieldname: field}, self.mesh_type)
+                else:
+                    return vtk_export(
+                        filename, self.pos, {fieldname: field}, self.mesh_type
+                    )
+            else:
+                print(
+                    "Field.to_vtk: No "
+                    + field_select
+                    + " stored in the class."
+                )
+        else:
+            raise ValueError(
+                "Unknown field value type: {}".format(self.value_type)
+            )
+
+    def to_pyvista(
+        self, field_select="field", fieldname="field"
+    ):  # pragma: no cover
+        """Create a VTK/PyVista grid of the stored field.
+
+        Parameters
+        ----------
+        field_select : :class:`str`, optional
+            Field that should be stored. Can be:
+            "field", "raw_field", "krige_field", "err_field" or "krige_var".
+            Default: "field"
+        fieldname : :class:`str`, optional
+            Name of the field in the VTK file. Default: "field"
+        """
+        grid = self._to_vtk_helper(
+            filename=None, field_select=field_select, fieldname=fieldname
+        )
+        return grid
+
+    def vtk_export(
+        self, filename, field_select="field", fieldname="field"
+    ):  # pragma: no cover
+        """Export the stored field to vtk.
+
+        Parameters
+        ----------
+        filename : :class:`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 : :class:`str`, optional
+            Field that should be stored. Can be:
+            "field", "raw_field", "krige_field", "err_field" or "krige_var".
+            Default: "field"
+        fieldname : :class:`str`, optional
+            Name of the field in the VTK file. Default: "field"
+        """
+        if not isinstance(filename, str):
+            raise TypeError("Please use a string filename.")
+        return self._to_vtk_helper(
+            filename=filename, field_select=field_select, fieldname=fieldname
+        )
+
     def _check_mesh_type(self, mesh_type):
         if mesh_type != "unstructured" and mesh_type != "structured":
             raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
@@ -438,116 +548,6 @@ def _pre_pos(self, pos, mesh_type="unstructured", make_unstruct=False):
             y, z = make_isotropic(self.model.dim, self.model.anis, y, z)
         return x, y, z, pos, mesh_type_gen, mesh_type_changed, axis_lens
 
-    def _to_vtk_helper(
-        self, filename=None, field_select="field", fieldname="field"
-    ):  # pragma: no cover
-        """Create a VTK/PyVista grid of the field or save it as a VTK file.
-
-        This is an internal helper that will handle saving or creating objects
-
-        Parameters
-        ----------
-        filename : :class:`str`
-            Filename of the file to be saved, including the path. Note that an
-            ending (.vtr or .vtu) will be added to the name. If ``None`` is
-            passed, a PyVista dataset of the appropriate type will be returned.
-        field_select : :class:`str`, optional
-            Field that should be stored. Can be:
-            "field", "raw_field", "krige_field", "err_field" or "krige_var".
-            Default: "field"
-        fieldname : :class:`str`, optional
-            Name of the field in the VTK file. Default: "field"
-        """
-        if self.value_type is None:
-            raise ValueError(
-                "Field value type not set! "
-                + "Specify 'scalar' or 'vector' before plotting."
-            )
-        elif self.value_type == "vector":
-            if hasattr(self, field_select):
-                field = getattr(self, field_select)
-            else:
-                field = None
-            if not (
-                self.pos is None or field is None or self.mesh_type is None
-            ):
-                suf = ["_X", "_Y", "_Z"]
-                fields = {}
-                for i in range(self.model.dim):
-                    fields[fieldname + suf[i]] = field[i]
-                if filename is None:
-                    return to_vtk(self.pos, fields, self.mesh_type)
-                else:
-                    return vtk_export(
-                        filename, self.pos, fields, self.mesh_type
-                    )
-        elif self.value_type == "scalar":
-            if hasattr(self, field_select):
-                field = getattr(self, field_select)
-            else:
-                field = None
-            if not (
-                self.pos is None or field is None or self.mesh_type is None
-            ):
-                if filename is None:
-                    return to_vtk(self.pos, {fieldname: field}, self.mesh_type)
-                else:
-                    return vtk_export(
-                        filename, self.pos, {fieldname: field}, self.mesh_type
-                    )
-            else:
-                print(
-                    "Field.to_vtk: No "
-                    + field_select
-                    + " stored in the class."
-                )
-        else:
-            raise ValueError(
-                "Unknown field value type: {}".format(self.value_type)
-            )
-
-    def to_pyvista(
-        self, field_select="field", fieldname="field"
-    ):  # pragma: no cover
-        """Create a VTK/PyVista grid of the stored field.
-
-        Parameters
-        ----------
-        field_select : :class:`str`, optional
-            Field that should be stored. Can be:
-            "field", "raw_field", "krige_field", "err_field" or "krige_var".
-            Default: "field"
-        fieldname : :class:`str`, optional
-            Name of the field in the VTK file. Default: "field"
-        """
-        grid = self._to_vtk_helper(
-            filename=None, field_select=field_select, fieldname=fieldname
-        )
-        return grid
-
-    def vtk_export(
-        self, filename, field_select="field", fieldname="field"
-    ):  # pragma: no cover
-        """Export the stored field to vtk.
-
-        Parameters
-        ----------
-        filename : :class:`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 : :class:`str`, optional
-            Field that should be stored. Can be:
-            "field", "raw_field", "krige_field", "err_field" or "krige_var".
-            Default: "field"
-        fieldname : :class:`str`, optional
-            Name of the field in the VTK file. Default: "field"
-        """
-        if not isinstance(filename, str):
-            raise TypeError("Please use a string filename.")
-        return self._to_vtk_helper(
-            filename=filename, field_select=field_select, fieldname=fieldname
-        )
-
     def plot(self, field="field", fig=None, ax=None):  # pragma: no cover
         """
         Plot the spatial random field.

From 870d892771a91ba1586e4a9d5c5b01ded4e26dfa Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 17 Apr 2020 18:56:42 +0200
Subject: [PATCH 33/49] Refactor Mesh export methods

---
 gstools/field/base.py | 44 ++++++++++++++++++-------------------------
 1 file changed, 18 insertions(+), 26 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index a4f0c6f3..6198166f 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -198,25 +198,24 @@ def mesh_type(self, value):
         self.point_data = {}
         self._mesh_type = value
 
-    def _to_vtk_helper(
-        self, filename=None, field_select="field", fieldname="field"
+    def _vtk_naming_helper(
+        self, field_select="field", fieldname="field"
     ):  # pragma: no cover
-        """Create a VTK/PyVista grid of the field or save it as a VTK file.
-
-        This is an internal helper that will handle saving or creating objects
+        """Prepare the field names for export.
 
         Parameters
         ----------
-        filename : :class:`str`
-            Filename of the file to be saved, including the path. Note that an
-            ending (.vtr or .vtu) will be added to the name. If ``None`` is
-            passed, a PyVista dataset of the appropriate type will be returned.
         field_select : :class:`str`, optional
             Field that should be stored. Can be:
             "field", "raw_field", "krige_field", "err_field" or "krige_var".
             Default: "field"
         fieldname : :class:`str`, optional
             Name of the field in the VTK file. Default: "field"
+
+        Returns
+        -------
+        fields : :class:`dict`
+            a dictionary containing the fields to be exported
         """
         if self.value_type is None:
             raise ValueError(
@@ -235,12 +234,7 @@ def _to_vtk_helper(
                 fields = {}
                 for i in range(self.model.dim):
                     fields[fieldname + suf[i]] = field[i]
-                if filename is None:
-                    return to_vtk(self.pos, fields, self.mesh_type)
-                else:
-                    return vtk_export(
-                        filename, self.pos, fields, self.mesh_type
-                    )
+                return fields
         elif self.value_type == "scalar":
             if hasattr(self, field_select):
                 field = getattr(self, field_select)
@@ -249,12 +243,7 @@ def _to_vtk_helper(
             if not (
                 self.pos is None or field is None or self.mesh_type is None
             ):
-                if filename is None:
-                    return to_vtk(self.pos, {fieldname: field}, self.mesh_type)
-                else:
-                    return vtk_export(
-                        filename, self.pos, {fieldname: field}, self.mesh_type
-                    )
+                return {fieldname: field}
             else:
                 print(
                     "Field.to_vtk: No "
@@ -280,9 +269,11 @@ def to_pyvista(
         fieldname : :class:`str`, optional
             Name of the field in the VTK file. Default: "field"
         """
-        grid = self._to_vtk_helper(
-            filename=None, field_select=field_select, fieldname=fieldname
+        field_names = self._vtk_naming_helper(
+            field_select=field_select, fieldname=fieldname
         )
+
+        grid = to_vtk(self.pos, field_names, self.mesh_type)
         return grid
 
     def vtk_export(
@@ -303,10 +294,11 @@ def vtk_export(
             Name of the field in the VTK file. Default: "field"
         """
         if not isinstance(filename, str):
-            raise TypeError("Please use a string filename.")
-        return self._to_vtk_helper(
-            filename=filename, field_select=field_select, fieldname=fieldname
+            raise TypeError("Please use a string as a filename.")
+        field_names = self._vtk_naming_helper(
+            field_select=field_select, fieldname=fieldname
         )
+        return vtk_export(filename, self.pos, field_names, self.mesh_type)
 
     def _check_mesh_type(self, mesh_type):
         if mesh_type != "unstructured" and mesh_type != "structured":

From 43afd7aef3727946566c4dd8f9691fcf8daaf151 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Sat, 18 Apr 2020 10:38:32 +0200
Subject: [PATCH 34/49] Move export tool fcts to Mesh

---
 gstools/field/base.py | 144 +++++++++++++++++++++++++++++++++++++++++-
 1 file changed, 141 insertions(+), 3 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 6198166f..9b317e5b 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -16,8 +16,9 @@
 
 import numpy as np
 
+from pyevtk.hl import gridToVTK, pointsToVTK
+
 from gstools.covmodel.base import CovModel
-from gstools.tools.export import to_vtk, vtk_export
 from gstools.field.tools import (
     _get_select,
     check_mesh,
@@ -85,6 +86,8 @@ class Mesh:
     def __init__(
         self, pos=None, name="field", values=None, mesh_type="unstructured",
     ):
+        self._vtk_export_fct = self._vtk_export_unstructured
+        self._to_vtk_fct = self._to_vtk_unstructured
         # mesh_type needs a special setter, therefore, `set_field_data` is not
         # used here
         self.mesh_type = mesh_type
@@ -195,6 +198,12 @@ def mesh_type(self, value):
         Warning: setting a new mesh type deletes all previously stored fields.
         """
         self._check_mesh_type(value)
+        if value == "structured":
+            self._vtk_export_fct = self._vtk_export_structured
+            self._to_vtk_fct = self._to_vtk_structured
+        else:
+            self._vtk_export_fct = self._vtk_export_unstructured
+            self._to_vtk_fct = self._to_vtk_unstructured
         self.point_data = {}
         self._mesh_type = value
 
@@ -273,7 +282,7 @@ def to_pyvista(
             field_select=field_select, fieldname=fieldname
         )
 
-        grid = to_vtk(self.pos, field_names, self.mesh_type)
+        grid = self._to_vtk_fct(field_names)
         return grid
 
     def vtk_export(
@@ -298,7 +307,136 @@ def vtk_export(
         field_names = self._vtk_naming_helper(
             field_select=field_select, fieldname=fieldname
         )
-        return vtk_export(filename, self.pos, field_names, self.mesh_type)
+        return self._vtk_export_fct(filename, field_names)
+
+    def _vtk_export_structured(self, filename, fields):
+        """Export a field to vtk structured rectilinear grid file.
+
+        Parameters
+        ----------
+        filename : :class:`str`
+            Filename of the file to be saved, including the path. Note that an
+            ending (.vtr) will be added to the name.
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Structured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+        """
+        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
+        return gridToVTK(filename, x, y, z, pointData=fields)
+
+    def _vtk_export_unstructured(self, filename, fields):
+        """Export a field to vtk unstructured grid file.
+
+        Parameters
+        ----------
+        filename : :class:`str`
+            Filename of the file to be saved, including the path. Note that an
+            ending (.vtu) will be added to the name.
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Unstructured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+        """
+        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
+        return pointsToVTK(filename, x, y, z, data=fields)
+
+    def _to_vtk_structured(self, fields):  # pragma: no cover
+        """Create a vtk structured rectilinear grid from a field.
+
+        Parameters
+        ----------
+        pos : :class:`list`
+            the position tuple, containing main direction and transversal
+            directions
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Structured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+
+        Returns
+        -------
+        :class:`pyvista.RectilinearGrid`
+            A PyVista rectilinear grid of the structured field data. Data arrays
+            live on the point data of this PyVista dataset.
+        """
+        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
+        try:
+            import pyvista as pv
+
+            grid = pv.RectilinearGrid(x, y, z)
+            grid.point_arrays.update(fields)
+        except ImportError:
+            raise ImportError("Please install PyVista to create VTK datasets.")
+        return grid
+
+    def _to_vtk_unstructured(self, fields):
+        """Export a field to vtk structured rectilinear grid file.
+
+        Parameters
+        ----------
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Unstructured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+
+        Returns
+        -------
+        :class:`pyvista.UnstructuredGrid`
+            A PyVista unstructured grid of the unstructured field data. Data arrays
+            live on the point data of this PyVista dataset. This is essentially
+            a point cloud with no topology.
+        """
+        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
+        try:
+            import pyvista as pv
+
+            grid = pv.PolyData(np.c_[x, y, z]).cast_to_unstructured_grid()
+            grid.point_arrays.update(fields)
+        except ImportError:
+            raise ImportError("Please install PyVista to create VTK datasets.")
+        return grid
+
+    def _vtk_structured_reshape(self, fields):
+        """An internal helper to extract what is needed for the vtk rectilinear grid
+        """
+        if not isinstance(fields, dict):
+            fields = {"field": fields}
+        x, y, z = pos2xyz(self.pos)
+        if y is None:
+            y = np.array([0])
+        if z is None:
+            z = np.array([0])
+        # need fortran order in VTK
+        for field in fields:
+            fields[field] = fields[field].reshape(-1, order="F")
+            if len(fields[field]) != len(x) * len(y) * len(z):
+                raise ValueError(
+                    "gstools.vtk_export_structured: "
+                    "field shape doesn't match the given mesh"
+                )
+        return x, y, z, fields
+
+    def _vtk_unstructured_reshape(self, fields):
+        if not isinstance(fields, dict):
+            fields = {"field": fields}
+        x, y, z = pos2xyz(self.pos)
+        if y is None:
+            y = np.zeros_like(x)
+        if z is None:
+            z = np.zeros_like(x)
+        for field in fields:
+            fields[field] = fields[field].reshape(-1)
+            if (
+                len(fields[field]) != len(x)
+                or len(fields[field]) != len(y)
+                or len(fields[field]) != len(z)
+            ):
+                raise ValueError(
+                    "gstools.vtk_export_unstructured: "
+                    "field shape doesn't match the given mesh"
+                )
+        return x, y, z, fields
 
     def _check_mesh_type(self, mesh_type):
         if mesh_type != "unstructured" and mesh_type != "structured":

From 89051e1cc4c03b38f89ff113f0cf017367cd5fd0 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Sat, 18 Apr 2020 10:40:50 +0200
Subject: [PATCH 35/49] Speed up export test

---
 tests/test_export.py | 7 +++----
 1 file changed, 3 insertions(+), 4 deletions(-)

diff --git a/tests/test_export.py b/tests/test_export.py
index 8c483d09..b3f6d18e 100644
--- a/tests/test_export.py
+++ b/tests/test_export.py
@@ -21,16 +21,15 @@
 class TestExport(unittest.TestCase):
     def setUp(self):
         self.test_dir = tempfile.mkdtemp()
-        # structured field with a size 100x100x100 and a grid-size of 1x1x1
-        x = y = z = range(50)
+        x = y = z = range(25)
         model = Gaussian(dim=3, var=0.6, len_scale=20)
         self.srf_structured = SRF(model)
         self.srf_structured((x, y, z), mesh_type="structured")
         # unstrucutred field
         seed = MasterRNG(19970221)
         rng = np.random.RandomState(seed())
-        x = rng.randint(0, 100, size=1000)
-        y = rng.randint(0, 100, size=1000)
+        x = rng.randint(0, 100, size=100)
+        y = rng.randint(0, 100, size=100)
         model = Exponential(
             dim=2, var=1, len_scale=[12.0, 3.0], angles=np.pi / 8.0
         )

From a01abbbb460f779134ee76491e7898c055f3b454 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Sat, 18 Apr 2020 10:41:03 +0200
Subject: [PATCH 36/49] Refactor export test

---
 tests/test_export.py | 50 ++++++++++++++++++++++++--------------------
 1 file changed, 27 insertions(+), 23 deletions(-)

diff --git a/tests/test_export.py b/tests/test_export.py
index b3f6d18e..838160db 100644
--- a/tests/test_export.py
+++ b/tests/test_export.py
@@ -6,8 +6,7 @@
 import tempfile
 import shutil
 
-from gstools import SRF, Gaussian, Exponential
-from gstools.random import MasterRNG
+import gstools as gs
 
 HAS_PYVISTA = False
 try:
@@ -21,41 +20,46 @@
 class TestExport(unittest.TestCase):
     def setUp(self):
         self.test_dir = tempfile.mkdtemp()
-        x = y = z = range(25)
-        model = Gaussian(dim=3, var=0.6, len_scale=20)
-        self.srf_structured = SRF(model)
-        self.srf_structured((x, y, z), mesh_type="structured")
-        # unstrucutred field
-        seed = MasterRNG(19970221)
+        self.x_grid = self.y_grid = self.z_grid = range(25)
+        model = gs.Gaussian(dim=3, var=0.6, len_scale=20)
+        self.srf_structured = gs.SRF(model)
+        self.srf_structured(
+            (self.x_grid, self.y_grid, self.z_grid), mesh_type="structured"
+        )
+        # unstructured field
+        seed = gs.random.MasterRNG(19970221)
         rng = np.random.RandomState(seed())
-        x = rng.randint(0, 100, size=100)
-        y = rng.randint(0, 100, size=100)
-        model = Exponential(
+        self.x_tuple = rng.randint(0, 100, size=100)
+        self.y_tuple = rng.randint(0, 100, size=100)
+        model = gs.Exponential(
             dim=2, var=1, len_scale=[12.0, 3.0], angles=np.pi / 8.0
         )
-        self.srf_unstructured = SRF(model, seed=20170519)
-        self.srf_unstructured([x, y])
+        self.srf_unstructured = gs.SRF(model, seed=20170519)
+        self.srf_unstructured([self.x_tuple, self.y_tuple])
 
     def tearDown(self):
         # Remove the test data directory after the test
         shutil.rmtree(self.test_dir)
 
     @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
-    def test_pyvista(self):
+    def test_pyvista_struct(self):
         mesh = self.srf_structured.to_pyvista()
         self.assertIsInstance(mesh, pv.RectilinearGrid)
+
+    @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
+    def test_pyvista_unstruct(self):
         mesh = self.srf_unstructured.to_pyvista()
         self.assertIsInstance(mesh, pv.UnstructuredGrid)
 
-    def test_pyevtk_export(self):
-        # Structured
-        sfilename = os.path.join(self.test_dir, "structured")
-        self.srf_structured.vtk_export(sfilename)
-        self.assertTrue(os.path.isfile(sfilename + ".vtr"))
-        # Unstructured
-        ufilename = os.path.join(self.test_dir, "unstructured")
-        self.srf_unstructured.vtk_export(ufilename)
-        self.assertTrue(os.path.isfile(ufilename + ".vtu"))
+    def test_pyevtk_export_struct(self):
+        filename = os.path.join(self.test_dir, "structured")
+        self.srf_structured.vtk_export(filename)
+        self.assertTrue(os.path.isfile(filename + ".vtr"))
+
+    def test_pyevtk_export_unstruct(self):
+        filename = os.path.join(self.test_dir, "unstructured")
+        self.srf_unstructured.vtk_export(filename)
+        self.assertTrue(os.path.isfile(filename + ".vtu"))
 
 
 if __name__ == "__main__":

From 68867afecdf8e93c703e20b2dc1d3628a1218975 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Sat, 18 Apr 2020 10:59:59 +0200
Subject: [PATCH 37/49] Add export tests for vector fields

---
 tests/test_export.py | 32 ++++++++++++++++++++++++++++++++
 1 file changed, 32 insertions(+)

diff --git a/tests/test_export.py b/tests/test_export.py
index 838160db..15bb8275 100644
--- a/tests/test_export.py
+++ b/tests/test_export.py
@@ -61,6 +61,38 @@ def test_pyevtk_export_unstruct(self):
         self.srf_unstructured.vtk_export(filename)
         self.assertTrue(os.path.isfile(filename + ".vtu"))
 
+    @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
+    def test_pyvista_vector_struct(self):
+        model = gs.Gaussian(dim=2, var=1, len_scale=10)
+        srf = gs.SRF(model, generator="VectorField")
+        srf((self.x_grid, self.y_grid), mesh_type="structured", seed=19841203)
+        mesh = srf.to_pyvista()
+        self.assertIsInstance(mesh, pv.RectilinearGrid)
+
+    @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
+    def test_pyvista_vector_unstruct(self):
+        model = gs.Gaussian(dim=2, var=1, len_scale=10)
+        srf = gs.SRF(model, generator="VectorField")
+        srf((self.x_tuple, self.y_tuple), mesh_type="unstructured", seed=19841203)
+        mesh = srf.to_pyvista()
+        self.assertIsInstance(mesh, pv.UnstructuredGrid)
+
+    def test_pyevtk_vector_export_struct(self):
+        filename = os.path.join(self.test_dir, "vector")
+        model = gs.Gaussian(dim=2, var=1, len_scale=10)
+        srf = gs.SRF(model, generator="VectorField")
+        srf((self.x_grid, self.y_grid), mesh_type="structured", seed=19841203)
+        srf.vtk_export(filename)
+        self.assertTrue(os.path.isfile(filename + ".vtr"))
+
+    def test_pyevtk_vector_export_unstruct(self):
+        filename = os.path.join(self.test_dir, "vector")
+        model = gs.Gaussian(dim=2, var=1, len_scale=10)
+        srf = gs.SRF(model, generator="VectorField")
+        srf((self.x_tuple, self.y_tuple), mesh_type="unstructured", seed=19841203)
+        srf.vtk_export(filename)
+        self.assertTrue(os.path.isfile(filename + ".vtu"))
+
 
 if __name__ == "__main__":
     unittest.main()

From e8328a966d8f486e82dd343dd5bd7249ac3b117f Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Wed, 22 Apr 2020 15:50:31 +0200
Subject: [PATCH 38/49] Remove argument from getter

---
 gstools/field/base.py | 8 +++++---
 1 file changed, 5 insertions(+), 3 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index 9b317e5b..5b366e54 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -179,10 +179,12 @@ def field(self, values):
         self.point_data[self.default_field] = values
 
     @property
-    def value_type(self, field="field"):
+    def value_type(self):
         """:any:`str`: The value type of the default field."""
-        if field in self.point_data:
-            r = value_type(self.mesh_type, self.point_data[field].shape)
+        if self.default_field in self.point_data:
+            r = value_type(
+                self.mesh_type, self.point_data[self.default_field].shape
+            )
         else:
             r = None
         return r

From e12abfc15a0b85eed34d23303744b502c2108600 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 23 Apr 2020 11:55:12 +0200
Subject: [PATCH 39/49] Put Mesh class into separate file

---
 gstools/__init__.py       |   3 +-
 gstools/field/__init__.py |   6 +-
 gstools/field/base.py     | 474 +-----------------------------------
 gstools/field/mesh.py     | 494 ++++++++++++++++++++++++++++++++++++++
 4 files changed, 502 insertions(+), 475 deletions(-)
 create mode 100644 gstools/field/mesh.py

diff --git a/gstools/__init__.py b/gstools/__init__.py
index 3e037124..6d94acc6 100644
--- a/gstools/__init__.py
+++ b/gstools/__init__.py
@@ -33,6 +33,7 @@
 
 .. autosummary::
    SRF
+   Field
    Mesh
 
 Covariance Base-Class
@@ -101,7 +102,7 @@
 """
 
 from gstools import field, variogram, random, covmodel, tools, krige, transform
-from gstools.field.base import Mesh
+from gstools.field.mesh import Mesh
 from gstools.field import SRF
 from gstools.tools import (
     vtk_export,
diff --git a/gstools/field/__init__.py b/gstools/field/__init__.py
index bba7f848..a90322a4 100644
--- a/gstools/field/__init__.py
+++ b/gstools/field/__init__.py
@@ -17,12 +17,14 @@
 
 .. autosummary::
    SRF
+   Field
    Mesh
 
 ----
 """
 
-from gstools.field.base import Mesh
+from gstools.field.mesh import Mesh
+from gstools.field.base import Field
 from gstools.field.srf import SRF
 
-__all__ = ["SRF", "Mesh"]
+__all__ = ["SRF", "Field", "Mesh"]
diff --git a/gstools/field/base.py b/gstools/field/base.py
index 5b366e54..8594179f 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -8,7 +8,6 @@
 
 .. autosummary::
    Field
-   Mesh
 """
 # pylint: disable=C0103
 
@@ -16,8 +15,7 @@
 
 import numpy as np
 
-from pyevtk.hl import gridToVTK, pointsToVTK
-
+from gstools.field import Mesh
 from gstools.covmodel.base import CovModel
 from gstools.field.tools import (
     _get_select,
@@ -28,475 +26,7 @@
 )
 from gstools.tools.geometric import pos2xyz, xyz2pos
 
-__all__ = ["Field", "Mesh"]
-
-
-def value_type(mesh_type, shape):
-    """Determine the value type ("scalar" or "vector")"""
-    r = "scalar"
-    if mesh_type == "unstructured":
-        # TODO is this the right place for doing these checks?!
-        if len(shape) == 2 and 2 <= shape[0] <= 3:
-            r = "vector"
-    else:
-        # for very small (2-3 points) meshes, this could break
-        # a 100% solution would require dim. information.
-        if len(shape) == shape[0] + 1:
-            r = "vector"
-    return r
-
-
-class Mesh:
-    """A base class encapsulating field data.
-
-    It holds a position array, which define the spatial locations of the
-    field values.
-    It can hold multiple fields in the :any:`self.point_data` dict. This assumes
-    that each field is defined on the same positions.
-    The mesh type must also be specified.
-
-    Parameters
-    ----------
-    pos : :class:`numpy.ndarray`, optional
-        positions of the field values
-    name : :any:`str`, optional
-        key of the field values
-    values : :any:`list`, optional
-        a list of the values of the fields
-    mesh_type : :any:`str`, optional
-        the type of mesh on which the field is defined on, can be
-        * unstructured
-        * structured
-
-    Examples
-    --------
-    >>> import numpy as np
-    >>> from gstools import Mesh
-    >>> pos = np.linspace(0., 100., 40)
-    >>> z = np.random.random(40)
-    >>> z2 = np.random.random(40)
-    >>> mesh = Mesh((pos,))
-    >>> mesh.add_field(z, "test_field1")
-    >>> mesh.add_field(z2, "test_field2")
-    >>> mesh.set_default_field("test_field2")
-    >>> print(mesh.field)
-
-    """
-
-    def __init__(
-        self, pos=None, name="field", values=None, mesh_type="unstructured",
-    ):
-        self._vtk_export_fct = self._vtk_export_unstructured
-        self._to_vtk_fct = self._to_vtk_unstructured
-        # mesh_type needs a special setter, therefore, `set_field_data` is not
-        # used here
-        self.mesh_type = mesh_type
-
-        # the pos/ points of the mesh
-        self._pos = pos
-
-        # data stored at each pos/ point, the "fields"
-        if values is not None:
-            self.point_data = {name: values}
-        else:
-            self.point_data = {}
-
-        # data valid for the global field
-        self.field_data = {}
-
-        self.set_field_data("default_field", name)
-
-        self.field_data["mesh_type"] = mesh_type
-
-    def set_field_data(self, name, value):
-        """Add an attribute to this instance and add it the `field_data`
-
-        This helper method is used to create attributes for easy access
-        while using this class, but it also adds an entry to the dictionary
-        `field_data`, which is used for exporting the data.
-        """
-        setattr(self, name, value)
-        self.field_data[name] = value
-
-    def add_field(
-        self, values, name="field", is_default_field=False,
-    ):
-        """Add a field (point data) to the mesh
-
-        .. note::
-            If no field has existed before calling this method,
-            the given field will be set to the default field.
-
-        .. warning::
-            If point data with the same `name` already exists, it will be
-            overwritten.
-
-        Parameters
-        ----------
-        values : :class:`numpy.ndarray`
-            the point data, has to be the same shape as the mesh
-        name : :class:`str`, optional
-            the name of the point data
-        is_default_field : :class:`bool`, optional
-            is this the default field of the mesh?
-
-        """
-        values = np.array(values)
-        self._check_point_data(values)
-        self.point_data[name] = values
-        # set the default_field to the first field added
-        if len(self.point_data) == 1 or is_default_field:
-            self.set_field_data("default_field", name)
-
-    def __getitem__(self, key):
-        """:any:`numpy.ndarray`: The values of the field."""
-        return self.point_data[key]
-
-    def __setitem__(self, key, value):
-        self.point_data[key] = value
-
-    @property
-    def pos(self):
-        """:any:`numpy.ndarray`: The pos. on which the field is defined."""
-        return self._pos
-
-    @pos.setter
-    def pos(self, value):
-        """
-        Warning: setting new positions deletes all previously stored fields.
-        """
-        self.point_data = {self.default_field: None}
-        self._pos = value
-
-    @property
-    def field(self):
-        """:class:`numpy.ndarray`: The point data of the default field."""
-        return self.point_data[self.default_field]
-
-    @field.setter
-    def field(self, values):
-        self._check_point_data(values)
-        self.point_data[self.default_field] = values
-
-    @property
-    def value_type(self):
-        """:any:`str`: The value type of the default field."""
-        if self.default_field in self.point_data:
-            r = value_type(
-                self.mesh_type, self.point_data[self.default_field].shape
-            )
-        else:
-            r = None
-        return r
-
-    @property
-    def mesh_type(self):
-        """:any:`str`: The mesh type of the fields."""
-        return self._mesh_type
-
-    @mesh_type.setter
-    def mesh_type(self, value):
-        """
-        Warning: setting a new mesh type deletes all previously stored fields.
-        """
-        self._check_mesh_type(value)
-        if value == "structured":
-            self._vtk_export_fct = self._vtk_export_structured
-            self._to_vtk_fct = self._to_vtk_structured
-        else:
-            self._vtk_export_fct = self._vtk_export_unstructured
-            self._to_vtk_fct = self._to_vtk_unstructured
-        self.point_data = {}
-        self._mesh_type = value
-
-    def _vtk_naming_helper(
-        self, field_select="field", fieldname="field"
-    ):  # pragma: no cover
-        """Prepare the field names for export.
-
-        Parameters
-        ----------
-        field_select : :class:`str`, optional
-            Field that should be stored. Can be:
-            "field", "raw_field", "krige_field", "err_field" or "krige_var".
-            Default: "field"
-        fieldname : :class:`str`, optional
-            Name of the field in the VTK file. Default: "field"
-
-        Returns
-        -------
-        fields : :class:`dict`
-            a dictionary containing the fields to be exported
-        """
-        if self.value_type is None:
-            raise ValueError(
-                "Field value type not set! "
-                + "Specify 'scalar' or 'vector' before plotting."
-            )
-        elif self.value_type == "vector":
-            if hasattr(self, field_select):
-                field = getattr(self, field_select)
-            else:
-                field = None
-            if not (
-                self.pos is None or field is None or self.mesh_type is None
-            ):
-                suf = ["_X", "_Y", "_Z"]
-                fields = {}
-                for i in range(self.model.dim):
-                    fields[fieldname + suf[i]] = field[i]
-                return fields
-        elif self.value_type == "scalar":
-            if hasattr(self, field_select):
-                field = getattr(self, field_select)
-            else:
-                field = None
-            if not (
-                self.pos is None or field is None or self.mesh_type is None
-            ):
-                return {fieldname: field}
-            else:
-                print(
-                    "Field.to_vtk: No "
-                    + field_select
-                    + " stored in the class."
-                )
-        else:
-            raise ValueError(
-                "Unknown field value type: {}".format(self.value_type)
-            )
-
-    def to_pyvista(
-        self, field_select="field", fieldname="field"
-    ):  # pragma: no cover
-        """Create a VTK/PyVista grid of the stored field.
-
-        Parameters
-        ----------
-        field_select : :class:`str`, optional
-            Field that should be stored. Can be:
-            "field", "raw_field", "krige_field", "err_field" or "krige_var".
-            Default: "field"
-        fieldname : :class:`str`, optional
-            Name of the field in the VTK file. Default: "field"
-        """
-        field_names = self._vtk_naming_helper(
-            field_select=field_select, fieldname=fieldname
-        )
-
-        grid = self._to_vtk_fct(field_names)
-        return grid
-
-    def vtk_export(
-        self, filename, field_select="field", fieldname="field"
-    ):  # pragma: no cover
-        """Export the stored field to vtk.
-
-        Parameters
-        ----------
-        filename : :class:`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 : :class:`str`, optional
-            Field that should be stored. Can be:
-            "field", "raw_field", "krige_field", "err_field" or "krige_var".
-            Default: "field"
-        fieldname : :class:`str`, optional
-            Name of the field in the VTK file. Default: "field"
-        """
-        if not isinstance(filename, str):
-            raise TypeError("Please use a string as a filename.")
-        field_names = self._vtk_naming_helper(
-            field_select=field_select, fieldname=fieldname
-        )
-        return self._vtk_export_fct(filename, field_names)
-
-    def _vtk_export_structured(self, filename, fields):
-        """Export a field to vtk structured rectilinear grid file.
-
-        Parameters
-        ----------
-        filename : :class:`str`
-            Filename of the file to be saved, including the path. Note that an
-            ending (.vtr) will be added to the name.
-        fields : :class:`dict` or :class:`numpy.ndarray`
-            Structured fields to be saved.
-            Either a single numpy array as returned by SRF,
-            or a dictionary of fields with theirs names as keys.
-        """
-        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
-        return gridToVTK(filename, x, y, z, pointData=fields)
-
-    def _vtk_export_unstructured(self, filename, fields):
-        """Export a field to vtk unstructured grid file.
-
-        Parameters
-        ----------
-        filename : :class:`str`
-            Filename of the file to be saved, including the path. Note that an
-            ending (.vtu) will be added to the name.
-        fields : :class:`dict` or :class:`numpy.ndarray`
-            Unstructured fields to be saved.
-            Either a single numpy array as returned by SRF,
-            or a dictionary of fields with theirs names as keys.
-        """
-        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
-        return pointsToVTK(filename, x, y, z, data=fields)
-
-    def _to_vtk_structured(self, fields):  # pragma: no cover
-        """Create a vtk structured rectilinear grid from a field.
-
-        Parameters
-        ----------
-        pos : :class:`list`
-            the position tuple, containing main direction and transversal
-            directions
-        fields : :class:`dict` or :class:`numpy.ndarray`
-            Structured fields to be saved.
-            Either a single numpy array as returned by SRF,
-            or a dictionary of fields with theirs names as keys.
-
-        Returns
-        -------
-        :class:`pyvista.RectilinearGrid`
-            A PyVista rectilinear grid of the structured field data. Data arrays
-            live on the point data of this PyVista dataset.
-        """
-        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
-        try:
-            import pyvista as pv
-
-            grid = pv.RectilinearGrid(x, y, z)
-            grid.point_arrays.update(fields)
-        except ImportError:
-            raise ImportError("Please install PyVista to create VTK datasets.")
-        return grid
-
-    def _to_vtk_unstructured(self, fields):
-        """Export a field to vtk structured rectilinear grid file.
-
-        Parameters
-        ----------
-        fields : :class:`dict` or :class:`numpy.ndarray`
-            Unstructured fields to be saved.
-            Either a single numpy array as returned by SRF,
-            or a dictionary of fields with theirs names as keys.
-
-        Returns
-        -------
-        :class:`pyvista.UnstructuredGrid`
-            A PyVista unstructured grid of the unstructured field data. Data arrays
-            live on the point data of this PyVista dataset. This is essentially
-            a point cloud with no topology.
-        """
-        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
-        try:
-            import pyvista as pv
-
-            grid = pv.PolyData(np.c_[x, y, z]).cast_to_unstructured_grid()
-            grid.point_arrays.update(fields)
-        except ImportError:
-            raise ImportError("Please install PyVista to create VTK datasets.")
-        return grid
-
-    def _vtk_structured_reshape(self, fields):
-        """An internal helper to extract what is needed for the vtk rectilinear grid
-        """
-        if not isinstance(fields, dict):
-            fields = {"field": fields}
-        x, y, z = pos2xyz(self.pos)
-        if y is None:
-            y = np.array([0])
-        if z is None:
-            z = np.array([0])
-        # need fortran order in VTK
-        for field in fields:
-            fields[field] = fields[field].reshape(-1, order="F")
-            if len(fields[field]) != len(x) * len(y) * len(z):
-                raise ValueError(
-                    "gstools.vtk_export_structured: "
-                    "field shape doesn't match the given mesh"
-                )
-        return x, y, z, fields
-
-    def _vtk_unstructured_reshape(self, fields):
-        if not isinstance(fields, dict):
-            fields = {"field": fields}
-        x, y, z = pos2xyz(self.pos)
-        if y is None:
-            y = np.zeros_like(x)
-        if z is None:
-            z = np.zeros_like(x)
-        for field in fields:
-            fields[field] = fields[field].reshape(-1)
-            if (
-                len(fields[field]) != len(x)
-                or len(fields[field]) != len(y)
-                or len(fields[field]) != len(z)
-            ):
-                raise ValueError(
-                    "gstools.vtk_export_unstructured: "
-                    "field shape doesn't match the given mesh"
-                )
-        return x, y, z, fields
-
-    def _check_mesh_type(self, mesh_type):
-        if mesh_type != "unstructured" and mesh_type != "structured":
-            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
-
-    def _check_point_data(self, values):
-        """Compare field shape to pos shape.
-
-        Parameters
-        ----------
-        values : :class:`numpy.ndarray`
-            the values of the field to be checked
-        """
-        err = True
-        if self.mesh_type == "unstructured":
-            # scalar
-            if len(values.shape) == 1:
-                if values.shape[0] == len(self.pos[0]):
-                    err = False
-            # vector
-            elif len(values.shape) == 2:
-                if (
-                    values.shape[1] == len(self.pos[0])
-                    and 2 <= values.shape[0] <= 3
-                ):
-                    err = False
-            if err:
-                raise ValueError(
-                    "Wrong field shape: {0} does not match mesh shape ({1},)".format(
-                        values.shape, len(self.pos[0])
-                    )
-                )
-        else:
-            # scalar
-            if len(values.shape) == len(self.pos):
-                if all(
-                    [
-                        values.shape[i] == len(self.pos[i])
-                        for i in range(len(self.pos))
-                    ]
-                ):
-                    err = False
-            # vector
-            elif len(values.shape) == len(self.pos) + 1:
-                if all(
-                    [
-                        values.shape[i + 1] == len(self.pos[i])
-                        for i in range(len(self.pos))
-                    ]
-                ) and values.shape[0] == len(self.pos):
-                    err = False
-            if err:
-                raise ValueError(
-                    "Wrong field shape: {0} does not match mesh shape [0/2/3]{1}".format(
-                        list(values.shape),
-                        [len(self.pos[i]) for i in range(len(self.pos))],
-                    )
-                )
+__all__ = ["Field"]
 
 
 class Field(Mesh):
diff --git a/gstools/field/mesh.py b/gstools/field/mesh.py
new file mode 100644
index 00000000..a7eb012f
--- /dev/null
+++ b/gstools/field/mesh.py
@@ -0,0 +1,494 @@
+# -*- coding: utf-8 -*-
+"""
+GStools subpackage providing a base class for spatial fields.
+
+.. currentmodule:: gstools.field.base
+
+The following classes are provided
+
+.. autosummary::
+   Mesh
+"""
+# pylint: disable=C0103
+
+import numpy as np
+
+from pyevtk.hl import gridToVTK, pointsToVTK
+
+from gstools.tools.geometric import pos2xyz
+
+__all__ = ["Mesh"]
+
+
+def value_type(mesh_type, shape):
+    """Determine the value type ("scalar" or "vector")"""
+    r = "scalar"
+    if mesh_type == "unstructured":
+        # TODO is this the right place for doing these checks?!
+        if len(shape) == 2 and 2 <= shape[0] <= 3:
+            r = "vector"
+    else:
+        # for very small (2-3 points) meshes, this could break
+        # a 100% solution would require dim. information.
+        if len(shape) == shape[0] + 1:
+            r = "vector"
+    return r
+
+
+class Mesh:
+    """A base class encapsulating field data.
+
+    It holds a position array, which define the spatial locations of the
+    field values.
+    It can hold multiple fields in the :any:`self.point_data` dict. This assumes
+    that each field is defined on the same positions.
+    The mesh type must also be specified.
+
+    Parameters
+    ----------
+    pos : :class:`numpy.ndarray`, optional
+        positions of the field values
+    name : :any:`str`, optional
+        key of the field values
+    values : :any:`list`, optional
+        a list of the values of the fields
+    mesh_type : :any:`str`, optional
+        the type of mesh on which the field is defined on, can be
+        * unstructured
+        * structured
+
+    Examples
+    --------
+    >>> import numpy as np
+    >>> from gstools import Mesh
+    >>> pos = np.linspace(0., 100., 40)
+    >>> z = np.random.random(40)
+    >>> z2 = np.random.random(40)
+    >>> mesh = Mesh((pos,))
+    >>> mesh.add_field(z, "test_field1")
+    >>> mesh.add_field(z2, "test_field2")
+    >>> mesh.set_default_field("test_field2")
+    >>> print(mesh.field)
+
+    """
+
+    def __init__(
+        self, pos=None, name="field", values=None, mesh_type="unstructured",
+    ):
+        self._vtk_export_fct = self._vtk_export_unstructured
+        self._to_vtk_fct = self._to_vtk_unstructured
+        # mesh_type needs a special setter, therefore, `set_field_data` is not
+        # used here
+        self.mesh_type = mesh_type
+
+        # the pos/ points of the mesh
+        self._pos = pos
+
+        # data stored at each pos/ point, the "fields"
+        if values is not None:
+            self.point_data = {name: values}
+        else:
+            self.point_data = {}
+
+        # data valid for the global field
+        self.field_data = {}
+
+        self.set_field_data("default_field", name)
+
+        self.field_data["mesh_type"] = mesh_type
+
+    def set_field_data(self, name, value):
+        """Add an attribute to this instance and add it the `field_data`
+
+        This helper method is used to create attributes for easy access
+        while using this class, but it also adds an entry to the dictionary
+        `field_data`, which is used for exporting the data.
+        """
+        setattr(self, name, value)
+        self.field_data[name] = value
+
+    def add_field(
+        self, values, name="field", is_default_field=False,
+    ):
+        """Add a field (point data) to the mesh
+
+        .. note::
+            If no field has existed before calling this method,
+            the given field will be set to the default field.
+
+        .. warning::
+            If point data with the same `name` already exists, it will be
+            overwritten.
+
+        Parameters
+        ----------
+        values : :class:`numpy.ndarray`
+            the point data, has to be the same shape as the mesh
+        name : :class:`str`, optional
+            the name of the point data
+        is_default_field : :class:`bool`, optional
+            is this the default field of the mesh?
+
+        """
+        values = np.array(values)
+        self._check_point_data(values)
+        self.point_data[name] = values
+        # set the default_field to the first field added
+        if len(self.point_data) == 1 or is_default_field:
+            self.set_field_data("default_field", name)
+
+    def __getitem__(self, key):
+        """:any:`numpy.ndarray`: The values of the field."""
+        return self.point_data[key]
+
+    def __setitem__(self, key, value):
+        self.point_data[key] = value
+
+    @property
+    def pos(self):
+        """:any:`numpy.ndarray`: The pos. on which the field is defined."""
+        return self._pos
+
+    @pos.setter
+    def pos(self, value):
+        """
+        Warning: setting new positions deletes all previously stored fields.
+        """
+        self.point_data = {self.default_field: None}
+        self._pos = value
+
+    @property
+    def field(self):
+        """:class:`numpy.ndarray`: The point data of the default field."""
+        return self.point_data[self.default_field]
+
+    @field.setter
+    def field(self, values):
+        self._check_point_data(values)
+        self.point_data[self.default_field] = values
+
+    @property
+    def value_type(self):
+        """:any:`str`: The value type of the default field."""
+        if self.default_field in self.point_data:
+            r = value_type(
+                self.mesh_type, self.point_data[self.default_field].shape
+            )
+        else:
+            r = None
+        return r
+
+    @property
+    def mesh_type(self):
+        """:any:`str`: The mesh type of the fields."""
+        return self._mesh_type
+
+    @mesh_type.setter
+    def mesh_type(self, value):
+        """
+        Warning: setting a new mesh type deletes all previously stored fields.
+        """
+        self._check_mesh_type(value)
+        if value == "structured":
+            self._vtk_export_fct = self._vtk_export_structured
+            self._to_vtk_fct = self._to_vtk_structured
+        else:
+            self._vtk_export_fct = self._vtk_export_unstructured
+            self._to_vtk_fct = self._to_vtk_unstructured
+        self.point_data = {}
+        self._mesh_type = value
+
+    def _vtk_naming_helper(
+        self, field_select="field", fieldname="field"
+    ):  # pragma: no cover
+        """Prepare the field names for export.
+
+        Parameters
+        ----------
+        field_select : :class:`str`, optional
+            Field that should be stored. Can be:
+            "field", "raw_field", "krige_field", "err_field" or "krige_var".
+            Default: "field"
+        fieldname : :class:`str`, optional
+            Name of the field in the VTK file. Default: "field"
+
+        Returns
+        -------
+        fields : :class:`dict`
+            a dictionary containing the fields to be exported
+        """
+        if self.value_type is None:
+            raise ValueError(
+                "Field value type not set! "
+                + "Specify 'scalar' or 'vector' before plotting."
+            )
+        elif self.value_type == "vector":
+            if hasattr(self, field_select):
+                field = getattr(self, field_select)
+            else:
+                field = None
+            if not (
+                self.pos is None or field is None or self.mesh_type is None
+            ):
+                suf = ["_X", "_Y", "_Z"]
+                fields = {}
+                for i in range(self.model.dim):
+                    fields[fieldname + suf[i]] = field[i]
+                return fields
+        elif self.value_type == "scalar":
+            if hasattr(self, field_select):
+                field = getattr(self, field_select)
+            else:
+                field = None
+            if not (
+                self.pos is None or field is None or self.mesh_type is None
+            ):
+                return {fieldname: field}
+            else:
+                print(
+                    "Field.to_vtk: No "
+                    + field_select
+                    + " stored in the class."
+                )
+        else:
+            raise ValueError(
+                "Unknown field value type: {}".format(self.value_type)
+            )
+
+    def to_pyvista(
+        self, field_select="field", fieldname="field"
+    ):  # pragma: no cover
+        """Create a VTK/PyVista grid of the stored field.
+
+        Parameters
+        ----------
+        field_select : :class:`str`, optional
+            Field that should be stored. Can be:
+            "field", "raw_field", "krige_field", "err_field" or "krige_var".
+            Default: "field"
+        fieldname : :class:`str`, optional
+            Name of the field in the VTK file. Default: "field"
+        """
+        field_names = self._vtk_naming_helper(
+            field_select=field_select, fieldname=fieldname
+        )
+
+        grid = self._to_vtk_fct(field_names)
+        return grid
+
+    def vtk_export(
+        self, filename, field_select="field", fieldname="field"
+    ):  # pragma: no cover
+        """Export the stored field to vtk.
+
+        Parameters
+        ----------
+        filename : :class:`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 : :class:`str`, optional
+            Field that should be stored. Can be:
+            "field", "raw_field", "krige_field", "err_field" or "krige_var".
+            Default: "field"
+        fieldname : :class:`str`, optional
+            Name of the field in the VTK file. Default: "field"
+        """
+        if not isinstance(filename, str):
+            raise TypeError("Please use a string as a filename.")
+        field_names = self._vtk_naming_helper(
+            field_select=field_select, fieldname=fieldname
+        )
+        return self._vtk_export_fct(filename, field_names)
+
+    def _vtk_export_structured(self, filename, fields):
+        """Export a field to vtk structured rectilinear grid file.
+
+        Parameters
+        ----------
+        filename : :class:`str`
+            Filename of the file to be saved, including the path. Note that an
+            ending (.vtr) will be added to the name.
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Structured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+        """
+        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
+        return gridToVTK(filename, x, y, z, pointData=fields)
+
+    def _vtk_export_unstructured(self, filename, fields):
+        """Export a field to vtk unstructured grid file.
+
+        Parameters
+        ----------
+        filename : :class:`str`
+            Filename of the file to be saved, including the path. Note that an
+            ending (.vtu) will be added to the name.
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Unstructured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+        """
+        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
+        return pointsToVTK(filename, x, y, z, data=fields)
+
+    def _to_vtk_structured(self, fields):  # pragma: no cover
+        """Create a vtk structured rectilinear grid from a field.
+
+        Parameters
+        ----------
+        pos : :class:`list`
+            the position tuple, containing main direction and transversal
+            directions
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Structured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+
+        Returns
+        -------
+        :class:`pyvista.RectilinearGrid`
+            A PyVista rectilinear grid of the structured field data. Data arrays
+            live on the point data of this PyVista dataset.
+        """
+        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
+        try:
+            import pyvista as pv
+
+            grid = pv.RectilinearGrid(x, y, z)
+            grid.point_arrays.update(fields)
+        except ImportError:
+            raise ImportError("Please install PyVista to create VTK datasets.")
+        return grid
+
+    def _to_vtk_unstructured(self, fields):
+        """Export a field to vtk structured rectilinear grid file.
+
+        Parameters
+        ----------
+        fields : :class:`dict` or :class:`numpy.ndarray`
+            Unstructured fields to be saved.
+            Either a single numpy array as returned by SRF,
+            or a dictionary of fields with theirs names as keys.
+
+        Returns
+        -------
+        :class:`pyvista.UnstructuredGrid`
+            A PyVista unstructured grid of the unstructured field data. Data arrays
+            live on the point data of this PyVista dataset. This is essentially
+            a point cloud with no topology.
+        """
+        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
+        try:
+            import pyvista as pv
+
+            grid = pv.PolyData(np.c_[x, y, z]).cast_to_unstructured_grid()
+            grid.point_arrays.update(fields)
+        except ImportError:
+            raise ImportError("Please install PyVista to create VTK datasets.")
+        return grid
+
+    def _vtk_structured_reshape(self, fields):
+        """An internal helper to extract what is needed for the vtk rectilinear grid
+        """
+        if not isinstance(fields, dict):
+            fields = {"field": fields}
+        x, y, z = pos2xyz(self.pos)
+        if y is None:
+            y = np.array([0])
+        if z is None:
+            z = np.array([0])
+        # need fortran order in VTK
+        for field in fields:
+            fields[field] = fields[field].reshape(-1, order="F")
+            if len(fields[field]) != len(x) * len(y) * len(z):
+                raise ValueError(
+                    "gstools.vtk_export_structured: "
+                    "field shape doesn't match the given mesh"
+                )
+        return x, y, z, fields
+
+    def _vtk_unstructured_reshape(self, fields):
+        if not isinstance(fields, dict):
+            fields = {"field": fields}
+        x, y, z = pos2xyz(self.pos)
+        if y is None:
+            y = np.zeros_like(x)
+        if z is None:
+            z = np.zeros_like(x)
+        for field in fields:
+            fields[field] = fields[field].reshape(-1)
+            if (
+                len(fields[field]) != len(x)
+                or len(fields[field]) != len(y)
+                or len(fields[field]) != len(z)
+            ):
+                raise ValueError(
+                    "gstools.vtk_export_unstructured: "
+                    "field shape doesn't match the given mesh"
+                )
+        return x, y, z, fields
+
+    def _check_mesh_type(self, mesh_type):
+        if mesh_type != "unstructured" and mesh_type != "structured":
+            raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
+
+    def _check_point_data(self, values):
+        """Compare field shape to pos shape.
+
+        Parameters
+        ----------
+        values : :class:`numpy.ndarray`
+            the values of the field to be checked
+        """
+        err = True
+        if self.mesh_type == "unstructured":
+            # scalar
+            if len(values.shape) == 1:
+                if values.shape[0] == len(self.pos[0]):
+                    err = False
+            # vector
+            elif len(values.shape) == 2:
+                if (
+                    values.shape[1] == len(self.pos[0])
+                    and 2 <= values.shape[0] <= 3
+                ):
+                    err = False
+            if err:
+                raise ValueError(
+                    "Wrong field shape: {0} does not match mesh shape ({1},)".format(
+                        values.shape, len(self.pos[0])
+                    )
+                )
+        else:
+            # scalar
+            if len(values.shape) == len(self.pos):
+                if all(
+                    [
+                        values.shape[i] == len(self.pos[i])
+                        for i in range(len(self.pos))
+                    ]
+                ):
+                    err = False
+            # vector
+            elif len(values.shape) == len(self.pos) + 1:
+                if all(
+                    [
+                        values.shape[i + 1] == len(self.pos[i])
+                        for i in range(len(self.pos))
+                    ]
+                ) and values.shape[0] == len(self.pos):
+                    err = False
+            if err:
+                raise ValueError(
+                    "Wrong field shape: {0} does not match mesh shape [0/2/3]{1}".format(
+                        list(values.shape),
+                        [len(self.pos[i]) for i in range(len(self.pos))],
+                    )
+                )
+
+
+if __name__ == "__main__":  # pragma: no cover
+    import doctest
+
+    doctest.testmod()

From a88641903c0b1e8b69029898d8ac995952710d30 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Apr 2020 15:15:20 +0200
Subject: [PATCH 40/49] Add an example for the Mesh class

---
 examples/03_variogram/README.rst  |  2 +
 examples/08_mesh/00_mesh_intro.py | 75 +++++++++++++++++++++++++++++++
 examples/08_mesh/README.rst       | 12 +++++
 3 files changed, 89 insertions(+)
 create mode 100644 examples/08_mesh/00_mesh_intro.py
 create mode 100644 examples/08_mesh/README.rst

diff --git a/examples/03_variogram/README.rst b/examples/03_variogram/README.rst
index d16c16c2..b0339f12 100644
--- a/examples/03_variogram/README.rst
+++ b/examples/03_variogram/README.rst
@@ -1,3 +1,5 @@
+.. _tutorial_03_variogram:
+
 Tutorial 3: Variogram Estimation
 ================================
 
diff --git a/examples/08_mesh/00_mesh_intro.py b/examples/08_mesh/00_mesh_intro.py
new file mode 100644
index 00000000..ba749a9e
--- /dev/null
+++ b/examples/08_mesh/00_mesh_intro.py
@@ -0,0 +1,75 @@
+"""
+Example: Using a mesh for GSTools
+---------------------------------
+
+# TODO write a summery once I'm finished with the example
+
+
+Pretending we have some Data
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+We will pretend that we have some external data, by creating some random data
+with GSTools, but we will delete the objects and only use the data, without the
+backend GSTools provides.
+"""
+from datetime import date
+import numpy as np
+import gstools as gs
+
+# create a circular mesh
+point_no = 10000
+rng = np.random.RandomState(20170521)
+r = 50.0 * np.sqrt(rng.uniform(size=point_no))
+theta = 2.0 * np.pi * rng.uniform(size=point_no)
+x = r * np.cos(theta)
+y = r * np.sin(theta)
+
+tmp_model = gs.Exponential(dim=2, var=1.5, len_scale=10.0)
+tmp_srf = gs.SRF(tmp_model)
+field = tmp_srf((x, y))
+tmp_srf.plot()
+
+# Now that we have our data, let's delete everything GSTools related and pretend
+# that this has never happend
+del(tmp_model)
+del(tmp_srf)
+
+
+# Creating the Mesh
+# ^^^^^^^^^^^^^^^^^
+#
+# Starting out fresh, we want to feed the mesh with our data
+mesh = gs.Mesh(pos=(x, y), values=field)
+
+# We can add meta data too
+mesh.set_field_data("location", "Süderbrarup")
+mesh.set_field_data("date", date(year=2020, month=2, day=28))
+
+# This can be conviniently accessed
+print(mesh.location)
+print(mesh.date)
+
+# But the meta data is also collected as a dictionary in case you want to export
+# it
+print(mesh.field_data)
+
+
+# Estimating the Variogram
+# ^^^^^^^^^^^^^^^^^^^^^^^^
+# Now, with our mesh, which was loaded from completely external sources, we can
+# estimate the variogram of the data.
+# To speed things up, we will only use a fraction of the available data
+
+bins = np.linspace(0, 50, 50)
+bin_centre, gamma = gs.vario_estimate_unstructured(
+    (x, y), field, bins, sampling_size=2000, sampling_seed=19900408)
+
+# As we are experts, we'll do an expert guess and say, that we will most likely
+# have data that has an exponential variogram. Non-experts can have a look at
+# the "Finding the best fitting variogram model" tutorial in
+# :ref:`tutorial_03_variogram`.
+fit_model = gs.Exponential(dim=2)
+fit_model.fit_variogram(bin_centre, gamma, nugget=False)
+
+ax = fit_model.plot(x_max=max(bin_centre))
+ax.plot(bin_centre, gamma)
diff --git a/examples/08_mesh/README.rst b/examples/08_mesh/README.rst
new file mode 100644
index 00000000..2b719627
--- /dev/null
+++ b/examples/08_mesh/README.rst
@@ -0,0 +1,12 @@
+.. _tutorial_08_mesh:
+
+Tutorial 8: Working With Meshes
+===============================
+
+In case you have external data, which does not come from the spatial random
+generation or similar methods of GSTools and you want to analyse this data, e.g.
+estimate the variogram, you can load the data into a mesh class and conveniently
+manipulate, analyse, visualise, and export the data with GSTools.
+
+Gallery
+-------

From b9e21958bebdc6c4d040998b70a3dbb6cad69995 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Apr 2020 15:17:58 +0200
Subject: [PATCH 41/49] Amend

---
 examples/08_mesh/00_mesh_intro.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/examples/08_mesh/00_mesh_intro.py b/examples/08_mesh/00_mesh_intro.py
index ba749a9e..bb7c3f9f 100644
--- a/examples/08_mesh/00_mesh_intro.py
+++ b/examples/08_mesh/00_mesh_intro.py
@@ -2,7 +2,7 @@
 Example: Using a mesh for GSTools
 ---------------------------------
 
-# TODO write a summery once I'm finished with the example
+This example shows how external data can be analysed with GSTools.
 
 
 Pretending we have some Data

From 0864e0cb64b29e419faa82e030d84c5fd95fcc92 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Apr 2020 16:18:15 +0200
Subject: [PATCH 42/49] Unify set_field_data args with add_field

---
 examples/08_mesh/00_mesh_intro.py | 4 ++--
 gstools/field/mesh.py             | 6 +++---
 tests/test_mesh.py                | 2 +-
 3 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/examples/08_mesh/00_mesh_intro.py b/examples/08_mesh/00_mesh_intro.py
index bb7c3f9f..e0d4d196 100644
--- a/examples/08_mesh/00_mesh_intro.py
+++ b/examples/08_mesh/00_mesh_intro.py
@@ -42,8 +42,8 @@
 mesh = gs.Mesh(pos=(x, y), values=field)
 
 # We can add meta data too
-mesh.set_field_data("location", "Süderbrarup")
-mesh.set_field_data("date", date(year=2020, month=2, day=28))
+mesh.set_field_data("Süderbrarup", "location")
+mesh.set_field_data(date(year=2020, month=2, day=28), "date")
 
 # This can be conviniently accessed
 print(mesh.location)
diff --git a/gstools/field/mesh.py b/gstools/field/mesh.py
index a7eb012f..e1387052 100644
--- a/gstools/field/mesh.py
+++ b/gstools/field/mesh.py
@@ -93,11 +93,11 @@ def __init__(
         # data valid for the global field
         self.field_data = {}
 
-        self.set_field_data("default_field", name)
+        self.set_field_data(name, "default_field")
 
         self.field_data["mesh_type"] = mesh_type
 
-    def set_field_data(self, name, value):
+    def set_field_data(self, value, name):
         """Add an attribute to this instance and add it the `field_data`
 
         This helper method is used to create attributes for easy access
@@ -135,7 +135,7 @@ def add_field(
         self.point_data[name] = values
         # set the default_field to the first field added
         if len(self.point_data) == 1 or is_default_field:
-            self.set_field_data("default_field", name)
+            self.set_field_data(name, "default_field")
 
     def __getitem__(self, key):
         """:any:`numpy.ndarray`: The values of the field."""
diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index be160ba7..8d9425ef 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -80,7 +80,7 @@ def test_default_value_type(self):
 
     def test_field_data_setter(self):
         # attribute creation by adding field_data
-        self.m2_tuple.set_field_data("mean", 3.14)
+        self.m2_tuple.set_field_data(3.14, "mean")
         self.assertEqual(self.m2_tuple.field_data["mean"], 3.14)
         self.assertEqual(self.m2_tuple.mean, 3.14)
 

From a76fb81ed8f04d6af7fc712a7b36e38dc7b66fa2 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Apr 2020 16:18:50 +0200
Subject: [PATCH 43/49] Meshanise GSTools

---
 gstools/field/condition.py |  4 ++--
 gstools/field/srf.py       | 21 ++++++++-------------
 2 files changed, 10 insertions(+), 15 deletions(-)

diff --git a/gstools/field/condition.py b/gstools/field/condition.py
index cbb0f076..55580d7d 100644
--- a/gstools/field/condition.py
+++ b/gstools/field/condition.py
@@ -53,7 +53,7 @@ def ordinary(srf):
         model=srf.model, cond_pos=srf.cond_pos, cond_val=err_data
     )
     err_field, __ = err_ok(srf.pos, srf.mesh_type)
-    cond_field = srf.raw_field + krige_field - err_field
+    cond_field = srf["raw_field"] + krige_field - err_field
     info = {"mean": krige_ok.mean}
     return cond_field, krige_field, err_field, krige_var, info
 
@@ -101,6 +101,6 @@ def simple(srf):
         cond_val=err_data,
     )
     err_field, __ = err_sk(srf.pos, srf.mesh_type)
-    cond_field = srf.raw_field + krige_field - err_field + srf.mean
+    cond_field = srf["raw_field"] + krige_field - err_field + srf.mean
     info = {}
     return cond_field, krige_field, err_field, krige_var, info
diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index 8c5a9a17..55f7db2a 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -91,11 +91,6 @@ def __init__(
         self._cond_pos = None
         self._cond_val = None
         self._krige_type = None
-        # initialize attributes
-        self.raw_field = None
-        self.krige_field = None
-        self.err_field = None
-        self.krige_var = None
         self.set_generator(generator, **generator_kwargs)
         self.upscaling = upscaling
         if self._value_type is None:
@@ -140,11 +135,11 @@ def __call__(
             pos, mesh_type
         )
         # generate the field
-        self.raw_field = self.generator.__call__(x, y, z, mt_gen)
+        self.add_field(self.generator.__call__(x, y, z, mt_gen), "raw_field")
         # reshape field if we got an unstructured mesh
         if mt_changed:
-            self.raw_field = reshape_field_from_unstruct_to_struct(
-                self.model.dim, self.raw_field, axis_lens
+            self["raw_field"] = reshape_field_from_unstruct_to_struct(
+                self.model.dim, self["raw_field"], axis_lens
             )
         # apply given conditions to the field
         if self.condition:
@@ -157,13 +152,13 @@ def __call__(
             ) = self.cond_func(self)
             # store everything in the class
             self.field = cond_field
-            self.krige_field = krige_field
-            self.err_field = err_field
-            self.krige_var = krigevar
+            self.add_field(krige_field, "krige_field")
+            self.add_field(err_field, "err_field")
+            self.set_field_data(krigevar, "krige_var")
             if "mean" in info:  # ordinary krging estimates mean
-                self.mean = info["mean"]
+                self.set_field_data(info["mean"], "mean")
         else:
-            self.field = self.raw_field + self.mean
+            self.field = self["raw_field"] + self.mean
         # upscaled variance
         if not np.isscalar(point_volumes) or not np.isclose(point_volumes, 0):
             scaled_var = self.upscaling_func(self.model, point_volumes)

From 5a33d33b61744babd74cc90021e9fb03f3f1bdd0 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Fri, 24 Apr 2020 16:30:11 +0200
Subject: [PATCH 44/49] Fix a bug when internally changing mesh_type

For rotating a structured field, it is converted to an unstructured grid
internally. This caused a problem, when saving the raw_field
prematurally.
---
 gstools/field/srf.py | 11 ++++++++---
 1 file changed, 8 insertions(+), 3 deletions(-)

diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index 55f7db2a..8d95621a 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -135,12 +135,17 @@ def __call__(
             pos, mesh_type
         )
         # generate the field
-        self.add_field(self.generator.__call__(x, y, z, mt_gen), "raw_field")
+        raw_field = self.generator.__call__(x, y, z, mt_gen)
         # reshape field if we got an unstructured mesh
         if mt_changed:
-            self["raw_field"] = reshape_field_from_unstruct_to_struct(
-                self.model.dim, self["raw_field"], axis_lens
+            self.add_field(
+                reshape_field_from_unstruct_to_struct(
+                    self.model.dim, raw_field, axis_lens
+                ),
+                "raw_field",
             )
+        else:
+            self.add_field(raw_field, "raw_field")
         # apply given conditions to the field
         if self.condition:
             (

From 1bb20e9f748c5b51dfc231dcbadd5c56400b5ca6 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Wed, 29 Apr 2020 14:42:09 +0200
Subject: [PATCH 45/49] Add method for del. field_data

---
 gstools/field/mesh.py | 15 +++++++++++++++
 tests/test_mesh.py    |  8 ++++++++
 2 files changed, 23 insertions(+)

diff --git a/gstools/field/mesh.py b/gstools/field/mesh.py
index e1387052..77b3a5f2 100644
--- a/gstools/field/mesh.py
+++ b/gstools/field/mesh.py
@@ -137,6 +137,21 @@ def add_field(
         if len(self.point_data) == 1 or is_default_field:
             self.set_field_data(name, "default_field")
 
+    def del_field_data(self):
+        """Delete the field data.
+
+        Deleting the field data resets the field data dictionary and deletes
+        the attributes too.
+        """
+        for attr in self.field_data.keys():
+            if attr == "default_field" or attr == "mesh_type":
+                continue
+            try:
+                delattr(self, attr)
+            except AttributeError:
+                pass
+        self.field_data = {"default_field": "field", "mesh_type": self.mesh_type}
+
     def __getitem__(self, key):
         """:any:`numpy.ndarray`: The values of the field."""
         return self.point_data[key]
diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index 8d9425ef..90d200e7 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -135,6 +135,14 @@ def test_default_field(self):
         self.assertEqual(self.m1_grid.default_field, "test_field2")
         self.assertEqual(self.m1_grid.field[5], self.m1_grid["test_field2"][5])
 
+    def test_reset(self):
+        self.m3_tuple.set_field_data("TestLoc", "location")
+        self.m3_tuple.del_field_data()
+        self.assertEqual(self.m3_tuple.field_data["default_field"], "field")
+        self.assertEqual(self.m3_tuple.field_data["mesh_type"], "unstructured")
+        self.assertEqual(self.m3_tuple.default_field, "field")
+        self.assertEqual(self.m3_tuple.mesh_type, "unstructured")
+
     def test_point_data_check(self):
         self.assertRaises(ValueError, self.m1_tuple.add_field, self.f1_grid)
         self.assertRaises(ValueError, self.m1_tuple.add_field, self.f2_grid)

From 974ec029e062e7fe37cf20fd07d13c5382a9d3e9 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Wed, 29 Apr 2020 18:19:18 +0200
Subject: [PATCH 46/49] 1d pos is possible again without tuple

---
 gstools/field/mesh.py | 18 +++++++++++-------
 1 file changed, 11 insertions(+), 7 deletions(-)

diff --git a/gstools/field/mesh.py b/gstools/field/mesh.py
index 77b3a5f2..6fb8fbf6 100644
--- a/gstools/field/mesh.py
+++ b/gstools/field/mesh.py
@@ -460,8 +460,12 @@ def _check_point_data(self, values):
         if self.mesh_type == "unstructured":
             # scalar
             if len(values.shape) == 1:
-                if values.shape[0] == len(self.pos[0]):
-                    err = False
+                try:
+                    if values.shape[0] == len(self.pos[0]):
+                        err = False
+                except TypeError:
+                    if values.shape[0] == len(self.pos):
+                        err = False
             # vector
             elif len(values.shape) == 2:
                 if (
@@ -477,7 +481,10 @@ def _check_point_data(self, values):
                 )
         else:
             # scalar
-            if len(values.shape) == len(self.pos):
+            # 1d case
+            if values.shape[0] == len(self.pos):
+                err = False
+            elif len(values.shape) == len(self.pos):
                 if all(
                     [
                         values.shape[i] == len(self.pos[i])
@@ -496,10 +503,7 @@ def _check_point_data(self, values):
                     err = False
             if err:
                 raise ValueError(
-                    "Wrong field shape: {0} does not match mesh shape [0/2/3]{1}".format(
-                        list(values.shape),
-                        [len(self.pos[i]) for i in range(len(self.pos))],
-                    )
+                    "Wrong field shape, it does not match mesh shape"
                 )
 
 

From e9e93343771eaaab16f1ed703b9b7fa7ca668ced Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Wed, 29 Apr 2020 18:20:11 +0200
Subject: [PATCH 47/49] Refactor SRF to be more mesh-centric

---
 gstools/field/condition.py |  4 +--
 gstools/field/srf.py       | 60 +++++++++++++++++++++++++-------------
 tests/test_condition.py    |  1 -
 tests/test_srf.py          | 13 +++++++++
 4 files changed, 54 insertions(+), 24 deletions(-)

diff --git a/gstools/field/condition.py b/gstools/field/condition.py
index 55580d7d..cbb0f076 100644
--- a/gstools/field/condition.py
+++ b/gstools/field/condition.py
@@ -53,7 +53,7 @@ def ordinary(srf):
         model=srf.model, cond_pos=srf.cond_pos, cond_val=err_data
     )
     err_field, __ = err_ok(srf.pos, srf.mesh_type)
-    cond_field = srf["raw_field"] + krige_field - err_field
+    cond_field = srf.raw_field + krige_field - err_field
     info = {"mean": krige_ok.mean}
     return cond_field, krige_field, err_field, krige_var, info
 
@@ -101,6 +101,6 @@ def simple(srf):
         cond_val=err_data,
     )
     err_field, __ = err_sk(srf.pos, srf.mesh_type)
-    cond_field = srf["raw_field"] + krige_field - err_field + srf.mean
+    cond_field = srf.raw_field + krige_field - err_field + srf.mean
     info = {}
     return cond_field, krige_field, err_field, krige_var, info
diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index 8d95621a..ecf856fb 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -11,6 +11,7 @@
 """
 # pylint: disable=C0103
 
+from warnings import warn
 import numpy as np
 from gstools.field.generator import RandMeth, IncomprRandMeth
 from gstools.field.tools import reshape_field_from_unstruct_to_struct
@@ -77,12 +78,14 @@ class SRF(Field):
     def __init__(
         self,
         model,
+        pos=None,
+        mesh_type="unstructured",
         mean=0.0,
         upscaling="no_scaling",
         generator="RandMeth",
         **generator_kwargs
     ):
-        super().__init__(model, mean=mean)
+        super().__init__(model, pos=pos, mesh_type=mesh_type, mean=mean)
         # initialize private attributes
         self._generator = None
         self._upscaling = None
@@ -100,7 +103,12 @@ def __init__(
             )
 
     def __call__(
-        self, pos, seed=np.nan, point_volumes=0.0, mesh_type="unstructured"
+        self,
+        pos=None,
+        name="field",
+        seed=np.nan,
+        point_volumes=0.0,
+        mesh_type=None,
     ):
         """Generate the spatial random field.
 
@@ -108,7 +116,7 @@ def __call__(
 
         Parameters
         ----------
-        pos : :class:`list`
+        pos : :class:`list`, optional
             the position tuple, containing main direction and transversal
             directions
         seed : :class:`int`, optional
@@ -127,25 +135,33 @@ def __call__(
         field : :class:`numpy.ndarray`
             the SRF
         """
-        self.mesh_type = mesh_type
+        if mesh_type is not None:
+            warn(
+                "'mesh_type' will be ignored in the srf.__call__, use " +
+                "__init__ or setter instead.",
+                DeprecationWarning,
+            )
+            self.mesh_type = mesh_type
+        if pos is not None:
+            warn(
+                "'pos' will be ignored in the srf.__call__, use "
+                "__init__ or setter instead.",
+                DeprecationWarning,
+            )
+            self.pos = pos
         # update the model/seed in the generator if any changes were made
         self.generator.update(self.model, seed)
         # internal conversation
-        x, y, z, self.pos, mt_gen, mt_changed, axis_lens = self._pre_pos(
-            pos, mesh_type
+        x, y, z, pos, mt_gen, mt_changed, axis_lens = self._pre_pos(
+            self.pos, self.mesh_type
         )
         # generate the field
-        raw_field = self.generator.__call__(x, y, z, mt_gen)
+        self.raw_field = self.generator.__call__(x, y, z, mt_gen)
         # reshape field if we got an unstructured mesh
         if mt_changed:
-            self.add_field(
-                reshape_field_from_unstruct_to_struct(
-                    self.model.dim, raw_field, axis_lens
-                ),
-                "raw_field",
-            )
-        else:
-            self.add_field(raw_field, "raw_field")
+            self.raw_field = reshape_field_from_unstruct_to_struct(
+                    self.model.dim, self.raw_field, axis_lens
+                )
         # apply given conditions to the field
         if self.condition:
             (
@@ -156,21 +172,23 @@ def __call__(
                 info,
             ) = self.cond_func(self)
             # store everything in the class
-            self.field = cond_field
+            self.add_field(cond_field, name)
             self.add_field(krige_field, "krige_field")
             self.add_field(err_field, "err_field")
             self.set_field_data(krigevar, "krige_var")
             if "mean" in info:  # ordinary krging estimates mean
                 self.set_field_data(info["mean"], "mean")
         else:
-            self.field = self["raw_field"] + self.mean
+            self.add_field(self.raw_field + self.mean, name)
         # upscaled variance
         if not np.isscalar(point_volumes) or not np.isclose(point_volumes, 0):
             scaled_var = self.upscaling_func(self.model, point_volumes)
-            self.field -= self.mean
-            self.field *= np.sqrt(scaled_var / self.model.sill)
-            self.field += self.mean
-        return self.field
+            self[name] -= self.mean
+            self[name] *= np.sqrt(scaled_var / self.model.sill)
+            self[name] += self.mean
+        self.add_field(self.raw_field, "raw_field")
+        del(self.raw_field)
+        return self[name]
 
     def set_condition(
         self, cond_pos=None, cond_val=None, krige_type="ordinary"
diff --git a/tests/test_condition.py b/tests/test_condition.py
index 3658df89..e7c58879 100644
--- a/tests/test_condition.py
+++ b/tests/test_condition.py
@@ -46,7 +46,6 @@ def test_simple(self):
             )
             srf = SRF(model, self.mean, seed=19970221)
             srf.set_condition(self.cond_pos[0], self.cond_val, "simple")
-            # TODO should this be possible?!?
             field_1 = srf.unstructured(self.pos[0])
             field_2 = srf.structured(self.pos[0])
             for i, val in enumerate(self.cond_val):
diff --git a/tests/test_srf.py b/tests/test_srf.py
index 082a7ec3..7a3d2440 100644
--- a/tests/test_srf.py
+++ b/tests/test_srf.py
@@ -73,6 +73,19 @@ def test_shape_3d(self):
         )
         self.assertEqual(field_unstr.shape, (len(self.x_tuple),))
 
+    def test_ens(self):
+        srf = SRF(
+            self.cov_model,
+            (self.x_tuple, self.y_tuple, self.z_tuple),
+            mean=self.mean,
+            mode_no=self.mode_no,
+        )
+        for i in range(3):
+            name = "ens_{}".format(i+1)
+            srf(name=name)
+            self.assertTrue(name in srf.point_data.keys())
+        self.assertEqual(srf.default_field, "ens_1")
+
     def test_anisotropy_2d(self):
         self.cov_model.dim = 2
         srf = SRF(self.cov_model, mean=self.mean, mode_no=self.mode_no)

From 5a6a6fe734b9a98b6d521f5d439887ed3dac897d Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Thu, 30 Apr 2020 12:25:36 +0200
Subject: [PATCH 48/49] Add dim to Mesh arg. list

---
 examples/08_mesh/00_mesh_intro.py |  2 +-
 gstools/field/base.py             |  6 +++++-
 gstools/field/mesh.py             | 10 ++++++++--
 tests/test_mesh.py                | 18 ++++++++++--------
 4 files changed, 24 insertions(+), 12 deletions(-)

diff --git a/examples/08_mesh/00_mesh_intro.py b/examples/08_mesh/00_mesh_intro.py
index e0d4d196..b45bb123 100644
--- a/examples/08_mesh/00_mesh_intro.py
+++ b/examples/08_mesh/00_mesh_intro.py
@@ -39,7 +39,7 @@
 # ^^^^^^^^^^^^^^^^^
 #
 # Starting out fresh, we want to feed the mesh with our data
-mesh = gs.Mesh(pos=(x, y), values=field)
+mesh = gs.Mesh(2, pos=(x, y), values=field)
 
 # We can add meta data too
 mesh.set_field_data("Süderbrarup", "location")
diff --git a/gstools/field/base.py b/gstools/field/base.py
index 8594179f..c31ff122 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -49,7 +49,11 @@ def __init__(
     ):
         # initialize attributes
         super().__init__(
-            pos=pos, name=name, values=values, mesh_type=mesh_type
+            dim=model.dim,
+            pos=pos,
+            name=name,
+            values=values,
+            mesh_type=mesh_type,
         )
         # initialize private attributes
         self._model = None
diff --git a/gstools/field/mesh.py b/gstools/field/mesh.py
index 6fb8fbf6..649c14ec 100644
--- a/gstools/field/mesh.py
+++ b/gstools/field/mesh.py
@@ -64,7 +64,7 @@ class Mesh:
     >>> pos = np.linspace(0., 100., 40)
     >>> z = np.random.random(40)
     >>> z2 = np.random.random(40)
-    >>> mesh = Mesh((pos,))
+    >>> mesh = Mesh(1, (pos,))
     >>> mesh.add_field(z, "test_field1")
     >>> mesh.add_field(z2, "test_field2")
     >>> mesh.set_default_field("test_field2")
@@ -73,8 +73,9 @@ class Mesh:
     """
 
     def __init__(
-        self, pos=None, name="field", values=None, mesh_type="unstructured",
+        self, dim, pos=None, name="field", values=None, mesh_type="unstructured",
     ):
+        self._dim = dim
         self._vtk_export_fct = self._vtk_export_unstructured
         self._to_vtk_fct = self._to_vtk_unstructured
         # mesh_type needs a special setter, therefore, `set_field_data` is not
@@ -159,6 +160,11 @@ def __getitem__(self, key):
     def __setitem__(self, key, value):
         self.point_data[key] = value
 
+    @property
+    def dim(self):
+        """:any:`int`: The dimension of the mesh."""
+        return self._dim
+
     @property
     def pos(self):
         """:any:`numpy.ndarray`: The pos. on which the field is defined."""
diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index 90d200e7..bed2916f 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -34,14 +34,16 @@ def setUp(self):
         self.f2_tuple = self.x_tuple * self.y_tuple
         self.f3_tuple = self.x_tuple * self.y_tuple * self.z_tuple
 
-        self.m1_grid = Mesh((self.x_grid,), mesh_type="structured")
-        self.m2_grid = Mesh((self.x_grid, self.y_grid), mesh_type="structured")
+        self.m1_grid = Mesh(1, (self.x_grid,), mesh_type="structured")
+        self.m2_grid = Mesh(
+            2, (self.x_grid, self.y_grid), mesh_type="structured"
+        )
         self.m3_grid = Mesh(
-            (self.x_grid, self.y_grid, self.z_grid), mesh_type="structured"
+            3, (self.x_grid, self.y_grid, self.z_grid), mesh_type="structured"
         )
-        self.m1_tuple = Mesh((self.x_tuple,))
-        self.m2_tuple = Mesh((self.x_tuple, self.y_tuple))
-        self.m3_tuple = Mesh((self.x_tuple, self.y_tuple, self.z_tuple))
+        self.m1_tuple = Mesh(1, (self.x_tuple,))
+        self.m2_tuple = Mesh(2, (self.x_tuple, self.y_tuple))
+        self.m3_tuple = Mesh(3, (self.x_tuple, self.y_tuple, self.z_tuple))
 
     def test_item_getter_setter(self):
         self.m3_grid.add_field(256.0 * self.f3_grid, name="2nd")
@@ -178,8 +180,8 @@ def test_point_data_check(self):
         z_tuple3 = self.rng.uniform(0.0, 10, (3, 100))
         f_tuple3 = np.vstack((x_tuple2, y_tuple2, z_tuple3))
 
-        m2_tuple = Mesh((x_tuple2, y_tuple2))
-        m3_tuple = Mesh((x_tuple3, y_tuple3, z_tuple3))
+        m2_tuple = Mesh(2, (x_tuple2, y_tuple2))
+        m3_tuple = Mesh(3, (x_tuple3, y_tuple3, z_tuple3))
 
         self.assertRaises(ValueError, m2_tuple.add_field, f_tuple3)
         self.assertRaises(ValueError, m3_tuple.add_field, f_tuple2)

From 344990f3c047903d9147ec8e4418dba7cd2777f8 Mon Sep 17 00:00:00 2001
From: LSchueler <mostem@posteo.de>
Date: Wed, 29 Jul 2020 15:35:44 +0200
Subject: [PATCH 49/49] [WIP] Try to get tests to run

---
 gstools/field/base.py      |  26 +--
 gstools/field/condition.py |   8 +-
 gstools/field/mesh.py      | 331 ++++++++++++++++++++++---------------
 gstools/field/srf.py       |  19 ++-
 gstools/field/tools.py     |  48 +++++-
 gstools/krige/base.py      |  12 +-
 gstools/krige/methods.py   |   2 +-
 tests/test_condition.py    |  13 +-
 tests/test_export.py       |  18 +-
 tests/test_krige.py        |   8 +-
 tests/test_mesh.py         |  46 ++++--
 tests/test_srf.py          |  21 ++-
 12 files changed, 345 insertions(+), 207 deletions(-)

diff --git a/gstools/field/base.py b/gstools/field/base.py
index c31ff122..47dd9db2 100755
--- a/gstools/field/base.py
+++ b/gstools/field/base.py
@@ -41,7 +41,7 @@ class Field(Mesh):
     def __init__(
         self,
         model,
-        pos=None,
+        points=None,
         name="field",
         values=None,
         mesh_type="unstructured",
@@ -50,7 +50,7 @@ def __init__(
         # initialize attributes
         super().__init__(
             dim=model.dim,
-            pos=pos,
+            points=points,
             name=name,
             values=values,
             mesh_type=mesh_type,
@@ -127,7 +127,7 @@ def mesh(
                 pnts = mesh.centroids_flat.T[select]
             else:
                 pnts = mesh.NODES.T[select]
-            out = self.unstructured(pos=pnts, **kwargs)
+            out = self.unstructured(points=pnts, **kwargs)
         else:
             if points == "centroids":
                 # define unique order of cells
@@ -140,9 +140,9 @@ def mesh(
                     offset.append(pnts.shape[0])
                     length.append(pnt.shape[0])
                     pnts = np.vstack((pnts, pnt))
-                # generate pos for __call__
+                # generate points for __call__
                 pnts = pnts.T[select]
-                out = self.unstructured(pos=pnts, **kwargs)
+                out = self.unstructured(points=pnts, **kwargs)
                 if isinstance(out, np.ndarray):
                     field = out
                 else:
@@ -153,7 +153,7 @@ def mesh(
                     field_dict[cell] = field[offset[i] : offset[i] + length[i]]
                 mesh.cell_data[name] = field_dict
             else:
-                out = self.unstructured(pos=mesh.points.T[select], **kwargs)
+                out = self.unstructured(points=mesh.points.T[select], **kwargs)
                 if isinstance(out, np.ndarray):
                     field = out
                 else:
@@ -162,13 +162,13 @@ def mesh(
                 mesh.point_data[name] = field
         return out
 
-    def _pre_pos(self, pos, mesh_type="unstructured", make_unstruct=False):
+    def _pre_points(self, points, mesh_type="unstructured", make_unstruct=False):
         """
         Preprocessing positions and mesh_type.
 
         Parameters
         ----------
-        pos : :any:`iterable`
+        points : :any:`iterable`
             the position tuple, containing main direction and transversal
             directions
         mesh_type : :class:`str`
@@ -184,7 +184,7 @@ def _pre_pos(self, pos, mesh_type="unstructured", make_unstruct=False):
             analog to x
         z : :class:`numpy.ndarray` or None
             analog to x
-        pos : :class:`tuple` of :class:`numpy.ndarray`
+        points : :class:`tuple` of :class:`numpy.ndarray`
             the normalized position tuple
         mesh_type_gen : :class:`str`
             'structured' / 'unstructured' for the generator
@@ -193,11 +193,11 @@ def _pre_pos(self, pos, mesh_type="unstructured", make_unstruct=False):
         axis_lens : :class:`tuple` or :any:`None`
             axis lengths of the structured mesh if mesh type was changed.
         """
-        x, y, z = pos2xyz(pos, max_dim=self.model.dim)
-        pos = xyz2pos(x, y, z)
+        x, y, z = pos2xyz(points, max_dim=self.model.dim)
+        points = xyz2pos(x, y, z)
         mesh_type_gen = mesh_type
         # format the positional arguments of the mesh
-        check_mesh(self.model.dim, x, y, z, mesh_type)
+        check_mesh(self.model.dim, points, mesh_type)
         mesh_type_changed = False
         axis_lens = None
         if (
@@ -212,7 +212,7 @@ def _pre_pos(self, pos, mesh_type="unstructured", make_unstruct=False):
             x, y, z = unrotate_mesh(self.model.dim, self.model.angles, x, y, z)
         if not self.model.is_isotropic:
             y, z = make_isotropic(self.model.dim, self.model.anis, y, z)
-        return x, y, z, pos, mesh_type_gen, mesh_type_changed, axis_lens
+        return x, y, z, points, mesh_type_gen, mesh_type_changed, axis_lens
 
     def plot(self, field="field", fig=None, ax=None):  # pragma: no cover
         """
diff --git a/gstools/field/condition.py b/gstools/field/condition.py
index cbb0f076..37aa857a 100644
--- a/gstools/field/condition.py
+++ b/gstools/field/condition.py
@@ -40,7 +40,7 @@ def ordinary(srf):
     krige_ok = Ordinary(
         model=srf.model, cond_pos=srf.cond_pos, cond_val=srf.cond_val
     )
-    krige_field, krige_var = krige_ok(srf.pos, srf.mesh_type)
+    krige_field, krige_var = krige_ok(srf.points, srf.mesh_type)
 
     # evaluate the field at the conditional points
     x, y, z = pos2xyz(srf.cond_pos, max_dim=srf.model.dim)
@@ -52,7 +52,7 @@ def ordinary(srf):
     err_ok = Ordinary(
         model=srf.model, cond_pos=srf.cond_pos, cond_val=err_data
     )
-    err_field, __ = err_ok(srf.pos, srf.mesh_type)
+    err_field, __ = err_ok(srf.points, srf.mesh_type)
     cond_field = srf.raw_field + krige_field - err_field
     info = {"mean": krige_ok.mean}
     return cond_field, krige_field, err_field, krige_var, info
@@ -85,7 +85,7 @@ def simple(srf):
         cond_pos=srf.cond_pos,
         cond_val=srf.cond_val,
     )
-    krige_field, krige_var = krige_sk(srf.pos, srf.mesh_type)
+    krige_field, krige_var = krige_sk(srf.points, srf.mesh_type)
 
     # evaluate the field at the conditional points
     x, y, z = pos2xyz(srf.cond_pos, max_dim=srf.model.dim)
@@ -100,7 +100,7 @@ def simple(srf):
         cond_pos=srf.cond_pos,
         cond_val=err_data,
     )
-    err_field, __ = err_sk(srf.pos, srf.mesh_type)
+    err_field, __ = err_sk(srf.points, srf.mesh_type)
     cond_field = srf.raw_field + krige_field - err_field + srf.mean
     info = {}
     return cond_field, krige_field, err_field, krige_var, info
diff --git a/gstools/field/mesh.py b/gstools/field/mesh.py
index 649c14ec..e322d73a 100644
--- a/gstools/field/mesh.py
+++ b/gstools/field/mesh.py
@@ -16,6 +16,7 @@
 from pyevtk.hl import gridToVTK, pointsToVTK
 
 from gstools.tools.geometric import pos2xyz
+from gstools.field.tools import check_mesh, check_point_data
 
 __all__ = ["Mesh"]
 
@@ -34,6 +35,21 @@ def value_type(mesh_type, shape):
             r = "vector"
     return r
 
+# TODO figure out what this is all about
+def convert_points(dim, points):
+    points = np.array(points)
+    shape = points.shape
+    if dim == 1:
+        #points = points.flatten()
+        if len(points.shape) > 1:
+            raise ValueError(
+                "points shape {} does not match dim {}".format(shape, dim)
+            )
+    else:
+        pass
+
+    return points
+
 
 class Mesh:
     """A base class encapsulating field data.
@@ -71,19 +87,15 @@ class Mesh:
     >>> print(mesh.field)
 
     """
-
     def __init__(
-        self, dim, pos=None, name="field", values=None, mesh_type="unstructured",
+        self, dim, points=None, name="field", values=None, mesh_type="unstructured",
     ):
         self._dim = dim
-        self._vtk_export_fct = self._vtk_export_unstructured
-        self._to_vtk_fct = self._to_vtk_unstructured
-        # mesh_type needs a special setter, therefore, `set_field_data` is not
-        # used here
-        self.mesh_type = mesh_type
 
         # the pos/ points of the mesh
-        self._pos = pos
+        if points is not None:
+            check_mesh(dim, points, mesh_type)
+        self._points = points
 
         # data stored at each pos/ point, the "fields"
         if values is not None:
@@ -94,9 +106,12 @@ def __init__(
         # data valid for the global field
         self.field_data = {}
 
-        self.set_field_data(name, "default_field")
+        # do following line manually in order to satisfy the linters
+        # self.set_field_data(name, "default_field")
+        self.default_field = name
+        self.field_data["default_field"] = name
 
-        self.field_data["mesh_type"] = mesh_type
+        self.mesh_type = mesh_type
 
     def set_field_data(self, value, name):
         """Add an attribute to this instance and add it the `field_data`
@@ -132,11 +147,22 @@ def add_field(
 
         """
         values = np.array(values)
-        self._check_point_data(values)
-        self.point_data[name] = values
+        check_point_data(
+            self.dim,
+            self.points,
+            values,
+            self.mesh_type,
+            value_type(self.mesh_type, values.shape),
+        )
         # set the default_field to the first field added
-        if len(self.point_data) == 1 or is_default_field:
+        if not self.point_data:
+            is_default_field = True
+        elif self.point_data.get(self.default_field, False) is None:
+            del(self.point_data[self.default_field])
+            is_default_field = True
+        if is_default_field:
             self.set_field_data(name, "default_field")
+        self.point_data[name] = values
 
     def del_field_data(self):
         """Delete the field data.
@@ -166,17 +192,18 @@ def dim(self):
         return self._dim
 
     @property
-    def pos(self):
+    def points(self):
         """:any:`numpy.ndarray`: The pos. on which the field is defined."""
-        return self._pos
+        return self._points
 
-    @pos.setter
-    def pos(self, value):
+    @points.setter
+    def points(self, value):
         """
         Warning: setting new positions deletes all previously stored fields.
         """
         self.point_data = {self.default_field: None}
-        self._pos = value
+        check_mesh(self.dim, value, self.mesh_type)
+        self._points = value
 
     @property
     def field(self):
@@ -185,13 +212,22 @@ def field(self):
 
     @field.setter
     def field(self, values):
-        self._check_point_data(values)
+        check_point_data(
+            self.dim,
+            self.points,
+            values,
+            self.mesh_type,
+            value_type(self.mesh_type, values.shape),
+        )
         self.point_data[self.default_field] = values
 
     @property
     def value_type(self):
         """:any:`str`: The value type of the default field."""
-        if self.default_field in self.point_data:
+        if (
+            self.default_field in self.point_data and
+            self.point_data[self.default_field] is not None
+        ):
             r = value_type(
                 self.mesh_type, self.point_data[self.default_field].shape
             )
@@ -218,6 +254,7 @@ def mesh_type(self, value):
             self._to_vtk_fct = self._to_vtk_unstructured
         self.point_data = {}
         self._mesh_type = value
+        self.field_data["mesh_type"] = value
 
     def _vtk_naming_helper(
         self, field_select="field", fieldname="field"
@@ -249,11 +286,11 @@ def _vtk_naming_helper(
             else:
                 field = None
             if not (
-                self.pos is None or field is None or self.mesh_type is None
+                self.points is None or field is None or self.mesh_type is None
             ):
                 suf = ["_X", "_Y", "_Z"]
                 fields = {}
-                for i in range(self.model.dim):
+                for i in range(self.dim):
                     fields[fieldname + suf[i]] = field[i]
                 return fields
         elif self.value_type == "scalar":
@@ -262,7 +299,7 @@ def _vtk_naming_helper(
             else:
                 field = None
             if not (
-                self.pos is None or field is None or self.mesh_type is None
+                self.points is None or field is None or self.mesh_type is None
             ):
                 return {fieldname: field}
             else:
@@ -276,13 +313,16 @@ def _vtk_naming_helper(
                 "Unknown field value type: {}".format(self.value_type)
             )
 
-    def to_pyvista(
-        self, field_select="field", fieldname="field"
+    def convert(
+        self, data_format="pyvista", field_select="field", fieldname="field"
     ):  # pragma: no cover
         """Create a VTK/PyVista grid of the stored field.
 
         Parameters
         ----------
+        data_format : :class:`str`
+            the format in which to convert the data, possible choices:
+            * 'pyvista'
         field_select : :class:`str`, optional
             Field that should be stored. Can be:
             "field", "raw_field", "krige_field", "err_field" or "krige_var".
@@ -290,6 +330,10 @@ def to_pyvista(
         fieldname : :class:`str`, optional
             Name of the field in the VTK file. Default: "field"
         """
+        if data_format.lower() != "pyvista":
+            raise NotImplementedError(
+                "Only pyvista format is supported at the moment."
+            )
         field_names = self._vtk_naming_helper(
             field_select=field_select, fieldname=fieldname
         )
@@ -297,8 +341,8 @@ def to_pyvista(
         grid = self._to_vtk_fct(field_names)
         return grid
 
-    def vtk_export(
-        self, filename, field_select="field", fieldname="field"
+    def export(
+        self, filename, data_format="vtk", field_select="field", fieldname="field"
     ):  # pragma: no cover
         """Export the stored field to vtk.
 
@@ -307,6 +351,9 @@ def vtk_export(
         filename : :class:`str`
             Filename of the file to be saved, including the path. Note that an
             ending (.vtr or .vtu) will be added to the name.
+        data_format : :class:`str`
+            the format in which to export the data, possible choices:
+            * 'vtk'
         field_select : :class:`str`, optional
             Field that should be stored. Can be:
             "field", "raw_field", "krige_field", "err_field" or "krige_var".
@@ -314,6 +361,10 @@ def vtk_export(
         fieldname : :class:`str`, optional
             Name of the field in the VTK file. Default: "field"
         """
+        if data_format.lower() != "vtk":
+            raise NotImplementedError(
+                "Only VTK format is supported at the moment."
+            )
         if not isinstance(filename, str):
             raise TypeError("Please use a string as a filename.")
         field_names = self._vtk_naming_helper(
@@ -321,118 +372,95 @@ def vtk_export(
         )
         return self._vtk_export_fct(filename, field_names)
 
-    def _vtk_export_structured(self, filename, fields):
-        """Export a field to vtk structured rectilinear grid file.
+    def _vtk_export_unstructured(self, filename, fields):
+        """Export a field to vtk unstructured grid file.
 
         Parameters
         ----------
         filename : :class:`str`
             Filename of the file to be saved, including the path. Note that an
-            ending (.vtr) will be added to the name.
+            ending (.vtu) will be added to the name.
         fields : :class:`dict` or :class:`numpy.ndarray`
-            Structured fields to be saved.
+            Unstructured fields to be saved.
             Either a single numpy array as returned by SRF,
             or a dictionary of fields with theirs names as keys.
         """
-        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
-        return gridToVTK(filename, x, y, z, pointData=fields)
+        x, y, z, fields = self._vtk_unstructured_reshape(fields)
+        return pointsToVTK(filename, x, y, z, data=fields)
 
-    def _vtk_export_unstructured(self, filename, fields):
-        """Export a field to vtk unstructured grid file.
+    def _vtk_export_structured(self, filename, fields):
+        """Export a field to vtk structured rectilinear grid file.
 
         Parameters
         ----------
         filename : :class:`str`
             Filename of the file to be saved, including the path. Note that an
-            ending (.vtu) will be added to the name.
+            ending (.vtr) will be added to the name.
         fields : :class:`dict` or :class:`numpy.ndarray`
-            Unstructured fields to be saved.
+            Structured fields to be saved.
             Either a single numpy array as returned by SRF,
             or a dictionary of fields with theirs names as keys.
         """
-        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
-        return pointsToVTK(filename, x, y, z, data=fields)
+        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
+        return gridToVTK(filename, x, y, z, pointData=fields)
 
-    def _to_vtk_structured(self, fields):  # pragma: no cover
-        """Create a vtk structured rectilinear grid from a field.
+    def _to_vtk_unstructured(self, fields):  # pragma: no cover
+        """Export a field to vtk structured rectilinear grid file.
 
         Parameters
         ----------
-        pos : :class:`list`
-            the position tuple, containing main direction and transversal
-            directions
         fields : :class:`dict` or :class:`numpy.ndarray`
-            Structured fields to be saved.
+            Unstructured fields to be saved.
             Either a single numpy array as returned by SRF,
             or a dictionary of fields with theirs names as keys.
 
         Returns
         -------
-        :class:`pyvista.RectilinearGrid`
-            A PyVista rectilinear grid of the structured field data. Data arrays
-            live on the point data of this PyVista dataset.
+        :class:`pyvista.UnstructuredGrid`
+            A PyVista unstructured grid of the unstructured field data. Data arrays
+            live on the point data of this PyVista dataset. This is essentially
+            a point cloud with no topology.
         """
-        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
+        x, y, z, fields = self._vtk_unstructured_reshape(fields)
         try:
             import pyvista as pv
 
-            grid = pv.RectilinearGrid(x, y, z)
+            grid = pv.PolyData(np.c_[x, y, z]).cast_to_unstructured_grid()
             grid.point_arrays.update(fields)
         except ImportError:
             raise ImportError("Please install PyVista to create VTK datasets.")
         return grid
 
-    def _to_vtk_unstructured(self, fields):
-        """Export a field to vtk structured rectilinear grid file.
+    def _to_vtk_structured(self, fields):  # pragma: no cover
+        """Create a vtk structured rectilinear grid from a field.
 
         Parameters
         ----------
         fields : :class:`dict` or :class:`numpy.ndarray`
-            Unstructured fields to be saved.
+            Structured fields to be saved.
             Either a single numpy array as returned by SRF,
             or a dictionary of fields with theirs names as keys.
 
         Returns
         -------
-        :class:`pyvista.UnstructuredGrid`
-            A PyVista unstructured grid of the unstructured field data. Data arrays
-            live on the point data of this PyVista dataset. This is essentially
-            a point cloud with no topology.
+        :class:`pyvista.RectilinearGrid`
+            A PyVista rectilinear grid of the structured field data. Data arrays
+            live on the point data of this PyVista dataset.
         """
-        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
+        x, y, z, fields = self._vtk_structured_reshape(fields=fields)
         try:
             import pyvista as pv
 
-            grid = pv.PolyData(np.c_[x, y, z]).cast_to_unstructured_grid()
+            grid = pv.RectilinearGrid(x, y, z)
             grid.point_arrays.update(fields)
         except ImportError:
             raise ImportError("Please install PyVista to create VTK datasets.")
         return grid
 
-    def _vtk_structured_reshape(self, fields):
-        """An internal helper to extract what is needed for the vtk rectilinear grid
-        """
-        if not isinstance(fields, dict):
-            fields = {"field": fields}
-        x, y, z = pos2xyz(self.pos)
-        if y is None:
-            y = np.array([0])
-        if z is None:
-            z = np.array([0])
-        # need fortran order in VTK
-        for field in fields:
-            fields[field] = fields[field].reshape(-1, order="F")
-            if len(fields[field]) != len(x) * len(y) * len(z):
-                raise ValueError(
-                    "gstools.vtk_export_structured: "
-                    "field shape doesn't match the given mesh"
-                )
-        return x, y, z, fields
-
     def _vtk_unstructured_reshape(self, fields):
         if not isinstance(fields, dict):
             fields = {"field": fields}
-        x, y, z = pos2xyz(self.pos)
+        x, y, z = pos2xyz(self.points)
         if y is None:
             y = np.zeros_like(x)
         if z is None:
@@ -450,67 +478,98 @@ def _vtk_unstructured_reshape(self, fields):
                 )
         return x, y, z, fields
 
+    def _vtk_structured_reshape(self, fields):
+        """An internal helper to extract what is needed for the vtk rectilinear grid
+        """
+        if not isinstance(fields, dict):
+            fields = {"field": fields}
+        x, y, z = pos2xyz(self.points)
+        if y is None:
+            y = np.array([0])
+        if z is None:
+            z = np.array([0])
+        # need fortran order in VTK
+        for field in fields:
+            fields[field] = fields[field].reshape(-1, order="F")
+            if len(fields[field]) != len(x) * len(y) * len(z):
+                raise ValueError(
+                    "gstools.vtk_export_structured: "
+                    "field shape doesn't match the given mesh"
+                )
+        return x, y, z, fields
+
     def _check_mesh_type(self, mesh_type):
         if mesh_type != "unstructured" and mesh_type != "structured":
             raise ValueError("Unknown 'mesh_type': {}".format(mesh_type))
 
-    def _check_point_data(self, values):
-        """Compare field shape to pos shape.
 
-        Parameters
-        ----------
-        values : :class:`numpy.ndarray`
-            the values of the field to be checked
-        """
-        err = True
-        if self.mesh_type == "unstructured":
-            # scalar
-            if len(values.shape) == 1:
-                try:
-                    if values.shape[0] == len(self.pos[0]):
-                        err = False
-                except TypeError:
-                    if values.shape[0] == len(self.pos):
-                        err = False
-            # vector
-            elif len(values.shape) == 2:
-                if (
-                    values.shape[1] == len(self.pos[0])
-                    and 2 <= values.shape[0] <= 3
-                ):
-                    err = False
-            if err:
-                raise ValueError(
-                    "Wrong field shape: {0} does not match mesh shape ({1},)".format(
-                        values.shape, len(self.pos[0])
-                    )
-                )
-        else:
-            # scalar
-            # 1d case
-            if values.shape[0] == len(self.pos):
-                err = False
-            elif len(values.shape) == len(self.pos):
-                if all(
-                    [
-                        values.shape[i] == len(self.pos[i])
-                        for i in range(len(self.pos))
-                    ]
-                ):
-                    err = False
-            # vector
-            elif len(values.shape) == len(self.pos) + 1:
-                if all(
-                    [
-                        values.shape[i + 1] == len(self.pos[i])
-                        for i in range(len(self.pos))
-                    ]
-                ) and values.shape[0] == len(self.pos):
-                    err = False
-            if err:
-                raise ValueError(
-                    "Wrong field shape, it does not match mesh shape"
-                )
+#class Mesh(MeshBase):
+#    def __init__(self):
+#        super().__init__()
+#
+#    def _vtk_reshape(self, fields):
+#        if not isinstance(fields, dict):
+#            fields = {"field": fields}
+#        x, y, z = pos2xyz(self.pos)
+#        if y is None:
+#            y = np.zeros_like(x)
+#        if z is None:
+#            z = np.zeros_like(x)
+#        for field in fields:
+#            fields[field] = fields[field].reshape(-1)
+#            if (
+#                len(fields[field]) != len(x)
+#                or len(fields[field]) != len(y)
+#                or len(fields[field]) != len(z)
+#            ):
+#                raise ValueError(
+#                    "gstools.vtk_export_unstructured: "
+#                    "field shape doesn't match the given mesh"
+#                )
+#        return x, y, z, fields
+#
+#    def _vtk_export(self, filename, fields):
+#        """Export a field to vtk unstructured grid file.
+#
+#        Parameters
+#        ----------
+#        filename : :class:`str`
+#            Filename of the file to be saved, including the path. Note that an
+#            ending (.vtu) will be added to the name.
+#        fields : :class:`dict` or :class:`numpy.ndarray`
+#            Unstructured fields to be saved.
+#            Either a single numpy array as returned by SRF,
+#            or a dictionary of fields with theirs names as keys.
+#        """
+#        x, y, z, fields = self._vtk_reshape(fields=fields)
+#        return pointsToVTK(filename, x, y, z, data=fields)
+#
+#    def _to_vtk(self, fields):
+#        """Export a field to vtk structured rectilinear grid file.
+#
+#        Parameters
+#        ----------
+#        fields : :class:`dict` or :class:`numpy.ndarray`
+#            Unstructured fields to be saved.
+#            Either a single numpy array as returned by SRF,
+#            or a dictionary of fields with theirs names as keys.
+#
+#        Returns
+#        -------
+#        :class:`pyvista.UnstructuredGrid`
+#            A PyVista unstructured grid of the unstructured field data. Data arrays
+#            live on the point data of this PyVista dataset. This is essentially
+#            a point cloud with no topology.
+#        """
+#        x, y, z, fields = self._vtk_unstructured_reshape(fields=fields)
+#        try:
+#            import pyvista as pv
+#
+#            grid = pv.PolyData(np.c_[x, y, z]).cast_to_unstructured_grid()
+#            grid.point_arrays.update(fields)
+#        except ImportError:
+#            raise ImportError("Please install PyVista to create VTK datasets.")
+#        return grid
 
 
 if __name__ == "__main__":  # pragma: no cover
diff --git a/gstools/field/srf.py b/gstools/field/srf.py
index ecf856fb..c5fda834 100644
--- a/gstools/field/srf.py
+++ b/gstools/field/srf.py
@@ -78,19 +78,20 @@ class SRF(Field):
     def __init__(
         self,
         model,
-        pos=None,
+        points=None,
         mesh_type="unstructured",
         mean=0.0,
         upscaling="no_scaling",
         generator="RandMeth",
         **generator_kwargs
     ):
-        super().__init__(model, pos=pos, mesh_type=mesh_type, mean=mean)
+        super().__init__(model, points=points, mesh_type=mesh_type, mean=mean)
         # initialize private attributes
         self._generator = None
         self._upscaling = None
         self._upscaling_func = None
         # condition related
+        # TODO rename to self._cond_points
         self._cond_pos = None
         self._cond_val = None
         self._krige_type = None
@@ -104,7 +105,7 @@ def __init__(
 
     def __call__(
         self,
-        pos=None,
+        points=None,
         name="field",
         seed=np.nan,
         point_volumes=0.0,
@@ -116,7 +117,7 @@ def __call__(
 
         Parameters
         ----------
-        pos : :class:`list`, optional
+        points : :class:`list`, optional
             the position tuple, containing main direction and transversal
             directions
         seed : :class:`int`, optional
@@ -142,18 +143,18 @@ def __call__(
                 DeprecationWarning,
             )
             self.mesh_type = mesh_type
-        if pos is not None:
+        if points is not None:
             warn(
-                "'pos' will be ignored in the srf.__call__, use "
+                "'points' will be ignored in the srf.__call__, use "
                 "__init__ or setter instead.",
                 DeprecationWarning,
             )
-            self.pos = pos
+            self.points = points
         # update the model/seed in the generator if any changes were made
         self.generator.update(self.model, seed)
         # internal conversation
-        x, y, z, pos, mt_gen, mt_changed, axis_lens = self._pre_pos(
-            self.pos, self.mesh_type
+        x, y, z, points, mt_gen, mt_changed, axis_lens = self._pre_points(
+            self.points, self.mesh_type
         )
         # generate the field
         self.raw_field = self.generator.__call__(x, y, z, mt_gen)
diff --git a/gstools/field/tools.py b/gstools/field/tools.py
index e9f289a2..4de6a7f1 100644
--- a/gstools/field/tools.py
+++ b/gstools/field/tools.py
@@ -74,34 +74,34 @@ def reshape_input_axis_from_struct(x, y=None, z=None):
 # SRF helpers #################################################################
 
 
-def check_mesh(dim, x, y, z, mesh_type):
+def check_mesh(dim, points, mesh_type):
     """Do a basic check of the shapes of the input arrays."""
     if dim >= 2:
-        if y is None:
+        if len(points) < 2:
             raise ValueError(
                 "The y-component is missing for " "{0} dimensions".format(dim)
             )
     if dim == 3:
-        if z is None:
+        if len(points) < 3:
             raise ValueError(
                 "The z-component is missing for " "{0} dimensions".format(dim)
             )
     if mesh_type == "unstructured":
         if dim >= 2:
             try:
-                if len(x) != len(y):
+                if len(points[0]) != len(points[1]):
                     raise ValueError(
                         "len(x) = {0} != len(y) = {1} "
-                        "for unstructured grids".format(len(x), len(y))
+                        "for unstructured grids".format(len(points[0]), len(points[1]))
                     )
             except TypeError:
                 pass
             if dim == 3:
                 try:
-                    if len(x) != len(z):
+                    if len(points[0]) != len(points[2]):
                         raise ValueError(
                             "len(x) = {0} != len(z) = {1} "
-                            "for unstructured grids".format(len(x), len(z))
+                            "for unstructured grids".format(len(points[0]), len(points[2]))
                         )
                 except TypeError:
                     pass
@@ -110,6 +110,40 @@ def check_mesh(dim, x, y, z, mesh_type):
     else:
         raise ValueError("Unknown mesh type {0}".format(mesh_type))
 
+def check_point_data(dim, points, data, mesh_type, value_type='scalar'):
+    """Do a basic check and compare shapes of the field data to mesh field."""
+    offset = 0
+    if value_type == 'vector':
+        if dim != data.shape[0]:
+            raise ValueError(
+                "Field is {} dim., but mesh is {}".
+                format(data.shape[0], dim)
+            )
+        offset = 1
+    if mesh_type == 'unstructured':
+        for d in range(dim):
+            if points[d].shape[0] != data.shape[offset]:
+                raise ValueError(
+                    "Field shape {} does not match existing mesh shape {}".
+                    format(data.shape[offset], points[d].shape[0])
+                )
+    elif mesh_type == 'structured':
+        try:
+            if dim != len(data.squeeze().shape)-offset:
+                raise ValueError(
+                    "Field is {} dim., but mesh {}".
+                    format(len(data.squeeze().shape), dim)
+                )
+        except AttributeError:
+            pass
+        for d in range(dim):
+            if points[d].shape[0] != data.shape[offset+d]:
+                raise ValueError(
+                    "Field shape {} does not match existing mesh shape {}".
+                    format(data.shape[offset+d], points[d].shape[0])
+                )
+    else:
+        raise ValueError("Unknown mesh type {0}".format(mesh_type))
 
 def make_isotropic(dim, anis, y, z):
     """Stretch given axes in order to implement anisotropy."""
diff --git a/gstools/krige/base.py b/gstools/krige/base.py
index 99c38a64..32e9d599 100755
--- a/gstools/krige/base.py
+++ b/gstools/krige/base.py
@@ -88,7 +88,7 @@ def __init__(
         self.set_condition(cond_pos, cond_val, ext_drift)
 
     def __call__(
-        self, pos, mesh_type="unstructured", ext_drift=None, chunk_size=None
+        self, points, mesh_type="unstructured", ext_drift=None, chunk_size=None
     ):
         """
         Generate the kriging field.
@@ -98,7 +98,7 @@ def __call__(
 
         Parameters
         ----------
-        pos : :class:`list`
+        points : :class:`list`
             the position tuple, containing main direction and transversal
             directions (x, [y, z])
         mesh_type : :class:`str`, optional
@@ -119,8 +119,8 @@ def __call__(
         """
         self.mesh_type = mesh_type
         # internal conversation
-        x, y, z, self.pos, __, mt_changed, axis_lens = self._pre_pos(
-            pos, mesh_type, make_unstruct=True
+        x, y, z, self.points, __, mt_changed, axis_lens = self._pre_points(
+            points, mesh_type, make_unstruct=True
         )
         point_no = len(x)
         # set chunk size
@@ -215,7 +215,7 @@ def _post_field(self, field, krige_var):
             self.field = field
         else:
             self.field = field + eval_func(
-                self.trend_function, self.pos, self.mesh_type
+                self.trend_function, self.points, self.mesh_type
             )
         self.krige_var = krige_var
 
@@ -309,7 +309,7 @@ def set_drift_functions(self, drift_functions=None):
 
     def update(self):
         """Update the kriging settings."""
-        x, y, z, __, __, __, __ = self._pre_pos(self.cond_pos)
+        x, y, z, __, __, __, __ = self._pre_points(self.cond_pos)
         # krige pos are the unrotated and isotropic condition positions
         self._krige_pos = (x, y, z)[: self.model.dim]
         self._krige_mat = self._get_krige_mat()
diff --git a/gstools/krige/methods.py b/gstools/krige/methods.py
index bad1b059..2be07b64 100644
--- a/gstools/krige/methods.py
+++ b/gstools/krige/methods.py
@@ -83,7 +83,7 @@ def _post_field(self, field, krige_var):
             self.field = (
                 field
                 + self.mean
-                + eval_func(self.trend_function, self.pos, self.mesh_type)
+                + eval_func(self.trend_function, self.points, self.mesh_type)
             )
         # add the given mean
         self.krige_var = self.model.sill - krige_var
diff --git a/tests/test_condition.py b/tests/test_condition.py
index e7c58879..5544366b 100644
--- a/tests/test_condition.py
+++ b/tests/test_condition.py
@@ -44,10 +44,11 @@ def test_simple(self):
             model = Model(
                 dim=1, var=0.5, len_scale=2, anis=[0.1, 1], angles=[0.5, 0, 0]
             )
+            # import pudb; pu.db
             srf = SRF(model, self.mean, seed=19970221)
             srf.set_condition(self.cond_pos[0], self.cond_val, "simple")
-            field_1 = srf.unstructured(self.pos[0])
-            field_2 = srf.structured(self.pos[0])
+            field_1 = srf.unstructured((self.pos[0],))
+            field_2 = srf.structured((self.pos[0],))
             for i, val in enumerate(self.cond_val):
                 self.assertAlmostEqual(val, field_1[i], places=2)
                 self.assertAlmostEqual(val, field_2[(i,)], places=2)
@@ -60,7 +61,7 @@ def test_simple(self):
                     anis=[0.1, 1],
                     angles=[0.5, 0, 0],
                 )
-                srf = SRF(model, self.mean, seed=19970221)
+                srf = SRF(model, mean=self.mean, seed=19970221)
                 srf.set_condition(self.cond_pos[:dim], self.cond_val, "simple")
                 field_1 = srf.unstructured(self.pos[:dim])
                 field_2 = srf.structured(self.pos[:dim])
@@ -73,10 +74,10 @@ def test_ordinary(self):
             model = Model(
                 dim=1, var=0.5, len_scale=2, anis=[0.1, 1], angles=[0.5, 0, 0]
             )
-            srf = SRF(model, seed=19970221)
+            srf = SRF(model, points=(self.pos[0],), seed=19970221)
             srf.set_condition(self.cond_pos[0], self.cond_val, "ordinary")
-            field_1 = srf.unstructured(self.pos[0])
-            field_2 = srf.structured(self.pos[0])
+            field_1 = srf.unstructured()
+            field_2 = srf.structured()
             for i, val in enumerate(self.cond_val):
                 self.assertAlmostEqual(val, field_1[i], places=2)
                 self.assertAlmostEqual(val, field_2[(i,)], places=2)
diff --git a/tests/test_export.py b/tests/test_export.py
index 15bb8275..c9f25231 100644
--- a/tests/test_export.py
+++ b/tests/test_export.py
@@ -20,7 +20,7 @@
 class TestExport(unittest.TestCase):
     def setUp(self):
         self.test_dir = tempfile.mkdtemp()
-        self.x_grid = self.y_grid = self.z_grid = range(25)
+        self.x_grid = self.y_grid = self.z_grid = np.arange(25)
         model = gs.Gaussian(dim=3, var=0.6, len_scale=20)
         self.srf_structured = gs.SRF(model)
         self.srf_structured(
@@ -43,22 +43,22 @@ def tearDown(self):
 
     @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
     def test_pyvista_struct(self):
-        mesh = self.srf_structured.to_pyvista()
+        mesh = self.srf_structured.convert()
         self.assertIsInstance(mesh, pv.RectilinearGrid)
 
     @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
     def test_pyvista_unstruct(self):
-        mesh = self.srf_unstructured.to_pyvista()
+        mesh = self.srf_unstructured.convert()
         self.assertIsInstance(mesh, pv.UnstructuredGrid)
 
     def test_pyevtk_export_struct(self):
         filename = os.path.join(self.test_dir, "structured")
-        self.srf_structured.vtk_export(filename)
+        self.srf_structured.export(filename, "vtk")
         self.assertTrue(os.path.isfile(filename + ".vtr"))
 
     def test_pyevtk_export_unstruct(self):
         filename = os.path.join(self.test_dir, "unstructured")
-        self.srf_unstructured.vtk_export(filename)
+        self.srf_unstructured.export(filename, "vtk")
         self.assertTrue(os.path.isfile(filename + ".vtu"))
 
     @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
@@ -66,7 +66,7 @@ def test_pyvista_vector_struct(self):
         model = gs.Gaussian(dim=2, var=1, len_scale=10)
         srf = gs.SRF(model, generator="VectorField")
         srf((self.x_grid, self.y_grid), mesh_type="structured", seed=19841203)
-        mesh = srf.to_pyvista()
+        mesh = srf.convert()
         self.assertIsInstance(mesh, pv.RectilinearGrid)
 
     @unittest.skipIf(not HAS_PYVISTA, "PyVista is not installed.")
@@ -74,7 +74,7 @@ def test_pyvista_vector_unstruct(self):
         model = gs.Gaussian(dim=2, var=1, len_scale=10)
         srf = gs.SRF(model, generator="VectorField")
         srf((self.x_tuple, self.y_tuple), mesh_type="unstructured", seed=19841203)
-        mesh = srf.to_pyvista()
+        mesh = srf.convert()
         self.assertIsInstance(mesh, pv.UnstructuredGrid)
 
     def test_pyevtk_vector_export_struct(self):
@@ -82,7 +82,7 @@ def test_pyevtk_vector_export_struct(self):
         model = gs.Gaussian(dim=2, var=1, len_scale=10)
         srf = gs.SRF(model, generator="VectorField")
         srf((self.x_grid, self.y_grid), mesh_type="structured", seed=19841203)
-        srf.vtk_export(filename)
+        srf.export(filename, "vtk")
         self.assertTrue(os.path.isfile(filename + ".vtr"))
 
     def test_pyevtk_vector_export_unstruct(self):
@@ -90,7 +90,7 @@ def test_pyevtk_vector_export_unstruct(self):
         model = gs.Gaussian(dim=2, var=1, len_scale=10)
         srf = gs.SRF(model, generator="VectorField")
         srf((self.x_tuple, self.y_tuple), mesh_type="unstructured", seed=19841203)
-        srf.vtk_export(filename)
+        srf.export(filename, "vtk")
         self.assertTrue(os.path.isfile(filename + ".vtu"))
 
 
diff --git a/tests/test_krige.py b/tests/test_krige.py
index 3afee2b0..75886c90 100644
--- a/tests/test_krige.py
+++ b/tests/test_krige.py
@@ -156,8 +156,12 @@ def test_extdrift(self):
                 anis=[0.9, 0.8],
                 angles=[2, 1, 0.5],
             )
-            srf = SRF(model)
-            field = srf(grid)
+            # TODO should meshgrids be accepted?!
+            if dim > 1:
+                grid = [g.flatten() for g in grid]
+            # TODO the pre_points shit has to be applied here too
+            srf = SRF(model, grid)
+            field = srf()
             ext_drift.append(field)
             field = field.reshape(self.grid_shape[:dim])
             cond_drift.append(field[self.data_idx[:dim]])
diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index bed2916f..28409745 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -7,6 +7,7 @@
 import numpy as np
 
 from gstools import Mesh
+from gstools.field.mesh import convert_points
 
 
 class TestMesh(unittest.TestCase):
@@ -45,6 +46,31 @@ def setUp(self):
         self.m2_tuple = Mesh(2, (self.x_tuple, self.y_tuple))
         self.m3_tuple = Mesh(3, (self.x_tuple, self.y_tuple, self.z_tuple))
 
+    def test_convert_points(self):
+        pass
+        # TODO figure out what this is all about
+        #1d
+        #p1_target = np.array((1, 2, 3, 4)).reshape(4, 1)
+        #print(p1_target.shape)
+        #p1_is = convert_points(1, p1_target)
+        #np.testing.assert_array_equal(p1_is, p1_target)
+
+        #p2 = [1, 2, 3, 4, 5]
+        #p2_target = np.array(p2)
+        #p2_is = convert_points(1, p2)
+        #np.testing.assert_array_equal(p2_is, p2_target)
+
+        #p3 = [[1], [2], [3], [4], [5]]
+        #p3_target = np.array(p3).reshape(-1)
+        #p3_is = convert_points(1, p3)
+        #np.testing.assert_array_equal(p3_is, p3_target)
+
+        ##2d
+        #p4 = [[1, 2], [3, 4], [5, 6]]
+        #p4_target = np.array(p4)
+        #print(p4_target.shape)
+
+
     def test_item_getter_setter(self):
         self.m3_grid.add_field(256.0 * self.f3_grid, name="2nd")
         self.m3_grid.add_field(512.0 * self.f3_grid, name="3rd")
@@ -60,15 +86,15 @@ def test_item_getter_setter(self):
             self.m3_tuple["tmp_field"].all(), (2.0 * self.f3_tuple).all()
         )
 
-    def test_pos_getter(self):
-        self.assertEqual(self.m1_grid.pos, (self.x_grid,))
-        self.assertEqual(self.m1_tuple.pos, (self.x_tuple,))
-        self.assertEqual(self.m2_grid.pos, (self.x_grid, self.y_grid))
+    def test_points_getter(self):
+        self.assertEqual(self.m1_grid.points, (self.x_grid,))
+        self.assertEqual(self.m1_tuple.points, (self.x_tuple,))
+        self.assertEqual(self.m2_grid.points, (self.x_grid, self.y_grid))
         self.assertEqual(
-            self.m3_grid.pos, (self.x_grid, self.y_grid, self.z_grid)
+            self.m3_grid.points, (self.x_grid, self.y_grid, self.z_grid)
         )
         self.assertEqual(
-            self.m3_tuple.pos, (self.x_tuple, self.y_tuple, self.z_tuple)
+            self.m3_tuple.points, (self.x_tuple, self.y_tuple, self.z_tuple)
         )
 
     def test_default_value_type(self):
@@ -86,15 +112,15 @@ def test_field_data_setter(self):
         self.assertEqual(self.m2_tuple.field_data["mean"], 3.14)
         self.assertEqual(self.m2_tuple.mean, 3.14)
 
-    def test_new_pos(self):
-        # set new pos. (which causes reset)
+    def test_new_points(self):
+        # set new points (which causes reset)
         x_tuple2 = self.rng.uniform(0.0, 10, 100)
         y_tuple2 = self.rng.uniform(0.0, 10, 100)
         self.m2_tuple.add_field(self.f2_tuple)
 
-        self.m2_tuple.pos = (x_tuple2, y_tuple2)
+        self.m2_tuple.points = (x_tuple2, y_tuple2)
 
-        self.assertEqual(self.m2_tuple.pos, (x_tuple2, y_tuple2))
+        self.assertEqual(self.m2_tuple.points, (x_tuple2, y_tuple2))
 
         # previous field has to be deleted
         self.assertEqual(self.m2_tuple.field, None)
diff --git a/tests/test_srf.py b/tests/test_srf.py
index 7a3d2440..53d4e2c0 100644
--- a/tests/test_srf.py
+++ b/tests/test_srf.py
@@ -32,8 +32,14 @@ def setUp(self):
 
     def test_shape_1d(self):
         self.cov_model.dim = 1
-        srf = SRF(self.cov_model, mean=self.mean, mode_no=self.mode_no)
-        field_str = srf([self.x_grid], seed=self.seed, mesh_type="structured")
+        srf = SRF(
+            self.cov_model,
+            (self.x_grid,),
+            mesh_type="unstructured",
+            mean=self.mean,
+            mode_no=self.mode_no,
+        )
+        field_str = srf(seed=self.seed)
         field_unstr = srf(
             [self.x_tuple], seed=self.seed, mesh_type="unstructured"
         )
@@ -85,6 +91,8 @@ def test_ens(self):
             srf(name=name)
             self.assertTrue(name in srf.point_data.keys())
         self.assertEqual(srf.default_field, "ens_1")
+        # 3 ens + 1 raw_field
+        self.assertEqual(len(srf.point_data), 4)
 
     def test_anisotropy_2d(self):
         self.cov_model.dim = 2
@@ -138,8 +146,13 @@ def test_rotation_unstruct_2d(self):
         field_str = np.reshape(field, (y_len, x_len))
 
         self.cov_model.angles = -np.pi / 2.0
-        srf = SRF(self.cov_model, mean=self.mean, mode_no=self.mode_no)
-        field_rot = srf((x_u, y_u), seed=self.seed, mesh_type="unstructured")
+        srf = SRF(
+            self.cov_model,
+            (x_u, y_u),
+            mean=self.mean,
+            mode_no=self.mode_no
+        )
+        field_rot = srf(seed=self.seed)
         field_rot_str = np.reshape(field_rot, (y_len, x_len))
 
         self.assertAlmostEqual(field_str[0, 0], field_rot_str[-1, 0])