diff --git a/docs/processing.md b/docs/processing.md
index 68535f9cdf20189d5944790e07a82a14677bd368..759c981c8bb76130ac8a850fa72dfa8ece9a535b 100644
--- a/docs/processing.md
+++ b/docs/processing.md
@@ -30,3 +30,10 @@ Here, we provide functionalities designed specifically for 3D image analysis and
             - watershed
             - fade_mask
             - overlay_rgb_images
+
+::: qim3d.processing.features
+    options:
+        members:
+            - area
+            - volume
+            - sphericity
\ No newline at end of file
diff --git a/mesh_from_segmentation.ipynb b/mesh_from_segmentation.ipynb
deleted file mode 100644
index 26a029ffdd3674cfa9e7bfbcdb6b6b01dfa89891..0000000000000000000000000000000000000000
Binary files a/mesh_from_segmentation.ipynb and /dev/null differ
diff --git a/qim3d/processing/__init__.py b/qim3d/processing/__init__.py
index 0b767fb1fe413699d3a56878ca216481fe1daae1..8f6edf04173e5781339f0ca33ff304a29851f0f9 100644
--- a/qim3d/processing/__init__.py
+++ b/qim3d/processing/__init__.py
@@ -6,3 +6,4 @@ from .operations import *
 from .cc import get_3d_cc
 from .layers2d import segment_layers, get_lines
 from .mesh import create_mesh
+from .features import volume, area, sphericity
diff --git a/qim3d/processing/features.py b/qim3d/processing/features.py
new file mode 100644
index 0000000000000000000000000000000000000000..9dd087028375eab5dfb949525cc79f5d20d924ba
--- /dev/null
+++ b/qim3d/processing/features.py
@@ -0,0 +1,151 @@
+import numpy as np
+import qim3d.processing
+from qim3d.utils.logger import log
+import trimesh
+import qim3d
+
+
+def volume(obj, **mesh_kwargs) -> float:
+    """
+    Compute the volume of a 3D volume or mesh.
+
+    Args:
+        obj: Either a np.ndarray volume or a mesh object of type trimesh.Trimesh.
+        **mesh_kwargs: Additional arguments for mesh creation if the input is a volume.
+
+    Returns:
+        volume: The volume of the object.
+
+    Example:
+        Compute volume from a mesh:
+        ```python
+        import qim3d
+
+        # Load a mesh from a file
+        mesh = qim3d.io.load_mesh('path/to/mesh.obj')
+
+        # Compute the volume of the mesh
+        volume = qim3d.processing.volume(mesh)
+        print('Volume:', volume)
+        ```
+
+        Compute volume from a np.ndarray:
+        ```python
+        import qim3d
+
+        # Generate a 3D blob
+        synthetic_blob = qim3d.generate.blob(noise_scale = 0.015)
+
+        # Compute the volume of the blob
+        volume = qim3d.processing.volume(synthetic_blob, level=0.5)
+        print('Volume:', volume)
+        ```
+
+    """
+    if isinstance(obj, np.ndarray):
+        log.info("Converting volume to mesh.")
+        obj = qim3d.processing.create_mesh(obj, **mesh_kwargs)
+
+    return obj.volume
+
+
+def area(obj, **mesh_kwargs) -> float:
+    """
+    Compute the surface area of a 3D volume or mesh.
+
+    Args:
+        obj: Either a np.ndarray volume or a mesh object of type trimesh.Trimesh.
+        **mesh_kwargs: Additional arguments for mesh creation if the input is a volume.
+
+    Returns:
+        area: The surface area of the object.
+
+    Example:
+        Compute area from a mesh:
+        ```python
+        import qim3d
+
+        # Load a mesh from a file
+        mesh = qim3d.io.load_mesh('path/to/mesh.obj')
+
+        # Compute the surface area of the mesh
+        area = qim3d.processing.area(mesh)
+        print(f"Area: {area}")
+        ```
+
+        Compute area from a np.ndarray:
+        ```python
+        import qim3d
+
+        # Generate a 3D blob
+        synthetic_blob = qim3d.generate.blob(noise_scale = 0.015)
+
+        # Compute the surface area of the blob
+        volume = qim3d.processing.area(synthetic_blob, level=0.5)
+        print('Area:', volume)
+        ```
+    """
+    if isinstance(obj, np.ndarray):
+        log.info("Converting volume to mesh.")
+        obj = qim3d.processing.create_mesh(obj, **mesh_kwargs)
+
+    return obj.area
+
+
+def sphericity(obj, **mesh_kwargs) -> float:
+    """
+    Compute the sphericity of a 3D volume or mesh.
+
+    Sphericity is a measure of how spherical an object is. It is defined as the ratio
+    of the surface area of a sphere with the same volume as the object to the object's
+    actual surface area.
+
+    Args:
+        obj: Either a np.ndarray volume or a mesh object of type trimesh.Trimesh.
+        **mesh_kwargs: Additional arguments for mesh creation if the input is a volume.
+
+    Returns:
+        sphericity: A float value representing the sphericity of the object.
+
+    Example:
+        Compute sphericity from a mesh:
+        ```python
+        import qim3d
+
+        # Load a mesh from a file
+        mesh = qim3d.io.load_mesh('path/to/mesh.obj')
+
+        # Compute the sphericity of the mesh
+        sphericity = qim3d.processing.sphericity(mesh)
+        ```
+
+        Compute sphericity from a np.ndarray:
+        ```python
+        import qim3d
+
+        # Generate a 3D blob
+        synthetic_blob = qim3d.generate.blob(noise_scale = 0.015)
+
+        # Compute the sphericity of the blob
+        sphericity = qim3d.processing.sphericity(synthetic_blob, level=0.5)
+        ```
+
+    !!! info "Limitations due to pixelation"
+        Sphericity is particularly sensitive to the resolution of the mesh, as it directly impacts the accuracy of surface area and volume calculations.
+        Since the mesh is generated from voxel-based 3D volume data, the discrete nature of the voxels leads to pixelation effects that reduce the precision of sphericity measurements.
+        Higher resolution meshes may mitigate these errors but often at the cost of increased computational demands.
+    """
+    if isinstance(obj, np.ndarray):
+        log.info("Converting volume to mesh.")
+        obj = qim3d.processing.create_mesh(obj, **mesh_kwargs)
+
+    volume = qim3d.processing.volume(obj)
+    area = qim3d.processing.area(obj)
+
+    if area == 0:
+        log.warning("Surface area is zero, sphericity is undefined.")
+        return np.nan
+
+    sphericity = (np.pi ** (1 / 3) * (6 * volume) ** (2 / 3)) / area
+    log.info(f"Sphericity: {sphericity}")
+    return sphericity
diff --git a/qim3d/processing/mesh.py b/qim3d/processing/mesh.py
index 88ff3de80351e3a345983b55c6e920b713967d3d..b4ff678398c765bebe496f7ee6d5d1015b7d7b1d 100644
--- a/qim3d/processing/mesh.py
+++ b/qim3d/processing/mesh.py
@@ -9,6 +9,7 @@ def create_mesh(
     volume: np.ndarray,
     level: float = None,
     step_size=1,
+    allow_degenerate=False,
     padding: Tuple[int, int, int] = (2, 2, 2),
     **kwargs: Any,
 ) -> trimesh.Trimesh:
@@ -18,6 +19,9 @@ def create_mesh(
     Args:
         volume (np.ndarray): The 3D numpy array representing the volume.
         level (float, optional): The threshold value for Marching Cubes. If None, Otsu's method is used.
+        step_size (int, optional): The step size for the Marching Cubes algorithm.
+        allow_degenerate (bool, optional): Whether to allow degenerate (i.e. zero-area) triangles in the end-result.
+        If False, degenerate triangles are removed, at the cost of making the algorithm slower. Default False.
         padding (tuple of int, optional): Padding to add around the volume.
         **kwargs: Additional keyword arguments to pass to `skimage.measure.marching_cubes`.
 
@@ -39,7 +43,7 @@ def create_mesh(
         mesh = qim3d.processing.create_mesh(vol, step_size=3)
         qim3d.viz.mesh(mesh.vertices, mesh.faces)
         ```
-
+        <iframe src="https://platform.qim.dk/k3d/mesh_visualization.html" width="100%" height="500" frameborder="0"></iframe>
     """
     if volume.ndim != 3:
         raise ValueError("The input volume must be a 3D numpy array.")
@@ -63,10 +67,13 @@ def create_mesh(
 
     # Call skimage.measure.marching_cubes with user-provided kwargs
     verts, faces, normals, values = measure.marching_cubes(
-        volume, level=level, step_size=step_size, **kwargs
+        volume, level=level, step_size=step_size, allow_degenerate=allow_degenerate, **kwargs
     )
 
     # Create the Trimesh object
     mesh = trimesh.Trimesh(vertices=verts, faces=faces)
 
+    # Fix face orientation to ensure normals point outwards
+    trimesh.repair.fix_inversion(mesh, multibody=True) 
+
     return mesh
diff --git a/qim3d/viz/k3d.py b/qim3d/viz/k3d.py
index e3f9b5119eab55a0d24b1ab04f233544cc7941c2..40ed5f9ac4b1439b4861706b997a0c5535272f41 100644
--- a/qim3d/viz/k3d.py
+++ b/qim3d/viz/k3d.py
@@ -191,6 +191,7 @@ def mesh(
         mesh = qim3d.processing.create_mesh(vol, step_size=3)
         qim3d.viz.mesh(mesh.vertices, mesh.faces)
         ```
+        <iframe src="https://platform.qim.dk/k3d/mesh_visualization.html" width="100%" height="500" frameborder="0"></iframe>
     """
     import k3d