Using PyVista meshes

PyVista is a helper module for the Visualization Toolkit (VTK) that takes a different approach on interfacing with VTK through NumPy and direct array access.

It provides mesh data structures and filtering methods for spatial datasets, makes 3D plotting simple and is built for large/complex data geometries.

The Field.mesh method enables easy field creation on PyVista meshes used by the SRF or Krige class.

import pyvista as pv

import gstools as gs

We create a structured grid with PyVista containing 50 segments on all three axes each with a length of 2 (whatever unit).

dim, spacing = (50, 50, 50), (2, 2, 2)
grid = pv.UniformGrid(dim, spacing)
/home/docs/checkouts/readthedocs.org/user_builds/gstools/envs/stable/lib/python3.8/site-packages/pyvista/core/grid.py:502: PyVistaDeprecationWarning: Behavior of pyvista.UniformGrid has changed. First argument must be either a ``vtk.vtkImageData`` or path.
  warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/gstools/envs/stable/lib/python3.8/site-packages/pyvista/core/grid.py:520: PyVistaDeprecationWarning: Behavior of pyvista.UniformGrid has changed. Use keyword arguments to specify dimensions, spacing, and origin. For example:

    >>> grid = pyvista.UniformGrid(
    ...     dimensions=(10, 10, 10),
    ...     spacing=(2, 1, 5),
    ...     origin=(10, 35, 50),
    ... )

  warnings.warn(

Now we set up the SRF class as always. We’ll use an anisotropic model.

model = gs.Gaussian(dim=3, len_scale=[16, 8, 4], angles=(0.8, 0.4, 0.2))
srf = gs.SRF(model, seed=19970221)

The PyVista mesh can now be directly passed to the SRF.mesh method. When dealing with meshes, one can choose if the field should be generated on the mesh-points (“points”) or the cell-centroids (“centroids”).

In addition we can set a name, under which the resulting field is stored in the mesh.

srf.mesh(grid, points="points", name="random-field")

Now we have access to PyVista’s abundancy of methods to explore the field.

Note

PyVista is not working on readthedocs, but you can try it out yourself by uncommenting the following line of code.

# grid.contour(isosurfaces=8).plot()

The result should look like this:

https://github.com/GeoStat-Framework/GeoStat-Framework.github.io/raw/master/img/GS_pyvista_cut.png

Total running time of the script: ( 0 minutes 7.411 seconds)

Gallery generated by Sphinx-Gallery