From 7c1c596bc7f04da3bebc2cae18e514b6ac275cc3 Mon Sep 17 00:00:00 2001
From: FZill <florian.zill@ufz.de>
Date: Mon, 5 Feb 2024 16:18:14 +0100
Subject: [PATCH] using pyvistas plot in meshlib examples

---
 .../plot_meshlib_pyvista_input.py             | 17 ++-----
 .../howto_meshlib/plot_meshlib_vtu_input.py   | 49 +++++++++----------
 2 files changed, 26 insertions(+), 40 deletions(-)

diff --git a/docs/examples/howto_meshlib/plot_meshlib_pyvista_input.py b/docs/examples/howto_meshlib/plot_meshlib_pyvista_input.py
index a457de915..6ad01daf3 100644
--- a/docs/examples/howto_meshlib/plot_meshlib_pyvista_input.py
+++ b/docs/examples/howto_meshlib/plot_meshlib_pyvista_input.py
@@ -8,15 +8,10 @@ For this example we create meshes from pyvista surfaces.
 """
 
 # %%
-import numpy as np  # For visualization only
-
-import ogstools.meshplotlib as mpl  # For visualization only
 from ogstools.meshlib.boundary import Layer
 from ogstools.meshlib.boundary_set import LayerSet
 from ogstools.meshlib.boundary_subset import Gaussian2D, Surface
-from ogstools.meshlib.region import (
-    to_region_tetraeder,
-)
+from ogstools.meshlib.region import to_region_tetraeder
 
 # See other examples for different meshing algorithms
 
@@ -30,13 +25,9 @@ surface3 = Surface(Gaussian2D(**args, height_offset=-200), material_id=2)
 
 ls = LayerSet([Layer(surface1, surface2), Layer(surface2, surface3)])
 mesh = to_region_tetraeder(ls, 40).mesh
+# interprets values as categories
+mesh["regions"] = [str(m) for m in mesh["MaterialIDs"]]
 
 # %%
 # Visualize the prism mesh
-
-slices = np.reshape(mesh.slice_along_axis(n=4, axis="y"), (1, -1))
-mpl.setup.aspect_limits = [0.2, 5.0]
-fig = mpl.plot(slices, "MaterialIDs")
-for ax, slice in zip(fig.axes, np.ravel(slices)):
-    ax.set_title(f"z = {slice.center[2]:.1f} {mpl.setup.length.data_unit}")
-_ = fig
+mesh.plot(scalars="regions", show_edges=True)
diff --git a/docs/examples/howto_meshlib/plot_meshlib_vtu_input.py b/docs/examples/howto_meshlib/plot_meshlib_vtu_input.py
index a2956693c..f8b60287e 100644
--- a/docs/examples/howto_meshlib/plot_meshlib_vtu_input.py
+++ b/docs/examples/howto_meshlib/plot_meshlib_vtu_input.py
@@ -10,9 +10,6 @@ For this example we create meshes from surface layers.
 # %%
 from pathlib import Path  # To get example vtu files
 
-import numpy as np  # For visualization only
-
-import ogstools.meshplotlib as mpl  # For visualization only
 from ogstools.definitions import TESTS_DIR  # To get example vtu files
 from ogstools.meshlib.boundary import Layer
 from ogstools.meshlib.boundary_set import LayerSet
@@ -24,10 +21,6 @@ from ogstools.meshlib.region import (
     to_region_voxel,
 )
 
-mpl.setup.reset()
-mpl.setup.length.output_unit = "km"
-mpl.setup.aspect_limits = [0.2, 5.0]
-
 # %%
 # The loaded surfaces are defined within VTU files and adhere to properties such as non-intersecting boundaries with consistent x and y bounds. Alternatively, surfaces can also be created using PyVista with the same properties.
 
@@ -54,28 +47,30 @@ pm = to_region_prism(layer_set1, resolution=200).mesh
 vm = to_region_voxel(layer_set1, resolution=[200, 200, 50]).mesh
 tm = to_region_tetraeder(layer_set1, resolution=200).mesh
 
-
+# %% [markdown]
+# Simplified mesh
+# ---------------
 # %%
-# Visualize the prism mesh
-
-mesh = pm
-slices = np.reshape(mesh.slice_along_axis(n=4, axis="y"), (-1, 1))
-fig = mpl.plot(slices, "MaterialIDs")
-for ax, slice in zip(fig.axes, np.ravel(slices)):
-    ax.set_title(f"z = {slice.center[2]:.1f} {mpl.setup.length.data_unit}")
+sm["regions"] = [str(m) for m in sm["MaterialIDs"]]
+sm.scale([1, 1, 5]).plot(scalars="regions", show_edges=True)
 
+# %% [markdown]
+# Voxel mesh
+# ---------------
 # %%
-# Visualize meshes with different meshing algorithm
+vm["regions"] = [str(m) for m in vm["MaterialIDs"]]
+vm.scale([1, 1, 5]).plot(scalars="regions", show_edges=True)
 
-meshes = [sm, vm, pm, tm]
-names = [
-    "to_region_simplified",
-    "to_region_voxel",
-    "to_region_prism",
-    "to_region_tetraeder",
-]
+# %% [markdown]
+# Prism mesh
+# ---------------
+# %%
+pm["regions"] = [str(m) for m in pm["MaterialIDs"]]
+pm.scale([1, 1, 5]).plot(scalars="regions", show_edges=True)
 
-x_slices = np.reshape([mesh.slice("x") for mesh in meshes], (-1, 1))
-fig = mpl.plot(x_slices, "MaterialIDs")
-for ax, name in zip(fig.axes, names):
-    ax.set_title(name)
+# %% [markdown]
+# Tetraeder mesh
+# ---------------
+# %%
+tm["regions"] = [str(m) for m in tm["MaterialIDs"]]
+tm.scale([1, 1, 5]).plot(scalars="regions", show_edges=True)
-- 
GitLab