add volume.slab()

docs/changes.md

@@ -13,7 +13,7 @@
 
 ### Other changes
 - add `core.apply_transform_from_actor()`
-
+- add `add volume.slab()`
 
 
 

examples/volumetric/earth_model.py

@@ -7,7 +7,7 @@
 conductor = tet.clone().threshold('cell_scalars', above=0, below=4)
 
 # Crop the initial mesh
-box = vedo.Box(size=[503500, 505000, 6414000, 6417000, -1830, 600])
+box = vedo.Box([503500, 505000, 6414000, 6417000, -1830, 600])
 tet.cut_with_mesh(box, whole_cells=True)
 
 # We need to build a look up table for our color bar

vedo/pyplot.py

@@ -2514,7 +2514,7 @@ def histogram(*args, **kwargs):
             data = args[0].pointdata[0]
         elif isinstance(args[0], vedo.Points):
             pd0 = args[0].pointdata[0]
-            if pd0:
+            if pd0 is not None:
                 data = pd0.ravel()
             else:
                 data = args[0].celldata[0].ravel()

vedo/shapes.py

@@ -3303,30 +3303,23 @@ def __init__(
         Build a box of dimensions `x=length, y=width and z=height`.
         Alternatively dimensions can be defined by setting `size` keyword with a tuple.
 
-        If `size` is a list of 6 numbers, this will be interpreted as the bounding box:
+        If `pos` is a list of 6 numbers, this will be interpreted as the bounding box:
         `[xmin,xmax, ymin,ymax, zmin,zmax]`
 
         Examples:
             - [aspring1.py](https://github.com/marcomusy/vedo/tree/master/examples/simulations/aspring1.py)
 
                 ![](https://vedo.embl.es/images/simulations/50738955-7e891800-11d9-11e9-85cd-02bd4f3f13ea.gif)
         """
-        if len(size) == 6:
-            bounds = size
-            length = bounds[1] - bounds[0]
-            width = bounds[3] - bounds[2]
-            height = bounds[5] - bounds[4]
-            xp = (bounds[1] + bounds[0]) / 2
-            yp = (bounds[3] + bounds[2]) / 2
-            zp = (bounds[5] + bounds[4]) / 2
-            pos = (xp, yp, zp)
+        src = vtk.new("CubeSource")
+        if len(pos) == 6:
+            src.SetBounds(pos)
+            pos = [0,0,0]
         elif len(size) == 3:
             length, width, height = size
-
-        src = vtk.new("CubeSource")
-        src.SetXLength(length)
-        src.SetYLength(width)
-        src.SetZLength(height)
+            src.SetXLength(length)
+            src.SetYLength(width)
+            src.SetZLength(height)
         src.Update()
         pd = src.GetOutput()
 

vedo/version.py

@@ -1 +1 @@
-_version = '2023.5.0+dev5'
+_version = '2023.5.0+dev6'

vedo/visual.py

@@ -48,13 +48,17 @@ def print(self):
     @property
     def LUT(self):
         """Return the lookup table of the object as a numpy object."""
-        _lut = self.mapper.GetLookupTable()
-        values = []
-        for i in range(_lut.GetTable().GetNumberOfTuples()):
-            # print("LUT i =", i, "value =", _lut.GetTableValue(i))
-            values.append(_lut.GetTableValue(i))
-        return np.array(values)
+        try:
+            _lut = self.mapper.GetLookupTable()
 
+            values = []
+            for i in range(_lut.GetTable().GetNumberOfTuples()):
+                # print("LUT i =", i, "value =", _lut.GetTableValue(i))
+                values.append(_lut.GetTableValue(i))
+            return np.array(values)
+        except AttributeError:
+            pass
+
     @LUT.setter
     def LUT(self, arr):
         """

vedo/volume.py

@@ -452,8 +452,118 @@ def slice_plane(self, origin=(0, 0, 0), normal=(1, 1, 1), autocrop=False):
         vslice.Update()
         msh = Mesh(vslice.GetOutput())
         msh.apply_transform(T)
-        msh.pipeline = utils.OperationNode("slice_plane", parents=[self], c="#4cc9f0:#e9c46a")
+        msh.pipeline = utils.OperationNode(
+            "slice_plane", parents=[self], c="#4cc9f0:#e9c46a")
         return msh
+
+    def slab(self, slice_range=(), axis='z', operation="mean"):
+        """
+        Extract a slab from a `Volume` by combining 
+        all of the slices of an image to create a single slice.
+
+        Returns a `Mesh` containing metadata which
+        can be accessed with e.g. `mesh.metadata["slab_range"]`.
+
+        Metadata:
+            slab_range : (list)
+                contains the range of slices extracted
+            slab_axis : (str)
+                contains the axis along which the slab was extracted
+            slab_operation : (str)
+                contains the operation performed on the slab
+            slab_bounding_box : (list)
+                contains the bounding box of the slab
+
+        Arguments:
+            slice_range : (list)
+                range of slices to extract
+            axis : (str)
+                axis along which to extract the slab
+            operation : (str)
+                operation to perform on the slab,
+                allowed values are: "sum", "min", "max", "mean".
+        """
+        if len(slice_range) != 2:
+            vedo.logger.error("in slab(): slice_range is empty or invalid")
+            raise ValueError()
+
+        islab = vtk.new("ImageSlab")
+        islab.SetInputData(self.dataset)
+        islab.SetSliceRange(slice_range)
+
+        if operation in ["+", "add", "sum"]:
+            islab.SetOperationToSum()
+        elif "min" in operation:
+            islab.SetOperationToMin()
+        elif "max" in operation:
+            islab.SetOperationToMax()
+        elif "mean" in operation:
+            islab.SetOperationToMean()
+        else:
+            vedo.logger.error(f"in slab(): unknown operation {operation}")
+            raise ValueError()
+
+        if axis == 'x':
+            islab.SetOrientationToX()
+        elif axis == 'y':
+            islab.SetOrientationToY()
+        elif axis == 'z':
+            islab.SetOrientationToZ()
+        else:
+            vedo.logger.error(f"Error in slab(): unknown axis {axis}")
+            raise RuntimeError()
+
+        islab.Update()
+
+        msh = Mesh(islab.GetOutput()).lighting('off')
+        msh.mapper.SetLookupTable(utils.ctf2lut(self, msh))
+        msh.mapper.SetScalarRange(self.scalar_range())
+
+        msh.metadata["slab_range"] = slice_range
+        msh.metadata["slab_axis"] = axis
+        msh.metadata["slab_operation"] = operation
+
+        # compute bounds of slices
+        dims = self.dimensions()
+        origin = self.origin()
+        spacing = self.spacing()
+        if axis == 'x':
+            msh.metadata["slab_bounding_box"] = [
+                origin[0] + slice_range[0]*spacing[0],
+                origin[0] + slice_range[1]*spacing[0],
+                origin[1],
+                origin[1] + dims[1]*spacing[1],
+                origin[2],
+                origin[2] + dims[2]*spacing[2],
+            ]
+        elif axis == 'y':
+            msh.metadata["slab_bounding_box"] = [
+                origin[0],
+                origin[0] + dims[0]*spacing[0],
+                origin[1] + slice_range[0]*spacing[1],
+                origin[1] + slice_range[1]*spacing[1],
+                origin[2],
+                origin[2] + dims[2]*spacing[2],
+            ]
+        elif axis == 'z':
+            msh.metadata["slab_bounding_box"] = [
+                origin[0],
+                origin[0] + dims[0]*spacing[0],
+                origin[1],
+                origin[1] + dims[1]*spacing[1],
+                origin[2] + slice_range[0]*spacing[2],
+                origin[2] + slice_range[1]*spacing[2],
+            ]
+
+        msh.pipeline = utils.OperationNode(
+            f"slab{slice_range}", 
+            comment=f"axis={axis}, operation={operation}",
+            parents=[self],
+            c="#4cc9f0:#e9c46a",
+        )
+        msh.name = "VolumeSlabMesh"
+        return msh
+
 
     def warp(self, source, target, sigma=1, mode="3d", fit=True):
         """
@@ -483,7 +593,7 @@ def warp(self, source, target, sigma=1, mode="3d", fit=True):
         self.pipeline = utils.OperationNode("warp", parents=[self], c="#4cc9f0")
         return self
 
-    def apply_transform(self, T, fit=True, interpolation="linear"):
+    def apply_transform(self, T, fit=True, interpolation="cubic"):
         """
         Apply a transform to the scalars in the volume.
 
@@ -494,16 +604,15 @@ def apply_transform(self, T, fit=True, interpolation="linear"):
                 fit/adapt the old bounding box to the modified geometry
             interpolation : (str)
                 one of the following: "linear", "nearest", "cubic"
-            
         """
         TI = T.compute_inverse()
         reslice = vtk.new("ImageReslice")
         reslice.SetInputData(self.dataset)
         reslice.SetResliceTransform(TI.T)
         reslice.SetOutputDimensionality(3)
-        if "linear" in interpolation.lower():
+        if "lin" in interpolation.lower():
             reslice.SetInterpolationModeToLinear()
-        elif "nearest" in interpolation.lower():
+        elif "near" in interpolation.lower():
             reslice.SetInterpolationModeToNearestNeighbor()
         elif "cubic" in interpolation.lower():
             reslice.SetInterpolationModeToCubic()

vedo/vtkclasses.py

@@ -631,6 +631,7 @@ def dump_hierarchy_to_file(fname=""):
     "vtkImageLaplacian",
     "vtkImageMedian3D",
     "vtkImageNormalize",
+    "vtkImageSlab",
 ]: location[name] = "vtkImagingGeneral"
 
 for name in ["vtkImageToPoints", "vtkSampleFunction"]: