Extract cubic volume from transformed imageData

Support
1

Hi,
I am new in VTK. I am using Python 3.7 for VTK 9.0.1 under anaconda in Ubuntu Linux.

I do have a 3D numpy array with a 3D scalar volume image and a cube source. I would like to split the image in two parts: the one inside the cube with opacity=1 and the one outside the cube with opacity=0.2. How can I split the image?
I must be able to move and rotate the image. The objective is to select a region of interest. I was able to show and move the image and also show the cube, but I have no idea how to split the image.
After that, I can obtain the transformation matrix and extract the region for another process.

Here is the code for showing and transforming the image and showing the cube. I do not have any clue to break the image:

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
@author: dani
"""

import numpy as np
import vtk
import vtk.util.numpy_support as numpy_support

def matrixExample():
    """
    Returns just a numpy 3D array as an example
    in real world this will become from a complex
    function
    """
    x=np.arange(-25,26)
    y=np.arange(-25,26)
    z=np.arange(-45,46)
    X,Y,Z=np.meshgrid(x,y,z)
    m=np.exp(-np.square(1-np.square(X)/250-np.square(Y)/250-np.square(Z)/700))
    max=m.max()
    min=m.min()
    m=255*(m-min)/(max-min)
    return m.astype(np.ubyte)

def boxCallback(obj, event):
    """
        Rotate the Prop with the same transformation than the widget
    """
    t = vtk.vtkTransform()
    obj.GetTransform(t)
    obj.GetProp3D().SetUserTransform(t)

def main():
    # Start from numpy 3D array and convert to imageData
    data=matrixExample()
    imdata = vtk.vtkImageData()
    depthArray = numpy_support.numpy_to_vtk( \
        data.ravel(), deep=True, \
        array_type=numpy_support.get_vtk_array_type(data.dtype) \
        )
    imdata.SetDimensions([data.shape[2],data.shape[1],data.shape[0]])
    imdata.SetSpacing(1,1,1)
    imdata.SetOrigin([-45,-25,-25])
    imdata.GetPointData().SetScalars(depthArray)
    
    # Colors and Opacity for the scalar range
    colorFunc = vtk.vtkColorTransferFunction()
    colorFunc.AddRGBPoint(0, 1, 0, 0.0) # Green
    colorFunc.AddRGBPoint(255, 0.0, 1, 0.0) # Green 
    opacity = vtk.vtkPiecewiseFunction()
    opacity.AddPoint(0, 0)
    opacity.AddPoint(10, 0)
    opacity.AddPoint(60, 0.5)
    opacity.AddPoint(255, 0.5)
    
    volumeProperty = vtk.vtkVolumeProperty()
    volumeProperty.SetColor(colorFunc)
    volumeProperty.SetScalarOpacity(opacity)
    volumeProperty.SetInterpolationTypeToLinear()
    volumeProperty.SetIndependentComponents(1)
    
    volumeMapper = vtk.vtkOpenGLGPUVolumeRayCastMapper()
    volumeMapper.SetInputData(imdata)
    volumeMapper.SetBlendModeToMaximumIntensity()
    
    volume = vtk.vtkVolume()
    volume.SetMapper(volumeMapper)
    volume.SetProperty(volumeProperty)
    
    # Cube Source
    cube = vtk.vtkCubeSource()
    cube.SetCenter(0,0,0)
    cube.SetBounds(-10,10,-10,10,-10,10)
    cubeMapper = vtk.vtkPolyDataMapper()
    cubeMapper.SetInputConnection(cube.GetOutputPort())
    cubeActor = vtk.vtkActor()
    cubeActor.GetProperty().SetOpacity(0.3)
    cubeActor.SetMapper(cubeMapper)

    # Renderer
    ren = vtk.vtkRenderer()
    ren.AddVolume(volume)
    ren.AddActor(cubeActor)
    ren.SetBackground(0,0,0)
    # Renderer Window
    renWin = vtk.vtkRenderWindow()
    renWin.AddRenderer(ren)
    renWin.SetSize(900, 900)
    # Interactor
    interactor = vtk.vtkRenderWindowInteractor()
    interactor.SetRenderWindow(renWin)
     
    #Add axes widget:
    axes = vtk.vtkAxesActor()
    widget = vtk.vtkOrientationMarkerWidget()
    widget.SetOrientationMarker(axes)
    widget.SetInteractor(interactor)
    widget.SetEnabled(1)
    widget.InteractiveOff()
    
    # BoxWidget to move the image
    boxWidget = vtk.vtkBoxWidget()
    boxWidget.SetInteractor(interactor)
    boxWidget.SetProp3D(volume)
    boxWidget.SetPlaceFactor(1)  # Make the box 1.25x larger than the actor
    boxWidget.PlaceWidget()
    boxWidget.HandlesOff() # to not be able to extend,  just to move and rotate
    boxWidget.On()
    # Connect the event to a function
    boxWidget.AddObserver('InteractionEvent', boxCallback)
    
    interactor.Initialize()
    ren.ResetCamera()
    renWin.Render()
    interactor.Start()

if __name__ == '__main__':
    main()

Thanks

2

The PolyDataToImageData example does some of what you need. The idea behind that example is that vtkPolyDataToImageStencil can be used to convert a polydata into a stencil (a form of mask) that can be applied to an image. There are a handful of image filters in VTK that can use a stencil to apply an operation to just part of an image: vtkImageStencil, vtkImageHistogram, vtkImageConnectivityFilter, vtkImageBlend, vtkImageReslice.

3

Thanks!
Since I am new to vtk, and conversion filters are a bit difficult for me at my stage, I have opted for a simple solution (or at least the one I can understand). I have added 6 clipping planes.