Commit af3ad6f5f9c784393dbf58dc76f264462fc26146

Authored by Thiago Franco de Moraes
2 parents 5be3e861 b59df909
Exists in master

Merge remote-tracking branch 'upstream/master'

invesalius/data/bases.py
1 1 import numpy as np
2 2 import invesalius.data.coordinates as dco
3 3 import invesalius.data.transformations as tr
4   -
  4 +import invesalius.data.coregistration as dcr
5 5  
6 6 def angle_calculation(ap_axis, coil_axis):
7 7 """
... ... @@ -100,78 +100,31 @@ def base_creation(fiducials):
100 100 return m, q
101 101  
102 102  
103   -# def calculate_fre(fiducials, minv, n, q, o):
104   -# """
105   -# Calculate the Fiducial Registration Error for neuronavigation.
106   -#
107   -# :param fiducials: array of 6 rows (image and tracker fiducials) and 3 columns (x, y, z) with coordinates
108   -# :param minv: inverse matrix given by base creation
109   -# :param n: base change matrix given by base creation
110   -# :param q: origin of first base
111   -# :param o: origin of second base
112   -# :return: float number of fiducial registration error
113   -# """
114   -#
115   -# img = np.zeros([3, 3])
116   -# dist = np.zeros([3, 1])
117   -#
118   -# q1 = np.mat(q).reshape(3, 1)
119   -# q2 = np.mat(o).reshape(3, 1)
120   -#
121   -# p1 = np.mat(fiducials[3, :]).reshape(3, 1)
122   -# p2 = np.mat(fiducials[4, :]).reshape(3, 1)
123   -# p3 = np.mat(fiducials[5, :]).reshape(3, 1)
124   -#
125   -# img[0, :] = np.asarray((q1 + (minv * n) * (p1 - q2)).reshape(1, 3))
126   -# img[1, :] = np.asarray((q1 + (minv * n) * (p2 - q2)).reshape(1, 3))
127   -# img[2, :] = np.asarray((q1 + (minv * n) * (p3 - q2)).reshape(1, 3))
128   -#
129   -# dist[0] = np.sqrt(np.sum(np.power((img[0, :] - fiducials[0, :]), 2)))
130   -# dist[1] = np.sqrt(np.sum(np.power((img[1, :] - fiducials[1, :]), 2)))
131   -# dist[2] = np.sqrt(np.sum(np.power((img[2, :] - fiducials[2, :]), 2)))
132   -#
133   -# return float(np.sqrt(np.sum(dist ** 2) / 3))
134   -
135   -
136   -def calculate_fre_m(fiducials):
  103 +def calculate_fre(fiducials_raw, fiducials, ref_mode_id, m_change, m_icp=None):
137 104 """
138 105 Calculate the Fiducial Registration Error for neuronavigation.
139 106  
  107 + :param fiducials_raw: array of 6 rows (tracker probe and reference) and 3 columns (x, y, z) with coordinates
  108 + :type fiducials_raw: numpy.ndarray
140 109 :param fiducials: array of 6 rows (image and tracker fiducials) and 3 columns (x, y, z) with coordinates
141   - :param minv: inverse matrix given by base creation
142   - :param n: base change matrix given by base creation
143   - :param q: origin of first base
144   - :param o: origin of second base
  110 + :type fiducials: numpy.ndarray
  111 + :param ref_mode_id: Reference mode ID
  112 + :type ref_mode_id: int
  113 + :param m_change: 3x3 array representing change of basis from head in tracking system to vtk head system
  114 + :type m_change: numpy.ndarray
  115 + :param m_icp: list with icp flag and 3x3 affine array
  116 + :type m_icp: list[int, numpy.ndarray]
145 117 :return: float number of fiducial registration error
146 118 """
  119 + if m_icp is not None:
  120 + icp = [True, m_icp]
  121 + else:
  122 + icp = [False, None]
147 123  
148   - m, q1, minv = base_creation_old(fiducials[:3, :])
149   - n, q2, ninv = base_creation_old(fiducials[3:, :])
150   -
151   - # TODO: replace the old by the new base creation
152   - # the values differ greatly if FRE is computed using the old or new base_creation
153   - # check the reason for the difference, because they should be the same
154   - # m, q1 = base_creation(fiducials[:3, :])
155   - # n, q2 = base_creation(fiducials[3:, :])
156   - # minv = np.linalg.inv(m)
157   -
158   - img = np.zeros([3, 3])
159 124 dist = np.zeros([3, 1])
160   -
161   - q1 = q1.reshape(3, 1)
162   - q2 = q2.reshape(3, 1)
163   -
164   - p1 = fiducials[3, :].reshape(3, 1)
165   - p2 = fiducials[4, :].reshape(3, 1)
166   - p3 = fiducials[5, :].reshape(3, 1)
167   -
168   - img[0, :] = (q1 + (minv @ n) * (p1 - q2)).reshape(1, 3)
169   - img[1, :] = (q1 + (minv @ n) * (p2 - q2)).reshape(1, 3)
170   - img[2, :] = (q1 + (minv @ n) * (p3 - q2)).reshape(1, 3)
171   -
172   - dist[0] = np.sqrt(np.sum(np.power((img[0, :] - fiducials[0, :]), 2)))
173   - dist[1] = np.sqrt(np.sum(np.power((img[1, :] - fiducials[1, :]), 2)))
174   - dist[2] = np.sqrt(np.sum(np.power((img[2, :] - fiducials[2, :]), 2)))
  125 + for i in range(0, 6, 2):
  126 + p_m, _ = dcr.corregistrate_dynamic((m_change, 0), fiducials_raw[i:i+2], ref_mode_id, icp)
  127 + dist[int(i/2)] = np.sqrt(np.sum(np.power((p_m[:3] - fiducials[int(i/2), :]), 2)))
175 128  
176 129 return float(np.sqrt(np.sum(dist ** 2) / 3))
177 130  
... ... @@ -203,6 +156,13 @@ def calculate_fre_m(fiducials):
203 156 #
204 157 # return point_rot
205 158  
  159 +def transform_icp(m_img, m_icp):
  160 + coord_img = [m_img[0, -1], -m_img[1, -1], m_img[2, -1], 1]
  161 + m_img[0, -1], m_img[1, -1], m_img[2, -1], _ = m_icp @ coord_img
  162 + m_img[0, -1], m_img[1, -1], m_img[2, -1] = m_img[0, -1], -m_img[1, -1], m_img[2, -1]
  163 +
  164 + return m_img
  165 +
206 166  
207 167 def object_registration(fiducials, orients, coord_raw, m_change):
208 168 """
... ...
invesalius/data/coordinates.py
... ... @@ -259,7 +259,7 @@ def DebugCoord(trk_init, trck_id, ref_mode):
259 259 coord4 = np.array([uniform(*dx), uniform(*dx), uniform(*dx),
260 260 uniform(*dt), uniform(*dt), uniform(*dt)])
261 261  
262   - sleep(0.05)
  262 + sleep(0.15)
263 263  
264 264 # coord1 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
265 265 # uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
... ...
invesalius/data/coregistration.py
... ... @@ -24,6 +24,7 @@ from time import sleep
24 24  
25 25 import invesalius.data.coordinates as dco
26 26 import invesalius.data.transformations as tr
  27 +import invesalius.data.bases as bases
27 28  
28 29  
29 30 # TODO: Replace the use of degrees by radians in every part of the navigation pipeline
... ... @@ -101,7 +102,7 @@ def tracker_to_image(m_change, m_probe_ref, r_obj_img, m_obj_raw, s0_dyn):
101 102 return m_img
102 103  
103 104  
104   -def corregistrate_object_dynamic(inp, coord_raw, ref_mode_id):
  105 +def corregistrate_object_dynamic(inp, coord_raw, ref_mode_id, icp):
105 106  
106 107 m_change, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img = inp
107 108  
... ... @@ -116,6 +117,8 @@ def corregistrate_object_dynamic(inp, coord_raw, ref_mode_id):
116 117 m_probe_ref[2, -1] = -m_probe_ref[2, -1]
117 118 # corregistrate from tracker to image space
118 119 m_img = tracker_to_image(m_change, m_probe_ref, r_obj_img, m_obj_raw, s0_dyn)
  120 + if icp[0]:
  121 + m_img = bases.transform_icp(m_img, icp[1])
119 122 # compute rotation angles
120 123 _, _, angles, _, _ = tr.decompose_matrix(m_img)
121 124 # create output coordiante list
... ... @@ -124,6 +127,9 @@ def corregistrate_object_dynamic(inp, coord_raw, ref_mode_id):
124 127  
125 128 return coord, m_img
126 129  
  130 +def UpdateICP(self, m_icp, flag):
  131 + self.m_icp = m_icp
  132 + self.icp = flag
127 133  
128 134 def compute_marker_transformation(coord_raw, obj_ref_mode):
129 135 psi, theta, phi = np.radians(coord_raw[obj_ref_mode, 3:])
... ... @@ -133,7 +139,7 @@ def compute_marker_transformation(coord_raw, obj_ref_mode):
133 139 return m_probe
134 140  
135 141  
136   -def corregistrate_dynamic(inp, coord_raw, ref_mode_id):
  142 +def corregistrate_dynamic(inp, coord_raw, ref_mode_id, icp):
137 143  
138 144 m_change, obj_ref_mode = inp
139 145  
... ... @@ -150,6 +156,10 @@ def corregistrate_dynamic(inp, coord_raw, ref_mode_id):
150 156 m_probe_ref[2, -1] = -m_probe_ref[2, -1]
151 157 # corregistrate from tracker to image space
152 158 m_img = m_change @ m_probe_ref
  159 +
  160 + if icp[0]:
  161 + m_img = bases.transform_icp(m_img, icp[1])
  162 +
153 163 # compute rotation angles
154 164 _, _, angles, _, _ = tr.decompose_matrix(m_img)
155 165 # create output coordiante list
... ... @@ -160,16 +170,19 @@ def corregistrate_dynamic(inp, coord_raw, ref_mode_id):
160 170  
161 171  
162 172 class CoordinateCorregistrate(threading.Thread):
163   - def __init__(self, ref_mode_id, trck_info, coreg_data, coord_queue, view_tracts, coord_tracts_queue, event, sle):
  173 + def __init__(self, ref_mode_id, trck_info, coreg_data, view_tracts, queues, event, sle):
164 174 threading.Thread.__init__(self, name='CoordCoregObject')
165 175 self.ref_mode_id = ref_mode_id
166 176 self.trck_info = trck_info
167 177 self.coreg_data = coreg_data
168   - self.coord_queue = coord_queue
  178 + self.coord_queue = queues[0]
169 179 self.view_tracts = view_tracts
170   - self.coord_tracts_queue = coord_tracts_queue
  180 + self.coord_tracts_queue = queues[1]
171 181 self.event = event
172 182 self.sle = sle
  183 + self.icp_queue = queues[2]
  184 + self.icp = None
  185 + self.m_icp = None
173 186  
174 187 def run(self):
175 188 trck_info = self.trck_info
... ... @@ -180,19 +193,20 @@ class CoordinateCorregistrate(threading.Thread):
180 193 # print('CoordCoreg: event {}'.format(self.event.is_set()))
181 194 while not self.event.is_set():
182 195 try:
  196 + if self.icp_queue.empty():
  197 + None
  198 + else:
  199 + self.icp, self.m_icp = self.icp_queue.get_nowait()
183 200 # print(f"Set the coordinate")
184 201 coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
185   - coord, m_img = corregistrate_object_dynamic(coreg_data, coord_raw, self.ref_mode_id)
186   - # m_img = np.array([[0.38, -0.8, -0.45, 40.17],
187   - # [0.82, 0.52, -0.24, 152.28],
188   - # [0.43, -0.28, 0.86, 235.78],
189   - # [0., 0., 0., 1.]])
190   - # angles = [-0.318, -0.441, 1.134]
191   - # coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], \
192   - # np.degrees(angles[0]), np.degrees(angles[1]), np.degrees(angles[2])
  202 + coord, m_img = corregistrate_object_dynamic(coreg_data, coord_raw, self.ref_mode_id, [self.icp, self.m_icp])
193 203 m_img_flip = m_img.copy()
194 204 m_img_flip[1, -1] = -m_img_flip[1, -1]
195 205 # self.pipeline.set_message(m_img_flip)
  206 +
  207 + if self.icp:
  208 + m_img = bases.transform_icp(m_img, self.m_icp)
  209 +
196 210 self.coord_queue.put_nowait([coord, m_img, view_obj])
197 211 # print('CoordCoreg: put {}'.format(count))
198 212 # count += 1
... ... @@ -200,6 +214,9 @@ class CoordinateCorregistrate(threading.Thread):
200 214 if self.view_tracts:
201 215 self.coord_tracts_queue.put_nowait(m_img_flip)
202 216  
  217 + if not self.icp_queue.empty():
  218 + self.icp_queue.task_done()
  219 +
203 220 # The sleep has to be in both threads
204 221 sleep(self.sle)
205 222 except queue.Full:
... ... @@ -207,16 +224,19 @@ class CoordinateCorregistrate(threading.Thread):
207 224  
208 225  
209 226 class CoordinateCorregistrateNoObject(threading.Thread):
210   - def __init__(self, ref_mode_id, trck_info, coreg_data, coord_queue, view_tracts, coord_tracts_queue, event, sle):
  227 + def __init__(self, ref_mode_id, trck_info, coreg_data, view_tracts, queues, event, sle):
211 228 threading.Thread.__init__(self, name='CoordCoregNoObject')
212 229 self.ref_mode_id = ref_mode_id
213 230 self.trck_info = trck_info
214 231 self.coreg_data = coreg_data
215   - self.coord_queue = coord_queue
  232 + self.coord_queue = queues[0]
216 233 self.view_tracts = view_tracts
217   - self.coord_tracts_queue = coord_tracts_queue
  234 + self.coord_tracts_queue = queues[1]
218 235 self.event = event
219 236 self.sle = sle
  237 + self.icp_queue = queues[2]
  238 + self.icp = None
  239 + self.m_icp = None
220 240  
221 241 def run(self):
222 242 trck_info = self.trck_info
... ... @@ -227,17 +247,28 @@ class CoordinateCorregistrateNoObject(threading.Thread):
227 247 # print('CoordCoreg: event {}'.format(self.event.is_set()))
228 248 while not self.event.is_set():
229 249 try:
  250 + if self.icp_queue.empty():
  251 + None
  252 + else:
  253 + self.icp, self.m_icp = self.icp_queue.get_nowait()
230 254 # print(f"Set the coordinate")
  255 + #print(self.icp, self.m_icp)
231 256 coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
232   - coord, m_img = corregistrate_dynamic(coreg_data, coord_raw, self.ref_mode_id)
  257 + coord, m_img = corregistrate_dynamic(coreg_data, coord_raw, self.ref_mode_id, [self.icp, self.m_icp])
233 258 # print("Coord: ", coord)
234 259 m_img_flip = m_img.copy()
235 260 m_img_flip[1, -1] = -m_img_flip[1, -1]
  261 +
  262 + if self.icp:
  263 + m_img = bases.transform_icp(m_img, self.m_icp)
  264 +
236 265 self.coord_queue.put_nowait([coord, m_img, view_obj])
237 266  
238 267 if self.view_tracts:
239 268 self.coord_tracts_queue.put_nowait(m_img_flip)
240 269  
  270 + if not self.icp_queue.empty():
  271 + self.icp_queue.task_done()
241 272 # The sleep has to be in both threads
242 273 sleep(self.sle)
243 274 except queue.Full:
... ...
invesalius/data/tractography.py
... ... @@ -230,7 +230,7 @@ def tracts_computation_branch(trk_list, alpha=255):
230 230  
231 231 class ComputeTractsThread(threading.Thread):
232 232  
233   - def __init__(self, inp, coord_tracts_queue, tracts_queue, event, sle):
  233 + def __init__(self, inp, queues, event, sle):
234 234 """Class (threading) to compute real time tractography data for visualization.
235 235  
236 236 Tracts are computed using the Trekker library by Baran Aydogan (https://dmritrekker.github.io/)
... ... @@ -243,10 +243,9 @@ class ComputeTractsThread(threading.Thread):
243 243  
244 244 :param inp: List of inputs: trekker instance, affine numpy array, seed_offset, seed_radius, n_threads
245 245 :type inp: list
246   - :param coord_tracts_queue: Queue instance that manage co-registered coordinates
247   - :type coord_tracts_queue: queue.Queue
248   - :param tracts_queue: Queue instance that manage the tracts to be visualized
249   - :type tracts_queue: queue.Queue
  246 + :param queues: Queue list with coord_tracts_queue (Queue instance that manage co-registered coordinates) and
  247 + tracts_queue (Queue instance that manage the tracts to be visualized)
  248 + :type queues: list[queue.Queue, queue.Queue]
250 249 :param event: Threading event to coordinate when tasks as done and allow UI release
251 250 :type event: threading.Event
252 251 :param sle: Sleep pause in seconds
... ... @@ -256,8 +255,8 @@ class ComputeTractsThread(threading.Thread):
256 255 threading.Thread.__init__(self, name='ComputeTractsThread')
257 256 self.inp = inp
258 257 # self.coord_queue = coord_queue
259   - self.coord_tracts_queue = coord_tracts_queue
260   - self.tracts_queue = tracts_queue
  258 + self.coord_tracts_queue = queues[0]
  259 + self.tracts_queue = queues[1]
261 260 # self.visualization_queue = visualization_queue
262 261 self.event = event
263 262 self.sle = sle
... ... @@ -362,7 +361,7 @@ class ComputeTractsThread(threading.Thread):
362 361  
363 362 class ComputeTractsACTThread(threading.Thread):
364 363  
365   - def __init__(self, inp, coord_tracts_queue, tracts_queue, event, sle):
  364 + def __init__(self, inp, queues, event, sle):
366 365 """Class (threading) to compute real time tractography data for visualization.
367 366  
368 367 Tracts are computed using the Trekker library by Baran Aydogan (https://dmritrekker.github.io/)
... ... @@ -375,10 +374,9 @@ class ComputeTractsACTThread(threading.Thread):
375 374  
376 375 :param inp: List of inputs: trekker instance, affine numpy array, seed_offset, seed_radius, n_threads
377 376 :type inp: list
378   - :param coord_tracts_queue: Queue instance that manage co-registered coordinates
379   - :type coord_tracts_queue: queue.Queue
380   - :param tracts_queue: Queue instance that manage the tracts to be visualized
381   - :type tracts_queue: queue.Queue
  377 + :param queues: Queue list with coord_tracts_queue (Queue instance that manage co-registered coordinates) and
  378 + tracts_queue (Queue instance that manage the tracts to be visualized)
  379 + :type queues: list[queue.Queue, queue.Queue]
382 380 :param event: Threading event to coordinate when tasks as done and allow UI release
383 381 :type event: threading.Event
384 382 :param sle: Sleep pause in seconds
... ... @@ -388,8 +386,8 @@ class ComputeTractsACTThread(threading.Thread):
388 386 threading.Thread.__init__(self, name='ComputeTractsThreadACT')
389 387 self.inp = inp
390 388 # self.coord_queue = coord_queue
391   - self.coord_tracts_queue = coord_tracts_queue
392   - self.tracts_queue = tracts_queue
  389 + self.coord_tracts_queue = queues[0]
  390 + self.tracts_queue = queues[1]
393 391 # on first pilots (january 12, 2021) used (-4, 4)
394 392 self.coord_list_w = img_utils.create_grid((-2, 2), (0, 20), inp[2]-5, 1)
395 393 # self.coord_list_sph = img_utils.create_spherical_grid(10, 1)
... ...
invesalius/data/viewer_slice.py
... ... @@ -837,8 +837,7 @@ class Viewer(wx.Panel):
837 837 # Publisher.subscribe(self.__update_cross_position,
838 838 # 'Update cross position %s' % self.orientation)
839 839 Publisher.subscribe(self.SetCrossFocalPoint, 'Set cross focal point')
840   - # Publisher.subscribe(self.UpdateSlicesPosition,
841   - # 'Co-registered points')
  840 + Publisher.subscribe(self.UpdateSlicesPosition, 'Update slices position')
842 841 ###
843 842 # Publisher.subscribe(self.ChangeBrushColour,
844 843 # 'Add mask')
... ... @@ -849,9 +848,9 @@ class Viewer(wx.Panel):
849 848 Publisher.subscribe(self.UpdateWindowLevelText,
850 849 'Update window level text')
851 850  
852   - #Publisher.subscribe(self._set_cross_visibility,\
853   - # 'Set cross visibility')
854   - ###
  851 + Publisher.subscribe(self._set_cross_visibility,
  852 + 'Set cross visibility')
  853 +
855 854 Publisher.subscribe(self.__set_layout,
856 855 'Set slice viewer layout')
857 856  
... ...
invesalius/gui/dialogs.py
... ... @@ -23,6 +23,7 @@ import os
23 23 import random
24 24 import sys
25 25 import time
  26 +from functools import partial
26 27  
27 28 from concurrent import futures
28 29  
... ... @@ -58,10 +59,12 @@ import invesalius.data.coordinates as dco
58 59 import invesalius.gui.widgets.gradient as grad
59 60 import invesalius.session as ses
60 61 import invesalius.utils as utils
  62 +import invesalius.data.bases as bases
61 63 from invesalius.gui.widgets.inv_spinctrl import InvSpinCtrl, InvFloatSpinCtrl
62 64 from invesalius.gui.widgets import clut_imagedata
63 65 from invesalius.gui.widgets.clut_imagedata import CLUTImageDataWidget, EVT_CLUT_NODE_CHANGED
64 66 import numpy as np
  67 +from numpy.core.umath_tests import inner1d
65 68  
66 69 from invesalius import inv_paths
67 70  
... ... @@ -887,6 +890,35 @@ def ShowNavigationTrackerWarning(trck_id, lib_mode):
887 890 dlg.ShowModal()
888 891 dlg.Destroy()
889 892  
  893 +def ICPcorregistration(fre):
  894 + msg = _("The fiducial registration error is: ") + str(round(fre, 2)) + '\n\n' + \
  895 + _("Would you like to improve accuracy?")
  896 + if sys.platform == 'darwin':
  897 + dlg = wx.MessageDialog(None, "", msg,
  898 + wx.YES_NO)
  899 + else:
  900 + dlg = wx.MessageDialog(None, msg, "InVesalius 3",
  901 + wx.YES_NO)
  902 +
  903 + if dlg.ShowModal() == wx.ID_YES:
  904 + flag = True
  905 + else:
  906 + flag = False
  907 +
  908 + dlg.Destroy()
  909 + return flag
  910 +
  911 +def ReportICPerror(prev_error, final_error):
  912 + msg = _("Error after refine: ") + str(round(final_error, 2)) + ' mm' + '\n\n' + \
  913 + _("Previous error: ") + str(round(prev_error, 2)) + ' mm'
  914 + if sys.platform == 'darwin':
  915 + dlg = wx.MessageDialog(None, "", msg,
  916 + wx.OK)
  917 + else:
  918 + dlg = wx.MessageDialog(None, msg, "InVesalius 3",
  919 + wx.OK)
  920 + dlg.ShowModal()
  921 + dlg.Destroy()
890 922  
891 923 def ShowEnterMarkerID(default):
892 924 msg = _("Edit marker ID")
... ... @@ -1155,7 +1187,7 @@ def ShowAboutDialog(parent):
1155 1187  
1156 1188 info.SetWebSite("https://www.cti.gov.br/invesalius")
1157 1189 info.SetIcon(icon)
1158   -
  1190 +
1159 1191 info.License = _("GNU GPL (General Public License) version 2")
1160 1192  
1161 1193 info.Developers = [u"Paulo Henrique Junqueira Amorim",
... ... @@ -1343,7 +1375,7 @@ class NewSurfaceDialog(wx.Dialog):
1343 1375 def ExportPicture(type_=""):
1344 1376 import invesalius.constants as const
1345 1377 import invesalius.project as proj
1346   -
  1378 +
1347 1379 INDEX_TO_EXTENSION = {0: "bmp", 1: "jpg", 2: "png", 3: "ps", 4:"povray", 5:"tiff"}
1348 1380 WILDCARD_SAVE_PICTURE = _("BMP image")+" (*.bmp)|*.bmp|"+\
1349 1381 _("JPG image")+" (*.jpg)|*.jpg|"+\
... ... @@ -1513,7 +1545,7 @@ class SurfaceCreationOptionsPanel(wx.Panel):
1513 1545 project = prj.Project()
1514 1546 index_list = project.mask_dict.keys()
1515 1547 self.mask_list = [project.mask_dict[index].name for index in sorted(index_list)]
1516   -
  1548 +
1517 1549 active_mask = slc.Slice().current_mask.index
1518 1550 #active_mask = len(self.mask_list)-1
1519 1551  
... ... @@ -2093,17 +2125,17 @@ class ImportBitmapParameters(wx.Dialog):
2093 2125  
2094 2126  
2095 2127 def _init_gui(self):
2096   -
  2128 +
2097 2129 import invesalius.project as prj
2098   -
  2130 +
2099 2131 p = wx.Panel(self, -1, style = wx.TAB_TRAVERSAL
2100 2132 | wx.CLIP_CHILDREN
2101 2133 | wx.FULL_REPAINT_ON_RESIZE)
2102   -
  2134 +
2103 2135 gbs_principal = self.gbs = wx.GridBagSizer(4,1)
2104 2136  
2105 2137 gbs = self.gbs = wx.GridBagSizer(5, 2)
2106   -
  2138 +
2107 2139 flag_labels = wx.ALIGN_RIGHT | wx.ALIGN_CENTER_VERTICAL
2108 2140  
2109 2141 stx_name = wx.StaticText(p, -1, _(u"Project name:"))
... ... @@ -2134,7 +2166,7 @@ class ImportBitmapParameters(wx.Dialog):
2134 2166  
2135 2167 #--- spacing --------------
2136 2168 gbs_spacing = wx.GridBagSizer(2, 6)
2137   -
  2169 +
2138 2170 stx_spacing_x = stx_spacing_x = wx.StaticText(p, -1, _(u"X:"))
2139 2171 fsp_spacing_x = self.fsp_spacing_x = InvFloatSpinCtrl(p, -1, min_value=0, max_value=1000000000,
2140 2172 increment=0.25, value=1.0, digits=8)
... ... @@ -2151,7 +2183,7 @@ class ImportBitmapParameters(wx.Dialog):
2151 2183  
2152 2184 try:
2153 2185 proj = prj.Project()
2154   -
  2186 +
2155 2187 sx = proj.spacing[0]
2156 2188 sy = proj.spacing[1]
2157 2189 sz = proj.spacing[2]
... ... @@ -2174,7 +2206,7 @@ class ImportBitmapParameters(wx.Dialog):
2174 2206  
2175 2207 #----- buttons ------------------------
2176 2208 gbs_button = wx.GridBagSizer(2, 4)
2177   -
  2209 +
2178 2210 btn_ok = self.btn_ok= wx.Button(p, wx.ID_OK)
2179 2211 btn_ok.SetDefault()
2180 2212  
... ... @@ -2197,14 +2229,14 @@ class ImportBitmapParameters(wx.Dialog):
2197 2229  
2198 2230 box = wx.BoxSizer()
2199 2231 box.Add(gbs_principal, 1, wx.ALL|wx.EXPAND, 10)
2200   -
  2232 +
2201 2233 p.SetSizer(box)
2202 2234 box.Fit(self)
2203 2235 self.Layout()
2204 2236  
2205 2237 def bind_evts(self):
2206 2238 self.btn_ok.Bind(wx.EVT_BUTTON, self.OnOk)
2207   -
  2239 +
2208 2240 def SetInterval(self, v):
2209 2241 self.interval = v
2210 2242  
... ... @@ -2228,10 +2260,10 @@ class ImportBitmapParameters(wx.Dialog):
2228 2260  
2229 2261  
2230 2262 def BitmapNotSameSize():
2231   -
  2263 +
2232 2264 dlg = wx.MessageDialog(None,_("All bitmaps files must be the same \n width and height size."), 'Error',\
2233 2265 wx.OK | wx.ICON_ERROR)
2234   -
  2266 +
2235 2267 dlg.ShowModal()
2236 2268 dlg.Destroy()
2237 2269  
... ... @@ -2908,7 +2940,7 @@ class FFillSegmentationOptionsDialog(wx.Dialog):
2908 2940  
2909 2941  
2910 2942 class CropOptionsDialog(wx.Dialog):
2911   -
  2943 +
2912 2944 def __init__(self, config, ID=-1, title=_(u"Crop mask"), style=wx.DEFAULT_DIALOG_STYLE|wx.FRAME_FLOAT_ON_PARENT|wx.STAY_ON_TOP):
2913 2945 self.config = config
2914 2946 wx.Dialog.__init__(self, wx.GetApp().GetTopWindow(), ID, title=title, style=style)
... ... @@ -3316,7 +3348,7 @@ class ObjectCalibrationDialog(wx.Dialog):
3316 3348 main_sizer = wx.BoxSizer(wx.VERTICAL)
3317 3349 main_sizer.Add(self.interactor, 0, wx.EXPAND)
3318 3350 main_sizer.Add(group_sizer, 0,
3319   - wx.EXPAND|wx.GROW|wx.LEFT|wx.TOP|wx.RIGHT|wx.BOTTOM|wx.ALIGN_CENTER_HORIZONTAL, 10)
  3351 + wx.EXPAND|wx.GROW|wx.LEFT|wx.TOP|wx.RIGHT|wx.BOTTOM, 10)
3320 3352  
3321 3353 self.SetSizer(main_sizer)
3322 3354 main_sizer.Fit(self)
... ... @@ -3493,6 +3525,363 @@ class ObjectCalibrationDialog(wx.Dialog):
3493 3525 def GetValue(self):
3494 3526 return self.obj_fiducials, self.obj_orients, self.obj_ref_id, self.obj_name, self.polydata
3495 3527  
  3528 +class ICPCorregistrationDialog(wx.Dialog):
  3529 +
  3530 + def __init__(self, nav_prop):
  3531 + import invesalius.project as prj
  3532 +
  3533 + self.__bind_events()
  3534 +
  3535 + self.tracker_id = nav_prop[0]
  3536 + self.trk_init = nav_prop[1]
  3537 + self.obj_ref_id = 2
  3538 + self.obj_name = None
  3539 + self.obj_actor = None
  3540 + self.polydata = None
  3541 + self.m_icp = None
  3542 + self.initial_focus = None
  3543 + self.prev_error = None
  3544 + self.final_error = None
  3545 + self.icp_mode = 0
  3546 + self.staticballs = []
  3547 + self.point_coord = []
  3548 + self.transformed_points = []
  3549 +
  3550 + self.obj_fiducials = np.full([5, 3], np.nan)
  3551 + self.obj_orients = np.full([5, 3], np.nan)
  3552 +
  3553 + wx.Dialog.__init__(self, wx.GetApp().GetTopWindow(), -1, _(u"Refine Corregistration"), size=(380, 440),
  3554 + style=wx.DEFAULT_DIALOG_STYLE | wx.FRAME_FLOAT_ON_PARENT|wx.STAY_ON_TOP)
  3555 +
  3556 + self.proj = prj.Project()
  3557 +
  3558 + self._init_gui()
  3559 +
  3560 + def __bind_events(self):
  3561 + Publisher.subscribe(self.UpdateCurrentCoord, 'Set cross focal point')
  3562 +
  3563 + def UpdateCurrentCoord(self, position):
  3564 + self.current_coord = position[:]
  3565 +
  3566 + def _init_gui(self):
  3567 + self.interactor = wxVTKRenderWindowInteractor(self, -1, size=self.GetSize())
  3568 + self.interactor.Enable(1)
  3569 + self.ren = vtk.vtkRenderer()
  3570 + self.interactor.GetRenderWindow().AddRenderer(self.ren)
  3571 +
  3572 + self.timer = wx.Timer(self)
  3573 + self.Bind(wx.EVT_TIMER, self.OnUpdate, self.timer)
  3574 +
  3575 + txt_surface = wx.StaticText(self, -1, _('Select the surface:'))
  3576 + txt_mode = wx.StaticText(self, -1, _('Registration mode:'))
  3577 +
  3578 + combo_surface_name = wx.ComboBox(self, -1, size=(210, 23),
  3579 + style=wx.CB_DROPDOWN | wx.CB_READONLY)
  3580 + # combo_surface_name.SetSelection(0)
  3581 + if sys.platform != 'win32':
  3582 + combo_surface_name.SetWindowVariant(wx.WINDOW_VARIANT_SMALL)
  3583 + combo_surface_name.Bind(wx.EVT_COMBOBOX, self.OnComboName)
  3584 + for n in range(len(self.proj.surface_dict)):
  3585 + combo_surface_name.Insert(str(self.proj.surface_dict[n].name), n)
  3586 +
  3587 + self.combo_surface_name = combo_surface_name
  3588 +
  3589 + init_surface = 0
  3590 + combo_surface_name.SetSelection(init_surface)
  3591 + self.surface = self.proj.surface_dict[init_surface].polydata
  3592 + self.LoadActor()
  3593 +
  3594 + tooltip = wx.ToolTip(_("Choose the registration mode:"))
  3595 + choice_icp_method = wx.ComboBox(self, -1, "", size=(100, 23),
  3596 + choices=([_("Affine"), _("Similarity"), _("RigidBody")]),
  3597 + style=wx.CB_DROPDOWN|wx.CB_READONLY)
  3598 + choice_icp_method.SetSelection(0)
  3599 + choice_icp_method.SetToolTip(tooltip)
  3600 + choice_icp_method.Bind(wx.EVT_COMBOBOX, self.OnChoiceICPMethod)
  3601 +
  3602 + # Buttons to acquire and remove points
  3603 + create_point = wx.Button(self, -1, label=_('Create point'))
  3604 + create_point.Bind(wx.EVT_BUTTON, self.OnCreatePoint)
  3605 +
  3606 + cont_point = wx.ToggleButton(self, -1, label=_('Continuous acquisition'))
  3607 + cont_point.Bind(wx.EVT_TOGGLEBUTTON, partial(self.OnContinuousAcquisition, btn=cont_point))
  3608 + self.cont_point = cont_point
  3609 +
  3610 + btn_reset = wx.Button(self, -1, label=_('Remove points'))
  3611 + btn_reset.Bind(wx.EVT_BUTTON, self.OnReset)
  3612 +
  3613 + btn_apply_icp = wx.Button(self, -1, label=_('Apply registration'))
  3614 + btn_apply_icp.Bind(wx.EVT_BUTTON, self.OnICP)
  3615 +
  3616 + # Buttons to finish or cancel object registration
  3617 + tooltip = wx.ToolTip(_(u"Refine done"))
  3618 + btn_ok = wx.Button(self, wx.ID_OK, _(u"Done"))
  3619 + btn_ok.SetToolTip(tooltip)
  3620 +
  3621 + btn_cancel = wx.Button(self, wx.ID_CANCEL)
  3622 + btn_cancel.SetHelpText("")
  3623 +
  3624 + top_sizer = wx.FlexGridSizer(rows=2, cols=2, hgap=50, vgap=5)
  3625 + top_sizer.AddMany([txt_surface, txt_mode,
  3626 + combo_surface_name, choice_icp_method])
  3627 +
  3628 + btn_acqui_sizer = wx.FlexGridSizer(rows=1, cols=3, hgap=15, vgap=15)
  3629 + btn_acqui_sizer.AddMany([create_point, cont_point, btn_reset])
  3630 +
  3631 + btn_ok_sizer = wx.FlexGridSizer(rows=1, cols=3, hgap=20, vgap=20)
  3632 + btn_ok_sizer.AddMany([btn_apply_icp, btn_ok, btn_cancel])
  3633 +
  3634 + btn_sizer = wx.FlexGridSizer(rows=2, cols=1, hgap=50, vgap=20)
  3635 + btn_sizer.AddMany([(btn_acqui_sizer, 1, wx.ALIGN_CENTER_HORIZONTAL),
  3636 + (btn_ok_sizer, 1, wx.ALIGN_RIGHT)])
  3637 +
  3638 + self.progress = wx.Gauge(self, -1)
  3639 +
  3640 + main_sizer = wx.BoxSizer(wx.VERTICAL)
  3641 + main_sizer.Add(top_sizer, 0, wx.LEFT|wx.TOP|wx.BOTTOM, 10)
  3642 + main_sizer.Add(self.interactor, 0, wx.EXPAND)
  3643 + main_sizer.Add(btn_sizer, 0,
  3644 + wx.EXPAND|wx.GROW|wx.LEFT|wx.TOP|wx.BOTTOM, 10)
  3645 + main_sizer.Add(self.progress, 0, wx.EXPAND | wx.ALL, 5)
  3646 +
  3647 + self.SetSizer(main_sizer)
  3648 + main_sizer.Fit(self)
  3649 +
  3650 + def LoadActor(self):
  3651 + '''
  3652 + Load the selected actor from the project (self.surface) into the scene
  3653 + :return:
  3654 + '''
  3655 + mapper = vtk.vtkPolyDataMapper()
  3656 + mapper.SetInputData(self.surface)
  3657 + mapper.ScalarVisibilityOff()
  3658 + #mapper.ImmediateModeRenderingOn()
  3659 +
  3660 + obj_actor = vtk.vtkActor()
  3661 + obj_actor.SetMapper(mapper)
  3662 + self.obj_actor = obj_actor
  3663 +
  3664 + self.ren.AddActor(obj_actor)
  3665 + self.ren.ResetCamera()
  3666 + self.interactor.Render()
  3667 +
  3668 + def RemoveActor(self):
  3669 + #self.ren.RemoveActor(self.obj_actor)
  3670 + self.ren.RemoveAllViewProps()
  3671 + self.point_coord = []
  3672 + self.transformed_points = []
  3673 + self.m_icp = None
  3674 + self.SetProgress(0)
  3675 + self.ren.ResetCamera()
  3676 + self.interactor.Render()
  3677 +
  3678 + def AddMarker(self, size, colour, coord):
  3679 + """
  3680 + Points are rendered into the scene. These points give visual information about the registration.
  3681 + :param size: value of the marker size
  3682 + :type size: int
  3683 + :param colour: RGB Color Code for the marker
  3684 + :type colour: tuple (int(R),int(G),int(B))
  3685 + :param coord: x, y, z of the marker
  3686 + :type coord: np.ndarray
  3687 + """
  3688 +
  3689 + x, y, z = coord[0], -coord[1], coord[2]
  3690 +
  3691 + ball_ref = vtk.vtkSphereSource()
  3692 + ball_ref.SetRadius(size)
  3693 + ball_ref.SetCenter(x, y, z)
  3694 +
  3695 + mapper = vtk.vtkPolyDataMapper()
  3696 + mapper.SetInputConnection(ball_ref.GetOutputPort())
  3697 +
  3698 + prop = vtk.vtkProperty()
  3699 + prop.SetColor(colour[0:3])
  3700 +
  3701 + #adding a new actor for the present ball
  3702 + sphere_actor = vtk.vtkActor()
  3703 +
  3704 + sphere_actor.SetMapper(mapper)
  3705 + sphere_actor.SetProperty(prop)
  3706 +
  3707 + self.ren.AddActor(sphere_actor)
  3708 + self.point_coord.append([x, y, z])
  3709 +
  3710 + self.Refresh()
  3711 +
  3712 + def SetProgress(self, progress):
  3713 + self.progress.SetValue(progress * 100)
  3714 + self.Refresh()
  3715 +
  3716 + def vtkmatrix_to_numpy(self, matrix):
  3717 + """
  3718 + Copies the elements of a vtkMatrix4x4 into a numpy array.
  3719 +
  3720 + :param matrix: The matrix to be copied into an array.
  3721 + :type matrix: vtk.vtkMatrix4x4
  3722 + :rtype: numpy.ndarray
  3723 + """
  3724 + m = np.ones((4, 4))
  3725 + for i in range(4):
  3726 + for j in range(4):
  3727 + m[i, j] = matrix.GetElement(i, j)
  3728 + return m
  3729 +
  3730 + def SetCameraVolume(self, position):
  3731 + """
  3732 + Positioning of the camera based on the acquired point
  3733 + :param position: x, y, z of the last acquired point
  3734 + :return:
  3735 + """
  3736 + cam_focus = np.array([position[0], -position[1], position[2]])
  3737 + cam = self.ren.GetActiveCamera()
  3738 +
  3739 + if self.initial_focus is None:
  3740 + self.initial_focus = np.array(cam.GetFocalPoint())
  3741 +
  3742 + cam_pos0 = np.array(cam.GetPosition())
  3743 + cam_focus0 = np.array(cam.GetFocalPoint())
  3744 + v0 = cam_pos0 - cam_focus0
  3745 + v0n = np.sqrt(inner1d(v0, v0))
  3746 +
  3747 + v1 = (cam_focus - self.initial_focus)
  3748 +
  3749 + v1n = np.sqrt(inner1d(v1, v1))
  3750 + if not v1n:
  3751 + v1n = 1.0
  3752 + cam_pos = (v1/v1n)*v0n + cam_focus
  3753 +
  3754 + cam.SetFocalPoint(cam_focus)
  3755 + cam.SetPosition(cam_pos)
  3756 + self.Refresh()
  3757 +
  3758 + def ErrorEstimation(self, surface, points):
  3759 + """
  3760 + Estimation of the average squared distance between the cloud of points to the closest mesh
  3761 + :param surface: Surface polydata of the scene
  3762 + :type surface: vtk.polydata
  3763 + :param points: Cloud of points
  3764 + :type points: np.ndarray
  3765 + :return: mean distance
  3766 + """
  3767 + cell_locator = vtk.vtkCellLocator()
  3768 + cell_locator.SetDataSet(surface)
  3769 + cell_locator.BuildLocator()
  3770 +
  3771 + cellId = vtk.mutable(0)
  3772 + c = [0.0, 0.0, 0.0]
  3773 + subId = vtk.mutable(0)
  3774 + d = vtk.mutable(0.0)
  3775 + error = []
  3776 + for i in range(len(points)):
  3777 + cell_locator.FindClosestPoint(points[i], c, cellId, subId, d)
  3778 + error.append(np.sqrt(float(d)))
  3779 +
  3780 + return np.mean(error)
  3781 +
  3782 + def OnComboName(self, evt):
  3783 + surface_name = evt.GetString()
  3784 + surface_index = evt.GetSelection()
  3785 + self.surface = self.proj.surface_dict[surface_index].polydata
  3786 + if self.obj_actor:
  3787 + self.RemoveActor()
  3788 + self.LoadActor()
  3789 +
  3790 + def OnChoiceICPMethod(self, evt):
  3791 + self.icp_mode = evt.GetSelection()
  3792 +
  3793 + def OnContinuousAcquisition(self, evt=None, btn=None):
  3794 + value = btn.GetValue()
  3795 + if value:
  3796 + self.timer.Start(500)
  3797 + else:
  3798 + self.timer.Stop()
  3799 +
  3800 + def OnUpdate(self, evt):
  3801 + self.AddMarker(3, (1, 0, 0), self.current_coord[:3])
  3802 + self.SetCameraVolume(self.current_coord[:3])
  3803 +
  3804 + def OnCreatePoint(self, evt):
  3805 + self.AddMarker(3,(1,0,0),self.current_coord[:3])
  3806 + self.SetCameraVolume(self.current_coord[:3])
  3807 +
  3808 + def OnReset(self, evt):
  3809 + self.RemoveActor()
  3810 + self.LoadActor()
  3811 +
  3812 + def OnICP(self, evt):
  3813 + self.SetProgress(0.3)
  3814 + if self.cont_point:
  3815 + self.cont_point.SetValue(False)
  3816 + self.OnContinuousAcquisition(evt=None, btn=self.cont_point)
  3817 + sourcePoints = np.array(self.point_coord)
  3818 + sourcePoints_vtk = vtk.vtkPoints()
  3819 +
  3820 + for i in range(len(sourcePoints)):
  3821 + id0 = sourcePoints_vtk.InsertNextPoint(sourcePoints[i])
  3822 +
  3823 + source = vtk.vtkPolyData()
  3824 + source.SetPoints(sourcePoints_vtk)
  3825 +
  3826 + icp = vtk.vtkIterativeClosestPointTransform()
  3827 + icp.SetSource(source)
  3828 + icp.SetTarget(self.surface)
  3829 +
  3830 + if self.icp_mode == 0:
  3831 + print("Affine mode")
  3832 + icp.GetLandmarkTransform().SetModeToAffine()
  3833 + elif self.icp_mode == 1:
  3834 + print("Similarity mode")
  3835 + icp.GetLandmarkTransform().SetModeToSimilarity()
  3836 + elif self.icp_mode == 2:
  3837 + print("Rigid mode")
  3838 + icp.GetLandmarkTransform().SetModeToRigidBody()
  3839 +
  3840 + #icp.DebugOn()
  3841 + icp.SetMaximumNumberOfIterations(1000)
  3842 +
  3843 + icp.Modified()
  3844 +
  3845 + icp.Update()
  3846 +
  3847 + self.m_icp = self.vtkmatrix_to_numpy(icp.GetMatrix())
  3848 +
  3849 + icpTransformFilter = vtk.vtkTransformPolyDataFilter()
  3850 + icpTransformFilter.SetInputData(source)
  3851 +
  3852 + icpTransformFilter.SetTransform(icp)
  3853 + icpTransformFilter.Update()
  3854 +
  3855 + transformedSource = icpTransformFilter.GetOutput()
  3856 +
  3857 + self.SetProgress(1)
  3858 +
  3859 + for i in range(transformedSource.GetNumberOfPoints()):
  3860 + p = [0, 0, 0]
  3861 + transformedSource.GetPoint(i, p)
  3862 + self.transformed_points.append(p)
  3863 + point = vtk.vtkSphereSource()
  3864 + point.SetCenter(p)
  3865 + point.SetRadius(3)
  3866 + point.SetPhiResolution(3)
  3867 + point.SetThetaResolution(3)
  3868 +
  3869 + mapper = vtk.vtkPolyDataMapper()
  3870 + mapper.SetInputConnection(point.GetOutputPort())
  3871 +
  3872 + actor = vtk.vtkActor()
  3873 + actor.SetMapper(mapper)
  3874 + actor.GetProperty().SetColor((0,1,0))
  3875 +
  3876 + self.ren.AddActor(actor)
  3877 +
  3878 + self.prev_error = self.ErrorEstimation(self.surface, sourcePoints)
  3879 + self.final_error = self.ErrorEstimation(self.surface, self.transformed_points)
  3880 +
  3881 + self.Refresh()
  3882 +
  3883 + def GetValue(self):
  3884 + return self.m_icp, self.point_coord, self.transformed_points, self.prev_error, self.final_error
3496 3885  
3497 3886 class SurfaceProgressWindow(object):
3498 3887 def __init__(self):
... ... @@ -3752,20 +4141,20 @@ class SetNDIconfigs(wx.Dialog):
3752 4141 wildcard="Rom files (*.rom)|*.rom", message="Select probe's rom file")
3753 4142 row_probe = wx.BoxSizer(wx.VERTICAL)
3754 4143 row_probe.Add(wx.StaticText(self, wx.ID_ANY, "Set probe's rom file"), 0, wx.TOP|wx.RIGHT, 5)
3755   - row_probe.Add(self.dir_probe, 0, wx.EXPAND|wx.ALIGN_CENTER)
  4144 + row_probe.Add(self.dir_probe, 0, wx.ALIGN_CENTER)
3756 4145  
3757 4146 self.dir_ref = wx.FilePickerCtrl(self, path=last_ndi_ref_marker, style=wx.FLP_USE_TEXTCTRL|wx.FLP_SMALL,
3758 4147 wildcard="Rom files (*.rom)|*.rom", message="Select reference's rom file")
3759 4148 row_ref = wx.BoxSizer(wx.VERTICAL)
3760 4149 row_ref.Add(wx.StaticText(self, wx.ID_ANY, "Set reference's rom file"), 0, wx.TOP | wx.RIGHT, 5)
3761   - row_ref.Add(self.dir_ref, 0, wx.EXPAND|wx.ALIGN_CENTER)
  4150 + row_ref.Add(self.dir_ref, 0, wx.ALIGN_CENTER)
3762 4151  
3763 4152 self.dir_obj = wx.FilePickerCtrl(self, path=last_ndi_obj_marker, style=wx.FLP_USE_TEXTCTRL|wx.FLP_SMALL,
3764 4153 wildcard="Rom files (*.rom)|*.rom", message="Select object's rom file")
3765 4154 #self.dir_probe.Bind(wx.EVT_FILEPICKER_CHANGED, self.Selected)
3766 4155 row_obj = wx.BoxSizer(wx.VERTICAL)
3767 4156 row_obj.Add(wx.StaticText(self, wx.ID_ANY, "Set object's rom file"), 0, wx.TOP|wx.RIGHT, 5)
3768   - row_obj.Add(self.dir_obj, 0, wx.EXPAND|wx.ALIGN_CENTER)
  4157 + row_obj.Add(self.dir_obj, 0, wx.ALIGN_CENTER)
3769 4158  
3770 4159 btn_ok = wx.Button(self, wx.ID_OK)
3771 4160 btn_ok.SetHelpText("")
... ...
invesalius/gui/task_navigator.py
... ... @@ -309,17 +309,24 @@ class NeuronavigationPanel(wx.Panel):
309 309  
310 310 # Initialize global variables
311 311 self.fiducials = np.full([6, 3], np.nan)
  312 + self.fiducials_raw = np.zeros((6, 6))
312 313 self.correg = None
313 314 self.current_coord = 0, 0, 0
314 315 self.trk_init = None
  316 + self.nav_status = False
315 317 self.trigger = None
316 318 self.trigger_state = False
317 319 self.obj_reg = None
318 320 self.obj_reg_status = False
319 321 self.track_obj = False
  322 + self.m_icp = None
  323 + self.fre = None
  324 + self.icp_fre = None
  325 + self.icp = False
320 326 self.event = threading.Event()
321 327  
322 328 self.coord_queue = QueueCustom(maxsize=1)
  329 + self.icp_queue = QueueCustom(maxsize=1)
323 330 # self.visualization_queue = QueueCustom(maxsize=1)
324 331 self.trigger_queue = QueueCustom(maxsize=1)
325 332 self.coord_tracts_queue = QueueCustom(maxsize=1)
... ... @@ -385,6 +392,7 @@ class NeuronavigationPanel(wx.Panel):
385 392  
386 393 # TODO: Find a better allignment between FRE, text and navigate button
387 394 txt_fre = wx.StaticText(self, -1, _('FRE:'))
  395 + txt_icp = wx.StaticText(self, -1, _('Refine:'))
388 396  
389 397 # Fiducial registration error text box
390 398 tooltip = wx.ToolTip(_("Fiducial registration error"))
... ... @@ -401,6 +409,14 @@ class NeuronavigationPanel(wx.Panel):
401 409 btn_nav.SetToolTip(tooltip)
402 410 btn_nav.Bind(wx.EVT_TOGGLEBUTTON, partial(self.OnNavigate, btn=(btn_nav, choice_trck, choice_ref)))
403 411  
  412 + tooltip = wx.ToolTip(_(u"Refine the coregistration"))
  413 + checkicp = wx.CheckBox(self, -1, _(' '))
  414 + checkicp.SetValue(False)
  415 + checkicp.Enable(False)
  416 + checkicp.Bind(wx.EVT_CHECKBOX, partial(self.Oncheckicp, ctrl=checkicp))
  417 + checkicp.SetToolTip(tooltip)
  418 + self.checkicp = checkicp
  419 +
404 420 # Image and tracker coordinates number controls
405 421 for m in range(len(self.btns_coord)):
406 422 for n in range(3):
... ... @@ -421,10 +437,12 @@ class NeuronavigationPanel(wx.Panel):
421 437 if m in range(1, 6):
422 438 self.numctrls_coord[m][n].SetEditable(False)
423 439  
424   - nav_sizer = wx.FlexGridSizer(rows=1, cols=3, hgap=5, vgap=5)
  440 + nav_sizer = wx.FlexGridSizer(rows=1, cols=5, hgap=5, vgap=5)
425 441 nav_sizer.AddMany([(txt_fre, 0, wx.ALIGN_CENTER_HORIZONTAL | wx.ALIGN_CENTER_VERTICAL),
426 442 (txtctrl_fre, 0, wx.ALIGN_CENTER_HORIZONTAL | wx.ALIGN_CENTER_VERTICAL),
427   - (btn_nav, 0, wx.ALIGN_CENTER_HORIZONTAL | wx.ALIGN_CENTER_VERTICAL)])
  443 + (btn_nav, 0, wx.ALIGN_CENTER_HORIZONTAL | wx.ALIGN_CENTER_VERTICAL),
  444 + (txt_icp, 0, wx.ALIGN_CENTER_HORIZONTAL | wx.ALIGN_CENTER_VERTICAL),
  445 + (checkicp, 0, wx.ALIGN_CENTER_HORIZONTAL | wx.ALIGN_CENTER_VERTICAL)])
428 446  
429 447 group_sizer = wx.FlexGridSizer(rows=9, cols=1, hgap=5, vgap=5)
430 448 group_sizer.AddGrowableCol(0, 1)
... ... @@ -459,6 +477,7 @@ class NeuronavigationPanel(wx.Panel):
459 477 Publisher.subscribe(self.UpdateTractsVisualization, 'Update tracts visualization')
460 478 Publisher.subscribe(self.EnableACT, 'Enable ACT')
461 479 Publisher.subscribe(self.UpdateACTData, 'Update ACT data')
  480 + Publisher.subscribe(self.UpdateNavigationStatus, 'Navigation status')
462 481  
463 482 def LoadImageFiducials(self, marker_id, coord):
464 483 for n in const.BTNS_IMG_MKS:
... ... @@ -470,6 +489,14 @@ class NeuronavigationPanel(wx.Panel):
470 489 for m in [0, 1, 2]:
471 490 self.numctrls_coord[btn_id][m].SetValue(coord[m])
472 491  
  492 + def UpdateNavigationStatus(self, nav_status, vis_status):
  493 + self.nav_status = nav_status
  494 + if nav_status and (self.m_icp is not None):
  495 + self.checkicp.Enable(True)
  496 + else:
  497 + self.checkicp.Enable(False)
  498 + #self.checkicp.SetValue(False)
  499 +
473 500 def UpdateFRE(self, fre):
474 501 # TODO: Exhibit FRE in a warning dialog and only starts navigation after user clicks ok
475 502 self.txtctrl_fre.SetValue(str(round(fre, 2)))
... ... @@ -661,9 +688,57 @@ class NeuronavigationPanel(wx.Panel):
661 688 # Update number controls with tracker coordinates
662 689 if coord is not None:
663 690 self.fiducials[btn_id, :] = coord[0:3]
  691 + if btn_id == 3:
  692 + self.fiducials_raw[0, :] = coord_raw[0, :]
  693 + self.fiducials_raw[1, :] = coord_raw[1, :]
  694 + elif btn_id == 4:
  695 + self.fiducials_raw[2, :] = coord_raw[0, :]
  696 + self.fiducials_raw[3, :] = coord_raw[1, :]
  697 + else:
  698 + self.fiducials_raw[4, :] = coord_raw[0, :]
  699 + self.fiducials_raw[5, :] = coord_raw[1, :]
  700 +
664 701 for n in [0, 1, 2]:
665 702 self.numctrls_coord[btn_id][n].SetValue(float(coord[n]))
666 703  
  704 + def OnICP(self):
  705 + dialog = dlg.ICPCorregistrationDialog(nav_prop=(self.tracker_id, self.trk_init, self.ref_mode_id))
  706 + if dialog.ShowModal() == wx.ID_OK:
  707 + self.m_icp, point_coord, transformed_points, prev_error, final_error = dialog.GetValue()
  708 + #TODO: checkbox in the dialog to transfer the icp points to 3D viewer
  709 + #create markers
  710 + # for i in range(len(point_coord)):
  711 + # img_coord = point_coord[i][0],-point_coord[i][1],point_coord[i][2], 0, 0, 0
  712 + # transf_coord = transformed_points[i][0],-transformed_points[i][1],transformed_points[i][2], 0, 0, 0
  713 + # Publisher.sendMessage('Create marker', coord=img_coord, marker_id=None, colour=(1,0,0))
  714 + # Publisher.sendMessage('Create marker', coord=transf_coord, marker_id=None, colour=(0,0,1))
  715 + if self.m_icp is not None:
  716 + dlg.ReportICPerror(prev_error, final_error)
  717 + self.checkicp.Enable(True)
  718 + self.checkicp.SetValue(True)
  719 + self.icp = True
  720 + else:
  721 + self.checkicp.Enable(False)
  722 + self.checkicp.SetValue(False)
  723 + self.icp = False
  724 +
  725 + return self.m_icp
  726 +
  727 + def Oncheckicp(self, evt, ctrl):
  728 + if ctrl.GetValue() and evt and (self.m_icp is not None):
  729 + self.icp = True
  730 + else:
  731 + self.icp = False
  732 + self.ctrl_icp()
  733 +
  734 + def ctrl_icp(self):
  735 + if self.icp:
  736 + self.UpdateFRE(self.icp_fre)
  737 + else:
  738 + self.UpdateFRE(self.fre)
  739 + self.icp_queue.put_nowait([self.icp, self.m_icp])
  740 + #print(self.icp, self.m_icp)
  741 +
667 742 def OnNavigate(self, evt, btn):
668 743 btn_nav = btn[0]
669 744 choice_trck = btn[1]
... ... @@ -673,7 +748,7 @@ class NeuronavigationPanel(wx.Panel):
673 748 # initialize jobs list
674 749 jobs_list = []
675 750 vis_components = [self.trigger_state, self.view_tracts]
676   - vis_queues = [self.coord_queue, self.trigger_queue, self.tracts_queue]
  751 + vis_queues = [self.coord_queue, self.trigger_queue, self.tracts_queue, self.icp_queue]
677 752  
678 753 nav_id = btn_nav.GetValue()
679 754 if not nav_id:
... ... @@ -752,10 +827,9 @@ class NeuronavigationPanel(wx.Panel):
752 827 tracker_mode = self.trk_init, self.tracker_id, self.ref_mode_id
753 828  
754 829 # compute fiducial registration error (FRE)
755   - # this is the old way to compute the fre, left here to recheck if new works fine.
756   - # fre = db.calculate_fre(self.fiducials, minv, n, q1, q2)
757   - fre = db.calculate_fre_m(self.fiducials)
758   - self.UpdateFRE(fre)
  830 + if not self.icp_fre:
  831 + self.fre = db.calculate_fre(self.fiducials_raw, self.fiducials, self.ref_mode_id, m_change)
  832 + self.UpdateFRE(self.fre)
759 833  
760 834 if self.track_obj:
761 835 # if object tracking is selected
... ... @@ -778,15 +852,15 @@ class NeuronavigationPanel(wx.Panel):
778 852 obj_data = db.object_registration(obj_fiducials, obj_orients, coord_raw, m_change)
779 853 coreg_data.extend(obj_data)
780 854  
781   - jobs_list.append(dcr.CoordinateCorregistrate(self.ref_mode_id, tracker_mode, coreg_data, self.coord_queue,
782   - self.view_tracts, self.coord_tracts_queue,
  855 + queues = [self.coord_queue, self.coord_tracts_queue, self.icp_queue]
  856 + jobs_list.append(dcr.CoordinateCorregistrate(self.ref_mode_id, tracker_mode, coreg_data,
  857 + self.view_tracts, queues,
783 858 self.event, self.sleep_nav))
784   -
785 859 else:
786 860 coreg_data = (m_change, 0)
  861 + queues = [self.coord_queue, self.coord_tracts_queue, self.icp_queue]
787 862 jobs_list.append(dcr.CoordinateCorregistrateNoObject(self.ref_mode_id, tracker_mode, coreg_data,
788   - self.coord_queue,
789   - self.view_tracts, self.coord_tracts_queue,
  863 + self.view_tracts, queues,
790 864 self.event, self.sleep_nav))
791 865  
792 866 if not errors:
... ... @@ -807,12 +881,11 @@ class NeuronavigationPanel(wx.Panel):
807 881 self.trk_inp = self.trekker, affine, self.seed_offset, self.n_tracts, self.seed_radius,\
808 882 self.n_threads, self.act_data, affine_vtk, matrix_shape[1]
809 883 # print("Appending the tract computation thread!")
  884 + queues = [self.coord_tracts_queue, self.tracts_queue]
810 885 if self.enable_act:
811   - jobs_list.append(dti.ComputeTractsACTThread(self.trk_inp, self.coord_tracts_queue,
812   - self.tracts_queue, self.event, self.sleep_nav))
  886 + jobs_list.append(dti.ComputeTractsACTThread(self.trk_inp, queues, self.event, self.sleep_nav))
813 887 else:
814   - jobs_list.append(dti.ComputeTractsThread(self.trk_inp, self.coord_tracts_queue,
815   - self.tracts_queue, self.event, self.sleep_nav))
  888 + jobs_list.append(dti.ComputeTractsThread(self.trk_inp, queues, self.event, self.sleep_nav))
816 889  
817 890 jobs_list.append(UpdateNavigationScene(vis_queues, vis_components,
818 891 self.event, self.sleep_nav))
... ... @@ -822,6 +895,13 @@ class NeuronavigationPanel(wx.Panel):
822 895 jobs.start()
823 896 # del jobs
824 897  
  898 + if not self.checkicp.GetValue():
  899 + if dlg.ICPcorregistration(self.fre):
  900 + m_icp = self.OnICP()
  901 + self.icp_fre = db.calculate_fre(self.fiducials_raw, self.fiducials, self.ref_mode_id,
  902 + m_change, m_icp)
  903 + self.ctrl_icp()
  904 +
825 905 def ResetImageFiducials(self):
826 906 for m in range(0, 3):
827 907 self.btns_coord[m].SetValue(False)
... ... @@ -838,15 +918,25 @@ class NeuronavigationPanel(wx.Panel):
838 918 self.txtctrl_fre.SetValue('')
839 919 self.txtctrl_fre.SetBackgroundColour('WHITE')
840 920  
  921 + def ResetIcp(self):
  922 + self.m_icp = None
  923 + self.fre = None
  924 + self.icp_fre = None
  925 + self.icp = False
  926 + self.checkicp.Enable(False)
  927 + self.checkicp.SetValue(False)
  928 +
841 929 def OnCloseProject(self):
842 930 self.ResetTrackerFiducials()
843 931 self.ResetImageFiducials()
  932 + self.ResetIcp()
844 933 self.OnChoiceTracker(False, self.choice_trck)
845 934 Publisher.sendMessage('Update object registration')
846 935 Publisher.sendMessage('Update track object state', flag=False, obj_name=False)
847 936 Publisher.sendMessage('Delete all markers')
848 937 Publisher.sendMessage("Update marker offset state", create=False)
849 938 Publisher.sendMessage("Remove tracts")
  939 + Publisher.sendMessage("Set cross visibility", visibility=0)
850 940 # TODO: Reset camera initial focus
851 941 Publisher.sendMessage('Reset cam clipping range')
852 942  
... ... @@ -1373,18 +1463,23 @@ class MarkersPanel(wx.Panel):
1373 1463 self.marker_ind -= 1
1374 1464 Publisher.sendMessage('Remove marker', index=index)
1375 1465  
1376   - def OnCreateMarker(self, evt=None, coord=None, marker_id=None):
  1466 + def OnCreateMarker(self, evt=None, coord=None, marker_id=None, colour=None):
1377 1467 # OnCreateMarker is used for both pubsub and button click events
1378 1468 # Pubsub is used for markers created with fiducial buttons, trigger and create marker button
  1469 + if not colour:
  1470 + colour = self.marker_colour
  1471 + if not coord:
  1472 + coord = self.current_coord
  1473 +
1379 1474 if evt is None:
1380 1475 if coord:
1381 1476 self.CreateMarker(coord=coord, colour=(0.0, 1.0, 0.0), size=self.marker_size,
1382 1477 marker_id=marker_id, seed=self.current_seed)
1383 1478 else:
1384   - self.CreateMarker(coord=self.current_coord, colour=self.marker_colour, size=self.marker_size,
  1479 + self.CreateMarker(coord=self.current_coord, colour=colour, size=self.marker_size,
1385 1480 seed=self.current_seed)
1386 1481 else:
1387   - self.CreateMarker(coord=self.current_coord, colour=self.marker_colour, size=self.marker_size,
  1482 + self.CreateMarker(coord=self.current_coord, colour=colour, size=self.marker_size,
1388 1483 seed=self.current_seed)
1389 1484  
1390 1485 def OnLoadMarkers(self, evt):
... ... @@ -1807,7 +1902,6 @@ class TractographyPanel(wx.Panel):
1807 1902 self.nav_status = nav_status
1808 1903  
1809 1904 def OnLinkBrain(self, event=None):
1810   -
1811 1905 Publisher.sendMessage('Update status text in GUI', label=_("Busy"))
1812 1906 Publisher.sendMessage('Begin busy cursor')
1813 1907 mask_path = dlg.ShowImportOtherFilesDialog(const.ID_NIFTI_IMPORT, _("Import brain mask"))
... ... @@ -2023,7 +2117,7 @@ class UpdateNavigationScene(threading.Thread):
2023 2117  
2024 2118 threading.Thread.__init__(self, name='UpdateScene')
2025 2119 self.trigger_state, self.view_tracts = vis_components
2026   - self.coord_queue, self.trigger_queue, self.tracts_queue = vis_queues
  2120 + self.coord_queue, self.trigger_queue, self.tracts_queue, self.icp_queue = vis_queues
2027 2121 self.sle = sle
2028 2122 self.event = event
2029 2123  
... ... @@ -2055,6 +2149,7 @@ class UpdateNavigationScene(threading.Thread):
2055 2149  
2056 2150 #TODO: If using the view_tracts substitute the raw coord from the offset coordinate, so the user
2057 2151 # see the red cross in the position of the offset marker
  2152 + wx.CallAfter(Publisher.sendMessage, 'Update slices position', position=coord[:3])
2058 2153 wx.CallAfter(Publisher.sendMessage, 'Set cross focal point', position=coord)
2059 2154 wx.CallAfter(Publisher.sendMessage, 'Update slice viewer')
2060 2155  
... ...