diff --git a/gustaf/edges.py b/gustaf/edges.py
index 28eea5e7..0f0c31f5 100644
--- a/gustaf/edges.py
+++ b/gustaf/edges.py
@@ -3,42 +3,55 @@
 Edges. Also known as lines.
 """
 
+from __future__ import annotations
+
 from copy import deepcopy
+from typing import TYPE_CHECKING
+
+import numpy as _np
 
-import numpy as np
+from gustaf import helpers as _helpers
+from gustaf import settings as _settings
+from gustaf import show as _show
+from gustaf import utils as _utils
+from gustaf.helpers.options import Option as _Option
+from gustaf.vertices import Vertices as _Vertices
 
-from gustaf import helpers, settings, show, utils
-from gustaf.helpers.options import Option
-from gustaf.vertices import Vertices
+if TYPE_CHECKING:
+    from typing import Any
 
+    from gustaf.faces import Faces
+    from gustaf.helpers.data import TrackedArray, Unique2DIntegers
+    from gustaf.volumes import Volumes
 
-class EdgesShowOption(helpers.options.ShowOption):
+
+class EdgesShowOption(_helpers.options.ShowOption):
     """
     Show options for vertices.
     """
 
-    _valid_options = helpers.options.make_valid_options(
-        *helpers.options.vedo_common_options,
-        Option(
+    _valid_options = _helpers.options.make_valid_options(
+        *_helpers.options.vedo_common_options,
+        _Option(
             "vedo",
             "lw",
             "Width of edges (lines) in pixel units.",
             (float, int),
         ),
-        Option("vedo", "as_arrows", "Show edges as arrows.", (bool,)),
-        Option(
+        _Option("vedo", "as_arrows", "Show edges as arrows.", (bool,)),
+        _Option(
             "vedo",
             "head_radius",
             "Radius of arrow head. Applicable if as_arrows is True",
             (float, int),
         ),
-        Option(
+        _Option(
             "vedo",
             "head_length",
             "Length of arrow head. Applicable if as_arrows is True",
             (float, int),
         ),
-        Option(
+        _Option(
             "vedo",
             "shaft_radius",
             "Radius of arrow shaft. Applicable if as_arrows is True",
@@ -62,20 +75,20 @@ def _initialize_showable(self):
         """
         if self.get("as_arrows", False):
             init_options = ("head_radius", "head_length", "shaft_radius")
-            return show.vedo.Arrows(
+            return _show.vedo.Arrows(
                 self._helpee.const_vertices[self._helpee.edges],
                 **self[init_options],
             )
 
         else:
             init_options = ("lw",)
-            return show.vedo.Lines(
+            return _show.vedo.Lines(
                 self._helpee.const_vertices[self._helpee.edges],
                 **self[init_options],
             )
 
 
-class Edges(Vertices):
+class Edges(_Vertices):
     kind = "edge"
 
     __slots__ = (
@@ -84,14 +97,14 @@ class Edges(Vertices):
     )
 
     __show_option__ = EdgesShowOption
-    __boundary_class__ = Vertices
+    __boundary_class__ = _Vertices
 
     def __init__(
         self,
-        vertices=None,
-        edges=None,
-        elements=None,
-    ):
+        vertices: list[list[float]] | TrackedArray | _np.ndarray = None,
+        edges: _np.ndarray | None = None,
+        elements: _np.ndarray | None = None,
+    ) -> None:
         """Edges. It has vertices and edges. Also known as lines.
 
         Parameters
@@ -108,7 +121,7 @@ def __init__(
             self.edges = elements
 
     @property
-    def edges(self):
+    def edges(self) -> TrackedArray:
         """Returns edges. If edges is not its original property.
 
         Parameters
@@ -123,7 +136,7 @@ def edges(self):
         return self._edges
 
     @edges.setter
-    def edges(self, es):
+    def edges(self, es: TrackedArray | _np.ndarray) -> None:
         """Edges setter. Similar to vertices, this is a tracked array.
 
         Parameters
@@ -136,20 +149,20 @@ def edges(self, es):
         """
         self._logd("setting edges")
 
-        self._edges = helpers.data.make_tracked_array(
-            es, settings.INT_DTYPE, copy=False
+        self._edges = _helpers.data.make_tracked_array(
+            es, _settings.INT_DTYPE, copy=False
         )
 
         # shape check
         if es is not None:
-            utils.arr.is_shape(es, (-1, 2), strict=True)
+            _utils.arr.is_shape(es, (-1, 2), strict=True)
 
         # same, but non-writeable view of tracked array
         self._const_edges = self._edges.view()
         self._const_edges.flags.writeable = False
 
     @property
-    def const_edges(self):
+    def const_edges(self) -> TrackedArray:
         """Returns non-writeable version of edges.
 
         Parameters
@@ -163,7 +176,7 @@ def const_edges(self):
         return self._const_edges
 
     @property
-    def whatami(self):
+    def whatami(self) -> str:
         """whatami?
 
         Parameters
@@ -176,8 +189,8 @@ def whatami(self):
         """
         return "edges"
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def sorted_edges(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def sorted_edges(self) -> _np.ndarray:
         """Sort edges along axis=1.
 
         Parameters
@@ -190,10 +203,10 @@ def sorted_edges(self):
         """
         edges = self._get_attr("edges")
 
-        return np.sort(edges, axis=1)
+        return _np.sort(edges, axis=1)
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def unique_edges(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def unique_edges(self) -> Unique2DIntegers:
         """Returns a named tuple of unique edge info. Info includes unique
         values, ids of unique edges, inverse ids, count of each unique values.
 
@@ -206,7 +219,7 @@ def unique_edges(self):
         unique_info: Unique2DIntegers
           valid attributes are {values, ids, inverse, counts}
         """
-        unique_info = utils.connec.sorted_unique(
+        unique_info = _utils.connec.sorted_unique(
             self.sorted_edges(), sorted_=True
         )
 
@@ -217,8 +230,8 @@ def unique_edges(self):
 
         return unique_info
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def single_edges(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def single_edges(self) -> _np.ndarray:
         """Returns indices of very unique edges: edges that appear only once.
         For well constructed faces, this can be considered as outlines.
 
@@ -235,7 +248,7 @@ def single_edges(self):
         return unique_info.ids[unique_info.counts == 1]
 
     @property
-    def elements(self):
+    def elements(self) -> TrackedArray:
         """Returns current connectivity. A short cut in FE friendly term.
         Elements mean different things for different classes: Vertices ->
         vertices Edges -> edges Faces -> faces Volumes -> volumes.
@@ -255,7 +268,7 @@ def elements(self):
         return getattr(self, elem_name)
 
     @elements.setter
-    def elements(self, elements):
+    def elements(self, elements: TrackedArray | _np.ndarray) -> Any | None:
         """Calls corresponding connectivity setter. A short cut in FEM friendly
         term. Vertices -> vertices Edges -> edges Faces -> faces Volumes ->
         volumes.
@@ -274,7 +287,7 @@ def elements(self, elements):
         return setattr(self, elem_name, elements)
 
     @property
-    def const_elements(self):
+    def const_elements(self) -> TrackedArray:
         """Returns non-mutable version of elements.
 
         Parameters
@@ -288,8 +301,8 @@ def const_elements(self):
         self._logd("returning const_elements")
         return getattr(self, "const_" + type(self).__qualname__.lower())
 
-    @helpers.data.ComputedMeshData.depends_on(["vertices", "elements"])
-    def centers(self):
+    @_helpers.data.ComputedMeshData.depends_on(["vertices", "elements"])
+    def centers(self) -> TrackedArray:
         """Center of elements.
 
         Parameters
@@ -304,10 +317,10 @@ def centers(self):
 
         return self.const_vertices[self.const_elements].mean(axis=1)
 
-    @helpers.data.ComputedMeshData.depends_on(["vertices", "elements"])
+    @_helpers.data.ComputedMeshData.depends_on(["vertices", "elements"])
     def referenced_vertices(
         self,
-    ):
+    ) -> _np.ndarray:
         """Returns mask of referenced vertices.
 
         Parameters
@@ -318,12 +331,12 @@ def referenced_vertices(
         --------
         referenced: (n,) np.ndarray
         """
-        referenced = np.zeros(len(self.const_vertices), dtype=bool)
+        referenced = _np.zeros(len(self.const_vertices), dtype=bool)
         referenced[self.const_elements] = True
 
         return referenced
 
-    def remove_unreferenced_vertices(self):
+    def remove_unreferenced_vertices(self) -> Edges | Faces | Volumes:
         """Remove unreferenced vertices. Adapted from
         `github.com/mikedh/trimesh`
 
@@ -337,15 +350,15 @@ def remove_unreferenced_vertices(self):
         """
         referenced = self.referenced_vertices()
 
-        inverse = np.zeros(len(self.vertices), dtype=settings.INT_DTYPE)
-        inverse[referenced] = np.arange(referenced.sum())
+        inverse = _np.zeros(len(self.vertices), dtype=_settings.INT_DTYPE)
+        inverse[referenced] = _np.arange(referenced.sum())
 
         return self.update_vertices(
             mask=referenced,
             inverse=inverse,
         )
 
-    def update_elements(self, mask):
+    def update_elements(self, mask: _np.ndarray) -> Edges | Faces | Volumes:
         """Similar to update_vertices, but for elements.
 
         Parameters
@@ -364,7 +377,7 @@ def update_edges(self, *args, **kwargs):
         """Alias to update_elements."""
         return self.update_elements(*args, **kwargs)
 
-    def dashed(self, spacing=None):
+    def dashed(self, spacing: Any | None = None) -> Edges:
         """Turn edges into dashed edges(=lines). Given spacing, it will try to
         chop edges as close to it as possible. Pattern should look:
 
@@ -395,30 +408,34 @@ def dashed(self, spacing=None):
         v0s = self.vertices[self.edges[:, 0]]
         v1s = self.vertices[self.edges[:, 1]]
 
-        distances = np.linalg.norm(v0s - v1s, axis=1)
-        linspaces = (((distances // (spacing * 1.5)) + 1) * 3).astype(np.int32)
+        distances = _np.linalg.norm(v0s - v1s, axis=1)
+        linspaces = (((distances // (spacing * 1.5)) + 1) * 3).astype(
+            _np.int32
+        )
 
         # chop vertices!
         new_vs = []
         for v0, v1, lins in zip(v0s, v1s, linspaces):
-            new_vs.append(np.linspace(v0, v1, lins))
+            new_vs.append(_np.linspace(v0, v1, lins))
 
         # we need all chopped vertices.
         # there might be duplicating vertices. you can use merge_vertices
-        new_vs = np.vstack(new_vs)
+        new_vs = _np.vstack(new_vs)
         # all mid points are explicitly defined, but they aren't required
         # so, rm.
-        mask = np.ones(len(new_vs), dtype=bool)
+        mask = _np.ones(len(new_vs), dtype=bool)
         mask[1::3] = False
         new_vs = new_vs[mask]
 
         # prepare edges
-        tmp_es = utils.connec.range_to_edges((0, len(new_vs)), closed=False)
+        tmp_es = _utils.connec.range_to_edges((0, len(new_vs)), closed=False)
         new_es = tmp_es[::2]
 
         return Edges(vertices=new_vs, edges=new_es)
 
-    def shrink(self, ratio=0.8, map_vertex_data=True):
+    def shrink(
+        self, ratio: float = 0.8, map_vertex_data: bool = True
+    ) -> Edges | Faces | Volumes:
         """Returns shrunk elements.
 
         Parameters
@@ -434,13 +451,13 @@ def shrink(self, ratio=0.8, map_vertex_data=True):
           shrunk elements
         """
         elements = self.const_elements
-        vs = np.vstack(self.vertices[elements])
-        es = np.arange(len(vs))
+        vs = _np.vstack(self.vertices[elements])
+        es = _np.arange(len(vs))
 
         nodes_per_element = elements.shape[1]
         es = es.reshape(-1, nodes_per_element)
 
-        mids = np.repeat(self.centers(), nodes_per_element, axis=0)
+        mids = _np.repeat(self.centers(), nodes_per_element, axis=0)
 
         vs -= mids
         vs *= ratio
@@ -460,7 +477,7 @@ def shrink(self, ratio=0.8, map_vertex_data=True):
 
         return s_elements
 
-    def to_vertices(self):
+    def to_vertices(self) -> _Vertices:
         """Returns Vertices obj.
 
         Parameters
@@ -471,9 +488,9 @@ def to_vertices(self):
         --------
         vertices: Vertices
         """
-        return Vertices(self.vertices)
+        return _Vertices(self.vertices)
 
-    def _get_attr(self, attr):
+    def _get_attr(self, attr: str) -> TrackedArray | _np.ndarray:
         """Internal function to get attribute that maybe property or callable.
         Some properties are replaced by callable in subclasses as it may depend
         on other properties of subclass.
diff --git a/gustaf/faces.py b/gustaf/faces.py
index f8013841..86c9e6a0 100644
--- a/gustaf/faces.py
+++ b/gustaf/faces.py
@@ -1,10 +1,21 @@
 """gustaf/gustaf/faces.py."""
 
-import numpy as np
+from __future__ import annotations
 
-from gustaf import helpers, settings, show, utils
-from gustaf.edges import Edges
-from gustaf.helpers.options import Option
+from typing import TYPE_CHECKING
+
+import numpy as _np
+
+from gustaf import helpers as _helpers
+from gustaf import settings as _settings
+from gustaf import show as _show
+from gustaf import utils as _utils
+from gustaf.edges import Edges as _Edges
+from gustaf.helpers.options import Option as _Option
+
+if TYPE_CHECKING:
+
+    from gustaf.helpers.data import TrackedArray, Unique2DIntegers
 
 # special types for face texture option
 try:
@@ -14,27 +25,29 @@
     # there are other ways to get here, but this is exact path for our use
     vtkTexture = vedo.vtkclasses.vtkTexture
 except ImportError as err:
-    vedoPicture = helpers.raise_if.ModuleImportRaiser("vedo", err)
-    vtkTexture = helpers.raise_if.ModuleImportRaiser("vedo", err)
+    vedoPicture = _helpers.raise_if.ModuleImportRaiser("vedo", err)
+    vtkTexture = _helpers.raise_if.ModuleImportRaiser("vedo", err)
 
 
-class FacesShowOption(helpers.options.ShowOption):
+class FacesShowOption(_helpers.options.ShowOption):
     """
     Show options for vertices.
     """
 
-    _valid_options = helpers.options.make_valid_options(
-        *helpers.options.vedo_common_options,
-        Option("vedo", "lw", "Width of edges (lines) in pixel units.", (int,)),
-        Option(
+    _valid_options = _helpers.options.make_valid_options(
+        *_helpers.options.vedo_common_options,
+        _Option(
+            "vedo", "lw", "Width of edges (lines) in pixel units.", (int,)
+        ),
+        _Option(
             "vedo", "lc", "Color of edges (lines).", (int, str, tuple, list)
         ),
-        Option(
+        _Option(
             "vedo",
             "texture",
             "Texture of faces in array, vedo.Picture, vtk.vtkTexture, "
             "or path to an image.",
-            (np.ndarray, tuple, list, str, vedoPicture, vtkTexture),
+            (_np.ndarray, tuple, list, str, vedoPicture, vtkTexture),
         ),
     )
 
@@ -53,7 +66,7 @@ def _initialize_showable(self):
         faces: vedo.Mesh
         """
 
-        faces = show.vedo.Mesh(
+        faces = _show.vedo.Mesh(
             [self._helpee.const_vertices, self._helpee.const_faces],
         )
 
@@ -65,20 +78,20 @@ def _initialize_showable(self):
         return faces
 
 
-class Faces(Edges):
+class Faces(_Edges):
     kind = "face"
 
-    const_edges = helpers.raise_if.invalid_inherited_attr(
+    const_edges = _helpers.raise_if.invalid_inherited_attr(
         "Edges.const_edges",
         __qualname__,
         property_=True,
     )
-    update_edges = helpers.raise_if.invalid_inherited_attr(
+    update_edges = _helpers.raise_if.invalid_inherited_attr(
         "Edges.update_edges",
         __qualname__,
         property_=False,
     )
-    dashed = helpers.raise_if.invalid_inherited_attr(
+    dashed = _helpers.raise_if.invalid_inherited_attr(
         "Edges.dashed",
         __qualname__,
         property_=False,
@@ -91,14 +104,14 @@ class Faces(Edges):
     )
 
     __show_option__ = FacesShowOption
-    __boundary_class__ = Edges
+    __boundary_class__ = _Edges
 
     def __init__(
         self,
-        vertices=None,
-        faces=None,
-        elements=None,
-    ):
+        vertices: list[list[float]] | TrackedArray | _np.ndarray = None,
+        faces: _np.ndarray | None = None,
+        elements: _np.ndarray | None = None,
+    ) -> None:
         """Faces. It has vertices and faces. Faces could be triangles or
         quadrilaterals.
 
@@ -116,8 +129,8 @@ def __init__(
 
         self.BC = {}
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def edges(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def edges(self) -> _np.ndarray:
         """Edges from here aren't main property. So this needs to be computed.
 
         Parameters
@@ -131,10 +144,10 @@ def edges(self):
         self._logd("computing edges")
         faces = self._get_attr("faces")
 
-        return utils.connec.faces_to_edges(faces)
+        return _utils.connec.faces_to_edges(faces)
 
     @property
-    def whatami(self):
+    def whatami(self) -> str:
         """Determines whatami.
 
         Parameters
@@ -148,7 +161,7 @@ def whatami(self):
         return type(self).whatareyou(self)
 
     @classmethod
-    def whatareyou(cls, face_obj):
+    def whatareyou(cls, face_obj: Faces) -> str:
         """classmethod that tells you if the Faces is tri or quad or invalid
         kind.
 
@@ -176,7 +189,7 @@ def whatareyou(cls, face_obj):
             )
 
     @property
-    def faces(self):
+    def faces(self) -> TrackedArray:
         """Returns faces.
 
         Parameters
@@ -191,7 +204,7 @@ def faces(self):
         return self._faces
 
     @faces.setter
-    def faces(self, fs):
+    def faces(self, fs: TrackedArray | _np.ndarray) -> None:
         """Faces setter. Similar to vertices, this will be a tracked array.
 
         Parameters
@@ -204,14 +217,14 @@ def faces(self, fs):
         """
         self._logd("setting faces")
 
-        self._faces = helpers.data.make_tracked_array(
+        self._faces = _helpers.data.make_tracked_array(
             fs,
-            settings.INT_DTYPE,
+            _settings.INT_DTYPE,
             copy=False,
         )
         # shape check
         if fs is not None:
-            utils.arr.is_one_of_shapes(
+            _utils.arr.is_one_of_shapes(
                 fs,
                 ((-1, 3), (-1, 4)),
                 strict=True,
@@ -222,7 +235,7 @@ def faces(self, fs):
         self._const_faces.flags.writeable = False
 
     @property
-    def const_faces(self):
+    def const_faces(self) -> TrackedArray:
         """Returns non-writeable view of faces.
 
         Parameters
@@ -235,8 +248,8 @@ def const_faces(self):
         """
         return self._const_faces
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def sorted_faces(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def sorted_faces(self) -> _np.ndarray:
         """Similar to edges_sorted but for faces.
 
         Parameters
@@ -249,10 +262,10 @@ def sorted_faces(self):
         """
         faces = self._get_attr("faces")
 
-        return np.sort(faces, axis=1)
+        return _np.sort(faces, axis=1)
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def unique_faces(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def unique_faces(self) -> Unique2DIntegers:
         """Returns a namedtuple of unique faces info. Similar to unique_edges.
 
         Parameters
@@ -264,7 +277,7 @@ def unique_faces(self):
         unique_info: Unique2DIntegers
           valid attributes are {values, ids, inverse, counts}
         """
-        unique_info = utils.connec.sorted_unique(
+        unique_info = _utils.connec.sorted_unique(
             self.sorted_faces(), sorted_=True
         )
 
@@ -274,8 +287,8 @@ def unique_faces(self):
 
         return unique_info
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def single_faces(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def single_faces(self) -> _np.ndarray:
         """Returns indices of very unique faces: faces that appear only once.
         For well constructed volumes, this can be considered as surfaces.
 
@@ -307,7 +320,7 @@ def to_edges(self, unique=True):
         --------
         edges: Edges
         """
-        return Edges(
+        return _Edges(
             self.vertices,
             edges=self.unique_edges().values if unique else self.edges(),
         )
diff --git a/gustaf/helpers/_base.py b/gustaf/helpers/_base.py
index cac37fcb..41641ee4 100644
--- a/gustaf/helpers/_base.py
+++ b/gustaf/helpers/_base.py
@@ -3,12 +3,13 @@
 Base class for helper
 """
 
+from typing import Any as _Any
 from weakref import ref
 
-from gustaf._base import GustafBase
+from gustaf._base import GustafBase as _GustafBase
 
 
-class HelperBase(GustafBase):
+class HelperBase(_GustafBase):
     """
     Minimal base layer for helper classes to avoid cyclic referencing.
     Instead of saving a pure reference to helpee object, this will create a
@@ -19,11 +20,11 @@ class HelperBase(GustafBase):
     __slots__ = ("_helpee_weak_ref",)
 
     @property
-    def _helpee(self):
+    def _helpee(self) -> _Any:
         """Returns dereferenced weak ref. Setter will create and save weakref
         of a helpee."""
         return self._helpee_weak_ref()
 
     @_helpee.setter
-    def _helpee(self, helpee):
+    def _helpee(self, helpee: _Any) -> None:
         self._helpee_weak_ref = ref(helpee)
diff --git a/gustaf/helpers/data.py b/gustaf/helpers/data.py
index 1d01f55b..697d3605 100644
--- a/gustaf/helpers/data.py
+++ b/gustaf/helpers/data.py
@@ -3,15 +3,18 @@
 Helps helpee to manage data. Some useful data structures.
 """
 
-from collections import namedtuple
-from functools import wraps
+from __future__ import annotations
 
-import numpy as np
+from collections import namedtuple as _namedtuple
+from functools import wraps as _wraps
+from typing import Any as _Any
 
-from gustaf.helpers._base import HelperBase
+import numpy as _np
 
+from gustaf.helpers._base import HelperBase as _HelperBase
 
-class TrackedArray(np.ndarray):
+
+class TrackedArray(_np.ndarray):
     """numpy array object that keeps mirroring inplace changes to the source.
     Meant to help control_points.
     """
@@ -21,7 +24,9 @@ class TrackedArray(np.ndarray):
         "_modified",
     )
 
-    def __array_finalize__(self, obj):
+    def __array_finalize__(
+        self, obj: TrackedArray | _np.ndarray
+    ) -> _Any | None:
         """Sets default flags for any arrays that maybe generated based on
         physical space array. For more information,
         see https://numpy.org/doc/stable/user/basics.subclassing.html"""
@@ -63,23 +68,23 @@ def modified(self):
         return self._modified
 
     @modified.setter
-    def modified(self, m):
+    def modified(self, m: bool) -> None:
         if self._super_arr is not None and self._super_arr is not True:
             self._super_arr._modified = m
         else:
             self._modified = m
 
-    def copy(self, *args, **kwargs):
+    def copy(self, *args: _Any, **kwargs: _Any) -> _np.ndarray:
         """copy creates regular numpy array"""
-        return np.array(self, *args, copy=True, **kwargs)
+        return _np.array(self, *args, copy=True, **kwargs)
 
-    def view(self, *args, **kwargs):
+    def view(self, *args: type, **kwargs: _Any) -> TrackedArray:
         """Set writeable flags to False for the view."""
         v = super(self.__class__, self).view(*args, **kwargs)
         v.flags.writeable = False
         return v
 
-    def __iadd__(self, *args, **kwargs):
+    def __iadd__(self, *args: TrackedArray, **kwargs: _Any) -> TrackedArray:
         sr = super(self.__class__, self).__iadd__(*args, **kwargs)
         self.modified = True
         return sr
@@ -149,14 +154,16 @@ def __ior__(self, *args, **kwargs):
         self.modified = True
         return sr
 
-    def __setitem__(self, key, value):
+    def __setitem__(self, key: _Any, value: TrackedArray) -> _Any | None:
         # set first. invalid setting will cause error
         sr = super(self.__class__, self).__setitem__(key, value)
         self.modified = True
         return sr
 
 
-def make_tracked_array(array, dtype=None, copy=True):
+def make_tracked_array(
+    array: _Any, dtype: str = None, copy: bool = True
+) -> TrackedArray:
     """Motivated by nice implementations of `trimesh` (see LICENSE.txt).
     `https://github.com/mikedh/trimesh/blob/main/trimesh/caching.py`.
 
@@ -182,9 +189,9 @@ def make_tracked_array(array, dtype=None, copy=True):
         array = []
 
     if copy:
-        array = np.array(array, dtype=dtype)
+        array = _np.array(array, dtype=dtype)
     else:
-        array = np.asanyarray(array, dtype=dtype)
+        array = _np.asanyarray(array, dtype=dtype)
 
     tracked = array.view(TrackedArray)
 
@@ -194,10 +201,10 @@ def make_tracked_array(array, dtype=None, copy=True):
     return tracked
 
 
-class DataHolder(HelperBase):
+class DataHolder(_HelperBase):
     __slots__ = ("_saved",)
 
-    def __init__(self, helpee):
+    def __init__(self, helpee: _Any) -> None:
         """Base class for any data holder. Behaves similar to dict.
 
         Parameters
@@ -309,7 +316,7 @@ def values(self):
         """
         return self._saved.values()
 
-    def items(self):
+    def items(self) -> _Any:
         """Returns items of data holding dict.
 
         Returns
@@ -325,7 +332,7 @@ class ComputedData(DataHolder):
 
     __slots__ = ()
 
-    def __init__(self, helpee, **_kwargs):
+    def __init__(self, helpee: _Any, **_kwargs: _Any) -> None:
         """Stores last computed values.
 
         Keys are expected to be the same as helpee's function that computes the
@@ -380,8 +387,10 @@ def inner(func):
 
                 cls._inv_depends[vn].append(func.__name__)
 
-            @wraps(func)
-            def compute_or_return_saved(*args, **kwargs):
+            @_wraps(func)
+            def compute_or_return_saved(
+                *args: _Any, **kwargs: _Any
+            ) -> Unique2DIntegers | _np.ndarray:
                 """Check if the key should be computed,"""
                 # extract some related info
                 self = args[0]  # the helpee itself
@@ -412,7 +421,7 @@ def compute_or_return_saved(*args, **kwargs):
 
                 # we've reached this point because we have to compute this
                 computed = func(*args, **kwargs)
-                if isinstance(computed, np.ndarray):
+                if isinstance(computed, _np.ndarray):
                     computed.flags.writeable = False  # configurable?
                 self._computed._saved[func.__name__] = computed
 
@@ -444,7 +453,7 @@ class VertexData(DataHolder):
 
     __slots__ = ()
 
-    def __init__(self, helpee):
+    def __init__(self, helpee: _Any) -> None:
         """Checks if helpee has vertices as attr beforehand.
 
         Parameters
@@ -457,7 +466,9 @@ def __init__(self, helpee):
 
         super().__init__(helpee)
 
-    def _validate_len(self, value=None, raise_=True):
+    def _validate_len(
+        self, value: _Any | None = None, raise_: bool = True
+    ) -> bool:
         """Checks if given value is a valid vertex_data based of its length.
 
         If raise_, throws error, else, deletes all incompatible values.
@@ -605,7 +616,7 @@ def as_scalar(self, key, default=None):
         if value.shape[1] == 1:
             value_norm = value
         else:
-            value_norm = np.linalg.norm(value, axis=1).reshape(-1, 1)
+            value_norm = _np.linalg.norm(value, axis=1).reshape(-1, 1)
 
         # save norm
         self[norm_key] = value_norm
@@ -642,7 +653,7 @@ def as_arrow(self, key, default=None, raise_=True):
         return value
 
 
-Unique2DFloats = namedtuple(
+Unique2DFloats = _namedtuple(
     "Unique2DFloats", ["values", "ids", "inverse", "intersection"]
 )
 Unique2DFloats.__doc__ = """
@@ -662,7 +673,7 @@ def as_arrow(self, key, default=None, raise_=True):
   Field number 3
 """
 
-Unique2DIntegers = namedtuple(
+Unique2DIntegers = _namedtuple(
     "Unique2DIntegers", ["values", "ids", "inverse", "counts"]
 )
 Unique2DIntegers.__doc__ = """
diff --git a/gustaf/helpers/options.py b/gustaf/helpers/options.py
index 4a497b10..b26de5dd 100644
--- a/gustaf/helpers/options.py
+++ b/gustaf/helpers/options.py
@@ -3,10 +3,13 @@
 Classes to help organize options.
 """
 
-from copy import deepcopy
+from __future__ import annotations
 
-from gustaf.helpers._base import HelperBase
-from gustaf.helpers.raise_if import ModuleImportRaiser
+from copy import deepcopy as _deepcopy
+from typing import Any as _Any
+
+from gustaf.helpers._base import HelperBase as _HelperBase
+from gustaf.helpers.raise_if import ModuleImportRaiser as _ModuleImportRaiser
 
 
 class Option:
@@ -21,8 +24,8 @@ class Option:
     description: str
     allowed_types: set
       set of types
-    default: one of allwed_types
-      Optional. Default is None
+    default: one of allowed_types
+      _Optional. Default is None
     """
 
     __slots__ = (
@@ -34,8 +37,13 @@ class Option:
     )
 
     def __init__(
-        self, backends, key, description, allowed_types, default=None
-    ):
+        self,
+        backends: int | str,
+        key: int | str,
+        description: int | str,
+        allowed_types: tuple[type, ...] | int,
+        default: int | str | None = None,
+    ) -> None:
         """
         Check types
         """
@@ -91,7 +99,7 @@ class SetDefault:
 
     __slots__ = ("key", "default")
 
-    def __init__(self, key, default):
+    def __init__(self, key: str, default: int | str) -> None:
         self.key = key
         self.default = default
 
@@ -192,7 +200,7 @@ def __init__(self, key, default):
 )
 
 
-def make_valid_options(*options):
+def make_valid_options(*options: _Any) -> dict[str, Option]:
     """
     Forms valid options. Should run only once during module loading.
 
@@ -212,12 +220,12 @@ def make_valid_options(*options):
             # and wrapped by ModuleImportRaiser
             allowed_types = []
             for at in opt.allowed_types:
-                if isinstance(at, ModuleImportRaiser):
+                if isinstance(at, _ModuleImportRaiser):
                     continue
                 allowed_types.append(at)
             opt.allowed_types = tuple(allowed_types)
 
-            valid_options[opt.key] = deepcopy(opt)
+            valid_options[opt.key] = _deepcopy(opt)
         elif isinstance(opt, SetDefault):
             # overwrite default of existing option.
             if opt.key not in valid_options:
@@ -232,7 +240,7 @@ def make_valid_options(*options):
     return valid_options
 
 
-class ShowOption(HelperBase):
+class ShowOption(_HelperBase):
     """
     Behaves similar to dict, but will only accept a set of options that's
     applicable to the helpee class. Intended use is to create a
@@ -248,7 +256,7 @@ class ShowOption(HelperBase):
 
     _helps = None
 
-    def __init__(self, helpee):
+    def __init__(self, helpee: _Any) -> None:
         """
         Parameters
         ----------
@@ -494,7 +502,7 @@ def copy_valid_options(self, copy_to, keys=None):
 
         for key, value in items:
             if key in valid_keys:
-                copy_to[key] = deepcopy(value)  # is deepcopy necessary?
+                copy_to[key] = _deepcopy(value)  # is deepcopy necessary?
 
     def _initialize_showable(self):
         """
diff --git a/gustaf/helpers/raise_if.py b/gustaf/helpers/raise_if.py
index bce64cf2..827c49f8 100644
--- a/gustaf/helpers/raise_if.py
+++ b/gustaf/helpers/raise_if.py
@@ -4,7 +4,9 @@
 behavior
 """
 
-from typing import Any, Optional
+from __future__ import annotations
+
+from typing import Any
 
 
 def invalid_inherited_attr(attr_name, qualname, property_=False):
@@ -48,7 +50,7 @@ class ModuleImportRaiser:
     them to function. Examples are `splinepy` and `vedo`.
     """
 
-    def __init__(self, lib_name: str, error_message: Optional[str] = None):
+    def __init__(self, lib_name: str, error_message: str | None = None):
         original_message = ""
         if error_message is not None:
             original_message = f"\nOriginal error message - {error_message}"
diff --git a/gustaf/settings.py b/gustaf/settings.py
index 556f5dbb..8a579f50 100644
--- a/gustaf/settings.py
+++ b/gustaf/settings.py
@@ -3,19 +3,21 @@
 Global variables/constants that's used throughout `gustaf`.
 """
 
-TOLERANCE = 1e-10
+from typing import Literal
 
-FLOAT_DTYPE = "float64"
-INT_DTYPE = "int32"
+TOLERANCE: float = 1e-10
+
+FLOAT_DTYPE: Literal["float64"] = "float64"
+INT_DTYPE: Literal["int32"] = "int32"
 
 # OPTIONS are <"vedo" | "trimesh" | "matplotlib">
-VISUALIZATION_BACKEND = "vedo"
+VISUALIZATION_BACKEND: Literal["vedo"] = "vedo"
 
-VEDO_DEFAULT_OPTIONS = {
+VEDO_DEFAULT_OPTIONS: dict[str, dict] = {
     "vertex": {},
     "edges": {},
     "faces": {},
     "volumes": {},
 }
 
-NTHREADS = 1
+NTHREADS: int = 1
diff --git a/gustaf/show.py b/gustaf/show.py
index 89a54fe9..95df312d 100644
--- a/gustaf/show.py
+++ b/gustaf/show.py
@@ -3,7 +3,10 @@
 Everything related to show/visualization.
 """
 
+from __future__ import annotations
+
 import sys
+from typing import TYPE_CHECKING
 
 import numpy as np
 
@@ -38,6 +41,11 @@
     vedo = ModuleImportRaiser("vedo", err)
     vedoUGrid = vedo
 
+if TYPE_CHECKING:
+    from typing import Any
+
+    from gustaf._base import GustafBase
+
 
 # enable `gus.show()`
 # taken from https://stackoverflow.com/questions/1060796/callable-modules
@@ -82,7 +90,9 @@ def show(*args, **kwargs):
     return_show_list = kwargs.get("return_showable_list", False)
     axes = kwargs.get("axes", None)
 
-    def clear_vedo_plotter(plotter, num_renderers, skip_cl=skip_clear):
+    def clear_vedo_plotter(
+        plotter: vedo.Plotter, num_renderers: int, skip_cl: bool = skip_clear
+    ):
         """enough said."""
         # for whatever reason it is desired
         if skip_cl:
@@ -98,7 +108,7 @@ def clear_vedo_plotter(plotter, num_renderers, skip_cl=skip_clear):
 
         return None
 
-    def cam_tuple_to_list(dict_cam):
+    def cam_tuple_to_list(dict_cam: None | dict[str, Any]):
         """if entity is tuple, turns it into list."""
         if dict_cam is None:
             return None
@@ -235,7 +245,9 @@ def cam_tuple_to_list(dict_cam):
         return plt
 
 
-def make_showable(obj, as_dict=False, **kwargs):
+# TODO the type of obj is not really easy to define since it in reality it
+# should be gustaf.Vertice but it might also be splinepy.Spline
+def make_showable(obj: GustafBase, as_dict: bool = False, **kwargs):
     """Generates a vedo obj based on `kind` attribute from given obj, as well
     as show_options.
 
@@ -411,7 +423,9 @@ def make_showable(obj, as_dict=False, **kwargs):
 
 # possibly relocate, is this actually used?
 # could not find any usage in this repo
-def interpolate_vedo_dictcam(cameras, resolutions, spline_degree=1):
+def interpolate_vedo_dictcam(
+    cameras: list[dict[str, Any]], resolutions: int, spline_degree: int = 1
+):
     """Interpolate between vedo dict cameras.
 
     Parameters
diff --git a/gustaf/utils/arr.py b/gustaf/utils/arr.py
index 3e8b6454..8f3a64d3 100644
--- a/gustaf/utils/arr.py
+++ b/gustaf/utils/arr.py
@@ -4,10 +4,16 @@
 `array` is python library and it sounds funny.
 """
 
-import numpy as np
+from __future__ import annotations
 
-from gustaf import settings
-from gustaf.helpers.raise_if import ModuleImportRaiser
+from typing import Any as _Any
+
+import numpy as _np
+from napf._napf import UIntVectorVector as _UIntVectorVector
+
+from gustaf import settings as _settings
+from gustaf.helpers.data import TrackedArray as _TrackedArray
+from gustaf.helpers.raise_if import ModuleImportRaiser as _ModuleImportRaiser
 
 has_funi = has_napf = has_scipy = False
 try:
@@ -15,22 +21,22 @@
 
     has_funi = True
 except ImportError:
-    funi = ModuleImportRaiser("funi")
+    funi = _ModuleImportRaiser("funi")
 try:
     import napf
 
     has_napf = True
 except ImportError:
-    napf = ModuleImportRaiser("napf")
+    napf = _ModuleImportRaiser("napf")
 try:
     import scipy
 
     has_scipy = True
 except ImportError:
-    scipy = ModuleImportRaiser("scipy")
+    scipy = _ModuleImportRaiser("scipy")
 
 
-def make_c_contiguous(array, dtype=None):
+def make_c_contiguous(array: _np.ndarray, dtype: str = None) -> _np.ndarray:
     """Make given array like object a c contiguous np.ndarray. dtype is
     optional. If None is given, just returns None.
 
@@ -38,7 +44,7 @@ def make_c_contiguous(array, dtype=None):
     -----------
     array: array-like
     dtype: type or str
-      (Optional) `numpy` interpretable type or str, describing type.
+      (_Optional) `numpy` interpretable type or str, describing type.
 
     Returns
     --------
@@ -47,17 +53,17 @@ def make_c_contiguous(array, dtype=None):
     if array is None:
         return None
 
-    if isinstance(array, np.ndarray) and array.flags.c_contiguous:
+    if isinstance(array, _np.ndarray) and array.flags.c_contiguous:
         if dtype is not None and array.dtype != dtype:
             return array.astype(dtype)
 
         return array
 
     if dtype:
-        return np.ascontiguousarray(array, dtype=dtype)
+        return _np.ascontiguousarray(array, dtype=dtype)
 
     else:
-        return np.ascontiguousarray(array)
+        return _np.ascontiguousarray(array)
 
 
 def enforce_len(value, n_len):
@@ -75,14 +81,14 @@ def enforce_len(value, n_len):
     len_n_array: (n_len,) np.ndarray
     """
     if isinstance(value, (int, float)):
-        return np.repeat(value, n_len)
-    elif isinstance(value, (np.ndarray, tuple, list)):
+        return _np.repeat(value, n_len)
+    elif isinstance(value, (_np.ndarray, tuple, list)):
         if len(value) != n_len:
             raise ValueError(
                 f"Invalid value length ({len(value)}). ",
                 f"Expected length is ({n_len})",
             )
-        return np.asarray(value)
+        return _np.asarray(value)
     else:
         raise TypeError(
             f"Invalid value type ({type(value)}). "
@@ -91,12 +97,12 @@ def enforce_len(value, n_len):
 
 
 def unique_rows(
-    in_arr,
-    return_index=True,
-    return_inverse=True,
-    return_counts=True,
-    dtype_name=None,
-):
+    in_arr: _np.ndarray,
+    return_index: bool = True,
+    return_inverse: bool = True,
+    return_counts: bool = True,
+    dtype_name: str = None,
+) -> list[_np.ndarray]:
     """
     Find unique rows using np.unique, but apply tricks. Adapted from
     `skimage.util.unique_rows`.
@@ -118,7 +124,7 @@ def unique_rows(
     unique_inv: (t,) np.ndarray
     """
     if dtype_name is None:
-        dtype_name = settings.INT_DTYPE
+        dtype_name = _settings.INT_DTYPE
 
     in_arr = make_c_contiguous(in_arr, dtype_name)
 
@@ -127,7 +133,7 @@ def unique_rows(
 
     in_arr_row_view = in_arr.view(f"|S{in_arr.itemsize * in_arr.shape[1]}")
 
-    unique_stuff = np.unique(
+    unique_stuff = _np.unique(
         in_arr_row_view,
         return_index=True,
         return_inverse=return_inverse,
@@ -145,8 +151,12 @@ def unique_rows(
 
 
 def close_rows(
-    arr, tolerance=None, return_intersection=False, nthreads=None, **_kwargs
-):
+    arr: _TrackedArray,
+    tolerance: float = None,
+    return_intersection: bool = False,
+    nthreads: _Any | None = None,
+    **_kwargs: _Any,
+) -> tuple[_np.ndarray, _np.ndarray, _np.ndarray, list | _UIntVectorVector]:
     """Similar to unique_rows, but if data type is floats, use this one.
     Performs radius search using KDTree. Currently uses
     `scipy.spatial.cKDTree`.
@@ -172,10 +182,10 @@ def close_rows(
         id of neighbors within the tolerance.
     """
     if tolerance is None:
-        tolerance = settings.TOLERANCE
+        tolerance = _settings.TOLERANCE
 
     if nthreads is None:
-        nthreads = settings.NTHREADS
+        nthreads = _settings.NTHREADS
 
     if has_funi and not return_intersection:
         return (
@@ -205,13 +215,13 @@ def close_rows(
         )
 
         # inverse based on original vertices.
-        o_inverse = np.array(
+        o_inverse = _np.array(
             [n[0] for n in neighbors],
-            dtype=settings.INT_DTYPE,
+            dtype=_settings.INT_DTYPE,
         )
 
         # unique of o_inverse, and inverse based on that
-        (_, uniq_id, inv) = np.unique(
+        (_, uniq_id, inv) = _np.unique(
             o_inverse,
             return_index=True,
             return_inverse=True,
@@ -228,7 +238,7 @@ def close_rows(
     )
 
 
-def bounds(arr):
+def bounds(arr: _TrackedArray) -> _np.ndarray:
     """Return bounds.
 
     Parameters
@@ -239,10 +249,10 @@ def bounds(arr):
     --------
     bounds: (2, d) np.ndarray
     """
-    return np.vstack(
+    return _np.vstack(
         (
-            np.min(arr, axis=0).ravel(),
-            np.max(arr, axis=0).ravel(),
+            _np.min(arr, axis=0).ravel(),
+            _np.max(arr, axis=0).ravel(),
         )
     )
 
@@ -275,7 +285,7 @@ def bounds_norm(arr):
     --------
     bounds_norm: float
     """
-    return np.linalg.norm(bounds_diagonal(arr))
+    return _np.linalg.norm(bounds_diagonal(arr))
 
 
 def bounds_mean(arr):
@@ -289,10 +299,12 @@ def bounds_mean(arr):
     --------
     bounds_mean: (n,) array-like
     """
-    return np.mean(bounds(arr), axis=0)
+    return _np.mean(bounds(arr), axis=0)
 
 
-def select_with_ranges(arr, ranges):
+def select_with_ranges(
+    arr: _TrackedArray, ranges: list[list[float]]
+) -> _np.ndarray:
     """Select array with ranges of each column. Always parsed as:
 
     [[greater_than, less_than], [....], ...]
@@ -315,21 +327,21 @@ def select_with_ranges(arr, ranges):
             lower = arr[:, i] > r[0]
             upper = arr[:, i] < r[1]
             if r[1] > r[0]:
-                masks.append(np.logical_and(lower, upper))
+                masks.append(_np.logical_and(lower, upper))
             else:
-                masks.append(np.logical_or(lower, upper))
+                masks.append(_np.logical_or(lower, upper))
 
     if len(masks) > 1:
-        mask = np.zeros(arr.shape[0], dtype=bool)
+        mask = _np.zeros(arr.shape[0], dtype=bool)
         for i, m in enumerate(masks):
             mask = (
-                np.logical_or(mask, m) if i == 0 else np.logical_and(mask, m)
+                _np.logical_or(mask, m) if i == 0 else _np.logical_and(mask, m)
             )
 
     else:
         mask = masks[0]
 
-    return np.arange(arr.shape[0])[mask]
+    return _np.arange(arr.shape[0])[mask]
 
 
 def rotation_matrix(rotation, degree=True):
@@ -343,7 +355,7 @@ def rotation_matrix(rotation, degree=True):
       Amount of rotation along [x,y,z] axis. Default is in degrees.
       In 2D, it can be float.
     degree: bool
-      (Optional) rotation given in degrees.
+      (_Optional) rotation given in degrees.
       Default is `True`. If `False`, in radian.
 
     Returns
@@ -352,10 +364,10 @@ def rotation_matrix(rotation, degree=True):
     """
     from scipy.spatial.transform import Rotation as R
 
-    rotation = np.asarray(rotation).ravel()
+    rotation = _np.asarray(rotation).ravel()
 
     if degree:
-        rotation = np.radians(rotation)
+        rotation = _np.radians(rotation)
 
     # 2D
     if len(rotation) == 1:
@@ -382,16 +394,16 @@ def rotate(arr, rotation, rotation_axis=None, degree=True):
     --------
     rotated_points: (n, d) np.ndarray
     """
-    arr = make_c_contiguous(arr, settings.FLOAT_DTYPE)
+    arr = make_c_contiguous(arr, _settings.FLOAT_DTYPE)
     if rotation_axis is not None:
-        rotation_axis = np.asanyarray(rotation_axis)
+        rotation_axis = _np.asanyarray(rotation_axis)
 
     if rotation_axis is None:
-        return np.matmul(arr, rotation_matrix(rotation, degree))
+        return _np.matmul(arr, rotation_matrix(rotation, degree))
 
     else:
         rotated_array = arr - rotation_axis
-        rotated_array = np.matmul(
+        rotated_array = _np.matmul(
             rotated_array, rotation_matrix(rotation, degree)
         )
         rotated_array += rotation_axis
@@ -412,7 +424,7 @@ def rotation_matrix_around_axis(axis=None, rotation=None, degree=True):
     rotation: float
       angle of rotation in either radiant or degrees
     degree: bool
-      (Optional) rotation given in degrees.
+      (_Optional) rotation given in degrees.
       Default is `True`. If `False`, in radian.
 
     Returns
@@ -423,28 +435,28 @@ def rotation_matrix_around_axis(axis=None, rotation=None, degree=True):
     if rotation is None:
         raise ValueError("No rotation angle specified.")
     elif degree:
-        rotation = np.radians(rotation)
+        rotation = _np.radians(rotation)
 
     # Check Problem dimensions
     if axis is None:
         problem_dimension = 2
     else:
-        axis = np.asarray(axis).ravel()
+        axis = _np.asarray(axis).ravel()
         if axis.shape[0] != 3:
             raise ValueError("Axis dimension must be 3D")
         problem_dimension = 3
 
     # Assemble rotation matrix
     if problem_dimension == 2:
-        rotation_matrix = np.array(
+        rotation_matrix = _np.array(
             [
-                [np.cos(rotation), -np.sin(rotation)],
-                [np.sin(rotation), np.cos(rotation)],
+                [_np.cos(rotation), -_np.sin(rotation)],
+                [_np.sin(rotation), _np.cos(rotation)],
             ]
         )
     else:
         # See Rodrigues' formula
-        rotation_matrix = np.array(
+        rotation_matrix = _np.array(
             [
                 [0, -axis[2], axis[1]],
                 [axis[2], 0, -axis[0]],
@@ -452,18 +464,22 @@ def rotation_matrix_around_axis(axis=None, rotation=None, degree=True):
             ]
         )
         rotation_matrix = (
-            np.eye(3)
-            + np.sin(rotation) * rotation_matrix
+            _np.eye(3)
+            + _np.sin(rotation) * rotation_matrix
             + (
-                (1 - np.cos(rotation))
-                * np.matmul(rotation_matrix, rotation_matrix)
+                (1 - _np.cos(rotation))
+                * _np.matmul(rotation_matrix, rotation_matrix)
             )
         )
 
     return rotation_matrix
 
 
-def is_shape(arr, shape, strict=False):
+def is_shape(
+    arr: _TrackedArray | _np.ndarray,
+    shape: tuple[int, int],
+    strict: bool = False,
+) -> bool:
     """Checks if arr matches given shape. shape can have negative numbers.
 
     Parameters
@@ -477,7 +493,7 @@ def is_shape(arr, shape, strict=False):
     --------
     matches: bool
     """
-    arr = np.asanyarray(arr)
+    arr = _np.asanyarray(arr)
 
     if arr.ndim != len(shape):
         if strict:
@@ -495,8 +511,12 @@ def is_shape(arr, shape, strict=False):
     return True
 
 
-def is_one_of_shapes(arr, shapes, strict=False):
-    """Tuple/list of given shapes, iterates and checks with is_shape. Useful if
+def is_one_of_shapes(
+    arr: list[list[int]] | _TrackedArray | _np.ndarray,
+    shapes: tuple[tuple[int, int], tuple[int, int]],
+    strict: bool = False,
+) -> bool:
+    """tuple/list of given shapes, iterates and checks with is_shape. Useful if
     you have multiple acceptable shapes.
 
     Parameters
@@ -510,7 +530,7 @@ def is_one_of_shapes(arr, shapes, strict=False):
     --------
     matches: bool
     """
-    arr = np.asanyarray(arr)
+    arr = _np.asanyarray(arr)
     matches = False
     for s in shapes:
         m = is_shape(arr, s, strict=False)
@@ -549,7 +569,7 @@ def derivatives_to_normals(derivatives, normalize=True):
     # 2D is simple index flip
     if shape[2] == 2:
         der = derivatives.reshape(-1, shape[2])
-        normals = np.empty_like(der)
+        normals = _np.empty_like(der)
         normals[:, 0] = der[:, 1]
         normals[:, 1] = -der[:, 0]
 
@@ -558,7 +578,7 @@ def derivatives_to_normals(derivatives, normalize=True):
         normals = cross3d(der[::2], der[1::2])
 
     if normalize:
-        normals /= np.linalg.norm(normals, axis=1).reshape(-1, 1)
+        normals /= _np.linalg.norm(normals, axis=1).reshape(-1, 1)
 
     return normals
 
@@ -580,28 +600,28 @@ def cross3d(a, b):
     # (1 5 - 2 4, 2 3 - 0 5, 0 4 - 1 3).
     # or from two arrays
     # (1 2 - 2 1, 2 0 - 0 2, 0 1 - 1 0).
-    o = np.empty_like(a)
+    o = _np.empty_like(a)
 
     # temporary aux arrays
     size = len(a)
-    t0 = np.empty(size)
-    t1 = np.empty(size)
+    t0 = _np.empty(size)
+    t1 = _np.empty(size)
 
     # short cuts
     a0, a1, a2 = a[..., 0], a[..., 1], a[..., 2]
     b0, b1, b2 = b[..., 0], b[..., 1], b[..., 2]
     o0, o1, o2 = o[..., 0], o[..., 1], o[..., 2]
 
-    np.multiply(a1, b2, out=t0)
-    np.multiply(a2, b1, out=t1)
-    np.subtract(t0, t1, out=o0)
+    _np.multiply(a1, b2, out=t0)
+    _np.multiply(a2, b1, out=t1)
+    _np.subtract(t0, t1, out=o0)
 
-    np.multiply(a2, b0, out=t0)
-    np.multiply(a0, b2, out=t1)
-    np.subtract(t0, t1, out=o1)
+    _np.multiply(a2, b0, out=t0)
+    _np.multiply(a0, b2, out=t1)
+    _np.subtract(t0, t1, out=o1)
 
-    np.multiply(a0, b1, out=t0)
-    np.multiply(a1, b0, out=t1)
-    np.subtract(t0, t1, out=o2)
+    _np.multiply(a0, b1, out=t0)
+    _np.multiply(a1, b0, out=t1)
+    _np.subtract(t0, t1, out=o2)
 
     return o
diff --git a/gustaf/utils/connec.py b/gustaf/utils/connec.py
index 99f337f5..127cf2e2 100644
--- a/gustaf/utils/connec.py
+++ b/gustaf/utils/connec.py
@@ -5,9 +5,15 @@
 even cooler, if it was palindrome.
 """
 
-import collections
+from __future__ import annotations
 
-import numpy as np
+import collections as _collections
+from typing import Any as _Any
+
+import numpy as _np
+
+from gustaf.helpers.data import TrackedArray as _TrackedArray
+from gustaf.helpers.data import Unique2DIntegers as _Unique2DIntegers
 
 try:
     import napf
@@ -20,7 +26,7 @@
 from gustaf.utils import arr
 
 
-def tet_to_tri(volumes):
+def tet_to_tri(volumes: _TrackedArray | _np.ndarray) -> _np.ndarray:
     """Computes tri faces based on following index scheme.
 
     ``Tetrahedron``
@@ -53,14 +59,14 @@ def tet_to_tri(volumes):
     --------
     faces: (n * 4, 3) np.ndarray
     """
-    volumes = np.asanyarray(volumes, settings.INT_DTYPE)
+    volumes = _np.asanyarray(volumes, settings.INT_DTYPE)
 
     if volumes.ndim != 2 or volumes.shape[1] != 4:
         raise ValueError("Given volumes are not `tet` volumes")
 
     fpe = 4  # faces per element
     faces = (
-        np.ones(((volumes.shape[0] * fpe), 3), dtype=settings.INT_DTYPE) * -1
+        _np.ones(((volumes.shape[0] * fpe), 3), dtype=settings.INT_DTYPE) * -1
     )  # -1 for safety check
 
     faces[:, 0] = volumes.ravel()
@@ -75,7 +81,7 @@ def tet_to_tri(volumes):
     return faces
 
 
-def hexa_to_quad(volumes):
+def hexa_to_quad(volumes: _TrackedArray | _np.ndarray) -> _np.ndarray:
     """Computes quad faces based on following index scheme.
 
     ``Hexahedron``
@@ -113,13 +119,13 @@ def hexa_to_quad(volumes):
     --------
     faces: (n * 8, 4) np.ndarray
     """
-    volumes = np.asanyarray(volumes, settings.INT_DTYPE)
+    volumes = _np.asanyarray(volumes, settings.INT_DTYPE)
 
     if volumes.ndim != 2 or volumes.shape[1] != 8:
         raise ValueError("Given volumes are not `hexa` volumes")
 
     fpe = 6  # faces per element
-    faces = np.empty(((volumes.shape[0] * fpe), 4), dtype=settings.INT_DTYPE)
+    faces = _np.empty(((volumes.shape[0] * fpe), 4), dtype=settings.INT_DTYPE)
 
     faces[::fpe] = volumes[:, [1, 0, 3, 2]]
     faces[1::fpe] = volumes[:, [0, 1, 5, 4]]
@@ -142,7 +148,7 @@ def volumes_to_faces(volumes):
     --------
     faces: (n*4, 3) or (m*6, 4) np.ndarray
     """
-    volumes = np.asanyarray(volumes, settings.INT_DTYPE)
+    volumes = _np.asanyarray(volumes, settings.INT_DTYPE)
     if volumes.shape[1] == 4:
         return tet_to_tri(volumes)
 
@@ -150,7 +156,7 @@ def volumes_to_faces(volumes):
         return hexa_to_quad(volumes)
 
 
-def faces_to_edges(faces):
+def faces_to_edges(faces: _TrackedArray | _np.ndarray) -> _np.ndarray:
     """Compute edges based on following edge scheme.
 
     .. code-block:: text
@@ -190,7 +196,7 @@ def faces_to_edges(faces):
     vertices_per_face = faces.shape[1]
 
     num_edges = int(num_faces * vertices_per_face)
-    edges = np.empty((num_edges, 2), dtype=settings.INT_DTYPE)
+    edges = _np.empty((num_edges, 2), dtype=settings.INT_DTYPE)
 
     edges[:, 0] = faces.ravel()
 
@@ -203,7 +209,9 @@ def faces_to_edges(faces):
     return edges
 
 
-def range_to_edges(range_, closed=False, continuous=True):
+def range_to_edges(
+    range_: tuple[int, int], closed: bool = False, continuous: bool = True
+) -> _np.ndarray:
     """Given range, for example (a, b), returns an edge sequence that
     sequentially connects indices. If int is given as range, it is considered
     as (0, value). Used to be called "closed/open_loop_index_train".
@@ -219,10 +227,10 @@ def range_to_edges(range_, closed=False, continuous=True):
     edges: (n, 2) np.ndarray
     """
     if isinstance(range_, int):
-        indices = np.arange(range_, dtype=settings.INT_DTYPE)
+        indices = _np.arange(range_, dtype=settings.INT_DTYPE)
     elif isinstance(range_, (list, tuple)):
         # pass range_ as is and check for valid output
-        indices = np.arange(*range_, dtype=settings.INT_DTYPE)
+        indices = _np.arange(*range_, dtype=settings.INT_DTYPE)
         if len(indices) < 2:
             raise ValueError(
                 f"{range_} is invalid range input. "
@@ -241,7 +249,7 @@ def range_to_edges(range_, closed=False, continuous=True):
     return sequence_to_edges(indices, closed)
 
 
-def sequence_to_edges(seq, closed=False):
+def sequence_to_edges(seq: _np.ndarray, closed: bool = False) -> _np.ndarray:
     """Given a sequence of int, "connect" to turn them into edges.
 
     Parameters
@@ -253,7 +261,7 @@ def sequence_to_edges(seq, closed=False):
     --------
     edges: (m, 2) np.ndarray
     """
-    edges = np.repeat(seq, 2)
+    edges = _np.repeat(seq, 2)
 
     if closed:
         first = int(edges[0])  # this is redundant copy to ensure detaching
@@ -269,7 +277,7 @@ def sequence_to_edges(seq, closed=False):
     return edges.reshape(-1, 2)
 
 
-def make_quad_faces(resolutions):
+def make_quad_faces(resolutions: list[int]) -> _np.ndarray:
     """Given number of nodes per each dimension, returns connectivity
     information of a structured mesh. Counter clock wise connectivity.
 
@@ -288,20 +296,20 @@ def make_quad_faces(resolutions):
     -------
     faces: (n, 4) np.ndarray
     """
-    nnpd = np.asarray(resolutions)  # number of nodes per dimension
+    nnpd = _np.asarray(resolutions)  # number of nodes per dimension
     if any(nnpd < 1):
         raise ValueError(f"The number of nodes per dimension is wrong: {nnpd}")
 
-    total_nodes = np.prod(nnpd)
+    total_nodes = _np.prod(nnpd)
     total_faces = (nnpd[0] - 1) * (nnpd[1] - 1)
     try:
-        node_indices = np.arange(
+        node_indices = _np.arange(
             total_nodes, dtype=settings.INT_DTYPE
         ).reshape(nnpd[1], nnpd[0])
     except ValueError as e:
         raise ValueError(f"Problem with generating node indices. {e}")
 
-    faces = np.empty((total_faces, 4), dtype=settings.INT_DTYPE)
+    faces = _np.empty((total_faces, 4), dtype=settings.INT_DTYPE)
 
     faces[:, 0] = node_indices[: (nnpd[1] - 1), : (nnpd[0] - 1)].ravel()
     faces[:, 1] = node_indices[: (nnpd[1] - 1), 1 : nnpd[0]].ravel()
@@ -311,7 +319,7 @@ def make_quad_faces(resolutions):
     return faces
 
 
-def make_hexa_volumes(resolutions):
+def make_hexa_volumes(resolutions: list[int]) -> _np.ndarray:
     """Given number of nodes per each dimension, returns connectivity
     information of structured hexahedron elements. Counter clock wise
     connectivity.
@@ -335,17 +343,17 @@ def make_hexa_volumes(resolutions):
     --------
     elements: (n, 8) np.ndarray
     """
-    nnpd = np.asarray(resolutions)  # number of nodes per dimension
+    nnpd = _np.asarray(resolutions)  # number of nodes per dimension
     if any(nnpd < 1):
         raise ValueError(f"The number of nodes per dimension is wrong: {nnpd}")
 
-    total_nodes = np.prod(nnpd)
-    total_volumes = np.prod(nnpd - 1)
-    node_indices = np.arange(total_nodes, dtype=settings.INT_DTYPE).reshape(
+    total_nodes = _np.prod(nnpd)
+    total_volumes = _np.prod(nnpd - 1)
+    node_indices = _np.arange(total_nodes, dtype=settings.INT_DTYPE).reshape(
         nnpd[::-1]
     )
 
-    volumes = np.empty((total_volumes, 8), dtype=settings.INT_DTYPE)
+    volumes = _np.empty((total_volumes, 8), dtype=settings.INT_DTYPE)
 
     volumes[:, 0] = node_indices[
         : (nnpd[2] - 1), : (nnpd[1] - 1), : (nnpd[0] - 1)
@@ -452,14 +460,14 @@ def subdivide_tri(mesh, return_dict=False):
 
     # Form new vertices
     edge_mid_v = mesh.vertices[mesh.unique_edges().values].mean(axis=1)
-    new_vertices = np.vstack((mesh.vertices, edge_mid_v))
+    new_vertices = _np.vstack((mesh.vertices, edge_mid_v))
 
-    subdivided_faces = np.empty(
+    subdivided_faces = _np.empty(
         (mesh.faces.shape[0] * 4, mesh.faces.shape[1]),
         dtype=settings.INT_DTYPE,
     )
 
-    mask = np.ones(subdivided_faces.shape[0], dtype=bool)
+    mask = _np.ones(subdivided_faces.shape[0], dtype=bool)
     mask[3::4] = False
 
     # 0th column minus (every 4th row, starting from 3rd row)
@@ -512,7 +520,7 @@ def subdivide_quad(
     # Form new vertices
     edge_mid_v = mesh.vertices[mesh.unique_edges().values].mean(axis=1)
     face_centers = mesh.centers()
-    new_vertices = np.vstack(
+    new_vertices = _np.vstack(
         (
             mesh.vertices,
             edge_mid_v,
@@ -520,15 +528,15 @@ def subdivide_quad(
         )
     )
 
-    subdivided_faces = np.empty(
+    subdivided_faces = _np.empty(
         (mesh.faces.shape[0] * 4, mesh.faces.shape[1]),
         dtype=settings.INT_DTYPE,
     )
 
     subdivided_faces[:, 0] = mesh.faces.ravel()
     subdivided_faces[:, 1] = mesh.unique_edges().inverse + len(mesh.vertices)
-    subdivided_faces[:, 2] = np.repeat(
-        np.arange(len(face_centers)) + (len(mesh.vertices) + len(edge_mid_v)),
+    subdivided_faces[:, 2] = _np.repeat(
+        _np.arange(len(face_centers)) + (len(mesh.vertices) + len(edge_mid_v)),
         4,
         dtype=settings.INT_DTYPE,
     )
@@ -546,7 +554,9 @@ def subdivide_quad(
         return new_vertices, subdivided_faces
 
 
-def sorted_unique(connectivity, sorted_=False):
+def sorted_unique(
+    connectivity: _np.ndarray, sorted_: bool = False
+) -> _Unique2DIntegers:
     """Given connectivity array, finds unique entries, based on its axis=1
     sorted values. Returned value will be sorted.
 
@@ -559,7 +569,7 @@ def sorted_unique(connectivity, sorted_=False):
     --------
     unique_info: Unique2DIntegers
     """
-    s_connec = connectivity if sorted_ else np.sort(connectivity, axis=1)
+    s_connec = connectivity if sorted_ else _np.sort(connectivity, axis=1)
 
     unique_stuff = arr.unique_rows(
         s_connec,
@@ -577,15 +587,23 @@ def sorted_unique(connectivity, sorted_=False):
     )
 
 
-def _sequentialize_directed_edges(edges, start=None, return_edges=False):
+def _sequentialize_directed_edges(
+    edges: _np.ndarray,
+    start: _Any | None = None,
+    return_edges: bool = False,
+) -> tuple[
+    list[list[int | _np.int32]]
+    | list[list[int | _np.int32] | list[_np.int32 | _np.int64]],
+    list[bool],
+]:
     """
     Sequentialize directed edges.
     """
     # we want to have an np array
-    edges = np.asanyarray(edges)
+    edges = _np.asanyarray(edges)
 
     # Build a lookup_array
-    lookup_array = np.full(edges.max() + 1, -1, dtype=settings.INT_DTYPE)
+    lookup_array = _np.full(edges.max() + 1, -1, dtype=settings.INT_DTYPE)
     lookup_array[edges[:, 0]] = edges[:, 1]
 
     # select starting point - lowest index
@@ -594,10 +612,10 @@ def _sequentialize_directed_edges(edges, start=None, return_edges=False):
     # initialize a set to keep track of processes vertices
     next_candidates = set(edges[:, 0])
     # we want to keep track of single occurrences, as they are line start
-    line_starts = set(np.where(np.bincount(edges.ravel()) == 1)[0])
+    line_starts = set(_np.where(_np.bincount(edges.ravel()) == 1)[0])
     # for this to work, we can't have a line that starts at column 1.
     # so, we remove those.
-    ls_col1 = set(np.where(np.bincount(edges[:, 1].ravel()) == 1)[0])
+    ls_col1 = set(_np.where(_np.bincount(edges[:, 1].ravel()) == 1)[0])
     line_starts.difference_update(ls_col1)
 
     polygons = []
@@ -650,15 +668,23 @@ def _sequentialize_directed_edges(edges, start=None, return_edges=False):
         return polygon_edges, is_polygon
 
 
-def _sequentialize_edges(edges, start=None, return_edges=False):
+def _sequentialize_edges(
+    edges: _np.ndarray,
+    start: _Any | None = None,
+    return_edges: bool = False,
+) -> tuple[
+    list[list[int | _np.int64]]
+    | list[list[int | _np.int64] | list[_np.int64]],
+    list[bool],
+]:
     """
     sequentialize undirected edges. No overlaps are allowed, for now.
     """
-    edges = np.asanyarray(edges)
+    edges = _np.asanyarray(edges)
 
     # only applicable to closed polygons and open lines
     # not for arbitrarily connected edges
-    bc = np.bincount(edges.ravel())
+    bc = _np.bincount(edges.ravel())
     if not all(bc < 3):
         raise ValueError(
             "This function currently supports individual lines/polygons "
@@ -666,25 +692,25 @@ def _sequentialize_edges(edges, start=None, return_edges=False):
         )
 
     # we want to keep track of single occurrences, as they are line start
-    line_starts = set(np.where(bc == 1)[0])
+    line_starts = set(_np.where(bc == 1)[0])
 
     # initialize a set to keep track of processes vertices
     next_candidates = set(edges.ravel())
 
     # create a look up to each edge column
-    edge_col = collections.namedtuple("a", "b")
+    edge_col = _collections.namedtuple("a", "b")
     edge_col.a = edges[:, 0]
     edge_col.b = edges[:, 1]
 
     # create trees for each edge column
-    tree = collections.namedtuple("a", "b")
+    tree = _collections.namedtuple("a", "b")
     tree.a = napf.KDT(edge_col.a.reshape(-1, 1))
     tree.b = napf.KDT(edge_col.b.reshape(-1, 1))
 
     # radius search size
     r = 0.1
 
-    current_id = np.argmin(edge_col.a) if start is None else start
+    current_id = _np.argmin(edge_col.a) if start is None else start
     start_value = int(edge_col.a[current_id])
     other_col = edge_col.b
 
@@ -788,7 +814,12 @@ def _sequentialize_edges(edges, start=None, return_edges=False):
         return polygon_edges
 
 
-def sequentialize_edges(edges, start=None, return_edges=False, directed=False):
+def sequentialize_edges(
+    edges: _np.ndarray,
+    start: _Any | None = None,
+    return_edges: bool = False,
+    directed: bool = False,
+) -> tuple[_Any, list[bool]]:
     """
     Organize edge connectivities to describe polygon or a line.
     This supports edges that describes separated/individual polygons and lines.
@@ -798,11 +829,11 @@ def sequentialize_edges(edges, start=None, return_edges=False, directed=False):
     -----------
     edges: (n, 2) list-like
     start: int
-      (Optional) Specify starting point. It will take minimum index otherwise.
+      (_Optional) Specify starting point. It will take minimum index otherwise.
     return_edges: bool
-      (Optional) Default is False. If set True, returns sequences as edges.
+      (_Optional) Default is False. If set True, returns sequences as edges.
     directed: bool
-      (Optional) Default is False. Set True, if given edges are directed.
+      (_Optional) Default is False. Set True, if given edges are directed.
       It should return the result faster.
 
     Returns
diff --git a/gustaf/utils/log.py b/gustaf/utils/log.py
index ea53c748..3c085c59 100644
--- a/gustaf/utils/log.py
+++ b/gustaf/utils/log.py
@@ -3,8 +3,8 @@
 Thin logging wrapper.
 """
 
-import logging
-from functools import partial
+import logging as _logging
+from functools import partial as _partial
 
 
 def configure(debug=False, logfile=None):
@@ -20,14 +20,14 @@ def configure(debug=False, logfile=None):
     None
     """
     # logger
-    logger = logging.getLogger("gustaf")
+    logger = _logging.getLogger("gustaf")
 
     # level
-    level = logging.DEBUG if debug else logging.INFO
+    level = _logging.DEBUG if debug else _logging.INFO
     logger.setLevel(level)
 
     # format
-    formatter = logging.Formatter(fmt="%(name)s [%(levelname)s] %(message)s")
+    formatter = _logging.Formatter(fmt="%(name)s [%(levelname)s] %(message)s")
 
     # apply format using stream handler
     # let's use only one stream handler so that calling configure multiple
@@ -35,7 +35,7 @@ def configure(debug=False, logfile=None):
     new_handlers = []
     for _i, h in enumerate(logger.handlers):
         # we skip all the stream handler.
-        if isinstance(h, logging.StreamHandler):
+        if isinstance(h, _logging.StreamHandler):
             continue
 
         # blindly keep other ones
@@ -43,7 +43,7 @@ def configure(debug=False, logfile=None):
             new_handlers.append(h)
 
     # add new stream handler
-    stream_handler = logging.StreamHandler()
+    stream_handler = _logging.StreamHandler()
     stream_handler.setLevel(level)
     stream_handler.setFormatter(formatter)
     new_handlers.append(stream_handler)
@@ -52,11 +52,11 @@ def configure(debug=False, logfile=None):
 
     # output logs
     if logfile is not None:
-        file_logger_handler = logging.FileHandler(logfile)
+        file_logger_handler = _logging.FileHandler(logfile)
         logger.addHandler(file_logger_handler)
 
 
-def debug(*log):
+def debug(*log: str) -> None:
     """Debug logger.
 
     Parameters
@@ -67,7 +67,7 @@ def debug(*log):
     --------
     None
     """
-    logger = logging.getLogger("gustaf")
+    logger = _logging.getLogger("gustaf")
     logger.debug(" ".join(map(str, log)))
 
 
@@ -82,7 +82,7 @@ def info(*log):
     --------
     None
     """
-    logger = logging.getLogger("gustaf")
+    logger = _logging.getLogger("gustaf")
     logger.info(" ".join(map(str, log)))
 
 
@@ -97,7 +97,7 @@ def warning(*log):
     --------
     None
     """
-    logger = logging.getLogger("gustaf")
+    logger = _logging.getLogger("gustaf")
     logger.warning(" ".join(map(str, log)))
 
 
@@ -115,4 +115,4 @@ def prepended_log(message, log_func):
     -------
     prepended: function
     """
-    return partial(log_func, message)
+    return _partial(log_func, message)
diff --git a/gustaf/vertices.py b/gustaf/vertices.py
index 762ea8d6..9e071415 100644
--- a/gustaf/vertices.py
+++ b/gustaf/vertices.py
@@ -3,35 +3,54 @@
 Vertices. Base of all "Mesh" geometries.
 """
 
+from __future__ import annotations
+
 import copy
+from typing import TYPE_CHECKING
+
+import numpy as _np
+
+from gustaf import helpers as _helpers
+from gustaf import settings as _settings
+from gustaf import show as _show
+from gustaf import utils as _utils
+from gustaf._base import GustafBase as _GustafBase
+from gustaf.helpers.data import TrackedArray as _TrackedArray
+from gustaf.helpers.data import VertexData as _VertexData
+from gustaf.helpers.options import Option as _Option
+
+if TYPE_CHECKING:
+    import sys
+    from typing import Any
 
-import numpy as np
+    if sys.version_info >= (3, 10):
+        from typing import Self
+    else:
+        from typing_extensions import Self
 
-from gustaf import helpers, settings, show, utils
-from gustaf._base import GustafBase
-from gustaf.helpers.options import Option
+    from gustaf.helpers.data import Unique2DFloats as _Unique2DFloats
 
 
-class VerticesShowOption(helpers.options.ShowOption):
+class VerticesShowOption(_helpers.options.ShowOption):
     """
     Show options for vertices.
     """
 
-    _valid_options = helpers.options.make_valid_options(
-        *helpers.options.vedo_common_options,
-        Option(
+    _valid_options = _helpers.options.make_valid_options(
+        *_helpers.options.vedo_common_options,
+        _Option(
             "vedo",
             "r",
             "Radius of vertices in units of pixels.",
             (float, int),
         ),
-        Option(
+        _Option(
             "vedo",
             "labels",
             "Places a label/description str at the place of vertices.",
-            (np.ndarray, tuple, list),
+            (_np.ndarray, tuple, list),
         ),
-        Option(
+        _Option(
             "vedo",
             "label_options",
             "Label kwargs to be passed during initialization."
@@ -61,7 +80,7 @@ def _initialize_showable(self):
         """
         init_options = ("r",)
 
-        vertices = show.vedo.Points(
+        vertices = _show.vedo.Points(
             self._helpee.const_vertices, **self[init_options]
         )
 
@@ -87,7 +106,7 @@ def _initialize_showable(self):
             return vertices
 
 
-class Vertices(GustafBase):
+class Vertices(_GustafBase):
     kind = "vertex"
 
     __slots__ = (
@@ -103,8 +122,8 @@ class Vertices(GustafBase):
 
     def __init__(
         self,
-        vertices=None,
-    ):
+        vertices: list[list[float]] | _TrackedArray | _np.ndarray = None,
+    ) -> None:
         """Vertices. It has vertices.
 
         Parameters
@@ -119,12 +138,12 @@ def __init__(
         self.vertices = vertices
 
         # init helpers
-        self._vertex_data = helpers.data.VertexData(self)
-        self._computed = helpers.data.ComputedMeshData(self)
+        self._vertex_data = _helpers.data.VertexData(self)
+        self._computed = _helpers.data.ComputedMeshData(self)
         self._show_options = self.__show_option__(self)
 
     @property
-    def vertices(self):
+    def vertices(self) -> _TrackedArray:
         """Returns vertices.
 
         Parameters
@@ -139,7 +158,9 @@ def vertices(self):
         return self._vertices
 
     @vertices.setter
-    def vertices(self, vs):
+    def vertices(
+        self, vs: list[list[float]] | _TrackedArray | _np.ndarray
+    ) -> None:
         """Vertices setter. This will saved as a tracked array. This tracked
         array is very sensitive and if we do anything with it that may hint an
         inplace operation, it will be marked as modified. This includes copying
@@ -157,13 +178,13 @@ def vertices(self, vs):
         self._logd("setting vertices")
 
         # we try not to make copy.
-        self._vertices = helpers.data.make_tracked_array(
-            vs, settings.FLOAT_DTYPE, copy=False
+        self._vertices = _helpers.data.make_tracked_array(
+            vs, _settings.FLOAT_DTYPE, copy=False
         )
 
         # shape check
         if self._vertices.size > 0:
-            utils.arr.is_shape(self._vertices, (-1, -1), strict=True)
+            _utils.arr.is_shape(self._vertices, (-1, -1), strict=True)
 
         # exact same, but not tracked.
         self._const_vertices = self._vertices.view()
@@ -175,7 +196,7 @@ def vertices(self, vs):
             self.vertex_data._validate_len(raise_=False)
 
     @property
-    def const_vertices(self):
+    def const_vertices(self) -> _TrackedArray:
         """Returns non-mutable view of `vertices`. Naming inspired by c/cpp
         sessions.
 
@@ -191,7 +212,7 @@ def const_vertices(self):
         return self._const_vertices
 
     @property
-    def vertex_data(self):
+    def vertex_data(self) -> _VertexData:
         """
         Returns vertex_data manager. Behaves similar to dict() and can be used
         to store values/data associated with each vertex.
@@ -208,7 +229,7 @@ def vertex_data(self):
         return self._vertex_data
 
     @property
-    def show_options(self):
+    def show_options(self) -> _helpers.options.ShowOption:
         """
         Returns a show option manager for this object. Behaves similar to
         dict.
@@ -226,7 +247,7 @@ def show_options(self):
         return self._show_options
 
     @property
-    def whatami(self):
+    def whatami(self) -> str:
         """Answers deep philosophical question: "what am i"?
 
         Parameters
@@ -240,8 +261,10 @@ def whatami(self):
         """
         return "vertices"
 
-    @helpers.data.ComputedMeshData.depends_on(["vertices"])
-    def unique_vertices(self, tolerance=None, **kwargs):
+    @_helpers.data.ComputedMeshData.depends_on(["vertices"])
+    def unique_vertices(
+        self, tolerance: float | None = None, **kwargs: Any
+    ) -> _Unique2DFloats:
         """Returns a namedtuple that holds unique vertices info. Unique here
         means "close-enough-within-tolerance".
 
@@ -259,21 +282,21 @@ def unique_vertices(self, tolerance=None, **kwargs):
         """
         self._logd("computing unique vertices")
         if tolerance is None:
-            tolerance = settings.TOLERANCE
+            tolerance = _settings.TOLERANCE
 
-        values, ids, inverse, intersection = utils.arr.close_rows(
+        values, ids, inverse, intersection = _utils.arr.close_rows(
             self.const_vertices, tolerance=tolerance, **kwargs
         )
 
-        return helpers.data.Unique2DFloats(
+        return _helpers.data.Unique2DFloats(
             values,
             ids,
             inverse,
             intersection,
         )
 
-    @helpers.data.ComputedMeshData.depends_on(["vertices"])
-    def bounds(self):
+    @_helpers.data.ComputedMeshData.depends_on(["vertices"])
+    def bounds(self) -> _np.ndarray:
         """Returns bounds of the vertices. Bounds means AABB of the geometry.
 
         Parameters
@@ -285,10 +308,10 @@ def bounds(self):
         bounds: (d,) np.ndarray
         """
         self._logd("computing bounds")
-        return utils.arr.bounds(self.vertices)
+        return _utils.arr.bounds(self.vertices)
 
-    @helpers.data.ComputedMeshData.depends_on(["vertices"])
-    def bounds_diagonal(self):
+    @_helpers.data.ComputedMeshData.depends_on(["vertices"])
+    def bounds_diagonal(self) -> _np.ndarray:
         """Returns diagonal vector of the bounding box.
 
         Parameters
@@ -304,8 +327,8 @@ def bounds_diagonal(self):
         bounds = self.bounds()
         return bounds[1] - bounds[0]
 
-    @helpers.data.ComputedMeshData.depends_on(["vertices"])
-    def bounds_diagonal_norm(self):
+    @_helpers.data.ComputedMeshData.depends_on(["vertices"])
+    def bounds_diagonal_norm(self) -> float:
         """Returns norm of bounds diagonal.
 
         Parameters
@@ -319,7 +342,9 @@ def bounds_diagonal_norm(self):
         self._logd("computing bounds_diagonal_norm")
         return float(sum(self.bounds_diagonal() ** 2) ** 0.5)
 
-    def update_vertices(self, mask, inverse=None):
+    def update_vertices(
+        self, mask: _np.ndarray, inverse: _np.ndarray | None = None
+    ) -> Vertices:
         """Update vertices with a mask. In other words, keeps only masked
         vertices. Adapted from `github.com/mikedh/trimesh`. Updates
         connectivity accordingly too.
@@ -336,7 +361,7 @@ def update_vertices(self, mask, inverse=None):
         vertices = self.const_vertices.copy()
 
         # make mask numpy array
-        mask = np.asarray(mask)
+        mask = _np.asarray(mask)
 
         if (mask.dtype.name == "bool" and mask.all()) or len(mask) == 0:
             return self
@@ -344,12 +369,12 @@ def update_vertices(self, mask, inverse=None):
         # create inverse mask if not passed
         check_neg = False
         if inverse is None and self.kind != "vertex":
-            inverse = np.full(len(vertices), -11, dtype=settings.INT_DTYPE)
+            inverse = _np.full(len(vertices), -11, dtype=_settings.INT_DTYPE)
             check_neg = True
             if mask.dtype.kind == "b":
-                inverse[mask] = np.arange(mask.sum())
+                inverse[mask] = _np.arange(mask.sum())
             elif mask.dtype.kind == "i":
-                inverse[mask] = np.arange(len(mask))
+                inverse[mask] = _np.arange(len(mask))
             else:
                 inverse = None
 
@@ -369,13 +394,15 @@ def update_vertices(self, mask, inverse=None):
         # apply mask
         vertices = vertices[mask]
 
-        def update_vertex_data(obj, m, vertex_data):
+        def update_vertex_data(
+            obj: Vertices, mask: _np.ndarray, vertex_data: dict
+        ) -> Any:
             """apply mask to vertex data if there's any."""
-            new_data = helpers.data.VertexData(obj)
+            new_data = _helpers.data.VertexData(obj)
 
             for key, values in vertex_data.items():
                 # should work, since this is called after updating vertices
-                new_data[key] = values[m]
+                new_data[key] = values[mask]
 
             obj._vertex_data = new_data
 
@@ -394,7 +421,7 @@ def update_vertex_data(obj, m, vertex_data):
 
         return self
 
-    def select_vertices(self, ranges):
+    def select_vertices(self, ranges: list[list[float]]) -> _np.ndarray:
         """Returns vertices inside the given range.
 
         Parameters
@@ -406,9 +433,9 @@ def select_vertices(self, ranges):
         --------
         ids: (n,) np.ndarray
         """
-        return utils.arr.select_with_ranges(self.vertices, ranges)
+        return _utils.arr.select_with_ranges(self.vertices, ranges)
 
-    def remove_vertices(self, ids):
+    def remove_vertices(self, ids: _np.ndarray) -> Vertices:
         """Removes vertices with given vertex ids.
 
         Parameters
@@ -419,12 +446,14 @@ def remove_vertices(self, ids):
         --------
         new_self: type(self)
         """
-        mask = np.ones(len(self.vertices), dtype=bool)
+        mask = _np.ones(len(self.vertices), dtype=bool)
         mask[ids] = False
 
         return self.update_vertices(mask)
 
-    def merge_vertices(self, tolerance=None, **kwargs):
+    def merge_vertices(
+        self, tolerance: float | None = None, **kwargs: Any
+    ) -> Vertices:
         """Based on unique vertices, merge vertices if it is mergeable.
 
         Parameters
@@ -459,7 +488,7 @@ def showable(self, **kwargs):
         showable: obj
           Obj of `gustaf.settings.VISUALIZATION_BACKEND`
         """
-        return show.make_showable(self, **kwargs)
+        return _show.make_showable(self, **kwargs)
 
     def show(self, **kwargs):
         """Show current object using visualization backend.
@@ -473,9 +502,9 @@ def show(self, **kwargs):
         --------
         None
         """
-        return show.show(self, **kwargs)
+        return _show.show(self, **kwargs)
 
-    def copy(self):
+    def copy(self) -> Vertices:
         """Returns deepcopy of self.
 
         Parameters
@@ -497,12 +526,12 @@ def copy(self):
         return copied
 
     @classmethod
-    def concat(cls, *instances):
+    def concat(cls, *instances: Vertices) -> Vertices:
         """Sequentially put them together to make one object.
 
         Parameters
         -----------
-        *instances: List[type(cls)]
+        *instances: list[type(cls)]
           Allows one iterable object also.
 
         Returns
@@ -510,7 +539,7 @@ def concat(cls, *instances):
         one_instance: type(cls)
         """
 
-        def is_concatable(inst):
+        def is_concatable(inst: Any) -> bool:
             """Return true, if it is same as type(cls)"""
             return bool(isinstance(inst, cls))
 
@@ -555,14 +584,14 @@ def is_concatable(inst):
 
         if has_elem:
             return cls(
-                vertices=np.vstack(vertices),
-                elements=np.vstack(elements),
+                vertices=_np.vstack(vertices),
+                elements=_np.vstack(elements),
             )
 
         else:
-            return Vertices(vertices=np.vstack(vertices))
+            return Vertices(vertices=_np.vstack(vertices))
 
-    def __add__(self, to_add):
+    def __add__(self, to_add: Self) -> Self:
         """Concat in form of +.
 
         Parameters
diff --git a/gustaf/volumes.py b/gustaf/volumes.py
index 87e87ec0..fee5d60b 100644
--- a/gustaf/volumes.py
+++ b/gustaf/volumes.py
@@ -1,21 +1,34 @@
 """gustaf/gustaf/volumes.py."""
 
-import numpy as np
+from __future__ import annotations
 
-from gustaf import helpers, settings, show, utils
-from gustaf.faces import Faces
-from gustaf.helpers.options import Option
+from typing import TYPE_CHECKING
 
+import numpy as _np
 
-class VolumesShowOption(helpers.options.ShowOption):
+from gustaf import helpers as _helpers
+from gustaf import settings as _settings
+from gustaf import show as _show
+from gustaf import utils as _utils
+from gustaf.faces import Faces as _Faces
+from gustaf.helpers.options import Option as _Option
+
+if TYPE_CHECKING:
+
+    from gustaf.helpers.data import TrackedArray, Unique2DIntegers
+
+
+class VolumesShowOption(_helpers.options.ShowOption):
     """
     Show options for vertices.
     """
 
-    _valid_options = helpers.options.make_valid_options(
-        *helpers.options.vedo_common_options,
-        Option("vedo", "lw", "Width of edges (lines) in pixel units.", (int,)),
-        Option(
+    _valid_options = _helpers.options.make_valid_options(
+        *_helpers.options.vedo_common_options,
+        _Option(
+            "vedo", "lw", "Width of edges (lines) in pixel units.", (int,)
+        ),
+        _Option(
             "vedo", "lc", "Color of edges (lines).", (int, str, tuple, list)
         ),
     )
@@ -40,14 +53,14 @@ def _initialize_showable(self):
             self.get("data", None) is None
             and not self.get("vertex_ids", False)
             and not self.get("arrow_data", False)
-            and not show.is_ipython
+            and not _show.is_ipython
         ):
             from vtk import VTK_HEXAHEDRON as herr_hexa
             from vtk import VTK_TETRA as frau_tetra
 
             to_vtktype = {"tet": frau_tetra, "hexa": herr_hexa}
             grid_type = to_vtktype[self._helpee.whatami]
-            u_grid = show.vedoUGrid(
+            u_grid = _show.vedoUGrid(
                 [
                     self._helpee.const_vertices,
                     self._helpee.const_volumes,
@@ -72,15 +85,15 @@ def _initialize_showable(self):
         return faces.show_options._initialize_showable()
 
 
-class Volumes(Faces):
+class Volumes(_Faces):
     kind = "volume"
 
-    const_faces = helpers.raise_if.invalid_inherited_attr(
+    const_faces = _helpers.raise_if.invalid_inherited_attr(
         "Faces.const_faces",
         __qualname__,
         property_=True,
     )
-    update_faces = helpers.raise_if.invalid_inherited_attr(
+    update_faces = _helpers.raise_if.invalid_inherited_attr(
         "Faces.update_edges",
         __qualname__,
         property_=False,
@@ -92,14 +105,14 @@ class Volumes(Faces):
     )
 
     __show_option__ = VolumesShowOption
-    __boundary_class__ = Faces
+    __boundary_class__ = _Faces
 
     def __init__(
         self,
-        vertices=None,
-        volumes=None,
-        elements=None,
-    ):
+        vertices: list[list[float]] | _np.ndarray = None,
+        volumes: list[list[int]] | None | _np.ndarray = None,
+        elements: _np.ndarray | None = None,
+    ) -> None:
         """Volumes. It has vertices and volumes. Volumes could be tetrahedrons
         or hexahedrons.
 
@@ -114,8 +127,8 @@ def __init__(
         elif elements is not None:
             self.volumes = elements
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def faces(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def faces(self) -> _np.ndarray:
         """Faces here aren't main property. So this needs to be computed.
 
         Parameters
@@ -129,14 +142,14 @@ def faces(self):
         whatami = self.whatami
         faces = None
         if whatami.startswith("tet"):
-            faces = utils.connec.tet_to_tri(self.volumes)
+            faces = _utils.connec.tet_to_tri(self.volumes)
         elif whatami.startswith("hexa"):
-            faces = utils.connec.hexa_to_quad(self.volumes)
+            faces = _utils.connec.hexa_to_quad(self.volumes)
 
         return faces
 
     @classmethod
-    def whatareyou(cls, volume_obj):
+    def whatareyou(cls, volume_obj: Volumes) -> str:
         """overwrites Faces.whatareyou to tell you is this volume is tet or
         hexa.
 
@@ -164,7 +177,7 @@ def whatareyou(cls, volume_obj):
             )
 
     @property
-    def volumes(self):
+    def volumes(self) -> TrackedArray:
         """Returns volumes.
 
         Parameters
@@ -178,7 +191,9 @@ def volumes(self):
         return self._volumes
 
     @volumes.setter
-    def volumes(self, vols):
+    def volumes(
+        self, vols: list[list[int]] | TrackedArray | _np.ndarray
+    ) -> None:
         """volumes setter. Similar to vertices, this will be a tracked array.
 
         Parameters
@@ -189,13 +204,13 @@ def volumes(self, vols):
         --------
         None
         """
-        self._volumes = helpers.data.make_tracked_array(
+        self._volumes = _helpers.data.make_tracked_array(
             vols,
-            settings.INT_DTYPE,
+            _settings.INT_DTYPE,
             copy=False,
         )
         if vols is not None:
-            utils.arr.is_one_of_shapes(
+            _utils.arr.is_one_of_shapes(
                 vols,
                 ((-1, 4), (-1, 8)),
                 strict=True,
@@ -205,7 +220,7 @@ def volumes(self, vols):
         self._const_volumes.flags.writeable = False
 
     @property
-    def const_volumes(self):
+    def const_volumes(self) -> TrackedArray:
         """Returns non-writeable view of volumes.
 
         Parameters
@@ -218,8 +233,8 @@ def const_volumes(self):
         """
         return self._const_volumes
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def sorted_volumes(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def sorted_volumes(self) -> _np.ndarray:
         """Sort volumes along axis=1.
 
         Parameters
@@ -232,10 +247,10 @@ def sorted_volumes(self):
         """
         volumes = self._get_attr("volumes")
 
-        return np.sort(volumes, axis=1)
+        return _np.sort(volumes, axis=1)
 
-    @helpers.data.ComputedMeshData.depends_on(["elements"])
-    def unique_volumes(self):
+    @_helpers.data.ComputedMeshData.depends_on(["elements"])
+    def unique_volumes(self) -> Unique2DIntegers:
         """Returns a namedtuple of unique volumes info. Similar to
         unique_edges.
 
@@ -248,7 +263,7 @@ def unique_volumes(self):
         unique_info: Unique2DIntegers
           valid attributes are {values, ids, inverse, counts}
         """
-        unique_info = utils.connec.sorted_unique(
+        unique_info = _utils.connec.sorted_unique(
             self.sorted_volumes(),
             sorted_=True,
         )
@@ -263,7 +278,7 @@ def update_volumes(self, *args, **kwargs):
         """Alias to update_elements."""
         self.update_elements(*args, **kwargs)
 
-    def to_faces(self, unique=True):
+    def to_faces(self, unique: bool = True) -> _Faces:
         """Returns Faces obj.
 
         Parameters
@@ -275,7 +290,7 @@ def to_faces(self, unique=True):
         --------
         faces: Faces
         """
-        return Faces(
+        return _Faces(
             self.vertices,
             faces=self.unique_faces().values if unique else self.faces(),
         )