FEAT: Depth 3D edit implemented

Implemented with a spin control from 0 to 1.

Function parameter improvements in mask_3d_edit.py

invesalius/data/styles_3d.py

@@ -680,7 +680,27 @@ def OnInsertPolygonPoint(self, obj, evt):
         self.poly.insert_point((mouse_x, mouse_y))
         self.viewer.UpdateCanvas()
 
-    def __get_model_to_screen_volume(self):
+    def __make_camera_matrix(self):
+        camera = self.viewer.ren.GetActiveCamera()
+
+        # Get camera position (translation)
+        position = camera.GetPosition()
+        print(f"Camera Position: {position}")
+
+        # Get camera focal point (where the camera is looking)
+        focal_point = camera.GetFocalPoint()
+        print(f"Camera Focal Point: {focal_point}")
+
+        # Get the view up vector
+        view_up = camera.GetViewUp()
+        print(f"Camera View Up: {view_up}")
+
+        m = create_model_view_matrix(position, focal_point, view_up)
+        print(m)
+
+        return m
+
+    def __get_cam_parameters(self):
         w, h = self.viewer.GetSize()
         self.viewer.ren.Render()
         cam = self.viewer.ren.GetActiveCamera()
@@ -693,7 +713,14 @@ def __get_model_to_screen_volume(self):
         MV = cam.GetViewTransformMatrix()
         MV = vtkarray_to_numpy(MV)
 
-        return M
+        params = {
+            "model_to_screen": M,
+            "model_view": MV,
+            "resolution": (w, h),
+            "clipping_range": (near, far),
+        }
+
+        return params
 
     def OnInsertPolygon(self, obj, evt):
         self.poly.complete_polygon()
@@ -702,14 +729,8 @@ def OnInsertPolygon(self, obj, evt):
         self.viewer.UpdateCanvas()
         self.viewer.ren.Render()
 
-        Publisher.sendMessage(
-            "M3E add polygon",
-            points=self.poly.polygon.points,
-            screen=self.viewer.ren.GetRenderWindow().GetSize(),
-        )
-        Publisher.sendMessage(
-            "M3E set model_to_screen", model_to_screen=self.__get_model_to_screen_volume()
-        )
+        Publisher.sendMessage("M3E add polygon", points=self.poly.polygon.points)
+        Publisher.sendMessage("M3E set camera", params=self.__get_cam_parameters())
 
 
 class Styles:

invesalius/gui/task_slice.py

@@ -828,7 +828,7 @@ def __bind_events_wx(self):
         self.btn_do_3d_edit.Bind(wx.EVT_BUTTON, self.OnDoMaskEdit3D)
         self.combo_mask_edit_3d_op.Bind(wx.EVT_COMBOBOX, self.OnComboMaskEdit3DMode)
         self.cbox_mask_edit_3d_depth.Bind(wx.EVT_CHECKBOX, self.OnUseDepthMaskEdit3D)
-        self.spin_mask_edit_3d_depth.Bind(wx.EVT_SPIN, self.OnSpinDepthMaskEdit3D)
+        self.spin_mask_edit_3d_depth.Bind(wx.EVT_SPINCTRLDOUBLE, self.OnSpinDepthMaskEdit3D)
         self.btn_clear_3d_poly.Bind(wx.EVT_BUTTON, self.OnClearPolyMaskEdit3D)
 
     def __bind_events(self):
@@ -950,8 +950,9 @@ def OnComboMaskEdit3DMode(self, evt):
         Publisher.sendMessage("M3E set edit mode", mode=op_id)
 
     def OnUseDepthMaskEdit3D(self, evt):
-        cb_val = self.cbox_mask_edit_3d.GetValue()
+        cb_val = self.cbox_mask_edit_3d_depth.GetValue()
         Publisher.sendMessage("M3E use depth", use=cb_val)
+        self.OnSpinDepthMaskEdit3D(evt)  # To pass the set value
 
     def OnSpinDepthMaskEdit3D(self, evt):
         spin_val = self.spin_mask_edit_3d_depth.GetValue()

invesalius/segmentation/mask_3d_edit.py

@@ -23,41 +23,59 @@
 import invesalius.constants as const
 import invesalius.data.slice_ as slc
 from invesalius.pubsub import pub as Publisher
-from invesalius_cy.mask_cut import mask_cut
+from invesalius_cy.mask_cut import mask_cut, mask_cut_with_depth
+
+
+class Mask3DEditException(Exception):
+    pass
 
 
 class Mask3DEditor:
     def __init__(self):
         self.__bind_events()
         self.polygons_to_operate = []
-        self.model_to_screen = None
         self.spacing = None
         self.edit_mode = const.MASK_3D_EDIT_INCLUDE
+        self.use_depth = False
+        self.depth_val = None
+        self.model_to_screen = None
+        self.model_view = None
+        self.resolution = None  # (w, h)
+        self.clipping_range = None  # (near, far)
 
     def __bind_events(self):
         Publisher.subscribe(self.AddPolygon, "M3E add polygon")
         Publisher.subscribe(self.DoMaskEdit, "M3E apply edit")
         Publisher.subscribe(self.ClearPolygons, "M3E clear polygons")
-        Publisher.subscribe(self.SetMVP, "M3E set model_to_screen")
+        Publisher.subscribe(self.SetCamParameters, "M3E set camera")
         Publisher.subscribe(self.SetEditMode, "M3E set edit mode")
+        Publisher.subscribe(self.SetUseDepthForEdit, "M3E use depth")
+        Publisher.subscribe(self.SetDepthValue, "M3E depth value")
 
-    def AddPolygon(self, points, screen):
+    def AddPolygon(self, points):
         """
         Adds polygon to be used in the edit
         """
-        self.polygons_to_operate.append((points, screen))
+        self.polygons_to_operate.append(points)
 
     def ClearPolygons(self):
         """
         Discards all added polygons
         """
         self.polygons_to_operate = []
 
-    def SetMVP(self, model_to_screen):
+    def SetCamParameters(self, params):
         """
         Sets the model-view-projection matrix used by the viewer
         """
-        self.model_to_screen = model_to_screen
+        if "model_to_screen" in params:
+            self.model_to_screen = params["model_to_screen"]
+        if "model_view" in params:
+            self.model_view = params["model_view"]
+        if "resolution" in params:
+            self.resolution = params["resolution"]
+        if "clipping_range" in params:
+            self.clipping_range = params["clipping_range"]
 
     def SetEditMode(self, mode):
         """
@@ -66,17 +84,25 @@ def SetEditMode(self, mode):
         """
         self.edit_mode = mode
 
+    def SetUseDepthForEdit(self, use):
+        """
+        Sets whether to perform a mask cut using depth or through all
+        """
+        self.use_depth = use
+
+    def SetDepthValue(self, value):
+        self.depth_val = value
+
     def __create_filter(self):
         """
         Based on the polygons and screen resolution,
         create a filter for the edit.
         """
-        w, h = tuple(self.polygons_to_operate[0][1])
-
+        w, h = self.resolution
         _filter = np.zeros((h, w), dtype="uint8")
 
         # Include all selected polygons to create the cut
-        for poly_points, _ in self.polygons_to_operate:
+        for poly_points in self.polygons_to_operate:
             print(f"Poly Points: {poly_points}")
             poly = np.array(poly_points)
             rr, cc = polygon(poly[:, 1], poly[:, 0], _filter.shape)
@@ -88,7 +114,7 @@ def __create_filter(self):
         return _filter
 
     # Unoptimized implementation
-    def __cut_mask_left_right(self, mask_data, spacing, _filter):
+    def __cut_mask(self, mask_data, spacing, _filter):
         """
         Cuts an Invesalius Mask with a filter (pixel-wise filter)
         mask_data: matrix data without metadata ([1:])
@@ -141,17 +167,29 @@ def DoMaskEdit(self):
         _cur_mask = s.current_mask
 
         if _cur_mask is None:
-            raise Exception("Attempted Slicing an empty mask")
+            raise Mask3DEditException("Attempted editing a non-existent mask")
 
         _filter = self.__create_filter()
 
         # Unoptimized implementation
-        # self.__cut_mask_left_right(_cur_mask.matrix[1:, 1:, 1:], s.spacing, _filter)
+        # self.__cut_mask(_cur_mask.matrix[1:, 1:, 1:], s.spacing, _filter)
 
         # Optimized implementation
         _mat = _cur_mask.matrix[1:, 1:, 1:].copy()
         sx, sy, sz = s.spacing
-        mask_cut(_mat, sx, sy, sz, _filter, self.model_to_screen, _mat)
+
+        print(s.spacing)
+        print(_mat.shape)
+
+        if self.use_depth:
+            near, far = self.clipping_range
+            depth = near + (far - near) * self.depth_val
+            mask_cut_with_depth(
+                _mat, sx, sy, sz, depth, _filter, self.model_to_screen, self.model_view, _mat
+            )
+        else:
+            mask_cut(_mat, sx, sy, sz, _filter, self.model_to_screen, _mat)
+
         _cur_mask.matrix[1:, 1:, 1:] = _mat
         _cur_mask.modified(all_volume=True)
 

invesalius_cy/mask_cut.pyx

@@ -59,3 +59,73 @@ def mask_cut(np.ndarray[image_t, ndim=3] mask_data,
                     if 0 <= px <= w and 0 <= py <= h:
                         if filter[<int>(py), <int>(px)]:
                             out[z, y, x] = 0
+
+@cython.boundscheck(False) # turn of bounds-checking for entire function
+@cython.cdivision(True)
+def mask_cut_with_depth(np.ndarray[image_t, ndim=3] image,
+                          float sx, float sy, float sz,
+                          float max_depth,
+                          np.ndarray[mask_t, ndim=2] mask,
+                          np.ndarray[DTYPEF64_t, ndim=2] M,
+                          np.ndarray[DTYPEF64_t, ndim=2] MV,
+                          np.ndarray[image_t, ndim=3] out):
+
+    cdef int dz = image.shape[0]
+    cdef int dy = image.shape[1]
+    cdef int dx = image.shape[2]
+
+    cdef int x
+    cdef int y
+    cdef int z
+
+    cdef int h = mask.shape[0]
+    cdef int w = mask.shape[1]
+
+    cdef DTYPEF64_t px
+    cdef DTYPEF64_t py
+
+    cdef DTYPEF64_t p0, p1, p2, p3
+    cdef DTYPEF64_t _q0, _q1, _q2, _q3
+    cdef DTYPEF64_t q0, q1, q2
+
+    cdef DTYPEF64_t _c0, _c1, _c2, _c3
+    cdef DTYPEF64_t c0, c1, c2
+
+    cdef DTYPEF64_t dist
+
+    for z in prange(dz, nogil=True):
+        for y in range(dy):
+            for x in range(dx):
+                p0 = <float>(x*sx)
+                p1 = <float>(y*sy)
+                p2 = <float>(z*sz)
+                p3 = 1.0
+
+                _q0 = p0 * M[0, 0] + p1 * M[0, 1] + p2 * M[0, 2] + p3 * M[0, 3]
+                _q1 = p0 * M[1, 0] + p1 * M[1, 1] + p2 * M[1, 2] + p3 * M[1, 3]
+                _q2 = p0 * M[2, 0] + p1 * M[2, 1] + p2 * M[2, 2] + p3 * M[2, 3]
+                _q3 = p0 * M[3, 0] + p1 * M[3, 1] + p2 * M[3, 2] + p3 * M[3, 3]
+
+                if _q3 > 0:
+                    q0 = _q0/_q3
+                    q1 = _q1/_q3
+                    q2 = _q2/_q3
+
+                    _c0 = p0 * MV[0, 0] + p1 * MV[0, 1] + p2 * MV[0, 2] + p3 * MV[0, 3]
+                    _c1 = p0 * MV[1, 0] + p1 * MV[1, 1] + p2 * MV[1, 2] + p3 * MV[1, 3]
+                    _c2 = p0 * MV[2, 0] + p1 * MV[2, 1] + p2 * MV[2, 2] + p3 * MV[2, 3]
+                    _c3 = p0 * MV[3, 0] + p1 * MV[3, 1] + p2 * MV[3, 2] + p3 * MV[3, 3]
+
+                    c0 = _c0/_c3
+                    c1 = _c1/_c3
+                    c2 = _c2/_c3
+
+                    dist = sqrt(c0*c0 + c1*c1 + c2*c2)
+
+                    if dist <= max_depth:
+                        px = (q0/2.0 + 0.5) * (w - 1)
+                        py = (q1/2.0 + 0.5) * (h - 1)
+
+                        if 0 <= px <= w and 0 <= py <= h:
+                            if mask[<int>round(py), <int>round(px)]:
+                                out[z, y, x] = 0