Commit 0673940c840118edf2b74d57dd93d5dd648a61a0

Authored by Victor Hugo Souza
Committed by GitHub
1 parent b7470270
Exists in master

Refactor of thread management in neuronavigation (#242)

* ADD: Initial support for real time tracts visualization during navigation

* ENH: Improved affine transformations for importing external polydata in world coordinates

* ENH: Tracts visualization with cross move started

* Create Trekker thread

* Tracts style working and save and load project with affine

* Tract computing incorporated to navigation pipeline

* Add offset to coil normal to get coordinate inside brain

* ADD: Pipeline for navigation with tract visualiztion and better debug tracker

* ENH: Refactor the cross update feature

* FIX: Fixed update cross during navigation and improved tract visualization

* FIX: Seed flip when applying offset and added list to offset user control

* ADD: Brain surface texture mapping with peeling depth

* Improved sleep for navigation with tractography

* ENH: Tract computation speed doubled and tracts colors fixed

* ENH: Print removed in tract coloring

* ENH: Improved visualization of tracts

* ENH: Divide tract computation in chuncks

* ADD: Threads using Queue for navigation

* ADD: Thread computing with class and pipeline

* ENH: Navigation working with separate threads and lock

* ENH: Navigation with tracts running in producer-consumer method with locks

* ENH: Refactoring tracts threads

* Update tractography data exchange

* ENH: Improved Trekker GUI and tractography navigation working

* Update optimal Trekker parameters

* ADD: Dialog to adjust Trekker parameters

* ENH: Dialog to load Trekker JSON configuration file

* ENH: Replaced Trekker default list by dictionary with values[C

* ENH: Remove tract code comments and refactor navigation pipeline

* ENH: Replaced dialog warnings for message box in navigation tab

* ENH: Split the removal and update of tracts

* ENH: Major changes in navigation array management and tracts cleaning
- Replaced all matrix to array, now using python matrix operator
- Refactored flip_x function
- Improved cleaning of tracts

* ENH: Tract visualization enabled in cross function and matrix operations replaced by arrays
- Additional code cleaning and commenting
- Refactoring of tip offset snippet

* ADD: Control of brain mask opacity

* FIX: Remove brain surface when closing project

* ENH: Update trekker parameters after loading files and major refactor in flip_x
- Additional refactor of FRE computation and use of base_creation

* ENH: Improved object load dialog, OnNavigate pipeline and better tracts cleaning

* FIX: Navigation with tractography working
- Improved offset computation for seed
- Improved flip_x. No need to use a function, only multiply Y coordinate by -1

* Refactored compute_seed, affine_to_vtk, and set_message for tracts computation
- Cleaned unecessary comments
- Providing m_img already flip_x while set_message

* FIX: Tract removal with update cross position was disabled

* ENH: Threading of navigation replaced to  LifoQueue, working with dev_tracker but not tested in real tracker

* ENH: Tractography running with new limited Queue and nowait
- Tracts actors are computed all at once, no optimization for computing blocks

* FIX: Update pubsub in tracker
- Navgation with tractography running

* ENH: Multiblock tract computation and fixed tract cleaning
- Tractography navigation working at 0.1-s sleep

* ENH: Documented tractography computation and cleaned code

* FIX: Navigation was not running if tracts were not found
- created thread in new format for trigger
- removed comments of unused pipelines and threading
- removed commented old FRE base creation

* FIX: Generalize new navigation to objects even without tracts
- Navigation with object working seamless with and without tracts
- Navigation in other modes (single ref, no object) not working
- Added one extra queue to generalize the pipeline in navigation

* FIX: Static object navigation working and refactor navigation pipeline
- Object corregistration split in functions
- Cleaned comments on reference computation

* FIX: Tractography navigation working with and without coil
- Applied the new threading to the non-object navigations
- Improved visualization of the offset marker for tractography
- Cleaned commented code

* ENH: Clean solved TODOs and commented unused navigation code

* FIX: Fixed open and save dialogs for navigation purposes

* FIX: Reading affine from old projects and replaced affine function form transformation

* UPD: Removed surface back culling prints

* DEL: Hided tractography panel and removed brain mesh mask
invesalius/constants.py
... ... @@ -19,6 +19,7 @@
19 19  
20 20 import os.path
21 21 import platform
  22 +import psutil
22 23 import sys
23 24 import wx
24 25 import itertools
... ... @@ -521,6 +522,11 @@ ID_GOTO_COORD = wx.NewId()
521 522  
522 523 ID_MANUAL_WWWL = wx.NewId()
523 524  
  525 +# Tractography with Trekker
  526 +ID_TREKKER_MASK = wx.NewId()
  527 +ID_TREKKER_IMG = wx.NewId()
  528 +ID_TREKKER_FOD = wx.NewId()
  529 +
524 530 #---------------------------------------------------------
525 531 STATE_DEFAULT = 1000
526 532 STATE_WL = 1001
... ... @@ -545,6 +551,7 @@ SLICE_STATE_REMOVE_MASK_PARTS = 3012
545 551 SLICE_STATE_SELECT_MASK_PARTS = 3013
546 552 SLICE_STATE_FFILL_SEGMENTATION = 3014
547 553 SLICE_STATE_CROP_MASK = 3015
  554 +SLICE_STATE_TRACTS = 3016
548 555  
549 556 VOLUME_STATE_SEED = 2001
550 557 # STATE_LINEAR_MEASURE = 3001
... ... @@ -559,7 +566,7 @@ TOOL_STATES = [STATE_WL, STATE_SPIN, STATE_ZOOM,
559 566  
560 567  
561 568 TOOL_SLICE_STATES = [SLICE_STATE_CROSS, SLICE_STATE_SCROLL,
562   - SLICE_STATE_REORIENT]
  569 + SLICE_STATE_REORIENT, SLICE_STATE_TRACTS]
563 570  
564 571  
565 572 SLICE_STYLES = TOOL_STATES + TOOL_SLICE_STATES
... ... @@ -651,6 +658,8 @@ BOOLEAN_XOR = 4
651 658  
652 659 #------------ Navigation defaults -------------------
653 660  
  661 +MARKER_COLOUR = (1.0, 1.0, 0.)
  662 +MARKER_SIZE = 2
654 663 SELECT = 0
655 664 MTC = 1
656 665 FASTRAK = 2
... ... @@ -738,3 +747,17 @@ COIL_COORD_THRESHOLD = 3
738 747 TIMESTAMP = 2.0
739 748  
740 749 CAM_MODE = True
  750 +
  751 +# Tractography visualization
  752 +N_TRACTS = 100
  753 +PEEL_DEPTH = 5
  754 +MAX_PEEL_DEPTH = 30
  755 +SEED_OFFSET = 15
  756 +SEED_RADIUS = 1.5
  757 +SLEEP_NAVIGATION = 0.1
  758 +BRAIN_OPACITY = 0.5
  759 +N_CPU = psutil.cpu_count()
  760 +TREKKER_CONFIG = {'seed_max': 1, 'step_size': 0.1, 'min_fod': 0.1, 'probe_quality': 3,
  761 + 'max_interval': 1, 'min_radius_curv': 0.8, 'probe_length': 0.4,
  762 + 'write_interval': 50, 'numb_threads': '', 'max_lenth': 200,
  763 + 'min_lenth': 20, 'max_sampling_step': 100}
... ...
invesalius/control.py
... ... @@ -32,6 +32,7 @@ import invesalius.data.mask as msk
32 32 import invesalius.data.measures as measures
33 33 import invesalius.data.slice_ as sl
34 34 import invesalius.data.surface as srf
  35 +import invesalius.data.transformations as tr
35 36 import invesalius.data.volume as volume
36 37 import invesalius.gui.dialogs as dialog
37 38 import invesalius.project as prj
... ... @@ -57,7 +58,6 @@ class Controller():
57 58 def __init__(self, frame):
58 59 self.surface_manager = srf.SurfaceManager()
59 60 self.volume = volume.Volume()
60   - self.plugin_manager = plugins.PluginManager()
61 61 self.__bind_events()
62 62 self.frame = frame
63 63 self.progress_dialog = None
... ... @@ -76,12 +76,9 @@ class Controller():
76 76  
77 77 Publisher.sendMessage('Load Preferences')
78 78  
79   - self.plugin_manager.find_plugins()
80   -
81 79 def __bind_events(self):
82 80 Publisher.subscribe(self.OnImportMedicalImages, 'Import directory')
83 81 Publisher.subscribe(self.OnImportGroup, 'Import group')
84   - Publisher.subscribe(self.OnImportFolder, 'Import folder')
85 82 Publisher.subscribe(self.OnShowDialogImportDirectory,
86 83 'Show import directory dialog')
87 84 Publisher.subscribe(self.OnShowDialogImportOtherFiles,
... ... @@ -123,8 +120,6 @@ class Controller():
123 120  
124 121 Publisher.subscribe(self.Send_affine, 'Get affine matrix')
125 122  
126   - Publisher.subscribe(self.create_project_from_matrix, 'Create project from matrix')
127   -
128 123 def SetBitmapSpacing(self, spacing):
129 124 proj = prj.Project()
130 125 proj.spacing = spacing
... ... @@ -336,6 +331,10 @@ class Controller():
336 331  
337 332 self.Slice.window_level = proj.level
338 333 self.Slice.window_width = proj.window
  334 + if proj.affine:
  335 + self.Slice.affine = np.asarray(proj.affine).reshape(4, 4)
  336 + else:
  337 + self.Slice.affine = None
339 338  
340 339 Publisher.sendMessage('Update threshold limits list',
341 340 threshold_range=proj.threshold_range)
... ... @@ -504,32 +503,7 @@ class Controller():
504 503 self.LoadProject()
505 504 Publisher.sendMessage("Enable state project", state=True)
506 505  
507   - def OnImportFolder(self, folder):
508   - Publisher.sendMessage('Begin busy cursor')
509   - folder = os.path.abspath(folder)
510   -
511   - proj = prj.Project()
512   - proj.load_from_folder(folder)
513   -
514   - self.Slice = sl.Slice()
515   - self.Slice._open_image_matrix(proj.matrix_filename,
516   - tuple(proj.matrix_shape),
517   - proj.matrix_dtype)
518   -
519   - self.Slice.window_level = proj.level
520   - self.Slice.window_width = proj.window
521   -
522   - Publisher.sendMessage('Update threshold limits list',
523   - threshold_range=proj.threshold_range)
524   -
525   - session = ses.Session()
526   - filename = proj.name+".inv3"
527   - filename = filename.replace("/", "") #Fix problem case other/Skull_DICOM
528   - dirpath = session.CreateProject(filename)
529   - self.LoadProject()
530   - Publisher.sendMessage("Enable state project", state=True)
531   -
532   - Publisher.sendMessage('End busy cursor')
  506 + #-------------------------------------------------------------------------------------
533 507  
534 508 def LoadProject(self):
535 509 proj = prj.Project()
... ... @@ -688,6 +662,9 @@ class Controller():
688 662 proj.matrix_dtype = matrix.dtype.name
689 663 proj.matrix_filename = matrix_filename
690 664  
  665 + if self.affine is not None:
  666 + proj.affine = self.affine.tolist()
  667 +
691 668 # Orientation must be CORONAL in order to as_closes_canonical and
692 669 # swap axis in img2memmap to work in a standardized way.
693 670 # TODO: Create standard import image for all acquisition orientations
... ... @@ -700,7 +677,6 @@ class Controller():
700 677 proj.level = self.Slice.window_level
701 678 proj.threshold_range = int(matrix.min()), int(matrix.max())
702 679 proj.spacing = self.Slice.spacing
703   - proj.affine = self.affine.tolist()
704 680  
705 681 ######
706 682 session = ses.Session()
... ... @@ -872,11 +848,13 @@ class Controller():
872 848 def OnOpenOtherFiles(self, filepath):
873 849 filepath = utils.decode(filepath, const.FS_ENCODE)
874 850 if not(filepath) == None:
875   - name = os.path.basename(filepath).split(".")[0]
  851 + name = filepath.rpartition('\\')[-1].split('.')
  852 +
876 853 group = oth.ReadOthers(filepath)
  854 +
877 855 if group:
878 856 matrix, matrix_filename = self.OpenOtherFiles(group)
879   - self.CreateOtherProject(name, matrix, matrix_filename)
  857 + self.CreateOtherProject(str(name[0]), matrix, matrix_filename)
880 858 self.LoadProject()
881 859 Publisher.sendMessage("Enable state project", state=True)
882 860 else:
... ... @@ -981,10 +959,10 @@ class Controller():
981 959 self.matrix, scalar_range, self.filename = image_utils.img2memmap(group)
982 960  
983 961 hdr = group.header
984   - if group.affine.any():
985   - self.affine = group.affine
986   - Publisher.sendMessage('Update affine matrix',
987   - affine=self.affine, status=True)
  962 + # if group.affine.any():
  963 + # self.affine = group.affine
  964 + # Publisher.sendMessage('Update affine matrix',
  965 + # affine=self.affine, status=True)
988 966 hdr.set_data_dtype('int16')
989 967 dims = hdr.get_zooms()
990 968 dimsf = tuple([float(s) for s in dims])
... ... @@ -1000,6 +978,18 @@ class Controller():
1000 978 self.Slice.window_level = wl
1001 979 self.Slice.window_width = ww
1002 980  
  981 + if group.affine.any():
  982 + # TODO: replace the inverse of the affine by the actual affine in the whole code
  983 + # remove scaling factor for non-unitary voxel dimensions
  984 + # self.affine = image_utils.world2invspace(affine=group.affine)
  985 + scale, shear, angs, trans, persp = tr.decompose_matrix(group.affine)
  986 + self.affine = np.linalg.inv(tr.compose_matrix(scale=None, shear=shear,
  987 + angles=angs, translate=trans, perspective=persp))
  988 + # print("repos_img: {}".format(repos_img))
  989 + self.Slice.affine = self.affine
  990 + Publisher.sendMessage('Update affine matrix',
  991 + affine=self.affine, status=True)
  992 +
1003 993 scalar_range = int(scalar_range[0]), int(scalar_range[1])
1004 994 Publisher.sendMessage('Update threshold limits list',
1005 995 threshold_range=scalar_range)
... ... @@ -1066,24 +1056,3 @@ class Controller():
1066 1056  
1067 1057 def ApplyReorientation(self):
1068 1058 self.Slice.apply_reorientation()
1069   -
1070   - def start_new_inv_instance(self, image, name, spacing, modality, orientation, window_width, window_level):
1071   - p = prj.Project()
1072   - project_folder = tempfile.mkdtemp()
1073   - p.create_project_file(name, spacing, modality, orientation, window_width, window_level, image, folder=project_folder)
1074   - err_msg = ''
1075   - try:
1076   - sp = subprocess.Popen([sys.executable, sys.argv[0], '--import-folder', project_folder],
1077   - stdout=subprocess.PIPE, stderr=subprocess.PIPE, cwd=os.getcwd())
1078   - except Exception as err:
1079   - err_msg = str(err)
1080   - else:
1081   - try:
1082   - if sp.wait(2):
1083   - err_msg = sp.stderr.read().decode('utf8')
1084   - sp.terminate()
1085   - except subprocess.TimeoutExpired:
1086   - pass
1087   -
1088   - if err_msg:
1089   - dialog.MessageBox(None, "It was not possible to launch new instance of InVesalius3 dsfa dfdsfa sdfas fdsaf asdfasf dsaa", err_msg)
... ...
invesalius/data/bases.py
1   -from math import sqrt, pi
2 1 import numpy as np
3 2 import invesalius.data.coordinates as dco
4 3 import invesalius.data.transformations as tr
... ... @@ -23,8 +22,7 @@ def angle_calculation(ap_axis, coil_axis):
23 22  
24 23 def base_creation_old(fiducials):
25 24 """
26   - Calculate the origin and matrix for coordinate system
27   - transformation.
  25 + Calculate the origin and matrix for coordinate system transformation.
28 26 q: origin of coordinate system
29 27 g1, g2, g3: orthogonal vectors of coordinate system
30 28  
... ... @@ -49,9 +47,9 @@ def base_creation_old(fiducials):
49 47  
50 48 g3 = np.cross(g2, g1)
51 49  
52   - g1 = g1/sqrt(np.dot(g1, g1))
53   - g2 = g2/sqrt(np.dot(g2, g2))
54   - g3 = g3/sqrt(np.dot(g3, g3))
  50 + g1 = g1/np.sqrt(np.dot(g1, g1))
  51 + g2 = g2/np.sqrt(np.dot(g2, g2))
  52 + g3 = g3/np.sqrt(np.dot(g3, g3))
55 53  
56 54 m = np.matrix([[g1[0], g1[1], g1[2]],
57 55 [g2[0], g2[1], g2[2]],
... ... @@ -79,7 +77,7 @@ def base_creation(fiducials):
79 77  
80 78 sub1 = p2 - p1
81 79 sub2 = p3 - p1
82   - lamb = (sub1[0]*sub2[0]+sub1[1]*sub2[1]+sub1[2]*sub2[2])/np.dot(sub1, sub1)
  80 + lamb = np.dot(sub1, sub2)/np.dot(sub1, sub1)
83 81  
84 82 q = p1 + lamb*sub1
85 83 g1 = p3 - q
... ... @@ -90,20 +88,52 @@ def base_creation(fiducials):
90 88  
91 89 g3 = np.cross(g1, g2)
92 90  
93   - g1 = g1/sqrt(np.dot(g1, g1))
94   - g2 = g2/sqrt(np.dot(g2, g2))
95   - g3 = g3/sqrt(np.dot(g3, g3))
96   -
97   - m = np.matrix([[g1[0], g2[0], g3[0]],
98   - [g1[1], g2[1], g3[1]],
99   - [g1[2], g2[2], g3[2]]])
100   -
101   - m_inv = m.I
102   -
103   - return m, q, m_inv
104   -
105   -
106   -def calculate_fre(fiducials, minv, n, q, o):
  91 + g1 = g1/np.sqrt(np.dot(g1, g1))
  92 + g2 = g2/np.sqrt(np.dot(g2, g2))
  93 + g3 = g3/np.sqrt(np.dot(g3, g3))
  94 +
  95 + m = np.zeros([3, 3])
  96 + m[:, 0] = g1/np.sqrt(np.dot(g1, g1))
  97 + m[:, 1] = g2/np.sqrt(np.dot(g2, g2))
  98 + m[:, 2] = g3/np.sqrt(np.dot(g3, g3))
  99 +
  100 + return m, q
  101 +
  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):
107 137 """
108 138 Calculate the Fiducial Registration Error for neuronavigation.
109 139  
... ... @@ -115,19 +145,29 @@ def calculate_fre(fiducials, minv, n, q, o):
115 145 :return: float number of fiducial registration error
116 146 """
117 147  
  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 +
118 158 img = np.zeros([3, 3])
119 159 dist = np.zeros([3, 1])
120 160  
121   - q1 = np.mat(q).reshape(3, 1)
122   - q2 = np.mat(o).reshape(3, 1)
  161 + q1 = q1.reshape(3, 1)
  162 + q2 = q2.reshape(3, 1)
123 163  
124   - p1 = np.mat(fiducials[3, :]).reshape(3, 1)
125   - p2 = np.mat(fiducials[4, :]).reshape(3, 1)
126   - p3 = np.mat(fiducials[5, :]).reshape(3, 1)
  164 + p1 = fiducials[3, :].reshape(3, 1)
  165 + p2 = fiducials[4, :].reshape(3, 1)
  166 + p3 = fiducials[5, :].reshape(3, 1)
127 167  
128   - img[0, :] = np.asarray((q1 + (minv * n) * (p1 - q2)).reshape(1, 3))
129   - img[1, :] = np.asarray((q1 + (minv * n) * (p2 - q2)).reshape(1, 3))
130   - img[2, :] = np.asarray((q1 + (minv * n) * (p3 - q2)).reshape(1, 3))
  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)
131 171  
132 172 dist[0] = np.sqrt(np.sum(np.power((img[0, :] - fiducials[0, :]), 2)))
133 173 dist[1] = np.sqrt(np.sum(np.power((img[1, :] - fiducials[1, :]), 2)))
... ... @@ -136,44 +176,68 @@ def calculate_fre(fiducials, minv, n, q, o):
136 176 return float(np.sqrt(np.sum(dist ** 2) / 3))
137 177  
138 178  
139   -def flip_x(point):
140   - """
141   - Flip coordinates of a vector according to X axis
142   - Coronal Images do not require this transformation - 1 tested
143   - and for this case, at navigation, the z axis is inverted
144   -
145   - It's necessary to multiply the z coordinate by (-1). Possibly
146   - because the origin of coordinate system of imagedata is
147   - located in superior left corner and the origin of VTK scene coordinate
148   - system (polygonal surface) is in the interior left corner. Second
149   - possibility is the order of slice stacking
150   -
151   - :param point: list of coordinates x, y and z
152   - :return: flipped coordinates
153   - """
154   -
155   - # TODO: check if the Flip function is related to the X or Y axis
156   -
157   - point = np.matrix(point + (0,))
158   - point[0, 2] = -point[0, 2]
159   -
160   - m_rot = np.matrix([[1.0, 0.0, 0.0, 0.0],
161   - [0.0, -1.0, 0.0, 0.0],
162   - [0.0, 0.0, -1.0, 0.0],
163   - [0.0, 0.0, 0.0, 1.0]])
164   - m_trans = np.matrix([[1.0, 0, 0, -point[0, 0]],
165   - [0.0, 1.0, 0, -point[0, 1]],
166   - [0.0, 0.0, 1.0, -point[0, 2]],
167   - [0.0, 0.0, 0.0, 1.0]])
168   - m_trans_return = np.matrix([[1.0, 0, 0, point[0, 0]],
169   - [0.0, 1.0, 0, point[0, 1]],
170   - [0.0, 0.0, 1.0, point[0, 2]],
171   - [0.0, 0.0, 0.0, 1.0]])
172   -
173   - point_rot = point*m_trans*m_rot*m_trans_return
174   - x, y, z = point_rot.tolist()[0][:3]
175   -
176   - return x, y, z
  179 +# def flip_x(point):
  180 +# """
  181 +# Flip coordinates of a vector according to X axis
  182 +# Coronal Images do not require this transformation - 1 tested
  183 +# and for this case, at navigation, the z axis is inverted
  184 +#
  185 +# It's necessary to multiply the z coordinate by (-1). Possibly
  186 +# because the origin of coordinate system of imagedata is
  187 +# located in superior left corner and the origin of VTK scene coordinate
  188 +# system (polygonal surface) is in the interior left corner. Second
  189 +# possibility is the order of slice stacking
  190 +#
  191 +# :param point: list of coordinates x, y and z
  192 +# :return: flipped coordinates
  193 +# """
  194 +#
  195 +# # TODO: check if the Flip function is related to the X or Y axis
  196 +#
  197 +# point = np.matrix(point + (0,))
  198 +# point[0, 2] = -point[0, 2]
  199 +#
  200 +# m_rot = np.matrix([[1.0, 0.0, 0.0, 0.0],
  201 +# [0.0, -1.0, 0.0, 0.0],
  202 +# [0.0, 0.0, -1.0, 0.0],
  203 +# [0.0, 0.0, 0.0, 1.0]])
  204 +# m_trans = np.matrix([[1.0, 0, 0, -point[0, 0]],
  205 +# [0.0, 1.0, 0, -point[0, 1]],
  206 +# [0.0, 0.0, 1.0, -point[0, 2]],
  207 +# [0.0, 0.0, 0.0, 1.0]])
  208 +# m_trans_return = np.matrix([[1.0, 0, 0, point[0, 0]],
  209 +# [0.0, 1.0, 0, point[0, 1]],
  210 +# [0.0, 0.0, 1.0, point[0, 2]],
  211 +# [0.0, 0.0, 0.0, 1.0]])
  212 +#
  213 +# point_rot = point*m_trans*m_rot*m_trans_return
  214 +# x, y, z = point_rot.tolist()[0][:3]
  215 +#
  216 +# return x, y, z
  217 +
  218 +
  219 +# def flip_x_m(point):
  220 +# """
  221 +# Rotate coordinates of a vector by pi around X axis in static reference frame.
  222 +#
  223 +# InVesalius also require to multiply the z coordinate by (-1). Possibly
  224 +# because the origin of coordinate system of imagedata is
  225 +# located in superior left corner and the origin of VTK scene coordinate
  226 +# system (polygonal surface) is in the interior left corner. Second
  227 +# possibility is the order of slice stacking
  228 +#
  229 +# :param point: list of coordinates x, y and z
  230 +# :return: rotated coordinates
  231 +# """
  232 +#
  233 +# point_4 = np.hstack((point, 1.)).reshape([4, 1])
  234 +# point_4[2, 0] = -point_4[2, 0]
  235 +#
  236 +# m_rot = tr.euler_matrix(pi, 0, 0)
  237 +#
  238 +# point_rot = m_rot @ point_4
  239 +#
  240 +# return point_rot[0, 0], point_rot[1, 0], point_rot[2, 0]
177 241  
178 242  
179 243 def flip_x_m(point):
... ... @@ -190,14 +254,14 @@ def flip_x_m(point):
190 254 :return: rotated coordinates
191 255 """
192 256  
193   - point_4 = np.hstack((point, 1.)).reshape([4, 1])
  257 + point_4 = np.hstack((point, 1.)).reshape(4, 1)
194 258 point_4[2, 0] = -point_4[2, 0]
195 259  
196   - m_rot = np.asmatrix(tr.euler_matrix(pi, 0, 0))
  260 + m_rot = tr.euler_matrix(np.pi, 0, 0)
197 261  
198   - point_rot = m_rot*point_4
  262 + point_rot = m_rot @ point_4
199 263  
200   - return point_rot[0, 0], point_rot[1, 0], point_rot[2, 0]
  264 + return point_rot
201 265  
202 266  
203 267 def object_registration(fiducials, orients, coord_raw, m_change):
... ... @@ -224,48 +288,48 @@ def object_registration(fiducials, orients, coord_raw, m_change):
224 288 fids_raw[ic, :] = dco.dynamic_reference_m2(coords[ic, :], coords[3, :])[:3]
225 289  
226 290 # compute initial alignment of probe fixed in the object in source frame
227   - t_s0_raw = np.asmatrix(tr.translation_matrix(coords[3, :3]))
228   - r_s0_raw = np.asmatrix(tr.euler_matrix(np.radians(coords[3, 3]), np.radians(coords[3, 4]),
229   - np.radians(coords[3, 5]), 'rzyx'))
230   - s0_raw = np.asmatrix(tr.concatenate_matrices(t_s0_raw, r_s0_raw))
  291 + t_s0_raw = tr.translation_matrix(coords[3, :3])
  292 + r_s0_raw = tr.euler_matrix(np.radians(coords[3, 3]), np.radians(coords[3, 4]),
  293 + np.radians(coords[3, 5]), 'rzyx')
  294 + s0_raw = tr.concatenate_matrices(t_s0_raw, r_s0_raw)
231 295  
232 296 # compute change of basis for object fiducials in source frame
233   - base_obj_raw, q_obj_raw, base_inv_obj_raw = base_creation(fids_raw[:3, :3])
234   - r_obj_raw = np.asmatrix(np.identity(4))
  297 + base_obj_raw, q_obj_raw = base_creation(fids_raw[:3, :3])
  298 + r_obj_raw = np.identity(4)
235 299 r_obj_raw[:3, :3] = base_obj_raw[:3, :3]
236   - t_obj_raw = np.asmatrix(tr.translation_matrix(q_obj_raw))
237   - m_obj_raw = np.asmatrix(tr.concatenate_matrices(t_obj_raw, r_obj_raw))
  300 + t_obj_raw = tr.translation_matrix(q_obj_raw)
  301 + m_obj_raw = tr.concatenate_matrices(t_obj_raw, r_obj_raw)
238 302  
239 303 for ic in range(0, 4):
240 304 if coord_raw.any():
241 305 # compute object fiducials in reference frame
242 306 fids_dyn[ic, :] = dco.dynamic_reference_m2(coords[ic, :], coord_raw[1, :])
243   - fids_dyn[ic, 2] = -fids_dyn[ic, 2]
244 307 else:
245 308 # compute object fiducials in source frame
246 309 fids_dyn[ic, :] = coords[ic, :]
  310 + fids_dyn[ic, 2] = -fids_dyn[ic, 2]
247 311  
248 312 # compute object fiducials in vtk head frame
249 313 a, b, g = np.radians(fids_dyn[ic, 3:])
250 314 T_p = tr.translation_matrix(fids_dyn[ic, :3])
251 315 R_p = tr.euler_matrix(a, b, g, 'rzyx')
252   - M_p = np.asmatrix(tr.concatenate_matrices(T_p, R_p))
253   - M_img = np.asmatrix(m_change) * M_p
  316 + M_p = tr.concatenate_matrices(T_p, R_p)
  317 + M_img = m_change @ M_p
254 318  
255 319 angles_img = np.degrees(np.asarray(tr.euler_from_matrix(M_img, 'rzyx')))
256   - coord_img = np.asarray(flip_x_m(tr.translation_from_matrix(M_img)))
  320 + coord_img = flip_x_m(tr.translation_from_matrix(M_img))
257 321  
258   - fids_img[ic, :] = np.hstack((coord_img, angles_img))
  322 + fids_img[ic, :] = np.hstack((coord_img[:3, 0], angles_img))
259 323  
260 324 # compute object base change in vtk head frame
261   - base_obj_img, q_obj_img, base_inv_obj_img = base_creation(fids_img[:3, :3])
262   - r_obj_img = np.asmatrix(np.identity(4))
  325 + base_obj_img, _ = base_creation(fids_img[:3, :3])
  326 + r_obj_img = np.identity(4)
263 327 r_obj_img[:3, :3] = base_obj_img[:3, :3]
264 328  
265 329 # compute initial alignment of probe fixed in the object in reference (or static) frame
266   - s0_trans_dyn = np.asmatrix(tr.translation_matrix(fids_dyn[3, :3]))
267   - s0_rot_dyn = np.asmatrix(tr.euler_matrix(np.radians(fids_dyn[3, 3]), np.radians(fids_dyn[3, 4]),
268   - np.radians(fids_dyn[3, 5]), 'rzyx'))
269   - s0_dyn = np.asmatrix(tr.concatenate_matrices(s0_trans_dyn, s0_rot_dyn))
  330 + s0_trans_dyn = tr.translation_matrix(fids_dyn[3, :3])
  331 + s0_rot_dyn = tr.euler_matrix(np.radians(fids_dyn[3, 3]), np.radians(fids_dyn[3, 4]),
  332 + np.radians(fids_dyn[3, 5]), 'rzyx')
  333 + s0_dyn = tr.concatenate_matrices(s0_trans_dyn, s0_rot_dyn)
270 334  
271 335 return t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img
... ...
invesalius/data/coordinates.py
... ... @@ -20,6 +20,7 @@
20 20 from math import sin, cos
21 21 import numpy as np
22 22  
  23 +import invesalius.data.bases as db
23 24 import invesalius.data.transformations as tr
24 25 import invesalius.constants as const
25 26  
... ... @@ -74,7 +75,7 @@ def PolarisCoord(trck_init, trck_id, ref_mode):
74 75 coord3 = np.hstack((trans_obj, angles_obj))
75 76  
76 77 coord = np.vstack([coord1, coord2, coord3])
77   - Publisher.sendMessage('Sensors ID', probe_id=trck.probeID, ref_id=trck.refID, obj_id=trck.objID)
  78 + # Publisher.sendMessage('Sensors ID', probe_id=trck.probeID, ref_id=trck.refID, obj_id=trck.objID)
78 79  
79 80 return coord
80 81  
... ... @@ -232,19 +233,45 @@ def DebugCoord(trk_init, trck_id, ref_mode):
232 233 :return: six coordinates x, y, z, alfa, beta and gama
233 234 """
234 235  
235   - sleep(0.05)
  236 + # Started to take a more reasonable, limited random coordinate generator based on
  237 + # the collected fiducials, but it is more complicated than this. It should account for the
  238 + # dynamic reference computation
  239 + # if trk_init:
  240 + # fiducials = trk_init[3:, :]
  241 + # fids_max = fiducials.max(axis=0)
  242 + # fids_min = fiducials.min(axis=0)
  243 + # fids_lim = np.hstack((fids_min[np.newaxis, :].T, fids_max[np.newaxis, :].T))
  244 + #
  245 + # dx = fids_max[]
  246 + # dt = [-180, 180]
  247 + #
  248 + # else:
236 249  
237   - coord1 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
238   - uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
  250 + dx = [-70, 70]
  251 + dt = [-180, 180]
239 252  
240   - coord2 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
241   - uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
  253 + coord1 = np.array([uniform(*dx), uniform(*dx), uniform(*dx),
  254 + uniform(*dt), uniform(*dt), uniform(*dt)])
  255 + coord2 = np.array([uniform(*dx), uniform(*dx), uniform(*dx),
  256 + uniform(*dt), uniform(*dt), uniform(*dt)])
  257 + coord3 = np.array([uniform(*dx), uniform(*dx), uniform(*dx),
  258 + uniform(*dt), uniform(*dt), uniform(*dt)])
  259 + coord4 = np.array([uniform(*dx), uniform(*dx), uniform(*dx),
  260 + uniform(*dt), uniform(*dt), uniform(*dt)])
242 261  
243   - coord3 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
244   - uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
  262 + sleep(0.05)
245 263  
246   - coord4 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
247   - uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
  264 + # coord1 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
  265 + # uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
  266 + #
  267 + # coord2 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
  268 + # uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
  269 + #
  270 + # coord3 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
  271 + # uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
  272 + #
  273 + # coord4 = np.array([uniform(1, 200), uniform(1, 200), uniform(1, 200),
  274 + # uniform(-180.0, 180.0), uniform(-180.0, 180.0), uniform(-180.0, 180.0)])
248 275  
249 276 Publisher.sendMessage('Sensors ID', probe_id=int(uniform(0, 5)), ref_id=int(uniform(0, 5)), obj_id=int(uniform(0, 5)))
250 277  
... ... @@ -297,19 +324,20 @@ def dynamic_reference_m(probe, reference):
297 324 :param reference: sensor two defined as reference
298 325 :return: rotated and translated coordinates
299 326 """
300   -
301 327 a, b, g = np.radians(reference[3:6])
302 328  
303   - T = tr.translation_matrix(reference[:3])
304   - R = tr.euler_matrix(a, b, g, 'rzyx')
305   - M = np.asmatrix(tr.concatenate_matrices(T, R))
306   - # M = tr.compose_matrix(angles=np.radians(reference[3:6]), translate=reference[:3])
307   - # print M
308   - probe_4 = np.vstack((np.asmatrix(probe[:3]).reshape([3, 1]), 1.))
309   - coord_rot = M.I * probe_4
310   - coord_rot = np.squeeze(np.asarray(coord_rot))
  329 + trans = tr.translation_matrix(reference[:3])
  330 + rot = tr.euler_matrix(a, b, g, 'rzyx')
  331 + affine = tr.concatenate_matrices(trans, rot)
  332 + probe_4 = np.vstack((probe[:3].reshape([3, 1]), 1.))
  333 + coord_rot = np.linalg.inv(affine) @ probe_4
  334 + # minus sign to the z coordinate
  335 + coord_rot[2, 0] = -coord_rot[2, 0]
  336 + coord_rot = coord_rot[:3, 0].tolist()
  337 + coord_rot.extend(probe[3:])
  338 +
  339 + return coord_rot
311 340  
312   - return coord_rot[0], coord_rot[1], -coord_rot[2], probe[3], probe[4], probe[5]
313 341  
314 342 def dynamic_reference_m2(probe, reference):
315 343 """
... ... @@ -331,70 +359,18 @@ def dynamic_reference_m2(probe, reference):
331 359 T_p = tr.translation_matrix(probe[:3])
332 360 R = tr.euler_matrix(a, b, g, 'rzyx')
333 361 R_p = tr.euler_matrix(a_p, b_p, g_p, 'rzyx')
334   - M = np.asmatrix(tr.concatenate_matrices(T, R))
335   - M_p = np.asmatrix(tr.concatenate_matrices(T_p, R_p))
336   - # M = tr.compose_matrix(angles=np.radians(reference[3:6]), translate=reference[:3])
337   - # print M
  362 + M = tr.concatenate_matrices(T, R)
  363 + M_p = tr.concatenate_matrices(T_p, R_p)
338 364  
339   - M_dyn = M.I * M_p
  365 + M_dyn = np.linalg.inv(M) @ M_p
340 366  
341 367 al, be, ga = tr.euler_from_matrix(M_dyn, 'rzyx')
342 368 coord_rot = tr.translation_from_matrix(M_dyn)
343 369  
344 370 coord_rot = np.squeeze(coord_rot)
345 371  
346   - # probe_4 = np.vstack((np.asmatrix(probe[:3]).reshape([3, 1]), 1.))
347   - # coord_rot_test = M.I * probe_4
348   - # coord_rot_test = np.squeeze(np.asarray(coord_rot_test))
349   - #
350   - # print "coord_rot: ", coord_rot
351   - # print "coord_rot_test: ", coord_rot_test
352   - # print "test: ", np.allclose(coord_rot, coord_rot_test[:3])
353   -
354 372 return coord_rot[0], coord_rot[1], coord_rot[2], np.degrees(al), np.degrees(be), np.degrees(ga)
355 373  
356   -# def dynamic_reference_m3(probe, reference):
357   -# """
358   -# Apply dynamic reference correction to probe coordinates. Uses the alpha, beta and gama
359   -# rotation angles of reference to rotate the probe coordinate and returns the x, y, z
360   -# difference between probe and reference. Angles sequences and equation was extracted from
361   -# Polhemus manual and Attitude matrix in Wikipedia.
362   -# General equation is:
363   -# coord = Mrot * (probe - reference)
364   -# :param probe: sensor one defined as probe
365   -# :param reference: sensor two defined as reference
366   -# :return: rotated and translated coordinates
367   -# """
368   -#
369   -# a, b, g = np.radians(reference[3:6])
370   -# a_p, b_p, g_p = np.radians(probe[3:6])
371   -#
372   -# T = tr.translation_matrix(reference[:3])
373   -# T_p = tr.translation_matrix(probe[:3])
374   -# R = tr.euler_matrix(a, b, g, 'rzyx')
375   -# R_p = tr.euler_matrix(a_p, b_p, g_p, 'rzyx')
376   -# M = np.asmatrix(tr.concatenate_matrices(T, R))
377   -# M_p = np.asmatrix(tr.concatenate_matrices(T_p, R_p))
378   -# # M = tr.compose_matrix(angles=np.radians(reference[3:6]), translate=reference[:3])
379   -# # print M
380   -#
381   -# M_dyn = M.I * M_p
382   -#
383   -# # al, be, ga = tr.euler_from_matrix(M_dyn, 'rzyx')
384   -# # coord_rot = tr.translation_from_matrix(M_dyn)
385   -# #
386   -# # coord_rot = np.squeeze(coord_rot)
387   -#
388   -# # probe_4 = np.vstack((np.asmatrix(probe[:3]).reshape([3, 1]), 1.))
389   -# # coord_rot_test = M.I * probe_4
390   -# # coord_rot_test = np.squeeze(np.asarray(coord_rot_test))
391   -# #
392   -# # print "coord_rot: ", coord_rot
393   -# # print "coord_rot_test: ", coord_rot_test
394   -# # print "test: ", np.allclose(coord_rot, coord_rot_test[:3])
395   -#
396   -# return M_dyn
397   -
398 374  
399 375 def str2float(data):
400 376 """
... ... @@ -414,3 +390,15 @@ def str2float(data):
414 390 data = [float(s) for s in data[1:len(data)]]
415 391  
416 392 return data
  393 +
  394 +
  395 +def offset_coordinate(p_old, norm_vec, offset):
  396 + """
  397 + Translate the coordinates of a point along a vector
  398 + :param p_old: (x, y, z) array with current point coordinates
  399 + :param norm_vec: (vx, vy, vz) array with normal vector coordinates
  400 + :param offset: double representing the magnitude of offset
  401 + :return: (x_new, y_new, z_new) array of offset coordinates
  402 + """
  403 + p_offset = p_old - offset * norm_vec
  404 + return p_offset
... ...
invesalius/data/coregistration.py
... ... @@ -17,347 +17,620 @@
17 17 # detalhes.
18 18 #--------------------------------------------------------------------------
19 19  
  20 +import numpy as np
  21 +import queue
20 22 import threading
21 23 from time import sleep
22 24  
23   -from numpy import asmatrix, mat, degrees, radians, identity
24   -import wx
25   -from pubsub import pub as Publisher
26   -
27 25 import invesalius.data.coordinates as dco
28 26 import invesalius.data.transformations as tr
29 27  
30   -# TODO: Optimize navigation thread. Remove the infinite loop and optimize sleep.
31   -
32 28  
33   -class CoregistrationStatic(threading.Thread):
34   - """
35   - Thread to update the coordinates with the fiducial points
36   - co-registration method while the Navigation Button is pressed.
37   - Sleep function in run method is used to avoid blocking GUI and
38   - for better real-time navigation
  29 +# TODO: Replace the use of degrees by radians in every part of the navigation pipeline
  30 +
  31 +def object_marker_to_center(coord_raw, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw):
  32 + """Translate and rotate the raw coordinate given by the tracking device to the reference system created during
  33 + the object registration.
  34 +
  35 + :param coord_raw: Coordinates returned by the tracking device
  36 + :type coord_raw: numpy.ndarray
  37 + :param obj_ref_mode:
  38 + :type obj_ref_mode: int
  39 + :param t_obj_raw:
  40 + :type t_obj_raw: numpy.ndarray
  41 + :param s0_raw:
  42 + :type s0_raw: numpy.ndarray
  43 + :param r_s0_raw: rotation transformation from marker to object basis
  44 + :type r_s0_raw: numpy.ndarray
  45 + :return: 4 x 4 numpy double array
  46 + :rtype: numpy.ndarray
39 47 """
40 48  
41   - def __init__(self, coreg_data, nav_id, trck_info):
42   - threading.Thread.__init__(self)
43   - self.coreg_data = coreg_data
44   - self.nav_id = nav_id
45   - self.trck_info = trck_info
46   - self._pause_ = False
47   - self.start()
48   -
49   - def stop(self):
50   - self._pause_ = True
51   -
52   - def run(self):
53   - # m_change = self.coreg_data[0]
54   - # obj_ref_mode = self.coreg_data[2]
55   - #
56   - # trck_init = self.trck_info[0]
57   - # trck_id = self.trck_info[1]
58   - # trck_mode = self.trck_info[2]
59   -
60   - m_change, obj_ref_mode = self.coreg_data
61   - trck_init, trck_id, trck_mode = self.trck_info
62   -
63   - while self.nav_id:
64   - coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
65   -
66   - psi, theta, phi = coord_raw[obj_ref_mode, 3:]
67   - t_probe_raw = asmatrix(tr.translation_matrix(coord_raw[obj_ref_mode, :3]))
  49 + as1, bs1, gs1 = np.radians(coord_raw[obj_ref_mode, 3:])
  50 + r_probe = tr.euler_matrix(as1, bs1, gs1, 'rzyx')
  51 + t_probe_raw = tr.translation_matrix(coord_raw[obj_ref_mode, :3])
  52 + t_offset_aux = np.linalg.inv(r_s0_raw) @ r_probe @ t_obj_raw
  53 + t_offset = np.identity(4)
  54 + t_offset[:, -1] = t_offset_aux[:, -1]
  55 + t_probe = s0_raw @ t_offset @ np.linalg.inv(s0_raw) @ t_probe_raw
  56 + m_probe = tr.concatenate_matrices(t_probe, r_probe)
68 57  
69   - t_probe_raw[2, -1] = -t_probe_raw[2, -1]
  58 + return m_probe
70 59  
71   - m_img = m_change * t_probe_raw
72 60  
73   - coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], psi, theta, phi
  61 +def object_to_reference(coord_raw, m_probe):
  62 + """Compute affine transformation matrix to the reference basis
74 63  
75   - wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
76   -
77   - # TODO: Optimize the value of sleep for each tracking device.
78   - sleep(0.175)
79   -
80   - if self._pause_:
81   - return
82   -
83   -
84   -class CoregistrationDynamic(threading.Thread):
85   - """
86   - Thread to update the coordinates with the fiducial points
87   - co-registration method while the Navigation Button is pressed.
88   - Sleep function in run method is used to avoid blocking GUI and
89   - for better real-time navigation
  64 + :param coord_raw: Coordinates returned by the tracking device
  65 + :type coord_raw: numpy.ndarray
  66 + :param m_probe: Probe coordinates
  67 + :type m_probe: numpy.ndarray
  68 + :return: 4 x 4 numpy double array
  69 + :rtype: numpy.ndarray
90 70 """
91 71  
92   - def __init__(self, coreg_data, nav_id, trck_info):
93   - threading.Thread.__init__(self)
94   - self.coreg_data = coreg_data
95   - self.nav_id = nav_id
96   - self.trck_info = trck_info
97   - self._pause_ = False
98   - self.start()
99   -
100   - def stop(self):
101   - self._pause_ = True
102   -
103   - def run(self):
104   - m_change, obj_ref_mode = self.coreg_data
105   - trck_init, trck_id, trck_mode = self.trck_info
106   -
107   - while self.nav_id:
108   - coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
109   -
110   - psi, theta, phi = radians(coord_raw[obj_ref_mode, 3:])
111   - r_probe = tr.euler_matrix(psi, theta, phi, 'rzyx')
112   - t_probe = tr.translation_matrix(coord_raw[obj_ref_mode, :3])
113   - m_probe = asmatrix(tr.concatenate_matrices(t_probe, r_probe))
114   -
115   - psi_ref, theta_ref, phi_ref = radians(coord_raw[1, 3:])
116   - r_ref = tr.euler_matrix(psi_ref, theta_ref, phi_ref, 'rzyx')
117   - t_ref = tr.translation_matrix(coord_raw[1, :3])
118   - m_ref = asmatrix(tr.concatenate_matrices(t_ref, r_ref))
119   -
120   - m_dyn = m_ref.I * m_probe
121   - m_dyn[2, -1] = -m_dyn[2, -1]
122   -
123   - m_img = m_change * m_dyn
124   -
125   - scale, shear, angles, trans, persp = tr.decompose_matrix(m_img)
126   -
127   - coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], \
128   - degrees(angles[0]), degrees(angles[1]), degrees(angles[2])
129   -
130   - wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
131   -
132   - # TODO: Optimize the value of sleep for each tracking device.
133   - sleep(0.175)
134   -
135   - if self._pause_:
136   - return
137   -
138   -
139   -class CoregistrationDynamic_old(threading.Thread):
140   - """
141   - Thread to update the coordinates with the fiducial points
142   - co-registration method while the Navigation Button is pressed.
143   - Sleep function in run method is used to avoid blocking GUI and
144   - for better real-time navigation
  72 + a, b, g = np.radians(coord_raw[1, 3:])
  73 + r_ref = tr.euler_matrix(a, b, g, 'rzyx')
  74 + t_ref = tr.translation_matrix(coord_raw[1, :3])
  75 + m_ref = tr.concatenate_matrices(t_ref, r_ref)
  76 +
  77 + m_dyn = np.linalg.inv(m_ref) @ m_probe
  78 + return m_dyn
  79 +
  80 +
  81 +def tracker_to_image(m_change, m_probe_ref, r_obj_img, m_obj_raw, s0_dyn):
  82 + """Compute affine transformation matrix to the reference basis
  83 +
  84 + :param m_change: Corregistration transformation obtained from fiducials
  85 + :type m_change: numpy.ndarray
  86 + :param m_probe_ref: Object or probe in reference coordinate system
  87 + :type m_probe_ref: numpy.ndarray
  88 + :param r_obj_img: Object coordinate system in image space (3d model)
  89 + :type r_obj_img: numpy.ndarray
  90 + :param m_obj_raw: Object basis in raw coordinates from tacker
  91 + :type m_obj_raw: numpy.ndarray
  92 + :param s0_dyn:
  93 + :type s0_dyn: numpy.ndarray
  94 + :return: 4 x 4 numpy double array
  95 + :rtype: numpy.ndarray
145 96 """
146 97  
147   - def __init__(self, bases, nav_id, trck_info):
148   - threading.Thread.__init__(self)
149   - self.bases = bases
150   - self.nav_id = nav_id
  98 + m_img = m_change @ m_probe_ref
  99 + r_obj = r_obj_img @ np.linalg.inv(m_obj_raw) @ np.linalg.inv(s0_dyn) @ m_probe_ref @ m_obj_raw
  100 + m_img[:3, :3] = r_obj[:3, :3]
  101 + return m_img
  102 +
  103 +
  104 +def corregistrate_object_dynamic(inp, coord_raw, ref_mode_id):
  105 +
  106 + m_change, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img = inp
  107 +
  108 + # transform raw marker coordinate to object center
  109 + m_probe = object_marker_to_center(coord_raw, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw)
  110 + # transform object center to reference marker if specified as dynamic reference
  111 + if ref_mode_id:
  112 + m_probe_ref = object_to_reference(coord_raw, m_probe)
  113 + else:
  114 + m_probe_ref = m_probe
  115 + # invert y coordinate
  116 + m_probe_ref[2, -1] = -m_probe_ref[2, -1]
  117 + # corregistrate from tracker to image space
  118 + m_img = tracker_to_image(m_change, m_probe_ref, r_obj_img, m_obj_raw, s0_dyn)
  119 + # compute rotation angles
  120 + _, _, angles, _, _ = tr.decompose_matrix(m_img)
  121 + # create output coordiante list
  122 + coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], \
  123 + np.degrees(angles[0]), np.degrees(angles[1]), np.degrees(angles[2])
  124 +
  125 + return coord, m_img
  126 +
  127 +
  128 +def compute_marker_transformation(coord_raw, obj_ref_mode):
  129 + psi, theta, phi = np.radians(coord_raw[obj_ref_mode, 3:])
  130 + r_probe = tr.euler_matrix(psi, theta, phi, 'rzyx')
  131 + t_probe = tr.translation_matrix(coord_raw[obj_ref_mode, :3])
  132 + m_probe = tr.concatenate_matrices(t_probe, r_probe)
  133 + return m_probe
  134 +
  135 +
  136 +def corregistrate_dynamic(inp, coord_raw, ref_mode_id):
  137 +
  138 + m_change, obj_ref_mode = inp
  139 +
  140 + # transform raw marker coordinate to object center
  141 + m_probe = compute_marker_transformation(coord_raw, obj_ref_mode)
  142 + # transform object center to reference marker if specified as dynamic reference
  143 + if ref_mode_id:
  144 + m_ref = compute_marker_transformation(coord_raw, 1)
  145 + m_probe_ref = np.linalg.inv(m_ref) @ m_probe
  146 + else:
  147 + m_probe_ref = m_probe
  148 +
  149 + # invert y coordinate
  150 + m_probe_ref[2, -1] = -m_probe_ref[2, -1]
  151 + # corregistrate from tracker to image space
  152 + m_img = m_change @ m_probe_ref
  153 + # compute rotation angles
  154 + _, _, angles, _, _ = tr.decompose_matrix(m_img)
  155 + # create output coordiante list
  156 + coord = m_img[0, -1], m_img[1, -1], m_img[2, -1],\
  157 + np.degrees(angles[0]), np.degrees(angles[1]), np.degrees(angles[2])
  158 +
  159 + return coord, m_img
  160 +
  161 +
  162 +class CoordinateCorregistrate(threading.Thread):
  163 + def __init__(self, ref_mode_id, trck_info, coreg_data, coord_queue, view_tracts, coord_tracts_queue, event, sle):
  164 + threading.Thread.__init__(self, name='CoordCoregObject')
  165 + self.ref_mode_id = ref_mode_id
151 166 self.trck_info = trck_info
152   - self._pause_ = False
153   - self.start()
154   -
155   - def stop(self):
156   - self._pause_ = True
157   -
158   - def run(self):
159   - m_inv = self.bases[0]
160   - n = self.bases[1]
161   - q1 = self.bases[2]
162   - q2 = self.bases[3]
163   - trck_init = self.trck_info[0]
164   - trck_id = self.trck_info[1]
165   - trck_mode = self.trck_info[2]
166   -
167   - while self.nav_id:
168   - # trck_coord, probe, reference = dco.GetCoordinates(trck_init, trck_id, trck_mode)
169   - coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
170   -
171   - trck_coord = dco.dynamic_reference(coord_raw[0, :], coord_raw[1, :])
172   -
173   - trck_xyz = mat([[trck_coord[0]], [trck_coord[1]], [trck_coord[2]]])
174   - img = q1 + (m_inv * n) * (trck_xyz - q2)
175   -
176   - coord = (float(img[0]), float(img[1]), float(img[2]), trck_coord[3],
177   - trck_coord[4], trck_coord[5])
178   - angles = coord_raw[0, 3:6]
179   -
180   - # Tried several combinations and different locations to send the messages,
181   - # however only this one does not block the GUI during navigation.
182   - wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=None, position=coord)
183   - wx.CallAfter(Publisher.sendMessage, 'Set camera in volume', coord)
184   - wx.CallAfter(Publisher.sendMessage, 'Update tracker angles', angles)
185   -
186   - # TODO: Optimize the value of sleep for each tracking device.
187   - # Debug tracker is not working with 0.175 so changed to 0.2
188   - # However, 0.2 is too low update frequency ~5 Hz. Need optimization URGENTLY.
189   - # sleep(.3)
190   - sleep(0.175)
191   -
192   - if self._pause_:
193   - return
194   -
195   -
196   -class CoregistrationObjectStatic(threading.Thread):
197   - """
198   - Thread to update the coordinates with the fiducial points
199   - co-registration method while the Navigation Button is pressed.
200   - Sleep function in run method is used to avoid blocking GUI and
201   - for better real-time navigation
202   - """
203   -
204   - def __init__(self, coreg_data, nav_id, trck_info):
205   - threading.Thread.__init__(self)
206 167 self.coreg_data = coreg_data
207   - self.nav_id = nav_id
208   - self.trck_info = trck_info
209   - self._pause_ = False
210   - self.start()
211   -
212   - def stop(self):
213   - self._pause_ = True
  168 + self.coord_queue = coord_queue
  169 + self.view_tracts = view_tracts
  170 + self.coord_tracts_queue = coord_tracts_queue
  171 + self.event = event
  172 + self.sle = sle
214 173  
215 174 def run(self):
216   - # m_change = self.coreg_data[0]
217   - # t_obj_raw = self.coreg_data[1]
218   - # s0_raw = self.coreg_data[2]
219   - # r_s0_raw = self.coreg_data[3]
220   - # s0_dyn = self.coreg_data[4]
221   - # m_obj_raw = self.coreg_data[5]
222   - # r_obj_img = self.coreg_data[6]
223   - # obj_ref_mode = self.coreg_data[7]
224   - #
225   - # trck_init = self.trck_info[0]
226   - # trck_id = self.trck_info[1]
227   - # trck_mode = self.trck_info[2]
228   -
229   - m_change, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img = self.coreg_data
230   - trck_init, trck_id, trck_mode = self.trck_info
231   -
232   - while self.nav_id:
233   - coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
234   -
235   - as1, bs1, gs1 = radians(coord_raw[obj_ref_mode, 3:])
236   - r_probe = asmatrix(tr.euler_matrix(as1, bs1, gs1, 'rzyx'))
237   - t_probe_raw = asmatrix(tr.translation_matrix(coord_raw[obj_ref_mode, :3]))
238   - t_offset_aux = r_s0_raw.I * r_probe * t_obj_raw
239   - t_offset = asmatrix(identity(4))
240   - t_offset[:, -1] = t_offset_aux[:, -1]
241   - t_probe = s0_raw * t_offset * s0_raw.I * t_probe_raw
242   - m_probe = asmatrix(tr.concatenate_matrices(t_probe, r_probe))
243   -
244   - m_probe[2, -1] = -m_probe[2, -1]
245   -
246   - m_img = m_change * m_probe
247   - r_obj = r_obj_img * m_obj_raw.I * s0_dyn.I * m_probe * m_obj_raw
248   -
249   - m_img[:3, :3] = r_obj[:3, :3]
250   -
251   - scale, shear, angles, trans, persp = tr.decompose_matrix(m_img)
252   -
253   - coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], \
254   - degrees(angles[0]), degrees(angles[1]), degrees(angles[2])
255   -
256   - wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
257   - wx.CallAfter(Publisher.sendMessage, 'Update object matrix', m_img=m_img, coord=coord)
258   -
259   - # TODO: Optimize the value of sleep for each tracking device.
260   - sleep(0.175)
261   -
262   - # Debug tracker is not working with 0.175 so changed to 0.2
263   - # However, 0.2 is too low update frequency ~5 Hz. Need optimization URGENTLY.
264   - # sleep(.3)
265   -
266   - # partially working for translate and offset,
267   - # but offset is kept always in same axis, have to fix for rotation
268   - # M_dyn = M_reference.I * T_stylus
269   - # M_dyn[2, -1] = -M_dyn[2, -1]
270   - # M_dyn_ch = M_change * M_dyn
271   - # ddd = M_dyn_ch[0, -1], M_dyn_ch[1, -1], M_dyn_ch[2, -1]
272   - # M_dyn_ch[:3, -1] = asmatrix(db.flip_x_m(ddd)).reshape([3, 1])
273   - # M_final = S0 * M_obj_trans_0 * S0.I * M_dyn_ch
274   -
275   - # this works for static reference object rotation
276   - # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_raw.I * R_stylus * M_obj_rot_raw
277   - # this works for dynamic reference in rotation but not in translation
278   - # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_dyn.I * R_reference.I * R_stylus * M_obj_rot_raw
279   -
280   - if self._pause_:
281   - return
282   -
283   -
284   -class CoregistrationObjectDynamic(threading.Thread):
285   - """
286   - Thread to update the coordinates with the fiducial points
287   - co-registration method while the Navigation Button is pressed.
288   - Sleep function in run method is used to avoid blocking GUI and
289   - for better real-time navigation
290   - """
291   -
292   - def __init__(self, coreg_data, nav_id, trck_info):
293   - threading.Thread.__init__(self)
294   - self.coreg_data = coreg_data
295   - self.nav_id = nav_id
  175 + trck_info = self.trck_info
  176 + coreg_data = self.coreg_data
  177 + view_obj = 1
  178 +
  179 + trck_init, trck_id, trck_mode = trck_info
  180 + # print('CoordCoreg: event {}'.format(self.event.is_set()))
  181 + while not self.event.is_set():
  182 + try:
  183 + # print(f"Set the coordinate")
  184 + 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_flip = m_img.copy()
  187 + m_img_flip[1, -1] = -m_img_flip[1, -1]
  188 + # self.pipeline.set_message(m_img_flip)
  189 + self.coord_queue.put_nowait([coord, m_img, view_obj])
  190 + # print('CoordCoreg: put {}'.format(count))
  191 + # count += 1
  192 +
  193 + if self.view_tracts:
  194 + self.coord_tracts_queue.put_nowait(m_img_flip)
  195 +
  196 + # The sleep has to be in both threads
  197 + sleep(self.sle)
  198 + except queue.Full:
  199 + pass
  200 +
  201 +
  202 +class CoordinateCorregistrateNoObject(threading.Thread):
  203 + def __init__(self, ref_mode_id, trck_info, coreg_data, coord_queue, view_tracts, coord_tracts_queue, event, sle):
  204 + threading.Thread.__init__(self, name='CoordCoregNoObject')
  205 + self.ref_mode_id = ref_mode_id
296 206 self.trck_info = trck_info
297   - self._pause_ = False
298   - self.start()
299   -
300   - def stop(self):
301   - self._pause_ = True
  207 + self.coreg_data = coreg_data
  208 + self.coord_queue = coord_queue
  209 + self.view_tracts = view_tracts
  210 + self.coord_tracts_queue = coord_tracts_queue
  211 + self.event = event
  212 + self.sle = sle
302 213  
303 214 def run(self):
  215 + trck_info = self.trck_info
  216 + coreg_data = self.coreg_data
  217 + view_obj = 0
  218 +
  219 + trck_init, trck_id, trck_mode = trck_info
  220 + # print('CoordCoreg: event {}'.format(self.event.is_set()))
  221 + while not self.event.is_set():
  222 + try:
  223 + # print(f"Set the coordinate")
  224 + coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
  225 + coord, m_img = corregistrate_dynamic(coreg_data, coord_raw, self.ref_mode_id)
  226 + # print("Coord: ", coord)
  227 + m_img_flip = m_img.copy()
  228 + m_img_flip[1, -1] = -m_img_flip[1, -1]
  229 + self.coord_queue.put_nowait([coord, m_img, view_obj])
  230 +
  231 + if self.view_tracts:
  232 + self.coord_tracts_queue.put_nowait(m_img_flip)
  233 +
  234 + # The sleep has to be in both threads
  235 + sleep(self.sle)
  236 + except queue.Full:
  237 + pass
  238 +
  239 +
  240 +# class CoregistrationStatic(threading.Thread):
  241 +# """
  242 +# Thread to update the coordinates with the fiducial points
  243 +# co-registration method while the Navigation Button is pressed.
  244 +# Sleep function in run method is used to avoid blocking GUI and
  245 +# for better real-time navigation
  246 +# """
  247 +#
  248 +# def __init__(self, coreg_data, nav_id, trck_info):
  249 +# threading.Thread.__init__(self)
  250 +# self.coreg_data = coreg_data
  251 +# self.nav_id = nav_id
  252 +# self.trck_info = trck_info
  253 +# self._pause_ = False
  254 +# self.start()
  255 +#
  256 +# def stop(self):
  257 +# self._pause_ = True
  258 +#
  259 +# def run(self):
  260 +# # m_change = self.coreg_data[0]
  261 +# # obj_ref_mode = self.coreg_data[2]
  262 +# #
  263 +# # trck_init = self.trck_info[0]
  264 +# # trck_id = self.trck_info[1]
  265 +# # trck_mode = self.trck_info[2]
  266 +#
  267 +# m_change, obj_ref_mode = self.coreg_data
  268 +# trck_init, trck_id, trck_mode = self.trck_info
  269 +#
  270 +# while self.nav_id:
  271 +# coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
  272 +#
  273 +# psi, theta, phi = coord_raw[obj_ref_mode, 3:]
  274 +# t_probe_raw = asmatrix(tr.translation_matrix(coord_raw[obj_ref_mode, :3]))
  275 +#
  276 +# t_probe_raw[2, -1] = -t_probe_raw[2, -1]
  277 +#
  278 +# m_img = m_change * t_probe_raw
  279 +#
  280 +# coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], psi, theta, phi
  281 +#
  282 +# wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
  283 +#
  284 +# # TODO: Optimize the value of sleep for each tracking device.
  285 +# sleep(0.175)
  286 +#
  287 +# if self._pause_:
  288 +# return
  289 +#
  290 +#
  291 +# class CoregistrationDynamic(threading.Thread):
  292 +# """
  293 +# Thread to update the coordinates with the fiducial points
  294 +# co-registration method while the Navigation Button is pressed.
  295 +# Sleep function in run method is used to avoid blocking GUI and
  296 +# for better real-time navigation
  297 +# """
  298 +#
  299 +# def __init__(self, coreg_data, nav_id, trck_info):
  300 +# threading.Thread.__init__(self)
  301 +# self.coreg_data = coreg_data
  302 +# self.nav_id = nav_id
  303 +# self.trck_info = trck_info
  304 +# # self.tracts_info = tracts_info
  305 +# # self.tracts = None
  306 +# self._pause_ = False
  307 +# # self.start()
  308 +#
  309 +# def stop(self):
  310 +# # self.tracts.stop()
  311 +# self._pause_ = True
  312 +#
  313 +# def run(self):
  314 +# m_change, obj_ref_mode = self.coreg_data
  315 +# trck_init, trck_id, trck_mode = self.trck_info
  316 +# # seed, tracker, affine, affine_vtk = self.tracts_info
  317 +#
  318 +# while self.nav_id:
  319 +# coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
  320 +#
  321 +# psi, theta, phi = radians(coord_raw[obj_ref_mode, 3:])
  322 +# r_probe = tr.euler_matrix(psi, theta, phi, 'rzyx')
  323 +# t_probe = tr.translation_matrix(coord_raw[obj_ref_mode, :3])
  324 +# m_probe = asmatrix(tr.concatenate_matrices(t_probe, r_probe))
  325 +#
  326 +# psi_ref, theta_ref, phi_ref = radians(coord_raw[1, 3:])
  327 +# r_ref = tr.euler_matrix(psi_ref, theta_ref, phi_ref, 'rzyx')
  328 +# t_ref = tr.translation_matrix(coord_raw[1, :3])
  329 +# m_ref = asmatrix(tr.concatenate_matrices(t_ref, r_ref))
  330 +#
  331 +# m_dyn = m_ref.I * m_probe
  332 +# m_dyn[2, -1] = -m_dyn[2, -1]
  333 +#
  334 +# m_img = m_change * m_dyn
  335 +#
  336 +# scale, shear, angles, trans, persp = tr.decompose_matrix(m_img)
  337 +#
  338 +# coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], \
  339 +# degrees(angles[0]), degrees(angles[1]), degrees(angles[2])
  340 +#
  341 +# # pos_world_aux = np.ones([4, 1])
  342 +# # pos_world_aux[:3, -1] = db.flip_x((m_img[0, -1], m_img[1, -1], m_img[2, -1]))[:3]
  343 +# # pos_world = np.linalg.inv(affine) @ pos_world_aux
  344 +# # seed_aux = pos_world.reshape([1, 4])[0, :3]
  345 +# # seed = seed_aux[np.newaxis, :]
  346 +# #
  347 +# # self.tracts = dtr.compute_tracts(tracker, seed, affine_vtk, True)
  348 +#
  349 +# # wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
  350 +# wx.CallAfter(Publisher.sendMessage, 'Update cross position', arg=m_img, position=coord)
  351 +#
  352 +# # TODO: Optimize the value of sleep for each tracking device.
  353 +# sleep(3.175)
  354 +#
  355 +# if self._pause_:
  356 +# return
  357 +#
  358 +#
  359 +# class CoregistrationDynamic_old(threading.Thread):
  360 +# """
  361 +# Thread to update the coordinates with the fiducial points
  362 +# co-registration method while the Navigation Button is pressed.
  363 +# Sleep function in run method is used to avoid blocking GUI and
  364 +# for better real-time navigation
  365 +# """
  366 +#
  367 +# def __init__(self, bases, nav_id, trck_info):
  368 +# threading.Thread.__init__(self)
  369 +# self.bases = bases
  370 +# self.nav_id = nav_id
  371 +# self.trck_info = trck_info
  372 +# self._pause_ = False
  373 +# self.start()
  374 +#
  375 +# def stop(self):
  376 +# self._pause_ = True
  377 +#
  378 +# def run(self):
  379 +# m_inv = self.bases[0]
  380 +# n = self.bases[1]
  381 +# q1 = self.bases[2]
  382 +# q2 = self.bases[3]
  383 +# trck_init = self.trck_info[0]
  384 +# trck_id = self.trck_info[1]
  385 +# trck_mode = self.trck_info[2]
  386 +#
  387 +# while self.nav_id:
  388 +# # trck_coord, probe, reference = dco.GetCoordinates(trck_init, trck_id, trck_mode)
  389 +# coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
  390 +#
  391 +# trck_coord = dco.dynamic_reference(coord_raw[0, :], coord_raw[1, :])
  392 +#
  393 +# trck_xyz = mat([[trck_coord[0]], [trck_coord[1]], [trck_coord[2]]])
  394 +# img = q1 + (m_inv * n) * (trck_xyz - q2)
  395 +#
  396 +# coord = (float(img[0]), float(img[1]), float(img[2]), trck_coord[3],
  397 +# trck_coord[4], trck_coord[5])
  398 +# angles = coord_raw[0, 3:6]
  399 +#
  400 +# # Tried several combinations and different locations to send the messages,
  401 +# # however only this one does not block the GUI during navigation.
  402 +# wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=None, position=coord)
  403 +# wx.CallAfter(Publisher.sendMessage, 'Set camera in volume', coord)
  404 +# wx.CallAfter(Publisher.sendMessage, 'Update tracker angles', angles)
  405 +#
  406 +# # TODO: Optimize the value of sleep for each tracking device.
  407 +# # Debug tracker is not working with 0.175 so changed to 0.2
  408 +# # However, 0.2 is too low update frequency ~5 Hz. Need optimization URGENTLY.
  409 +# # sleep(.3)
  410 +# sleep(0.175)
  411 +#
  412 +# if self._pause_:
  413 +# return
  414 +#
  415 +#
  416 +# class CoregistrationObjectStatic(threading.Thread):
  417 +# """
  418 +# Thread to update the coordinates with the fiducial points
  419 +# co-registration method while the Navigation Button is pressed.
  420 +# Sleep function in run method is used to avoid blocking GUI and
  421 +# for better real-time navigation
  422 +# """
  423 +#
  424 +# def __init__(self, coreg_data, nav_id, trck_info):
  425 +# threading.Thread.__init__(self)
  426 +# self.coreg_data = coreg_data
  427 +# self.nav_id = nav_id
  428 +# self.trck_info = trck_info
  429 +# self._pause_ = False
  430 +# self.start()
  431 +#
  432 +# def stop(self):
  433 +# self._pause_ = True
  434 +#
  435 +# def run(self):
  436 +# # m_change = self.coreg_data[0]
  437 +# # t_obj_raw = self.coreg_data[1]
  438 +# # s0_raw = self.coreg_data[2]
  439 +# # r_s0_raw = self.coreg_data[3]
  440 +# # s0_dyn = self.coreg_data[4]
  441 +# # m_obj_raw = self.coreg_data[5]
  442 +# # r_obj_img = self.coreg_data[6]
  443 +# # obj_ref_mode = self.coreg_data[7]
  444 +# #
  445 +# # trck_init = self.trck_info[0]
  446 +# # trck_id = self.trck_info[1]
  447 +# # trck_mode = self.trck_info[2]
  448 +#
  449 +# m_change, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img = self.coreg_data
  450 +# trck_init, trck_id, trck_mode = self.trck_info
  451 +#
  452 +# while self.nav_id:
  453 +# coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
  454 +#
  455 +# as1, bs1, gs1 = radians(coord_raw[obj_ref_mode, 3:])
  456 +# r_probe = asmatrix(tr.euler_matrix(as1, bs1, gs1, 'rzyx'))
  457 +# t_probe_raw = asmatrix(tr.translation_matrix(coord_raw[obj_ref_mode, :3]))
  458 +# t_offset_aux = r_s0_raw.I * r_probe * t_obj_raw
  459 +# t_offset = asmatrix(identity(4))
  460 +# t_offset[:, -1] = t_offset_aux[:, -1]
  461 +# t_probe = s0_raw * t_offset * s0_raw.I * t_probe_raw
  462 +# m_probe = asmatrix(tr.concatenate_matrices(t_probe, r_probe))
  463 +#
  464 +# m_probe[2, -1] = -m_probe[2, -1]
  465 +#
  466 +# m_img = m_change * m_probe
  467 +# r_obj = r_obj_img * m_obj_raw.I * s0_dyn.I * m_probe * m_obj_raw
  468 +#
  469 +# m_img[:3, :3] = r_obj[:3, :3]
  470 +#
  471 +# scale, shear, angles, trans, persp = tr.decompose_matrix(m_img)
  472 +#
  473 +# coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], \
  474 +# degrees(angles[0]), degrees(angles[1]), degrees(angles[2])
  475 +#
  476 +# wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
  477 +# wx.CallAfter(Publisher.sendMessage, 'Update object matrix', m_img=m_img, coord=coord)
  478 +#
  479 +# # TODO: Optimize the value of sleep for each tracking device.
  480 +# sleep(0.175)
  481 +#
  482 +# # Debug tracker is not working with 0.175 so changed to 0.2
  483 +# # However, 0.2 is too low update frequency ~5 Hz. Need optimization URGENTLY.
  484 +# # sleep(.3)
  485 +#
  486 +# # partially working for translate and offset,
  487 +# # but offset is kept always in same axis, have to fix for rotation
  488 +# # M_dyn = M_reference.I * T_stylus
  489 +# # M_dyn[2, -1] = -M_dyn[2, -1]
  490 +# # M_dyn_ch = M_change * M_dyn
  491 +# # ddd = M_dyn_ch[0, -1], M_dyn_ch[1, -1], M_dyn_ch[2, -1]
  492 +# # M_dyn_ch[:3, -1] = asmatrix(db.flip_x_m(ddd)).reshape([3, 1])
  493 +# # M_final = S0 * M_obj_trans_0 * S0.I * M_dyn_ch
  494 +#
  495 +# # this works for static reference object rotation
  496 +# # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_raw.I * R_stylus * M_obj_rot_raw
  497 +# # this works for dynamic reference in rotation but not in translation
  498 +# # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_dyn.I * R_reference.I * R_stylus * M_obj_rot_raw
  499 +#
  500 +# if self._pause_:
  501 +# return
  502 +#
  503 +#
  504 +# class CoregistrationObjectDynamic(threading.Thread):
  505 +# """
  506 +# Thread to update the coordinates with the fiducial points
  507 +# co-registration method while the Navigation Button is pressed.
  508 +# Sleep function in run method is used to avoid blocking GUI and
  509 +# for better real-time navigation
  510 +# """
  511 +#
  512 +# def __init__(self, coreg_data, nav_id, trck_info, tracts_info):
  513 +# threading.Thread.__init__(self)
  514 +# self.coreg_data = coreg_data
  515 +# self.nav_id = nav_id
  516 +# self.trck_info = trck_info
  517 +# # self.tracts_info = tracts_info
  518 +# # self.tracts = None
  519 +# self._pause_ = False
  520 +# self.start()
  521 +#
  522 +# def stop(self):
  523 +# # self.tracts.stop()
  524 +# self._pause_ = True
  525 +#
  526 +# def run(self):
  527 +#
  528 +# m_change, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img = self.coreg_data
  529 +# trck_init, trck_id, trck_mode = self.trck_info
  530 +# # seed, tracker, affine, affine_vtk = self.tracts_info
  531 +#
  532 +# while self.nav_id:
  533 +# coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
  534 +#
  535 +# as1, bs1, gs1 = radians(coord_raw[obj_ref_mode, 3:])
  536 +# r_probe = asmatrix(tr.euler_matrix(as1, bs1, gs1, 'rzyx'))
  537 +# t_probe_raw = asmatrix(tr.translation_matrix(coord_raw[obj_ref_mode, :3]))
  538 +# t_offset_aux = r_s0_raw.I * r_probe * t_obj_raw
  539 +# t_offset = asmatrix(identity(4))
  540 +# t_offset[:, -1] = t_offset_aux[:, -1]
  541 +# t_probe = s0_raw * t_offset * s0_raw.I * t_probe_raw
  542 +# m_probe = asmatrix(tr.concatenate_matrices(t_probe, r_probe))
  543 +#
  544 +# a, b, g = radians(coord_raw[1, 3:])
  545 +# r_ref = tr.euler_matrix(a, b, g, 'rzyx')
  546 +# t_ref = tr.translation_matrix(coord_raw[1, :3])
  547 +# m_ref = asmatrix(tr.concatenate_matrices(t_ref, r_ref))
  548 +#
  549 +# m_dyn = m_ref.I * m_probe
  550 +# m_dyn[2, -1] = -m_dyn[2, -1]
  551 +#
  552 +# m_img = m_change * m_dyn
  553 +# r_obj = r_obj_img * m_obj_raw.I * s0_dyn.I * m_dyn * m_obj_raw
  554 +#
  555 +# m_img[:3, :3] = r_obj[:3, :3]
  556 +#
  557 +# scale, shear, angles, trans, persp = tr.decompose_matrix(m_img)
  558 +#
  559 +# coord = m_img[0, -1], m_img[1, -1], m_img[2, -1],\
  560 +# degrees(angles[0]), degrees(angles[1]), degrees(angles[2])
  561 +#
  562 +# # norm_vec = m_img[:3, 2].reshape([1, 3]).tolist()
  563 +# # p0 = m_img[:3, -1].reshape([1, 3]).tolist()
  564 +# # p2 = [x - 30 * y for x, y in zip(p0[0], norm_vec[0])]
  565 +# # m_tract = m_img.copy()
  566 +# # m_tract[:3, -1] = np.reshape(np.asarray(p2)[np.newaxis, :], [3, 1])
  567 +#
  568 +# # pos_world_aux = np.ones([4, 1])
  569 +# # pos_world_aux[:3, -1] = db.flip_x((p2[0], p2[1], p2[2]))[:3]
  570 +# # pos_world = np.linalg.inv(affine) @ pos_world_aux
  571 +# # seed_aux = pos_world.reshape([1, 4])[0, :3]
  572 +# # seed = seed_aux[np.newaxis, :]
  573 +#
  574 +# # self.tracts = dtr.compute_tracts(tracker, seed, affine_vtk, True)
  575 +#
  576 +# # wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
  577 +# wx.CallAfter(Publisher.sendMessage, 'Update cross position', arg=m_img, position=coord)
  578 +# wx.CallAfter(Publisher.sendMessage, 'Update object matrix', m_img=m_img, coord=coord)
  579 +#
  580 +# # TODO: Optimize the value of sleep for each tracking device.
  581 +# #sleep(2.175)
  582 +# sleep(0.175)
  583 +#
  584 +# # Debug tracker is not working with 0.175 so changed to 0.2
  585 +# # However, 0.2 is too low update frequency ~5 Hz. Need optimization URGENTLY.
  586 +# # sleep(.3)
  587 +#
  588 +# # partially working for translate and offset,
  589 +# # but offset is kept always in same axis, have to fix for rotation
  590 +# # M_dyn = M_reference.I * T_stylus
  591 +# # M_dyn[2, -1] = -M_dyn[2, -1]
  592 +# # M_dyn_ch = M_change * M_dyn
  593 +# # ddd = M_dyn_ch[0, -1], M_dyn_ch[1, -1], M_dyn_ch[2, -1]
  594 +# # M_dyn_ch[:3, -1] = asmatrix(db.flip_x_m(ddd)).reshape([3, 1])
  595 +# # M_final = S0 * M_obj_trans_0 * S0.I * M_dyn_ch
  596 +#
  597 +# # this works for static reference object rotation
  598 +# # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_raw.I * R_stylus * M_obj_rot_raw
  599 +# # this works for dynamic reference in rotation but not in translation
  600 +# # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_dyn.I * R_reference.I * R_stylus * M_obj_rot_raw
  601 +#
  602 +# if self._pause_:
  603 +# return
  604 +#
  605 +#
  606 +# def corregistrate_object(inp, coord_raw):
  607 +# m_change, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img = inp
  608 +# as1, bs1, gs1 = radians(coord_raw[obj_ref_mode, 3:])
  609 +# r_probe = tr.euler_matrix(as1, bs1, gs1, 'rzyx')
  610 +# t_probe_raw = tr.translation_matrix(coord_raw[obj_ref_mode, :3])
  611 +# t_offset_aux = np.linalg.inv(r_s0_raw) @ r_probe @ t_obj_raw
  612 +# t_offset = identity(4)
  613 +# t_offset[:, -1] = t_offset_aux[:, -1]
  614 +# t_probe = s0_raw @ t_offset @ np.linalg.inv(s0_raw) @ t_probe_raw
  615 +# m_probe = tr.concatenate_matrices(t_probe, r_probe)
  616 +#
  617 +# a, b, g = radians(coord_raw[1, 3:])
  618 +# r_ref = tr.euler_matrix(a, b, g, 'rzyx')
  619 +# t_ref = tr.translation_matrix(coord_raw[1, :3])
  620 +# m_ref = tr.concatenate_matrices(t_ref, r_ref)
  621 +#
  622 +# m_dyn = np.linalg.inv(m_ref) @ m_probe
  623 +# m_dyn[2, -1] = -m_dyn[2, -1]
  624 +#
  625 +# m_img = m_change @ m_dyn
  626 +# r_obj = r_obj_img @ np.linalg.inv(m_obj_raw) @ np.linalg.inv(s0_dyn) @ m_dyn @ m_obj_raw
  627 +#
  628 +# m_img[:3, :3] = r_obj[:3, :3]
  629 +#
  630 +# scale, shear, angles, trans, persp = tr.decompose_matrix(m_img)
  631 +#
  632 +# coord = m_img[0, -1], m_img[1, -1], m_img[2, -1], \
  633 +# degrees(angles[0]), degrees(angles[1]), degrees(angles[2])
  634 +#
  635 +# return coord, m_img
304 636  
305   - m_change, obj_ref_mode, t_obj_raw, s0_raw, r_s0_raw, s0_dyn, m_obj_raw, r_obj_img = self.coreg_data
306   - trck_init, trck_id, trck_mode = self.trck_info
307   -
308   - while self.nav_id:
309   - coord_raw = dco.GetCoordinates(trck_init, trck_id, trck_mode)
310   -
311   - as1, bs1, gs1 = radians(coord_raw[obj_ref_mode, 3:])
312   - r_probe = asmatrix(tr.euler_matrix(as1, bs1, gs1, 'rzyx'))
313   - t_probe_raw = asmatrix(tr.translation_matrix(coord_raw[obj_ref_mode, :3]))
314   - t_offset_aux = r_s0_raw.I * r_probe * t_obj_raw
315   - t_offset = asmatrix(identity(4))
316   - t_offset[:, -1] = t_offset_aux[:, -1]
317   - t_probe = s0_raw * t_offset * s0_raw.I * t_probe_raw
318   - m_probe = asmatrix(tr.concatenate_matrices(t_probe, r_probe))
319   -
320   - a, b, g = radians(coord_raw[1, 3:])
321   - r_ref = tr.euler_matrix(a, b, g, 'rzyx')
322   - t_ref = tr.translation_matrix(coord_raw[1, :3])
323   - m_ref = asmatrix(tr.concatenate_matrices(t_ref, r_ref))
324   -
325   - m_dyn = m_ref.I * m_probe
326   - m_dyn[2, -1] = -m_dyn[2, -1]
327   -
328   - m_img = m_change * m_dyn
329   - r_obj = r_obj_img * m_obj_raw.I * s0_dyn.I * m_dyn * m_obj_raw
330   -
331   - m_img[:3, :3] = r_obj[:3, :3]
332   -
333   - scale, shear, angles, trans, persp = tr.decompose_matrix(m_img)
334   -
335   - coord = m_img[0, -1], m_img[1, -1], m_img[2, -1],\
336   - degrees(angles[0]), degrees(angles[1]), degrees(angles[2])
337   -
338   - wx.CallAfter(Publisher.sendMessage, 'Co-registered points', arg=m_img, position=coord)
339   - wx.CallAfter(Publisher.sendMessage, 'Update object matrix', m_img=m_img, coord=coord)
340   -
341   - # TODO: Optimize the value of sleep for each tracking device.
342   - sleep(0.175)
343   -
344   - # Debug tracker is not working with 0.175 so changed to 0.2
345   - # However, 0.2 is too low update frequency ~5 Hz. Need optimization URGENTLY.
346   - # sleep(.3)
347   -
348   - # partially working for translate and offset,
349   - # but offset is kept always in same axis, have to fix for rotation
350   - # M_dyn = M_reference.I * T_stylus
351   - # M_dyn[2, -1] = -M_dyn[2, -1]
352   - # M_dyn_ch = M_change * M_dyn
353   - # ddd = M_dyn_ch[0, -1], M_dyn_ch[1, -1], M_dyn_ch[2, -1]
354   - # M_dyn_ch[:3, -1] = asmatrix(db.flip_x_m(ddd)).reshape([3, 1])
355   - # M_final = S0 * M_obj_trans_0 * S0.I * M_dyn_ch
356   -
357   - # this works for static reference object rotation
358   - # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_raw.I * R_stylus * M_obj_rot_raw
359   - # this works for dynamic reference in rotation but not in translation
360   - # R_dyn = M_vtk * M_obj_rot_raw.I * S0_rot_dyn.I * R_reference.I * R_stylus * M_obj_rot_raw
361   -
362   - if self._pause_:
363   - return
... ...
invesalius/data/imagedata_utils.py
... ... @@ -37,6 +37,8 @@ from invesalius.data import vtk_utils as vtk_utils
37 37 import invesalius.reader.bitmap_reader as bitmap_reader
38 38 import invesalius.utils as utils
39 39 import invesalius.data.converters as converters
  40 +import invesalius.data.slice_ as sl
  41 +import invesalius.data.transformations as tr
40 42  
41 43 from invesalius import inv_paths
42 44  
... ... @@ -516,3 +518,67 @@ def image_normalize(image, min_=0.0, max_=1.0, output_dtype=np.int16):
516 518 imin, imax = image.min(), image.max()
517 519 output[:] = (image - imin) * ((max_ - min_) / (imax - imin)) + min_
518 520 return output
  521 +
  522 +
  523 +def world2invspace(affine=None):
  524 + """
  525 + Normalize image pixel intensity for int16 gray scale values.
  526 +
  527 + :param repos: list of translation and rotation [trans_x, trans_y, trans_z, rot_x, rot_y, rot_z] to reposition the
  528 + vtk object prior to applying the affine matrix transformation. Note: rotation given in degrees
  529 + :param user_matrix: affine matrix from image header, prefered QForm matrix
  530 + :return: vtk transform filter for repositioning the polydata and affine matrix to be used as SetUserMatrix in actor
  531 + """
  532 +
  533 + # remove scaling factor for non-unitary voxel dimensions
  534 + scale, shear, angs, trans, persp = tr.decompose_matrix(affine)
  535 + affine_noscale = tr.compose_matrix(scale=None, shear=shear, angles=angs, translate=trans, perspective=persp)
  536 + # repos_img = [0.] * 6
  537 + # repos_img[1] = -float(shape[1])
  538 + #
  539 + # repos_mat = np.identity(4)
  540 + # # translation
  541 + # repos_mat[:3, -1] = repos_img[:3]
  542 + # # rotation (in principle for invesalius space no rotation is needed)
  543 + # repos_mat[:3, :3] = tr.euler_matrix(*np.deg2rad(repos_img[3:]), axes='sxyz')[:3, :3]
  544 +
  545 + # if repos:
  546 + # transx, transy, transz, rotx, roty, rotz = repos
  547 + # # create a transform that rotates the stl source
  548 + # transform = vtk.vtkTransform()
  549 + # transform.PostMultiply()
  550 + # transform.RotateX(rotx)
  551 + # transform.RotateY(roty)
  552 + # transform.RotateZ(rotz)
  553 + # transform.Translate(transx, transy, transz)
  554 + #
  555 + # transform_filt = vtk.vtkTransformPolyDataFilter()
  556 + # transform_filt.SetTransform(transform)
  557 + # transform_filt.Update()
  558 +
  559 + # assuming vtk default transformation order is PreMultiply, the user matrix is set so:
  560 + # 1. Replace the object -> 2. Transform the object to desired position/orientation
  561 + # PreMultiplty: M = M*A where M is current transformation and A is applied transformation
  562 + # user_matrix = np.linalg.inv(user_matrix) @ repos_mat
  563 +
  564 + return np.linalg.inv(affine_noscale)
  565 +
  566 +
  567 +def convert_world_to_voxel(xyz, affine):
  568 + """
  569 + Convert a coordinate from the world space ((x, y, z); scanner space; millimeters) to the
  570 + voxel space ((i, j, k)). This is achieved by multiplying a coordinate by the inverse
  571 + of the affine transformation.
  572 + More information: https://nipy.org/nibabel/coordinate_systems.html
  573 + :param xyz: a list or array of 3 coordinates (x, y, z) in the world coordinates
  574 + :param affine: a 4x4 array containing the image affine transformation in homogeneous coordinates
  575 + :return: a 1x3 array with the point coordinates in image space (i, j, k)
  576 + """
  577 +
  578 + # print("xyz: ", xyz, "\naffine", affine)
  579 + # convert xyz coordinate to 1x4 homogeneous coordinates array
  580 + xyz_homo = np.hstack((xyz, 1.)).reshape([4, 1])
  581 + ijk_homo = np.linalg.inv(affine) @ xyz_homo
  582 + ijk = ijk_homo.T[np.newaxis, 0, :3]
  583 +
  584 + return ijk
... ...
invesalius/data/record_coords.py
... ... @@ -42,7 +42,8 @@ class Record(threading.Thread):
42 42 self.start()
43 43  
44 44 def __bind_events(self):
45   - Publisher.subscribe(self.UpdateCurrentCoords, 'Co-registered points')
  45 + # Publisher.subscribe(self.UpdateCurrentCoords, 'Co-registered points')
  46 + Publisher.subscribe(self.UpdateCurrentCoords, 'Update cross position')
46 47  
47 48 def UpdateCurrentCoords(self, arg, position):
48 49 self.coord = asarray(position)
... ... @@ -50,7 +51,10 @@ class Record(threading.Thread):
50 51 def stop(self):
51 52 self._pause_ = True
52 53 #save coords dialog
53   - filename = dlg.ShowSaveCoordsDialog("coords.csv")
  54 + filename = dlg.ShowLoadSaveDialog(message=_(u"Save coords as..."),
  55 + wildcard=_("Coordinates files (*.csv)|*.csv"),
  56 + style=wx.FD_SAVE | wx.FD_OVERWRITE_PROMPT,
  57 + default_filename="coords.csv", save_ext="csv")
54 58 if filename:
55 59 savetxt(filename, self.coord_list, delimiter=',', fmt='%.4f', header="time, x, y, z, a, b, g", comments="")
56 60  
... ...
invesalius/data/slice_.py
... ... @@ -82,6 +82,9 @@ class Slice(metaclass=utils.Singleton):
82 82 self.blend_filter = None
83 83 self.histogram = None
84 84 self._matrix = None
  85 + self._affine = np.identity(4)
  86 + self._n_tracts = 0
  87 + self._tracker = None
85 88 self.aux_matrices = {}
86 89 self.aux_matrices_colours = {}
87 90 self.state = const.STATE_DEFAULT
... ... @@ -144,6 +147,30 @@ class Slice(metaclass=utils.Singleton):
144 147 (s * d / 2.0) for (d, s) in zip(self.matrix.shape[::-1], self.spacing)
145 148 ]
146 149  
  150 + @property
  151 + def affine(self):
  152 + return self._affine
  153 +
  154 + @affine.setter
  155 + def affine(self, value):
  156 + self._affine = value
  157 +
  158 + @property
  159 + def n_tracts(self):
  160 + return self._n_tracts
  161 +
  162 + @n_tracts.setter
  163 + def n_tracts(self, value):
  164 + self._n_tracts = value
  165 +
  166 + @property
  167 + def tracker(self):
  168 + return self._tracker
  169 +
  170 + @tracker.setter
  171 + def tracker(self, value):
  172 + self._tracker = value
  173 +
147 174 def __bind_events(self):
148 175 # General slice control
149 176 Publisher.subscribe(self.CreateSurfaceFromIndex, "Create surface from index")
... ...
invesalius/data/styles.py
... ... @@ -47,6 +47,14 @@ from invesalius.data.measures import (CircleDensityMeasure, MeasureData,
47 47  
48 48 from invesalius_cy import floodfill
49 49  
  50 +# For tracts
  51 +import invesalius.data.tractography as dtr
  52 +# import invesalius.project as prj
  53 +import invesalius.data.slice_ as sl
  54 +import invesalius.data.bases as bases
  55 +# from time import sleep
  56 +# ---
  57 +
50 58 ORIENTATIONS = {
51 59 "AXIAL": const.AXIAL,
52 60 "CORONAL": const.CORONAL,
... ... @@ -465,6 +473,18 @@ class CrossInteractorStyle(DefaultInteractorStyle):
465 473 self.slice_actor = viewer.slice_data.actor
466 474 self.slice_data = viewer.slice_data
467 475  
  476 + # tracts
  477 + # self.seed = [0., 0., 0.]
  478 + # slic = sl.Slice()
  479 + # self.affine = slic.affine
  480 + # self.tracker = slic.tracker
  481 + #
  482 + # self.affine_vtk = vtk.vtkMatrix4x4()
  483 + # for row in range(0, 4):
  484 + # for col in range(0, 4):
  485 + # self.affine_vtk.SetElement(row, col, self.affine[row, col])
  486 + # ---
  487 +
468 488 self.picker = vtk.vtkWorldPointPicker()
469 489  
470 490 self.AddObserver("MouseMoveEvent", self.OnCrossMove)
... ... @@ -478,6 +498,7 @@ class CrossInteractorStyle(DefaultInteractorStyle):
478 498  
479 499 def CleanUp(self):
480 500 self.viewer._set_cross_visibility(0)
  501 + Publisher.sendMessage('Remove tracts')
481 502 Publisher.sendMessage('Toggle toolbar item',
482 503 _id=self.state_code, value=False)
483 504  
... ... @@ -497,12 +518,20 @@ class CrossInteractorStyle(DefaultInteractorStyle):
497 518 wx, wy, wz = self.viewer.get_coordinate_cursor(mouse_x, mouse_y, self.picker)
498 519 px, py = self.viewer.get_slice_pixel_coord_by_world_pos(wx, wy, wz)
499 520 coord = self.viewer.calcultate_scroll_position(px, py)
500   - Publisher.sendMessage('Update cross position', position=(wx, wy, wz))
  521 +
  522 + # Tracts
  523 + # pos_world_aux = np.ones([4, 1])
  524 + # pos_world_aux[:3, -1] = bases.flip_x((wx, wy, wz))[:3]
  525 + # pos_world = np.linalg.inv(self.affine) @ pos_world_aux
  526 + # seed_aux = pos_world.reshape([1, 4])[0, :3]
  527 + # self.seed = seed_aux[np.newaxis, :]
  528 + # print("Check the seed: ", self.seed)
  529 + #
  530 +
  531 + Publisher.sendMessage('Update cross position', arg=None, position=(wx, wy, wz, 0., 0., 0.))
501 532 self.ScrollSlice(coord)
502   - Publisher.sendMessage('Set ball reference position', position=(wx, wy, wz))
503   - Publisher.sendMessage('Co-registered points', arg=None, position=(wx, wy, wz, 0., 0., 0.))
504 533  
505   - iren.Render()
  534 + # iren.Render()
506 535  
507 536 def ScrollSlice(self, coord):
508 537 if self.orientation == "AXIAL":
... ... @@ -522,6 +551,92 @@ class CrossInteractorStyle(DefaultInteractorStyle):
522 551 index=coord[0])
523 552  
524 553  
  554 +class TractsInteractorStyle(CrossInteractorStyle):
  555 + """
  556 + Interactor style responsible for tracts visualization.
  557 + """
  558 + def __init__(self, viewer):
  559 + CrossInteractorStyle.__init__(self, viewer)
  560 +
  561 + # self.state_code = const.SLICE_STATE_TRACTS
  562 +
  563 + self.viewer = viewer
  564 + # print("Im fucking brilliant!")
  565 + self.tracts = None
  566 +
  567 + # data_dir = b'C:\Users\deoliv1\OneDrive\data\dti'
  568 + # FOD_path = b"sub-P0_dwi_FOD.nii"
  569 + # full_path = os.path.join(data_dir, FOD_path)
  570 + # self.tracker = Trekker.tracker(full_path)
  571 + # self.orientation = viewer.orientation
  572 + # self.slice_actor = viewer.slice_data.actor
  573 + # self.slice_data = viewer.slice_data
  574 +
  575 + # self.picker = vtk.vtkWorldPointPicker()
  576 +
  577 + self.AddObserver("MouseMoveEvent", self.OnTractsMove)
  578 + self.AddObserver("LeftButtonPressEvent", self.OnTractsMouseClick)
  579 + self.AddObserver("LeftButtonReleaseEvent", self.OnTractsReleaseLeftButton)
  580 +
  581 + # def SetUp(self):
  582 + # self.viewer._set_cross_visibility(1)
  583 + # Publisher.sendMessage('Toggle toolbar item',
  584 + # _id=self.state_code, value=True)
  585 +
  586 + # def CleanUp(self):
  587 + # self.viewer._set_cross_visibility(0)
  588 + # Publisher.sendMessage('Toggle toolbar item',
  589 + # _id=self.state_code, value=False)
  590 +
  591 + def OnTractsMove(self, obj, evt):
  592 + # The user moved the mouse with left button pressed
  593 + if self.left_pressed:
  594 + # print("OnTractsMove interactor style")
  595 + # iren = obj.GetInteractor()
  596 + self.ChangeTracts(True)
  597 +
  598 + def OnTractsMouseClick(self, obj, evt):
  599 + # print("Single mouse click")
  600 + # self.tracts = dtr.compute_tracts(self.tracker, self.seed, self.left_pressed)
  601 + self.ChangeTracts(True)
  602 +
  603 + def OnTractsReleaseLeftButton(self, obj, evt):
  604 + # time.sleep(3.)
  605 + self.tracts.stop()
  606 + # self.ChangeCrossPosition(iren)
  607 +
  608 + def ChangeTracts(self, pressed):
  609 + # print("Trying to compute tracts")
  610 + self.tracts = dtr.compute_tracts(self.tracker, self.seed, self.affine_vtk, pressed)
  611 + # mouse_x, mouse_y = iren.GetEventPosition()
  612 + # wx, wy, wz = self.viewer.get_coordinate_cursor(mouse_x, mouse_y, self.picker)
  613 + # px, py = self.viewer.get_slice_pixel_coord_by_world_pos(wx, wy, wz)
  614 + # coord = self.viewer.calcultate_scroll_position(px, py)
  615 + # Publisher.sendMessage('Update cross position', position=(wx, wy, wz))
  616 + # # self.ScrollSlice(coord)
  617 + # Publisher.sendMessage('Set ball reference position', position=(wx, wy, wz))
  618 + # Publisher.sendMessage('Co-registered points', arg=None, position=(wx, wy, wz, 0., 0., 0.))
  619 +
  620 + # iren.Render()
  621 +
  622 + # def ScrollSlice(self, coord):
  623 + # if self.orientation == "AXIAL":
  624 + # Publisher.sendMessage(('Set scroll position', 'SAGITAL'),
  625 + # index=coord[0])
  626 + # Publisher.sendMessage(('Set scroll position', 'CORONAL'),
  627 + # index=coord[1])
  628 + # elif self.orientation == "SAGITAL":
  629 + # Publisher.sendMessage(('Set scroll position', 'AXIAL'),
  630 + # index=coord[2])
  631 + # Publisher.sendMessage(('Set scroll position', 'CORONAL'),
  632 + # index=coord[1])
  633 + # elif self.orientation == "CORONAL":
  634 + # Publisher.sendMessage(('Set scroll position', 'AXIAL'),
  635 + # index=coord[2])
  636 + # Publisher.sendMessage(('Set scroll position', 'SAGITAL'),
  637 + # index=coord[0])
  638 +
  639 +
525 640 class WWWLInteractorStyle(DefaultInteractorStyle):
526 641 """
527 642 Interactor style responsible for Window Level & Width functionality.
... ... @@ -2784,6 +2899,7 @@ class Styles:
2784 2899 const.SLICE_STATE_SELECT_MASK_PARTS: SelectMaskPartsInteractorStyle,
2785 2900 const.SLICE_STATE_FFILL_SEGMENTATION: FloodFillSegmentInteractorStyle,
2786 2901 const.SLICE_STATE_CROP_MASK: CropMaskInteractorStyle,
  2902 + const.SLICE_STATE_TRACTS: TractsInteractorStyle,
2787 2903 }
2788 2904  
2789 2905 @classmethod
... ...
invesalius/data/surface.py
... ... @@ -383,6 +383,7 @@ class SurfaceManager():
383 383  
384 384 actor = vtk.vtkActor()
385 385 actor.SetMapper(mapper)
  386 + actor.GetProperty().SetBackfaceCulling(1)
386 387  
387 388 print("BOunds", actor.GetBounds())
388 389  
... ... @@ -494,6 +495,7 @@ class SurfaceManager():
494 495  
495 496 # Represent an object (geometry & properties) in the rendered scene
496 497 actor = vtk.vtkActor()
  498 + actor.GetProperty().SetBackfaceCulling(1)
497 499 actor.SetMapper(mapper)
498 500  
499 501 # Set actor colour and transparency
... ... @@ -536,6 +538,7 @@ class SurfaceManager():
536 538  
537 539 # Represent an object (geometry & properties) in the rendered scene
538 540 actor = vtk.vtkActor()
  541 + actor.GetProperty().SetBackfaceCulling(1)
539 542 actor.SetMapper(mapper)
540 543 del mapper
541 544 #Create Surface instance
... ...
invesalius/data/trackers.py
... ... @@ -62,6 +62,7 @@ def DefaultTracker(tracker_id):
62 62 # return tracker initialization variable and type of connection
63 63 return trck_init, 'wrapper'
64 64  
  65 +
65 66 def PolarisTracker(tracker_id):
66 67 from wx import ID_OK
67 68 trck_init = None
... ... @@ -91,6 +92,7 @@ def PolarisTracker(tracker_id):
91 92 # return tracker initialization variable and type of connection
92 93 return trck_init, lib_mode
93 94  
  95 +
94 96 def CameraTracker(tracker_id):
95 97 trck_init = None
96 98 try:
... ... @@ -105,6 +107,7 @@ def CameraTracker(tracker_id):
105 107 # return tracker initialization variable and type of connection
106 108 return trck_init, 'wrapper'
107 109  
  110 +
108 111 def ClaronTracker(tracker_id):
109 112 import invesalius.constants as const
110 113 from invesalius import inv_paths
... ...
invesalius/data/tractography.py 0 → 100644
... ... @@ -0,0 +1,471 @@
  1 +# -*- coding: utf-8 -*-
  2 +
  3 +#--------------------------------------------------------------------------
  4 +# Software: InVesalius - Software de Reconstrucao 3D de Imagens Medicas
  5 +# Copyright: (C) 2001 Centro de Pesquisas Renato Archer
  6 +# Homepage: http://www.softwarepublico.gov.br
  7 +# Contact: invesalius@cti.gov.br
  8 +# License: GNU - GPL 2 (LICENSE.txt/LICENCA.txt)
  9 +#--------------------------------------------------------------------------
  10 +# Este programa e software livre; voce pode redistribui-lo e/ou
  11 +# modifica-lo sob os termos da Licenca Publica Geral GNU, conforme
  12 +# publicada pela Free Software Foundation; de acordo com a versao 2
  13 +# da Licenca.
  14 +#
  15 +# Este programa eh distribuido na expectativa de ser util, mas SEM
  16 +# QUALQUER GARANTIA; sem mesmo a garantia implicita de
  17 +# COMERCIALIZACAO ou de ADEQUACAO A QUALQUER PROPOSITO EM
  18 +# PARTICULAR. Consulte a Licenca Publica Geral GNU para obter mais
  19 +# detalhes.
  20 +#--------------------------------------------------------------------------
  21 +
  22 +# Author: Victor Hugo Souza (victorhos-at-hotmail.com)
  23 +# Contributions: Dogu Baran Aydogan
  24 +# Initial date: 8 May 2020
  25 +
  26 +import threading
  27 +import time
  28 +
  29 +import numpy as np
  30 +import queue
  31 +from pubsub import pub as Publisher
  32 +import vtk
  33 +
  34 +import invesalius.constants as const
  35 +import invesalius.data.imagedata_utils as img_utils
  36 +
  37 +# Nice print for arrays
  38 +# np.set_printoptions(precision=2)
  39 +# np.set_printoptions(suppress=True)
  40 +
  41 +
  42 +def compute_directions(trk_n):
  43 + """Compute direction of a single tract in each point and return as an RGB color
  44 +
  45 + :param trk_n: nx3 array of doubles (x, y, z) point coordinates composing the tract
  46 + :type trk_n: numpy.ndarray
  47 + :return: nx3 array of int (x, y, z) RGB colors in the range 0 - 255
  48 + :rtype: numpy.ndarray
  49 + """
  50 +
  51 + # trk_d = np.diff(trk_n, axis=0, append=2*trk_n[np.newaxis, -1, :])
  52 + trk_d = np.diff(trk_n, axis=0, append=trk_n[np.newaxis, -2, :])
  53 + trk_d[-1, :] *= -1
  54 + # check that linalg norm makes second norm
  55 + # https://stackoverflow.com/questions/21030391/how-to-normalize-an-array-in-numpy
  56 + direction = 255 * np.absolute((trk_d / np.linalg.norm(trk_d, axis=1)[:, None]))
  57 + return direction.astype(int)
  58 +
  59 +
  60 +def compute_tubes(trk, direction):
  61 + """Compute and assign colors to a vtkTube for visualization of a single tract
  62 +
  63 + :param trk: nx3 array of doubles (x, y, z) point coordinates composing the tract
  64 + :type trk: numpy.ndarray
  65 + :param direction: nx3 array of int (x, y, z) RGB colors in the range 0 - 255
  66 + :type direction: numpy.ndarray
  67 + :return: a vtkTubeFilter instance
  68 + :rtype: vtkTubeFilter
  69 + """
  70 +
  71 + numb_points = trk.shape[0]
  72 + points = vtk.vtkPoints()
  73 + lines = vtk.vtkCellArray()
  74 +
  75 + colors = vtk.vtkUnsignedCharArray()
  76 + colors.SetNumberOfComponents(3)
  77 +
  78 + k = 0
  79 + lines.InsertNextCell(numb_points)
  80 + for j in range(numb_points):
  81 + points.InsertNextPoint(trk[j, :])
  82 + colors.InsertNextTuple(direction[j, :])
  83 + lines.InsertCellPoint(k)
  84 + k += 1
  85 +
  86 + trk_data = vtk.vtkPolyData()
  87 + trk_data.SetPoints(points)
  88 + trk_data.SetLines(lines)
  89 + trk_data.GetPointData().SetScalars(colors)
  90 +
  91 + # make it a tube
  92 + trk_tube = vtk.vtkTubeFilter()
  93 + trk_tube.SetRadius(0.5)
  94 + trk_tube.SetNumberOfSides(4)
  95 + trk_tube.SetInputData(trk_data)
  96 + trk_tube.Update()
  97 +
  98 + return trk_tube
  99 +
  100 +
  101 +def combine_tracts_root(out_list, root, n_block):
  102 + """Adds a set of tracts to given position in a given vtkMultiBlockDataSet
  103 +
  104 + :param out_list: List of vtkTubeFilters representing the tracts
  105 + :type out_list: list
  106 + :param root: A collection of tracts as a vtkMultiBlockDataSet
  107 + :type root: vtkMultiBlockDataSet
  108 + :param n_block: The location in the given vtkMultiBlockDataSet to insert the new tracts
  109 + :type n_block: int
  110 + :return: The updated collection of tracts as a vtkMultiBlockDataSet
  111 + :rtype: vtkMultiBlockDataSet
  112 + """
  113 +
  114 + # create tracts only when at least one was computed
  115 + # print("Len outlist in root: ", len(out_list))
  116 + if not out_list.count(None) == len(out_list):
  117 + for n, tube in enumerate(out_list):
  118 + root.SetBlock(n_block + n, tube.GetOutput())
  119 +
  120 + return root
  121 +
  122 +
  123 +def combine_tracts_branch(out_list):
  124 + """Combines a set of tracts in vtkMultiBlockDataSet
  125 +
  126 + :param out_list: List of vtkTubeFilters representing the tracts
  127 + :type out_list: list
  128 + :return: A collection of tracts as a vtkMultiBlockDataSet
  129 + :rtype: vtkMultiBlockDataSet
  130 + """
  131 +
  132 + branch = vtk.vtkMultiBlockDataSet()
  133 + # create tracts only when at least one was computed
  134 + # print("Len outlist in root: ", len(out_list))
  135 + if not out_list.count(None) == len(out_list):
  136 + for n, tube in enumerate(out_list):
  137 + branch.SetBlock(n, tube.GetOutput())
  138 +
  139 + return branch
  140 +
  141 +
  142 +def tracts_computation(trk_list, root, n_tracts):
  143 + """Convert the list of all computed tracts given by Trekker run and returns a vtkMultiBlockDataSet
  144 +
  145 + :param trk_list: List of lists containing the computed tracts and corresponding coordinates
  146 + :type trk_list: list
  147 + :param root: A collection of tracts as a vtkMultiBlockDataSet
  148 + :type root: vtkMultiBlockDataSet
  149 + :param n_tracts:
  150 + :type n_tracts: int
  151 + :return: The updated collection of tracts as a vtkMultiBlockDataSet
  152 + :rtype: vtkMultiBlockDataSet
  153 + """
  154 +
  155 + # Transform tracts to array
  156 + trk_arr = [np.asarray(trk_n).T if trk_n else None for trk_n in trk_list]
  157 +
  158 + # Compute the directions
  159 + trk_dir = [compute_directions(trk_n) for trk_n in trk_arr]
  160 +
  161 + # Compute the vtk tubes
  162 + out_list = [compute_tubes(trk_arr_n, trk_dir_n) for trk_arr_n, trk_dir_n in zip(trk_arr, trk_dir)]
  163 +
  164 + root = combine_tracts_root(out_list, root, n_tracts)
  165 +
  166 + return root
  167 +
  168 +
  169 +def compute_tracts(trekker, position, affine, affine_vtk, n_tracts):
  170 + """ Compute tractograms using the Trekker library.
  171 +
  172 + :param trekker: Trekker library instance
  173 + :type trekker: Trekker.T
  174 + :param position: 3 double coordinates (x, y, z) in list or array
  175 + :type position: list
  176 + :param affine: 4 x 4 numpy double array
  177 + :type affine: numpy.ndarray
  178 + :param affine_vtk: vtkMatrix4x4 isntance with affine transformation matrix
  179 + :type affine_vtk: vtkMatrix4x4
  180 + :param n_tracts: number of tracts to compute
  181 + :type n_tracts: int
  182 + """
  183 +
  184 + # during neuronavigation, root needs to be initialized outside the while loop so the new tracts
  185 + # can be appended to the root block set
  186 + root = vtk.vtkMultiBlockDataSet()
  187 + # Juuso's
  188 + # seed = np.array([[-8.49, -8.39, 2.5]])
  189 + # Baran M1
  190 + # seed = np.array([[27.53, -77.37, 46.42]])
  191 + seed_trk = img_utils.convert_world_to_voxel(position, affine)
  192 + # print("seed example: {}".format(seed_trk))
  193 + trekker.seed_coordinates(np.repeat(seed_trk, n_tracts, axis=0))
  194 + # print("trk list len: ", len(trekker.run()))
  195 + trk_list = trekker.run()
  196 + if trk_list:
  197 + root = tracts_computation(trk_list, root, 0)
  198 + Publisher.sendMessage('Remove tracts')
  199 + Publisher.sendMessage('Update tracts', flag=True, root=root, affine_vtk=affine_vtk)
  200 + else:
  201 + Publisher.sendMessage('Remove tracts')
  202 +
  203 +
  204 +def tracts_computation_branch(trk_list):
  205 + """Convert the list of all computed tracts given by Trekker run and returns a vtkMultiBlockDataSet
  206 +
  207 + :param trk_list: List of lists containing the computed tracts and corresponding coordinates
  208 + :type trk_list: list
  209 + :return: The collection of tracts as a vtkMultiBlockDataSet
  210 + :rtype: vtkMultiBlockDataSet
  211 + """
  212 + # Transform tracts to array
  213 + trk_arr = [np.asarray(trk_n).T if trk_n else None for trk_n in trk_list]
  214 + # Compute the directions
  215 + trk_dir = [compute_directions(trk_n) for trk_n in trk_arr]
  216 + # Compute the vtk tubes
  217 + tube_list = [compute_tubes(trk_arr_n, trk_dir_n) for trk_arr_n, trk_dir_n in zip(trk_arr, trk_dir)]
  218 + branch = combine_tracts_branch(tube_list)
  219 +
  220 + return branch
  221 +
  222 +
  223 +class ComputeTractsThread(threading.Thread):
  224 +
  225 + def __init__(self, inp, affine_vtk, coord_tracts_queue, tracts_queue, event, sle):
  226 + """Class (threading) to compute real time tractography data for visualization.
  227 +
  228 + Tracts are computed using the Trekker library by Baran Aydogan (https://dmritrekker.github.io/)
  229 + For VTK visualization, each tract (fiber) is a constructed as a tube and many tubes combined in one
  230 + vtkMultiBlockDataSet named as a branch. Several branches are combined in another vtkMultiBlockDataSet named as
  231 + bundle, to obtain fast computation and visualization. The bundle dataset is mapped to a single vtkActor.
  232 + Mapper and Actor are computer in the data/viewer_volume.py module for easier handling in the invesalius 3D scene.
  233 +
  234 + Sleep function in run method is used to avoid blocking GUI and more fluent, real-time navigation
  235 +
  236 + :param inp: List of inputs: trekker instance, affine numpy array, seed_offset, seed_radius, n_threads
  237 + :type inp: list
  238 + :param affine_vtk: Affine matrix in vtkMatrix4x4 instance to update objects position in 3D scene
  239 + :type affine_vtk: vtkMatrix4x4
  240 + :param coord_queue: Queue instance that manage coordinates read from tracking device and coregistered
  241 + :type coord_queue: queue.Queue
  242 + :param visualization_queue: Queue instance that manage coordinates to be visualized
  243 + :type visualization_queue: queue.Queue
  244 + :param event: Threading event to coordinate when tasks as done and allow UI release
  245 + :type event: threading.Event
  246 + :param sle: Sleep pause in seconds
  247 + :type sle: float
  248 + """
  249 +
  250 + threading.Thread.__init__(self, name='ComputeTractsThread')
  251 + self.inp = inp
  252 + self.affine_vtk = affine_vtk
  253 + # self.coord_queue = coord_queue
  254 + self.coord_tracts_queue = coord_tracts_queue
  255 + self.tracts_queue = tracts_queue
  256 + # self.visualization_queue = visualization_queue
  257 + self.event = event
  258 + self.sle = sle
  259 +
  260 + def run(self):
  261 +
  262 + trekker, affine, offset, n_tracts_total, seed_radius, n_threads = self.inp
  263 + # n_threads = n_tracts_total
  264 + p_old = np.array([[0., 0., 0.]])
  265 + n_tracts = 0
  266 +
  267 + # Compute the tracts
  268 + # print('ComputeTractsThread: event {}'.format(self.event.is_set()))
  269 + while not self.event.is_set():
  270 + try:
  271 + # print("Computing tracts")
  272 + # get from the queue the coordinates, coregistration transformation matrix, and flipped matrix
  273 + m_img_flip = self.coord_tracts_queue.get_nowait()
  274 + # coord, m_img, m_img_flip = self.coord_queue.get_nowait()
  275 + # print('ComputeTractsThread: get {}'.format(count))
  276 +
  277 + # 20200402: in this new refactored version the m_img comes different than the position
  278 + # the new version m_img is already flixped in y, which means that Y is negative
  279 + # if only the Y is negative maybe no need for the flip_x funtcion at all in the navigation
  280 + # but check all coord_queue before why now the m_img comes different than position
  281 + # 20200403: indeed flip_x is just a -1 multiplication to the Y coordinate, remove function flip_x
  282 + # m_img_flip = m_img.copy()
  283 + # m_img_flip[1, -1] = -m_img_flip[1, -1]
  284 +
  285 + # translate the coordinate along the normal vector of the object/coil
  286 + coord_offset = m_img_flip[:3, -1] - offset * m_img_flip[:3, 2]
  287 + # coord_offset = np.array([[27.53, -77.37, 46.42]])
  288 + dist = abs(np.linalg.norm(p_old - np.asarray(coord_offset)))
  289 + p_old = coord_offset.copy()
  290 +
  291 + # print("p_new_shape", coord_offset.shape)
  292 + # print("m_img_flip_shape", m_img_flip.shape)
  293 + seed_trk = img_utils.convert_world_to_voxel(coord_offset, affine)
  294 + # Juuso's
  295 + # seed_trk = np.array([[-8.49, -8.39, 2.5]])
  296 + # Baran M1
  297 + # seed_trk = np.array([[27.53, -77.37, 46.42]])
  298 + # print("Seed: {}".format(seed))
  299 +
  300 + # set the seeds for trekker, one seed is repeated n_threads times
  301 + # trekker has internal multiprocessing approach done in C. Here the number of available threads is give,
  302 + # but in case a large number of tracts is requested, it will compute all in parallel automatically
  303 + # for a more fluent navigation, better to compute the maximum number the computer handles
  304 + trekker.seed_coordinates(np.repeat(seed_trk, n_threads, axis=0))
  305 +
  306 + # run the trekker, this is the slowest line of code, be careful to just use once!
  307 + trk_list = trekker.run()
  308 +
  309 + if trk_list:
  310 + # print("dist: {}".format(dist))
  311 + if dist >= seed_radius:
  312 + # when moving the coil further than the seed_radius restart the bundle computation
  313 + bundle = vtk.vtkMultiBlockDataSet()
  314 + n_branches = 0
  315 + branch = tracts_computation_branch(trk_list)
  316 + bundle.SetBlock(n_branches, branch)
  317 + n_branches += 1
  318 + n_tracts = branch.GetNumberOfBlocks()
  319 +
  320 + # TODO: maybe keep computing even if reaches the maximum
  321 + elif dist < seed_radius and n_tracts < n_tracts_total:
  322 + # compute tracts blocks and add to bungle until reaches the maximum number of tracts
  323 + branch = tracts_computation_branch(trk_list)
  324 + if bundle:
  325 + bundle.SetBlock(n_branches, branch)
  326 + n_tracts += branch.GetNumberOfBlocks()
  327 + n_branches += 1
  328 +
  329 + else:
  330 + bundle = None
  331 +
  332 + # rethink if this should be inside the if condition, it may lock the thread if no tracts are found
  333 + # use no wait to ensure maximum speed and avoid visualizing old tracts in the queue, this might
  334 + # be more evident in slow computer or for heavier tract computations, it is better slow update
  335 + # than visualizing old data
  336 + # self.visualization_queue.put_nowait([coord, m_img, bundle])
  337 + self.tracts_queue.put_nowait((bundle, self.affine_vtk, coord_offset))
  338 + # print('ComputeTractsThread: put {}'.format(count))
  339 +
  340 + self.coord_tracts_queue.task_done()
  341 + # self.coord_queue.task_done()
  342 + # print('ComputeTractsThread: done {}'.format(count))
  343 +
  344 + # sleep required to prevent user interface from being unresponsive
  345 + time.sleep(self.sle)
  346 + # if no coordinates pass
  347 + except queue.Empty:
  348 + # print("Empty queue in tractography")
  349 + pass
  350 + # if queue is full mark as done (may not be needed in this new "nowait" method)
  351 + except queue.Full:
  352 + # self.coord_queue.task_done()
  353 + self.coord_tracts_queue.task_done()
  354 +
  355 +
  356 +class ComputeTractsThreadSingleBlock(threading.Thread):
  357 +
  358 + def __init__(self, inp, affine_vtk, coord_queue, visualization_queue, event, sle):
  359 + """Class (threading) to compute real time tractography data for visualization in a single loop.
  360 +
  361 + Different than ComputeTractsThread because it does not keep adding tracts to the bundle until maximum,
  362 + is reached. It actually compute all requested tracts at once. (Might be deleted in the future)!
  363 + Tracts are computed using the Trekker library by Baran Aydogan (https://dmritrekker.github.io/)
  364 + For VTK visualization, each tract (fiber) is a constructed as a tube and many tubes combined in one
  365 + vtkMultiBlockDataSet named as a branch. Several branches are combined in another vtkMultiBlockDataSet named as
  366 + bundle, to obtain fast computation and visualization. The bundle dataset is mapped to a single vtkActor.
  367 + Mapper and Actor are computer in the data/viewer_volume.py module for easier handling in the invesalius 3D scene.
  368 +
  369 + Sleep function in run method is used to avoid blocking GUI and more fluent, real-time navigation
  370 +
  371 + :param inp: List of inputs: trekker instance, affine numpy array, seed_offset, seed_radius, n_threads
  372 + :type inp: list
  373 + :param affine_vtk: Affine matrix in vtkMatrix4x4 instance to update objects position in 3D scene
  374 + :type affine_vtk: vtkMatrix4x4
  375 + :param coord_queue: Queue instance that manage coordinates read from tracking device and coregistered
  376 + :type coord_queue: queue.Queue
  377 + :param visualization_queue: Queue instance that manage coordinates to be visualized
  378 + :type visualization_queue: queue.Queue
  379 + :param event: Threading event to coordinate when tasks as done and allow UI release
  380 + :type event: threading.Event
  381 + :param sle: Sleep pause in seconds
  382 + :type sle: float
  383 + """
  384 +
  385 + threading.Thread.__init__(self, name='ComputeTractsThread')
  386 + self.inp = inp
  387 + self.affine_vtk = affine_vtk
  388 + self.coord_queue = coord_queue
  389 + self.visualization_queue = visualization_queue
  390 + self.event = event
  391 + self.sle = sle
  392 +
  393 + def run(self):
  394 +
  395 + trekker, affine, offset, n_tracts_total, seed_radius, n_threads = self.inp
  396 + # as a single block, computes all the maximum number of tracts at once, not optimal for navigation
  397 + n_threads = n_tracts_total
  398 + p_old = np.array([[0., 0., 0.]])
  399 + root = vtk.vtkMultiBlockDataSet()
  400 +
  401 + # Compute the tracts
  402 + # print('ComputeTractsThread: event {}'.format(self.event.is_set()))
  403 + while not self.event.is_set():
  404 + try:
  405 + coord, m_img, m_img_flip = self.coord_queue.get_nowait()
  406 +
  407 + # translate the coordinate along the normal vector of the object/coil
  408 + coord_offset = m_img_flip[:3, -1] - offset * m_img_flip[:3, 2]
  409 + # coord_offset = np.array([[27.53, -77.37, 46.42]])
  410 + dist = abs(np.linalg.norm(p_old - np.asarray(coord_offset)))
  411 + p_old = coord_offset.copy()
  412 + seed_trk = img_utils.convert_world_to_voxel(coord_offset, affine)
  413 + # Juuso's
  414 + # seed_trk = np.array([[-8.49, -8.39, 2.5]])
  415 + # Baran M1
  416 + # seed_trk = np.array([[27.53, -77.37, 46.42]])
  417 + # print("Seed: {}".format(seed))
  418 +
  419 + # set the seeds for trekker, one seed is repeated n_threads times
  420 + # trekker has internal multiprocessing approach done in C. Here the number of available threads is give,
  421 + # but in case a large number of tracts is requested, it will compute all in parallel automatically
  422 + # for a more fluent navigation, better to compute the maximum number the computer handles
  423 + trekker.seed_coordinates(np.repeat(seed_trk, n_threads, axis=0))
  424 + # run the trekker, this is the slowest line of code, be careful to just use once!
  425 + trk_list = trekker.run()
  426 +
  427 + if trk_list:
  428 + # if the seed is outside the defined radius, restart the bundle computation
  429 + if dist >= seed_radius:
  430 + root = tracts_computation(trk_list, root, 0)
  431 + self.visualization_queue.put_nowait((coord, m_img, root))
  432 +
  433 + self.coord_queue.task_done()
  434 + time.sleep(self.sle)
  435 + except queue.Empty:
  436 + pass
  437 + except queue.Full:
  438 + self.coord_queue.task_done()
  439 +
  440 +
  441 +def set_trekker_parameters(trekker, params):
  442 + """Set all user-defined parameters for tractography computation using the Trekker library
  443 +
  444 + :param trekker: Trekker instance
  445 + :type trekker: Trekker.T
  446 + :param params: Dictionary containing the parameters values to set in Trekker. Initial values are in constants.py
  447 + :type params: dict
  448 + :return: List containing the Trekker instance and number of threads for parallel processing in the computer
  449 + :rtype: list
  450 + """
  451 + trekker.seed_maxTrials(params['seed_max'])
  452 + # trekker.stepSize(params['step_size'])
  453 + trekker.minFODamp(params['min_fod'])
  454 + # trekker.probeQuality(params['probe_quality'])
  455 + # trekker.maxEstInterval(params['max_interval'])
  456 + # trekker.minRadiusOfCurvature(params['min_radius_curv'])
  457 + # trekker.probeLength(params['probe_length'])
  458 + trekker.writeInterval(params['write_interval'])
  459 + trekker.maxLength(params['max_lenth'])
  460 + trekker.minLength(params['min_lenth'])
  461 + trekker.maxSamplingPerStep(params['max_sampling_step'])
  462 +
  463 + # check number if number of cores is valid in configuration file,
  464 + # otherwise use the maximum number of threads which is usually 2*N_CPUS
  465 + n_threads = 2 * const.N_CPU
  466 + if isinstance((params['numb_threads']), int) and params['numb_threads'] <= 2*const.N_CPU:
  467 + n_threads = params['numb_threads']
  468 +
  469 + trekker.numberOfThreads(n_threads)
  470 + # print("Trekker config updated: n_threads, {}; seed_max, {}".format(n_threads, params['seed_max']))
  471 + return trekker, n_threads
... ...
invesalius/data/trigger.py
... ... @@ -39,7 +39,7 @@ class Trigger(threading.Thread):
39 39 try:
40 40 import serial
41 41  
42   - self.trigger_init = serial.Serial('COM1', baudrate=9600, timeout=0)
  42 + self.trigger_init = serial.Serial('COM5', baudrate=115200, timeout=0)
43 43 self.COM = True
44 44  
45 45 except:
... ... @@ -82,3 +82,91 @@ class Trigger(threading.Thread):
82 82 if self.trigger_init:
83 83 self.trigger_init.close()
84 84 return
  85 +
  86 +
  87 +class TriggerNew(threading.Thread):
  88 +
  89 + def __init__(self, trigger_queue, event, sle):
  90 + """Class (threading) to compute real time tractography data for visualization in a single loop.
  91 +
  92 + Different than ComputeTractsThread because it does not keep adding tracts to the bundle until maximum,
  93 + is reached. It actually compute all requested tracts at once. (Might be deleted in the future)!
  94 + Tracts are computed using the Trekker library by Baran Aydogan (https://dmritrekker.github.io/)
  95 + For VTK visualization, each tract (fiber) is a constructed as a tube and many tubes combined in one
  96 + vtkMultiBlockDataSet named as a branch. Several branches are combined in another vtkMultiBlockDataSet named as
  97 + bundle, to obtain fast computation and visualization. The bundle dataset is mapped to a single vtkActor.
  98 + Mapper and Actor are computer in the data/viewer_volume.py module for easier handling in the invesalius 3D scene.
  99 +
  100 + Sleep function in run method is used to avoid blocking GUI and more fluent, real-time navigation
  101 +
  102 + :param inp: List of inputs: trekker instance, affine numpy array, seed_offset, seed_radius, n_threads
  103 + :type inp: list
  104 + :param affine_vtk: Affine matrix in vtkMatrix4x4 instance to update objects position in 3D scene
  105 + :type affine_vtk: vtkMatrix4x4
  106 + :param coord_queue: Queue instance that manage coordinates read from tracking device and coregistered
  107 + :type coord_queue: queue.Queue
  108 + :param visualization_queue: Queue instance that manage coordinates to be visualized
  109 + :type visualization_queue: queue.Queue
  110 + :param event: Threading event to coordinate when tasks as done and allow UI release
  111 + :type event: threading.Event
  112 + :param sle: Sleep pause in seconds
  113 + :type sle: float
  114 + """
  115 +
  116 + threading.Thread.__init__(self, name='Trigger')
  117 +
  118 + self.trigger_init = None
  119 + self.stylusplh = False
  120 + # self.COM = False
  121 + self.__bind_events()
  122 + try:
  123 + import serial
  124 +
  125 + self.trigger_init = serial.Serial('COM5', baudrate=115200, timeout=0)
  126 + # self.COM = True
  127 +
  128 + except:
  129 + #wx.MessageBox(_('Connection with port COM1 failed'), _('Communication error'), wx.OK | wx.ICON_ERROR)
  130 + print("Trigger init error: Connection with port COM failed")
  131 + # self.COM = False
  132 + pass
  133 +
  134 + # self.coord_queue = coord_queue
  135 + self.trigger_queue = trigger_queue
  136 + self.event = event
  137 + self.sle = sle
  138 +
  139 + def run(self):
  140 +
  141 + while not self.event.is_set():
  142 + trigger_on = False
  143 + try:
  144 + self.trigger_init.write(b'0')
  145 + sleep(0.3)
  146 + lines = self.trigger_init.readlines()
  147 + # Following lines can simulate a trigger in 3 sec repetitions
  148 + # sleep(3)
  149 + # lines = True
  150 + if lines:
  151 + trigger_on = True
  152 + # wx.CallAfter(Publisher.sendMessage, 'Create marker')
  153 +
  154 + if self.stylusplh:
  155 + trigger_on = True
  156 + # wx.CallAfter(Publisher.sendMessage, 'Create marker')
  157 + self.stylusplh = False
  158 +
  159 + self.trigger_queue.put_nowait(trigger_on)
  160 + sleep(self.sle)
  161 +
  162 + except:
  163 + print("Trigger not read, error")
  164 + pass
  165 +
  166 + # except queue.Empty:
  167 + # pass
  168 + # except queue.Full:
  169 + # self.coord_queue.task_done()
  170 + else:
  171 + if self.trigger_init:
  172 + self.trigger_init.close()
... ...
invesalius/data/viewer_slice.py
... ... @@ -573,7 +573,7 @@ class Viewer(wx.Panel):
573 573  
574 574 self.cross.SetFocalPoint((ux, uy, uz))
575 575 self.ScrollSlice(coord)
576   - Publisher.sendMessage('Set ball reference position', position=(ux, uy, uz))
  576 + self.interactor.Render()
577 577  
578 578 def ScrollSlice(self, coord):
579 579 if self.orientation == "AXIAL":
... ... @@ -817,9 +817,9 @@ class Viewer(wx.Panel):
817 817 ('Set scroll position',
818 818 self.orientation))
819 819 Publisher.subscribe(self.__update_cross_position,
820   - 'Update cross position')
821   - Publisher.subscribe(self.UpdateSlicesNavigation,
822   - 'Co-registered points')
  820 + 'Update cross position')
  821 + # Publisher.subscribe(self.UpdateSlicesNavigation,
  822 + # 'Co-registered points')
823 823 ###
824 824 # Publisher.subscribe(self.ChangeBrushColour,
825 825 # 'Add mask')
... ... @@ -1136,8 +1136,9 @@ class Viewer(wx.Panel):
1136 1136  
1137 1137 renderer.AddActor(cross_actor)
1138 1138  
1139   - def __update_cross_position(self, position):
1140   - self.cross.SetFocalPoint(position)
  1139 + def __update_cross_position(self, arg, position):
  1140 + # self.cross.SetFocalPoint(position[:3])
  1141 + self.UpdateSlicesNavigation(None, position)
1141 1142  
1142 1143 def _set_cross_visibility(self, visibility):
1143 1144 self.cross_actor.SetVisibility(visibility)
... ... @@ -1296,8 +1297,8 @@ class Viewer(wx.Panel):
1296 1297 if self.state == const.SLICE_STATE_CROSS:
1297 1298 # Update other slice's cross according to the new focal point from
1298 1299 # the actual orientation.
1299   - focal_point = self.cross.GetFocalPoint()
1300   - Publisher.sendMessage('Update cross position', position=focal_point)
  1300 + x, y, z = self.cross.GetFocalPoint()
  1301 + Publisher.sendMessage('Update cross position', arg=None, position=(x, y, z, 0., 0., 0.))
1301 1302 Publisher.sendMessage('Update slice viewer')
1302 1303 else:
1303 1304 self.interactor.Render()
... ... @@ -1509,7 +1510,6 @@ class Viewer(wx.Panel):
1509 1510  
1510 1511 def UpdateCross(self, coord):
1511 1512 self.cross.SetFocalPoint(coord)
1512   - Publisher.sendMessage('Set ball reference position', position=self.cross.GetFocalPoint())
1513 1513 Publisher.sendMessage('Co-registered points', arg=None, position=(coord[0], coord[1], coord[2], 0., 0., 0.))
1514 1514 self.OnScrollBar()
1515 1515 self.interactor.Render()
... ...
invesalius/data/viewer_volume.py
... ... @@ -19,9 +19,10 @@
19 19 # detalhes.
20 20 #--------------------------------------------------------------------------
21 21  
22   -from math import cos, sin
  22 +# from math import cos, sin
23 23 import os
24 24 import sys
  25 +import time
25 26  
26 27 import numpy as np
27 28 from numpy.core.umath_tests import inner1d
... ... @@ -36,6 +37,7 @@ from imageio import imsave
36 37  
37 38 import invesalius.constants as const
38 39 import invesalius.data.bases as bases
  40 +import invesalius.data.slice_ as sl
39 41 import invesalius.data.transformations as tr
40 42 import invesalius.data.vtk_utils as vtku
41 43 import invesalius.project as prj
... ... @@ -168,6 +170,10 @@ class Viewer(wx.Panel):
168 170 #self.pTarget = [0., 0., 0.]
169 171  
170 172 # self.obj_axes = None
  173 + self.x_actor = None
  174 + self.y_actor = None
  175 + self.z_actor = None
  176 + self.mark_actor = None
171 177  
172 178 self._mode_cross = False
173 179 self._to_show_ball = 0
... ... @@ -188,12 +194,31 @@ class Viewer(wx.Panel):
188 194 self.anglethreshold = const.COIL_ANGLES_THRESHOLD
189 195 self.distthreshold = const.COIL_COORD_THRESHOLD
190 196  
  197 + # for DTI support tests
  198 + # self.ntimes = False
  199 + # self._to_show_stream = True
  200 + # data_dir = b'C:\Users\deoliv1\OneDrive\data\dti'
  201 + # nii_path = b'sub-P0_dwi_FOD.nii'
  202 + # trk_path = os.path.join(data_dir, nii_path)
  203 + # self.tracker_FOD = Trekker.tracker(trk_path)
  204 + # proj = prj.Project()
  205 + # self.affine = np.identity(4)
  206 + self.actor_tracts = None
  207 + self.actor_peel = None
  208 + self.seed_offset = const.SEED_OFFSET
  209 + # Publisher.sendMessage('Get affine matrix')
  210 +
  211 + # initialize Trekker parameters
  212 + slic = sl.Slice()
  213 + affine = slic.affine
  214 + self.affine_vtk = vtku.numpy_to_vtkMatrix4x4(affine)
  215 +
191 216 def __bind_events(self):
192 217 Publisher.subscribe(self.LoadActor,
193 218 'Load surface actor into viewer')
194 219 Publisher.subscribe(self.RemoveActor,
195 220 'Remove surface actor from viewer')
196   - Publisher.subscribe(self.OnShowSurface, 'Show surface')
  221 + # Publisher.subscribe(self.OnShowSurface, 'Show surface')
197 222 Publisher.subscribe(self.UpdateRender,
198 223 'Render volume viewer')
199 224 Publisher.subscribe(self.ChangeBackgroundColour,
... ... @@ -225,7 +250,8 @@ class Viewer(wx.Panel):
225 250 Publisher.subscribe(self.LoadSlicePlane, 'Load slice plane')
226 251  
227 252 Publisher.subscribe(self.ResetCamClippingRange, 'Reset cam clipping range')
228   - Publisher.subscribe(self.SetVolumeCamera, 'Co-registered points')
  253 + # Publisher.subscribe(self.SetVolumeCamera, 'Update cross position')
  254 + # Publisher.subscribe(self.SetVolumeCamera, 'Co-registered points')
229 255 # Publisher.subscribe(self.SetVolumeCamera, 'Set camera in volume')
230 256 Publisher.subscribe(self.SetVolumeCameraState, 'Update volume camera state')
231 257  
... ... @@ -255,8 +281,8 @@ class Viewer(wx.Panel):
255 281  
256 282 Publisher.subscribe(self.RemoveVolume, 'Remove Volume')
257 283  
258   - Publisher.subscribe(self.SetBallReferencePosition,
259   - 'Set ball reference position')
  284 + Publisher.subscribe(self.UpdateCameraBallPosition,
  285 + 'Update cross position')
260 286 Publisher.subscribe(self._check_ball_reference, 'Enable style')
261 287 Publisher.subscribe(self._uncheck_ball_reference, 'Disable style')
262 288  
... ... @@ -283,11 +309,17 @@ class Viewer(wx.Panel):
283 309 Publisher.subscribe(self.OnUpdateObjectTargetGuide, 'Update object matrix')
284 310 Publisher.subscribe(self.OnUpdateTargetCoordinates, 'Update target')
285 311 Publisher.subscribe(self.OnRemoveTarget, 'Disable or enable coil tracker')
286   - # Publisher.subscribe(self.UpdateObjectTargetView, 'Co-registered points')
287 312 Publisher.subscribe(self.OnTargetMarkerTransparency, 'Set target transparency')
288 313 Publisher.subscribe(self.OnUpdateAngleThreshold, 'Update angle threshold')
289 314 Publisher.subscribe(self.OnUpdateDistThreshold, 'Update dist threshold')
290 315  
  316 + Publisher.subscribe(self.OnUpdateTracts, 'Update tracts')
  317 + Publisher.subscribe(self.OnRemoveTracts, 'Remove tracts')
  318 + Publisher.subscribe(self.UpdateSeedOffset, 'Update seed offset')
  319 + Publisher.subscribe(self.UpdateMarkerOffsetState, 'Update marker offset state')
  320 + Publisher.subscribe(self.UpdateMarkerOffsetPosition, 'Update marker offset')
  321 + Publisher.subscribe(self.AddPeeledSurface, 'Update peel')
  322 +
291 323 def SetStereoMode(self, mode):
292 324 ren_win = self.interactor.GetRenderWindow()
293 325  
... ... @@ -319,13 +351,23 @@ class Viewer(wx.Panel):
319 351 def _check_ball_reference(self, style):
320 352 if style == const.SLICE_STATE_CROSS:
321 353 self._mode_cross = True
322   - self._check_and_set_ball_visibility()
  354 + # self._check_and_set_ball_visibility()
  355 + self._ball_ref_visibility = True
  356 + # if self._to_show_ball:
  357 + if not self.ball_actor:
  358 + self.CreateBallReference()
  359 +
323 360 self.interactor.Render()
324 361  
325 362 def _uncheck_ball_reference(self, style):
326 363 if style == const.SLICE_STATE_CROSS:
327 364 self._mode_cross = False
328   - self.RemoveBallReference()
  365 + # self.RemoveBallReference()
  366 + self._ball_ref_visibility = False
  367 + if self.ball_actor:
  368 + self.ren.RemoveActor(self.ball_actor)
  369 + self.ball_actor = None
  370 +
329 371 self.interactor.Render()
330 372  
331 373 def OnSensors(self, probe_id, ref_id, obj_id=0):
... ... @@ -385,12 +427,12 @@ class Viewer(wx.Panel):
385 427 self.probe = self.ref = self.obj = False
386 428 self.interactor.Render()
387 429  
388   - def OnShowSurface(self, index, visibility):
389   - if visibility:
390   - self._to_show_ball += 1
391   - else:
392   - self._to_show_ball -= 1
393   - self._check_and_set_ball_visibility()
  430 + # def OnShowSurface(self, index, visibility):
  431 + # if visibility:
  432 + # self._to_show_ball += 1
  433 + # else:
  434 + # self._to_show_ball -= 1
  435 + # self._check_and_set_ball_visibility()
394 436  
395 437 def OnStartSeed(self):
396 438 self.seed_points = []
... ... @@ -485,8 +527,8 @@ class Viewer(wx.Panel):
485 527 if (volumes.GetNumberOfItems()):
486 528 self.ren.RemoveVolume(volumes.GetLastProp())
487 529 self.interactor.Render()
488   - self._to_show_ball -= 1
489   - self._check_and_set_ball_visibility()
  530 + # self._to_show_ball -= 1
  531 + # self._check_and_set_ball_visibility()
490 532  
491 533 def RemoveActors(self, actors):
492 534 "Remove a list of actors"
... ... @@ -534,11 +576,11 @@ class Viewer(wx.Panel):
534 576 Markers created by navigation tools and rendered in volume viewer.
535 577 """
536 578 self.ball_id = ball_id
537   - x, y, z = bases.flip_x(coord)
  579 + coord_flip = bases.flip_x_m(coord)[:3, 0]
538 580  
539 581 ball_ref = vtk.vtkSphereSource()
540 582 ball_ref.SetRadius(size)
541   - ball_ref.SetCenter(x, y, z)
  583 + ball_ref.SetCenter(coord_flip)
542 584  
543 585 mapper = vtk.vtkPolyDataMapper()
544 586 mapper.SetInputConnection(ball_ref.GetOutputPort())
... ... @@ -557,6 +599,35 @@ class Viewer(wx.Panel):
557 599 #self.UpdateRender()
558 600 self.Refresh()
559 601  
  602 + def add_marker(self, coord, color):
  603 + """Simplified version for creating a spherical marker in the 3D scene
  604 +
  605 + :param coord:
  606 + :param color:
  607 + :return: vtkActor
  608 + """
  609 +
  610 + # x, y, z = coord
  611 +
  612 + ball_ref = vtk.vtkSphereSource()
  613 + ball_ref.SetRadius(2)
  614 + ball_ref.SetCenter(coord)
  615 +
  616 + mapper = vtk.vtkPolyDataMapper()
  617 + mapper.SetInputConnection(ball_ref.GetOutputPort())
  618 +
  619 + prop = vtk.vtkProperty()
  620 + prop.SetColor(color)
  621 +
  622 + actor = vtk.vtkActor()
  623 + actor.SetMapper(mapper)
  624 + actor.SetProperty(prop)
  625 + actor.GetProperty().SetOpacity(.5)
  626 +
  627 + # ren.AddActor(actor)
  628 +
  629 + return actor
  630 +
560 631 def HideAllMarkers(self, indexes):
561 632 ballid = indexes
562 633 for i in range(0, ballid):
... ... @@ -610,7 +681,6 @@ class Viewer(wx.Panel):
610 681 self.staticballs[index].GetProperty().SetColor(color)
611 682 self.Refresh()
612 683  
613   -
614 684 def OnTargetMarkerTransparency(self, status, index):
615 685 if status:
616 686 self.staticballs[index].GetProperty().SetOpacity(1)
... ... @@ -1104,7 +1174,6 @@ class Viewer(wx.Panel):
1104 1174 cam.SetFocalPoint(cam_focus)
1105 1175 cam.SetPosition(cam_pos)
1106 1176  
1107   -
1108 1177 def CreateBallReference(self):
1109 1178 """
1110 1179 Red sphere on volume visualization to reference center of
... ... @@ -1129,31 +1198,14 @@ class Viewer(wx.Panel):
1129 1198  
1130 1199 self.ren.AddActor(self.ball_actor)
1131 1200  
1132   - def ActivateBallReference(self):
1133   - self._mode_cross = True
1134   - self._ball_ref_visibility = True
1135   - if self._to_show_ball:
1136   - if not self.ball_actor:
1137   - self.CreateBallReference()
1138   -
1139   - def RemoveBallReference(self):
1140   - self._mode_cross = False
1141   - self._ball_ref_visibility = False
1142   - if self.ball_actor:
1143   - self.ren.RemoveActor(self.ball_actor)
1144   - self.ball_actor = None
  1201 + def UpdateCameraBallPosition(self, arg, position):
  1202 + coord_flip = bases.flip_x_m(position[:3])[:3, 0]
  1203 + self.ball_actor.SetPosition(coord_flip)
  1204 + self.SetVolumeCamera(coord_flip)
1145 1205  
1146 1206 def SetBallReferencePosition(self, position):
1147   - if self._to_show_ball:
1148   - if not self.ball_actor:
1149   - self.ActivateBallReference()
1150   -
1151   - coord = position
1152   - x, y, z = bases.flip_x(coord)
1153   - self.ball_actor.SetPosition(x, y, z)
1154   -
1155   - else:
1156   - self.RemoveBallReference()
  1207 + coord_flip = bases.flip_x_m(position[:3])[:3, 0]
  1208 + self.ball_actor.SetPosition(coord_flip)
1157 1209  
1158 1210 def CreateObjectPolyData(self, filename):
1159 1211 """
... ... @@ -1225,7 +1277,14 @@ class Viewer(wx.Panel):
1225 1277 self.obj_actor.GetProperty().SetOpacity(.4)
1226 1278 self.obj_actor.SetVisibility(0)
1227 1279  
  1280 + self.x_actor = self.add_line([0., 0., 0.], [1., 0., 0.], color=[.0, .0, 1.0])
  1281 + self.y_actor = self.add_line([0., 0., 0.], [0., 1., 0.], color=[.0, 1.0, .0])
  1282 + self.z_actor = self.add_line([0., 0., 0.], [0., 0., 1.], color=[1.0, .0, .0])
  1283 +
1228 1284 self.ren.AddActor(self.obj_actor)
  1285 + self.ren.AddActor(self.x_actor)
  1286 + self.ren.AddActor(self.y_actor)
  1287 + self.ren.AddActor(self.z_actor)
1229 1288  
1230 1289 # self.obj_axes = vtk.vtkAxesActor()
1231 1290 # self.obj_axes.SetShaftTypeToCylinder()
... ... @@ -1236,25 +1295,88 @@ class Viewer(wx.Panel):
1236 1295  
1237 1296 # self.ren.AddActor(self.obj_axes)
1238 1297  
1239   - def OnNavigationStatus(self, status):
1240   - self.nav_status = status
  1298 + def add_line(self, p1, p2, color=[0.0, 0.0, 1.0]):
  1299 + line = vtk.vtkLineSource()
  1300 + line.SetPoint1(p1)
  1301 + line.SetPoint2(p2)
  1302 +
  1303 + mapper = vtk.vtkPolyDataMapper()
  1304 + mapper.SetInputConnection(line.GetOutputPort())
  1305 +
  1306 + actor = vtk.vtkActor()
  1307 + actor.SetMapper(mapper)
  1308 + actor.GetProperty().SetColor(color)
  1309 +
  1310 + return actor
  1311 +
  1312 + def AddPeeledSurface(self, flag, actor):
  1313 + if self.actor_peel:
  1314 + self.ren.RemoveActor(self.actor_peel)
  1315 + self.actor_peel = None
  1316 + if flag and actor:
  1317 + self.ren.AddActor(actor)
  1318 + self.actor_peel = actor
  1319 + self.Refresh()
  1320 +
  1321 + def OnNavigationStatus(self, nav_status, vis_status):
  1322 + self.nav_status = nav_status
  1323 + self.tracts_status = vis_status[1]
1241 1324 self.pTarget = self.CenterOfMass()
1242   - if self.obj_actor and self.nav_status:
1243   - self.obj_actor.SetVisibility(self.obj_state)
1244   - if not self.obj_state:
1245   - self.Refresh()
  1325 +
  1326 + if self.nav_status:
  1327 + if self.obj_actor:
  1328 + self.obj_actor.SetVisibility(self.obj_state)
  1329 + self.x_actor.SetVisibility(self.obj_state)
  1330 + self.y_actor.SetVisibility(self.obj_state)
  1331 + self.z_actor.SetVisibility(self.obj_state)
  1332 +
  1333 + self.Refresh()
  1334 +
  1335 + def UpdateSeedOffset(self, data):
  1336 + self.seed_offset = data
  1337 +
  1338 + def UpdateMarkerOffsetState(self, create=False):
  1339 + if create:
  1340 + if not self.mark_actor:
  1341 + self.mark_actor = self.add_marker([0., 0., 0.], color=[0., 1., 1.])
  1342 + self.ren.AddActor(self.mark_actor)
  1343 + else:
  1344 + if self.mark_actor:
  1345 + self.ren.RemoveActor(self.mark_actor)
  1346 + self.mark_actor = None
  1347 + self.Refresh()
  1348 +
  1349 + def CreateMarkerOffset(self):
  1350 + self.mark_actor = self.add_marker([0., 0., 0.], color=[0., 1., 1.])
  1351 + self.ren.AddActor(self.mark_actor)
  1352 + self.Refresh()
  1353 +
  1354 + def UpdateMarkerOffsetPosition(self, coord_offset):
  1355 + self.mark_actor.SetPosition(coord_offset)
  1356 + self.Refresh()
1246 1357  
1247 1358 def UpdateObjectOrientation(self, m_img, coord):
1248   - m_img[:3, -1] = np.asmatrix(bases.flip_x_m((m_img[0, -1], m_img[1, -1], m_img[2, -1]))).reshape([3, 1])
  1359 + # print("Update object orientation")
1249 1360  
1250   - m_img_vtk = vtk.vtkMatrix4x4()
  1361 + m_img_flip = m_img.copy()
  1362 + m_img_flip[1, -1] = -m_img_flip[1, -1]
1251 1363  
1252   - for row in range(0, 4):
1253   - for col in range(0, 4):
1254   - m_img_vtk.SetElement(row, col, m_img[row, col])
  1364 + # translate coregistered coordinate to display a marker where Trekker seed is computed
  1365 + # coord_offset = m_img_flip[:3, -1] - self.seed_offset * m_img_flip[:3, 2]
  1366 +
  1367 + # print("m_img copy in viewer_vol: {}".format(m_img_copy))
  1368 +
  1369 + # m_img[:3, 0] is from posterior to anterior direction of the coil
  1370 + # m_img[:3, 1] is from left to right direction of the coil
  1371 + # m_img[:3, 2] is from bottom to up direction of the coil
  1372 +
  1373 + m_img_vtk = vtku.numpy_to_vtkMatrix4x4(m_img_flip)
1255 1374  
1256 1375 self.obj_actor.SetUserMatrix(m_img_vtk)
1257 1376 # self.obj_axes.SetUserMatrix(m_rot_vtk)
  1377 + self.x_actor.SetUserMatrix(m_img_vtk)
  1378 + self.y_actor.SetUserMatrix(m_img_vtk)
  1379 + self.z_actor.SetUserMatrix(m_img_vtk)
1258 1380  
1259 1381 self.Refresh()
1260 1382  
... ... @@ -1267,13 +1389,42 @@ class Viewer(wx.Panel):
1267 1389 else:
1268 1390 if self.obj_actor:
1269 1391 self.ren.RemoveActor(self.obj_actor)
  1392 + self.ren.RemoveActor(self.x_actor)
  1393 + self.ren.RemoveActor(self.y_actor)
  1394 + self.ren.RemoveActor(self.z_actor)
  1395 + self.ren.RemoveActor(self.mark_actor)
1270 1396 self.obj_actor = None
  1397 + self.x_actor = None
  1398 + self.y_actor = None
  1399 + self.z_actor = None
  1400 + self.mark_actor = None
1271 1401 self.Refresh()
1272 1402  
1273 1403 def UpdateShowObjectState(self, state):
1274 1404 self.obj_state = state
1275 1405 if self.obj_actor and not self.obj_state:
1276 1406 self.obj_actor.SetVisibility(self.obj_state)
  1407 + self.x_actor.SetVisibility(self.obj_state)
  1408 + self.y_actor.SetVisibility(self.obj_state)
  1409 + self.z_actor.SetVisibility(self.obj_state)
  1410 +
  1411 + self.Refresh()
  1412 +
  1413 + def OnUpdateTracts(self, evt=None, flag=None, actor=None, root=None, affine_vtk=None, count=0):
  1414 + mapper = vtk.vtkCompositePolyDataMapper2()
  1415 + mapper.SetInputDataObject(root)
  1416 +
  1417 + self.actor_tracts = vtk.vtkActor()
  1418 + self.actor_tracts.SetMapper(mapper)
  1419 + self.actor_tracts.SetUserMatrix(affine_vtk)
  1420 +
  1421 + self.ren.AddActor(self.actor_tracts)
  1422 + self.Refresh()
  1423 +
  1424 + def OnRemoveTracts(self):
  1425 + if self.actor_tracts:
  1426 + self.ren.RemoveActor(self.actor_tracts)
  1427 + self.actor_tracts = None
1277 1428 self.Refresh()
1278 1429  
1279 1430 def __bind_events_wx(self):
... ... @@ -1448,10 +1599,12 @@ class Viewer(wx.Panel):
1448 1599 def SetVolumeCameraState(self, camera_state):
1449 1600 self.camera_state = camera_state
1450 1601  
1451   - def SetVolumeCamera(self, arg, position):
  1602 + # def SetVolumeCamera(self, arg, position):
  1603 + def SetVolumeCamera(self, cam_focus):
1452 1604 if self.camera_state:
1453 1605 # TODO: exclude dependency on initial focus
1454   - cam_focus = np.array(bases.flip_x(position[:3]))
  1606 + # cam_focus = np.array(bases.flip_x(position[:3]))
  1607 + # cam_focus = np.array(bases.flip_x(position))
1455 1608 cam = self.ren.GetActiveCamera()
1456 1609  
1457 1610 if self.initial_focus is None:
... ... @@ -1545,16 +1698,16 @@ class Viewer(wx.Panel):
1545 1698 def OnShowRaycasting(self):
1546 1699 if not self.raycasting_volume:
1547 1700 self.raycasting_volume = True
1548   - self._to_show_ball += 1
1549   - self._check_and_set_ball_visibility()
  1701 + # self._to_show_ball += 1
  1702 + # self._check_and_set_ball_visibility()
1550 1703 if self.on_wl:
1551 1704 self.text.Show()
1552 1705  
1553 1706 def OnHideRaycasting(self):
1554 1707 self.raycasting_volume = False
1555 1708 self.text.Hide()
1556   - self._to_show_ball -= 1
1557   - self._check_and_set_ball_visibility()
  1709 + # self._to_show_ball -= 1
  1710 + # self._check_and_set_ball_visibility()
1558 1711  
1559 1712 def OnSize(self, evt):
1560 1713 self.UpdateRender()
... ... @@ -1581,16 +1734,16 @@ class Viewer(wx.Panel):
1581 1734  
1582 1735 #self.ShowOrientationCube()
1583 1736 self.interactor.Render()
1584   - self._to_show_ball += 1
1585   - self._check_and_set_ball_visibility()
  1737 + # self._to_show_ball += 1
  1738 + # self._check_and_set_ball_visibility()
1586 1739  
1587 1740 def RemoveActor(self, actor):
1588 1741 utils.debug("RemoveActor")
1589 1742 ren = self.ren
1590 1743 ren.RemoveActor(actor)
1591 1744 self.interactor.Render()
1592   - self._to_show_ball -= 1
1593   - self._check_and_set_ball_visibility()
  1745 + # self._to_show_ball -= 1
  1746 + # self._check_and_set_ball_visibility()
1594 1747  
1595 1748 def RemoveAllActor(self):
1596 1749 utils.debug("RemoveAllActor")
... ... @@ -1602,7 +1755,8 @@ class Viewer(wx.Panel):
1602 1755  
1603 1756 def LoadVolume(self, volume, colour, ww, wl):
1604 1757 self.raycasting_volume = True
1605   - self._to_show_ball += 1
  1758 + # self._to_show_ball += 1
  1759 + # self._check_and_set_ball_visibility()
1606 1760  
1607 1761 self.light = self.ren.GetLights().GetNextItem()
1608 1762  
... ... @@ -1622,15 +1776,14 @@ class Viewer(wx.Panel):
1622 1776 self.ren.ResetCamera()
1623 1777 self.ren.ResetCameraClippingRange()
1624 1778  
1625   - self._check_and_set_ball_visibility()
1626 1779 self.UpdateRender()
1627 1780  
1628 1781 def UnloadVolume(self, volume):
1629 1782 self.ren.RemoveVolume(volume)
1630 1783 del volume
1631 1784 self.raycasting_volume = False
1632   - self._to_show_ball -= 1
1633   - self._check_and_set_ball_visibility()
  1785 + # self._to_show_ball -= 1
  1786 + # self._check_and_set_ball_visibility()
1634 1787  
1635 1788 def OnSetViewAngle(self, view):
1636 1789 self.SetViewAngle(view)
... ... @@ -1777,16 +1930,17 @@ class Viewer(wx.Panel):
1777 1930 self.SetViewAngle(const.VOL_ISO)
1778 1931 self.repositioned_coronal_plan = 1
1779 1932  
1780   - def _check_and_set_ball_visibility(self):
1781   - #TODO: When creating Raycasting volume and cross is pressed, it is not
1782   - # automatically creating the ball reference.
1783   - if self._mode_cross:
1784   - if self._to_show_ball > 0 and not self._ball_ref_visibility:
1785   - self.ActivateBallReference()
1786   - self.interactor.Render()
1787   - elif not self._to_show_ball and self._ball_ref_visibility:
1788   - self.RemoveBallReference()
1789   - self.interactor.Render()
  1933 + # def _check_and_set_ball_visibility(self):
  1934 + # #TODO: When creating Raycasting volume and cross is pressed, it is not
  1935 + # # automatically creating the ball reference.
  1936 + # # print("mode_cross, show_ball, ball_vis ", self._mode_cross, self._to_show_ball, self._ball_ref_visibility)
  1937 + # if self._mode_cross:
  1938 + # if self._to_show_ball > 0 and not self._ball_ref_visibility:
  1939 + # self.ActivateBallReference()
  1940 + # self.interactor.Render()
  1941 + # elif not self._to_show_ball and self._ball_ref_visibility:
  1942 + # self.RemoveBallReference()
  1943 + # self.interactor.Render()
1790 1944  
1791 1945 class SlicePlane:
1792 1946 def __init__(self):
... ...
invesalius/data/vtk_utils.py
... ... @@ -287,3 +287,21 @@ class TextZero(object):
287 287 # font.SetWeight(wx.FONTWEIGHT_BOLD)
288 288 font.SetSymbolicSize(self.symbolic_syze)
289 289 canvas.draw_text(self.text, (x, y), font=font)
  290 +
  291 +
  292 +def numpy_to_vtkMatrix4x4(affine):
  293 + """
  294 + Convert a numpy 4x4 array to a vtk 4x4 matrix
  295 + :param affine: 4x4 array
  296 + :return: vtkMatrix4x4 object representing the affine
  297 + """
  298 + # test for type and shape of affine matrix
  299 + # assert isinstance(affine, np.ndarray)
  300 + assert affine.shape == (4, 4)
  301 +
  302 + affine_vtk = vtk.vtkMatrix4x4()
  303 + for row in range(0, 4):
  304 + for col in range(0, 4):
  305 + affine_vtk.SetElement(row, col, affine[row, col])
  306 +
  307 + return affine_vtk
... ...
invesalius/gui/dialogs.py
... ... @@ -387,6 +387,15 @@ def ShowImportOtherFilesDialog(id_type):
387 387 if id_type == const.ID_NIFTI_IMPORT:
388 388 dlg.SetMessage(_("Import NIFTi 1 file"))
389 389 dlg.SetWildcard(WILDCARD_NIFTI)
  390 + elif id_type == const.ID_TREKKER_MASK:
  391 + dlg.SetMessage(_("Import Trekker mask"))
  392 + dlg.SetWildcard(WILDCARD_NIFTI)
  393 + elif id_type == const.ID_TREKKER_IMG:
  394 + dlg.SetMessage(_("Import Trekker anatomical image"))
  395 + dlg.SetWildcard(WILDCARD_NIFTI)
  396 + elif id_type == const.ID_TREKKER_FOD:
  397 + dlg.SetMessage(_("Import Trekker FOD"))
  398 + dlg.SetWildcard(WILDCARD_NIFTI)
390 399 elif id_type == const.ID_PARREC_IMPORT:
391 400 dlg.SetMessage(_("Import PAR/REC file"))
392 401 dlg.SetWildcard(WILDCARD_PARREC)
... ... @@ -508,79 +517,11 @@ def ShowSaveAsProjectDialog(default_filename=None):
508 517 return filename, wildcard == INV_COMPRESSED
509 518  
510 519  
511   -# Dialog for neuronavigation markers
512   -def ShowSaveMarkersDialog(default_filename=None):
513   - current_dir = os.path.abspath(".")
514   - dlg = wx.FileDialog(None,
515   - _("Save markers as..."), # title
516   - "", # last used directory
517   - default_filename,
518   - _("Markers files (*.mks)|*.mks"),
519   - wx.FD_SAVE | wx.FD_OVERWRITE_PROMPT)
520   - # dlg.SetFilterIndex(0) # default is VTI
521   -
522   - filename = None
523   - try:
524   - if dlg.ShowModal() == wx.ID_OK:
525   - filename = dlg.GetPath()
526   - ok = 1
527   - else:
528   - ok = 0
529   - except(wx._core.PyAssertionError): # TODO: fix win64
530   - filename = dlg.GetPath()
531   - ok = 1
532   -
533   - if (ok):
534   - extension = "mks"
535   - if sys.platform != 'win32':
536   - if filename.split(".")[-1] != extension:
537   - filename = filename + "." + extension
538   -
539   - os.chdir(current_dir)
540   - return filename
541   -
542   -def ShowSaveCoordsDialog(default_filename=None):
543   - current_dir = os.path.abspath(".")
544   - dlg = wx.FileDialog(None,
545   - _("Save coords as..."), # title
546   - "", # last used directory
547   - default_filename,
548   - _("Coordinates files (*.csv)|*.csv"),
549   - wx.FD_SAVE | wx.FD_OVERWRITE_PROMPT)
550   - # dlg.SetFilterIndex(0) # default is VTI
551   -
552   - filename = None
553   - try:
554   - if dlg.ShowModal() == wx.ID_OK:
555   - filename = dlg.GetPath()
556   - ok = 1
557   - else:
558   - ok = 0
559   - except(wx._core.PyAssertionError): # TODO: fix win64
560   - filename = dlg.GetPath()
561   - ok = 1
  520 +def ShowLoadSaveDialog(message=_(u"Load File"), current_dir=os.path.abspath("."), style=wx.FD_OPEN | wx.FD_CHANGE_DIR,
  521 + wildcard=_("Registration files (*.obr)|*.obr"), default_filename="", save_ext=None):
562 522  
563   - if (ok):
564   - extension = "csv"
565   - if sys.platform != 'win32':
566   - if filename.split(".")[-1] != extension:
567   - filename = filename + "." + extension
568   -
569   - os.chdir(current_dir)
570   - return filename
571   -
572   -
573   -def ShowLoadMarkersDialog():
574   - current_dir = os.path.abspath(".")
575   -
576   - dlg = wx.FileDialog(None, message=_("Load markers"),
577   - defaultDir="",
578   - defaultFile="",
579   - wildcard=_("Markers files (*.mks)|*.mks"),
580   - style=wx.FD_OPEN|wx.FD_CHANGE_DIR)
581   -
582   - # inv3 filter is default
583   - dlg.SetFilterIndex(0)
  523 + dlg = wx.FileDialog(None, message=message, defaultDir="", defaultFile=default_filename,
  524 + wildcard=wildcard, style=style)
584 525  
585 526 # Show the dialog and retrieve the user response. If it is the OK response,
586 527 # process the data.
... ... @@ -589,73 +530,25 @@ def ShowLoadMarkersDialog():
589 530 if dlg.ShowModal() == wx.ID_OK:
590 531 # This returns a Python list of files that were selected.
591 532 filepath = dlg.GetPath()
  533 + ok_press = 1
  534 + else:
  535 + ok_press = 0
592 536 except(wx._core.PyAssertionError): # FIX: win64
593 537 filepath = dlg.GetPath()
  538 + ok_press = 1
594 539  
595   - # Destroy the dialog. Don't do this until you are done with it!
596   - # BAD things can happen otherwise!
597   - dlg.Destroy()
598   - os.chdir(current_dir)
599   - return filepath
600   -
601   -
602   -def ShowSaveRegistrationDialog(default_filename=None):
603   - current_dir = os.path.abspath(".")
604   - dlg = wx.FileDialog(None,
605   - _("Save object registration as..."), # title
606   - "", # last used directory
607   - default_filename,
608   - _("Registration files (*.obr)|*.obr"),
609   - wx.FD_SAVE | wx.FD_OVERWRITE_PROMPT)
610   - # dlg.SetFilterIndex(0) # default is VTI
611   -
612   - filename = None
613   - try:
614   - if dlg.ShowModal() == wx.ID_OK:
615   - filename = dlg.GetPath()
616   - ok = 1
617   - else:
618   - ok = 0
619   - except(wx._core.PyAssertionError): # TODO: fix win64
620   - filename = dlg.GetPath()
621   - ok = 1
622   -
623   - if (ok):
624   - extension = "obr"
  540 + # fix the extension if set different than expected
  541 + if save_ext and ok_press:
  542 + extension = save_ext
625 543 if sys.platform != 'win32':
626   - if filename.split(".")[-1] != extension:
627   - filename = filename + "." + extension
628   -
629   - os.chdir(current_dir)
630   - return filename
631   -
632   -
633   -def ShowLoadRegistrationDialog():
634   - current_dir = os.path.abspath(".")
635   -
636   - dlg = wx.FileDialog(None, message=_("Load object registration"),
637   - defaultDir="",
638   - defaultFile="",
639   - wildcard=_("Registration files (*.obr)|*.obr"),
640   - style=wx.FD_OPEN|wx.FD_CHANGE_DIR)
641   -
642   - # inv3 filter is default
643   - dlg.SetFilterIndex(0)
644   -
645   - # Show the dialog and retrieve the user response. If it is the OK response,
646   - # process the data.
647   - filepath = None
648   - try:
649   - if dlg.ShowModal() == wx.ID_OK:
650   - # This returns a Python list of files that were selected.
651   - filepath = dlg.GetPath()
652   - except(wx._core.PyAssertionError): # FIX: win64
653   - filepath = dlg.GetPath()
  544 + if filepath.split(".")[-1] != extension:
  545 + filepath = filepath + "." + extension
654 546  
655 547 # Destroy the dialog. Don't do this until you are done with it!
656 548 # BAD things can happen otherwise!
657 549 dlg.Destroy()
658 550 os.chdir(current_dir)
  551 +
659 552 return filepath
660 553  
661 554  
... ... @@ -920,31 +813,9 @@ def SurfaceSelectionRequiredForDuplication():
920 813  
921 814  
922 815 # Dialogs for neuronavigation mode
923   -def InvalidFiducials():
924   - msg = _("Fiducials are invalid. Select all coordinates.")
925   - if sys.platform == 'darwin':
926   - dlg = wx.MessageDialog(None, "", msg,
927   - wx.ICON_INFORMATION | wx.OK)
928   - else:
929   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
930   - wx.ICON_INFORMATION | wx.OK)
931   - dlg.ShowModal()
932   - dlg.Destroy()
933   -
  816 +# ----------------------------------
934 817  
935   -def InvalidObjectRegistration():
936   - msg = _("Perform coil registration before navigation.")
937   - if sys.platform == 'darwin':
938   - dlg = wx.MessageDialog(None, "", msg,
939   - wx.ICON_INFORMATION | wx.OK)
940   - else:
941   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
942   - wx.ICON_INFORMATION | wx.OK)
943   - dlg.ShowModal()
944   - dlg.Destroy()
945   -
946   -
947   -def NavigationTrackerWarning(trck_id, lib_mode):
  818 +def ShowNavigationTrackerWarning(trck_id, lib_mode):
948 819 """
949 820 Spatial Tracker connection error
950 821 """
... ... @@ -976,85 +847,46 @@ def NavigationTrackerWarning(trck_id, lib_mode):
976 847 dlg.Destroy()
977 848  
978 849  
979   -def InvalidMarkersFile():
980   - msg = _("The TXT file is invalid.")
981   - if sys.platform == 'darwin':
982   - dlg = wx.MessageDialog(None, "", msg,
983   - wx.ICON_INFORMATION | wx.OK)
984   - else:
985   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
986   - wx.ICON_INFORMATION | wx.OK)
987   - dlg.ShowModal()
988   - dlg.Destroy()
989   -
990   -
991   -def NoMarkerSelected():
992   - msg = _("No data selected")
  850 +def ShowEnterMarkerID(default):
  851 + msg = _("Edit marker ID")
993 852 if sys.platform == 'darwin':
994   - dlg = wx.MessageDialog(None, "", msg,
995   - wx.ICON_INFORMATION | wx.OK)
  853 + dlg = wx.TextEntryDialog(None, "", msg, defaultValue=default)
996 854 else:
997   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
998   - wx.ICON_INFORMATION | wx.OK)
  855 + dlg = wx.TextEntryDialog(None, msg, "InVesalius 3", value=default)
999 856 dlg.ShowModal()
  857 + result = dlg.GetValue()
1000 858 dlg.Destroy()
  859 + return result
1001 860  
1002 861  
1003   -def DeleteAllMarkers():
1004   - msg = _("Do you really want to delete all markers?")
  862 +def ShowConfirmationDialog(msg=_('Proceed?')):
  863 + # msg = _("Do you want to delete all markers?")
1005 864 if sys.platform == 'darwin':
1006 865 dlg = wx.MessageDialog(None, "", msg,
1007 866 wx.OK | wx.CANCEL | wx.ICON_QUESTION)
1008 867 else:
1009   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
  868 + dlg = wx.MessageDialog(None, msg, "InVesalius 3",
1010 869 wx.OK | wx.CANCEL | wx.ICON_QUESTION)
1011 870 result = dlg.ShowModal()
1012 871 dlg.Destroy()
1013 872 return result
1014 873  
1015   -def DeleteTarget():
1016   - msg = _("Target deleted")
1017   - if sys.platform == 'darwin':
1018   - dlg = wx.MessageDialog(None, "", msg,
1019   - wx.ICON_INFORMATION | wx.OK)
1020   - else:
1021   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
1022   - wx.ICON_INFORMATION | wx.OK)
1023   - dlg.ShowModal()
1024   - dlg.Destroy()
1025 874  
1026   -def NewTarget():
1027   - msg = _("New target selected")
1028   - if sys.platform == 'darwin':
1029   - dlg = wx.MessageDialog(None, "", msg,
1030   - wx.ICON_INFORMATION | wx.OK)
1031   - else:
1032   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
1033   - wx.ICON_INFORMATION | wx.OK)
1034   - dlg.ShowModal()
1035   - dlg.Destroy()
  875 +def ShowColorDialog(color_current):
  876 + cdata = wx.ColourData()
  877 + cdata.SetColour(wx.Colour(color_current))
  878 + dlg = wx.ColourDialog(None, data=cdata)
  879 + dlg.GetColourData().SetChooseFull(True)
1036 880  
1037   -def InvalidTargetID():
1038   - msg = _("Sorry, you cannot use 'TARGET' ID")
1039   - if sys.platform == 'darwin':
1040   - dlg = wx.MessageDialog(None, "", msg,
1041   - wx.ICON_INFORMATION | wx.OK)
  881 + if dlg.ShowModal() == wx.ID_OK:
  882 + color_new = dlg.GetColourData().GetColour().Get(includeAlpha=False)
1042 883 else:
1043   - dlg = wx.MessageDialog(None, msg, "InVesalius 3 - Neuronavigator",
1044   - wx.ICON_INFORMATION | wx.OK)
1045   - dlg.ShowModal()
1046   - dlg.Destroy()
  884 + color_new = None
1047 885  
1048   -def EnterMarkerID(default):
1049   - msg = _("Edit marker ID")
1050   - if sys.platform == 'darwin':
1051   - dlg = wx.TextEntryDialog(None, "", msg, defaultValue=default)
1052   - else:
1053   - dlg = wx.TextEntryDialog(None, msg, "InVesalius 3", value=default)
1054   - dlg.ShowModal()
1055   - result = dlg.GetValue()
1056 886 dlg.Destroy()
1057   - return result
  887 + return color_new
  888 +
  889 +# ----------------------------------
1058 890  
1059 891  
1060 892 class NewMask(wx.Dialog):
... ... @@ -3569,7 +3401,10 @@ class ObjectCalibrationDialog(wx.Dialog):
3569 3401 else:
3570 3402 coord = coord_raw[0, :]
3571 3403 else:
3572   - NavigationTrackerWarning(0, 'choose')
  3404 + ShowNavigationTrackerWarning(0, 'choose')
  3405 +
  3406 + if btn_id == 3:
  3407 + coord = np.zeros([6,])
3573 3408  
3574 3409 # Update text controls with tracker coordinates
3575 3410 if coord is not None or np.sum(coord) != 0.0:
... ... @@ -3582,7 +3417,7 @@ class ObjectCalibrationDialog(wx.Dialog):
3582 3417 self.ball_actors[btn_id].GetProperty().SetColor(0.0, 1.0, 0.0)
3583 3418 self.Refresh()
3584 3419 else:
3585   - NavigationTrackerWarning(0, 'choose')
  3420 + ShowNavigationTrackerWarning(0, 'choose')
3586 3421  
3587 3422 def OnChoiceRefMode(self, evt):
3588 3423 # When ref mode is changed the tracker coordinates are set to nan
... ...
invesalius/gui/frame.py
... ... @@ -1185,12 +1185,13 @@ class MenuBar(wx.MenuBar):
1185 1185 else:
1186 1186 self.FindItemById(const.ID_GOTO_COORD).Enable(False)
1187 1187  
1188   - def OnEnableNavigation(self, status):
  1188 + def OnEnableNavigation(self, nav_status, vis_status):
1189 1189 """
1190 1190 Disable mode menu when navigation is on.
1191   - :param status: Navigation status
  1191 + :param nav_status: Navigation status
  1192 + :param vis_status: Status of items visualization during navigation
1192 1193 """
1193   - value = status
  1194 + value = nav_status
1194 1195 if value:
1195 1196 self.FindItemById(const.ID_MODE_NAVIGATION).Enable(False)
1196 1197 else:
... ...
invesalius/gui/task_navigator.py
... ... @@ -18,10 +18,13 @@
18 18 #--------------------------------------------------------------------------
19 19  
20 20 from functools import partial
  21 +# import os
  22 +import queue
21 23 import sys
22   -import os
  24 +import threading
23 25  
24 26 import numpy as np
  27 +# import Trekker
25 28 import wx
26 29  
27 30 try:
... ... @@ -31,22 +34,23 @@ except ImportError:
31 34 import wx.lib.hyperlink as hl
32 35 import wx.lib.foldpanelbar as fpb
33 36  
  37 +import wx.lib.colourselect as csel
34 38 import wx.lib.masked.numctrl
35 39 from pubsub import pub as Publisher
36   -import wx.lib.colourselect as csel
37   -import wx.lib.platebtn as pbtn
38   -
39   -from math import cos, sin, pi
40 40 from time import sleep
41 41  
42   -import invesalius.data.transformations as tr
43 42 import invesalius.constants as const
44 43 import invesalius.data.bases as db
  44 +# import invesalius.data.brainmesh_handler as brain
45 45 import invesalius.data.coordinates as dco
46 46 import invesalius.data.coregistration as dcr
  47 +import invesalius.data.slice_ as sl
47 48 import invesalius.data.trackers as dt
  49 +import invesalius.data.tractography as dti
  50 +import invesalius.data.transformations as tr
48 51 import invesalius.data.trigger as trig
49 52 import invesalius.data.record_coords as rec
  53 +import invesalius.data.vtk_utils as vtk_utils
50 54 import invesalius.gui.dialogs as dlg
51 55 from invesalius import utils
52 56  
... ... @@ -136,7 +140,7 @@ class InnerFoldPanel(wx.Panel):
136 140 # Study this.
137 141  
138 142 fold_panel = fpb.FoldPanelBar(self, -1, wx.DefaultPosition,
139   - (10, 290), 0, fpb.FPB_SINGLE_FOLD)
  143 + (10, 310), 0, fpb.FPB_SINGLE_FOLD)
140 144 # Fold panel style
141 145 style = fpb.CaptionBarStyle()
142 146 style.SetCaptionStyle(fpb.CAPTIONBAR_GRADIENT_V)
... ... @@ -161,15 +165,23 @@ class InnerFoldPanel(wx.Panel):
161 165 leftSpacing=0, rightSpacing=0)
162 166  
163 167 # Fold 3 - Markers panel
164   - item = fold_panel.AddFoldPanel(_("Extra tools"), collapsed=True)
  168 + item = fold_panel.AddFoldPanel(_("Markers"), collapsed=True)
165 169 mtw = MarkersPanel(item)
166 170  
167 171 fold_panel.ApplyCaptionStyle(item, style)
168 172 fold_panel.AddFoldPanelWindow(item, mtw, spacing= 0,
169 173 leftSpacing=0, rightSpacing=0)
170 174  
171   - # Fold 4 - DBS
  175 + # Fold 4 - Tractography panel
  176 + item = fold_panel.AddFoldPanel(_("Tractography"), collapsed=True)
  177 + otw = TractographyPanel(item)
172 178  
  179 + fold_panel.ApplyCaptionStyle(item, style)
  180 + fold_panel.AddFoldPanelWindow(item, otw, spacing=0,
  181 + leftSpacing=0, rightSpacing=0)
  182 + item.Hide()
  183 +
  184 + # Fold 5 - DBS
173 185 self.dbs_item = fold_panel.AddFoldPanel(_("Deep Brain Stimulation"), collapsed=True)
174 186 dtw = DbsPanel(self.dbs_item) #Atribuir nova var, criar panel
175 187  
... ... @@ -239,8 +251,8 @@ class InnerFoldPanel(wx.Panel):
239 251 def OnHideDbs(self):
240 252 self.dbs_item.Hide()
241 253  
242   - def OnCheckStatus(self, status):
243   - if status:
  254 + def OnCheckStatus(self, nav_status, vis_status):
  255 + if nav_status:
244 256 self.checktrigger.Enable(False)
245 257 self.checkobj.Enable(False)
246 258 else:
... ... @@ -295,6 +307,23 @@ class NeuronavigationPanel(wx.Panel):
295 307 self.obj_reg = None
296 308 self.obj_reg_status = False
297 309 self.track_obj = False
  310 + self.event = threading.Event()
  311 +
  312 + self.coord_queue = QueueCustom(maxsize=1)
  313 + # self.visualization_queue = QueueCustom(maxsize=1)
  314 + self.trigger_queue = QueueCustom(maxsize=1)
  315 + self.coord_tracts_queue = QueueCustom(maxsize=1)
  316 + self.tracts_queue = QueueCustom(maxsize=1)
  317 +
  318 + # Tractography parameters
  319 + self.trk_inp = None
  320 + self.trekker = None
  321 + self.n_threads = None
  322 + self.view_tracts = False
  323 + self.n_tracts = const.N_TRACTS
  324 + self.seed_offset = const.SEED_OFFSET
  325 + self.seed_radius = const.SEED_RADIUS
  326 + self.sleep_nav = const.SLEEP_NAVIGATION
298 327  
299 328 self.tracker_id = const.DEFAULT_TRACKER
300 329 self.ref_mode_id = const.DEFAULT_REF_MODE
... ... @@ -405,10 +434,17 @@ class NeuronavigationPanel(wx.Panel):
405 434 Publisher.subscribe(self.LoadImageFiducials, 'Load image fiducials')
406 435 Publisher.subscribe(self.UpdateTriggerState, 'Update trigger state')
407 436 Publisher.subscribe(self.UpdateTrackObjectState, 'Update track object state')
408   - Publisher.subscribe(self.UpdateImageCoordinates, 'Set ball reference position')
  437 + Publisher.subscribe(self.UpdateImageCoordinates, 'Update cross position')
409 438 Publisher.subscribe(self.OnDisconnectTracker, 'Disconnect tracker')
410 439 Publisher.subscribe(self.UpdateObjectRegistration, 'Update object registration')
411 440 Publisher.subscribe(self.OnCloseProject, 'Close project data')
  441 + Publisher.subscribe(self.UpdateTrekkerObject, 'Update Trekker object')
  442 + Publisher.subscribe(self.UpdateNumTracts, 'Update number of tracts')
  443 + Publisher.subscribe(self.UpdateSeedOffset, 'Update seed offset')
  444 + Publisher.subscribe(self.UpdateSeedRadius, 'Update seed radius')
  445 + Publisher.subscribe(self.UpdateSleep, 'Update sleep')
  446 + Publisher.subscribe(self.UpdateNumberThreads, 'Update number of threads')
  447 + Publisher.subscribe(self.UpdateTractsVisualization, 'Update tracts visualization')
412 448  
413 449 def LoadImageFiducials(self, marker_id, coord):
414 450 for n in const.BTNS_IMG_MKS:
... ... @@ -420,7 +456,37 @@ class NeuronavigationPanel(wx.Panel):
420 456 for m in [0, 1, 2]:
421 457 self.numctrls_coord[btn_id][m].SetValue(coord[m])
422 458  
423   - def UpdateImageCoordinates(self, position):
  459 + def UpdateFRE(self, fre):
  460 + # TODO: Exhibit FRE in a warning dialog and only starts navigation after user clicks ok
  461 + self.txtctrl_fre.SetValue(str(round(fre, 2)))
  462 + if fre <= 3:
  463 + self.txtctrl_fre.SetBackgroundColour('GREEN')
  464 + else:
  465 + self.txtctrl_fre.SetBackgroundColour('RED')
  466 +
  467 + def UpdateTrekkerObject(self, data):
  468 + # self.trk_inp = data
  469 + self.trekker = data
  470 +
  471 + def UpdateNumTracts(self, data):
  472 + self.n_tracts = data
  473 +
  474 + def UpdateSeedOffset(self, data):
  475 + self.seed_offset = data
  476 +
  477 + def UpdateSeedRadius(self, data):
  478 + self.seed_radius = data
  479 +
  480 + def UpdateSleep(self, data):
  481 + self.sleep_nav = data
  482 +
  483 + def UpdateNumberThreads(self, data):
  484 + self.n_threads = data
  485 +
  486 + def UpdateTractsVisualization(self, data):
  487 + self.view_tracts = data
  488 +
  489 + def UpdateImageCoordinates(self, arg, position):
424 490 # TODO: Change from world coordinates to matrix coordinates. They are better for multi software communication.
425 491 self.current_coord = position
426 492 for m in [0, 1, 2]:
... ... @@ -473,7 +539,7 @@ class NeuronavigationPanel(wx.Panel):
473 539 label=_("Tracker disconnected successfully"))
474 540 self.trk_init = dt.TrackerConnection(self.tracker_id, None, 'connect')
475 541 if not self.trk_init[0]:
476   - dlg.NavigationTrackerWarning(self.tracker_id, self.trk_init[1])
  542 + dlg.ShowNavigationTrackerWarning(self.tracker_id, self.trk_init[1])
477 543 ctrl.SetSelection(0)
478 544 print("Tracker not connected!")
479 545 else:
... ... @@ -488,7 +554,7 @@ class NeuronavigationPanel(wx.Panel):
488 554 self.trk_init = dt.TrackerConnection(self.tracker_id, trck, 'disconnect')
489 555 if not self.trk_init[0]:
490 556 if evt is not False:
491   - dlg.NavigationTrackerWarning(self.tracker_id, 'disconnect')
  557 + dlg.ShowNavigationTrackerWarning(self.tracker_id, 'disconnect')
492 558 self.tracker_id = 0
493 559 ctrl.SetSelection(self.tracker_id)
494 560 Publisher.sendMessage('Update status text in GUI',
... ... @@ -507,7 +573,7 @@ class NeuronavigationPanel(wx.Panel):
507 573 self.tracker_id = choice
508 574 self.trk_init = dt.TrackerConnection(self.tracker_id, None, 'connect')
509 575 if not self.trk_init[0]:
510   - dlg.NavigationTrackerWarning(self.tracker_id, self.trk_init[1])
  576 + dlg.ShowNavigationTrackerWarning(self.tracker_id, self.trk_init[1])
511 577 self.tracker_id = 0
512 578 ctrl.SetSelection(self.tracker_id)
513 579  
... ... @@ -549,7 +615,18 @@ class NeuronavigationPanel(wx.Panel):
549 615 coord = None
550 616  
551 617 if self.trk_init and self.tracker_id:
  618 + # if self.tracker_id == const.DEBUGTRACK:
  619 + # if btn_id == 3:
  620 + # coord1 = np.array([-120., 0., 0., 0., 0., 0.])
  621 + # elif btn_id == 4:
  622 + # coord1 = np.array([120., 0., 0., 0., 0., 0.])
  623 + # elif btn_id == 5:
  624 + # coord1 = np.array([0., 120., 0., 0., 0., 0.])
  625 + # coord2 = np.zeros([3, 6])
  626 + # coord_raw = np.vstack([coord1, coord2])
  627 + # else:
552 628 coord_raw = dco.GetCoordinates(self.trk_init, self.tracker_id, self.ref_mode_id)
  629 +
553 630 if self.ref_mode_id:
554 631 coord = dco.dynamic_reference_m(coord_raw[0, :], coord_raw[1, :])
555 632 else:
... ... @@ -557,7 +634,7 @@ class NeuronavigationPanel(wx.Panel):
557 634 coord[2] = -coord[2]
558 635  
559 636 else:
560   - dlg.NavigationTrackerWarning(0, 'choose')
  637 + dlg.ShowNavigationTrackerWarning(0, 'choose')
561 638  
562 639 # Update number controls with tracker coordinates
563 640 if coord is not None:
... ... @@ -569,112 +646,152 @@ class NeuronavigationPanel(wx.Panel):
569 646 btn_nav = btn[0]
570 647 choice_trck = btn[1]
571 648 choice_ref = btn[2]
  649 + errors = False
  650 +
  651 + # initialize jobs list
  652 + jobs_list = []
  653 + vis_components = [self.trigger_state, self.view_tracts]
  654 + vis_queues = [self.coord_queue, self.trigger_queue, self.tracts_queue]
572 655  
573 656 nav_id = btn_nav.GetValue()
574   - if nav_id:
  657 + if not nav_id:
  658 + self.event.set()
  659 +
  660 + # print("coord unfinished: {}, queue {}", self.coord_queue.unfinished_tasks, self.coord_queue.qsize())
  661 + # print("coord_tracts unfinished: {}, queue {}", self.coord_tracts_queue.unfinished_tasks, self.coord_tracts_queue.qsize())
  662 + # print("tracts unfinished: {}, queue {}", self.tracts_queue.unfinished_tasks, self.tracts_queue.qsize())
  663 + self.coord_queue.clear()
  664 + # self.visualization_queue.clear()
  665 + if self.trigger_state:
  666 + self.trigger_queue.clear()
  667 + if self.view_tracts:
  668 + self.coord_tracts_queue.clear()
  669 + self.tracts_queue.clear()
  670 +
  671 + # print("coord after unfinished: {}, queue {}", self.coord_queue.unfinished_tasks, self.coord_queue.qsize())
  672 + # print("coord_tracts after unfinished: {}, queue {}", self.coord_tracts_queue.unfinished_tasks, self.coord_tracts_queue.qsize())
  673 + # print("tracts after unfinished: {}, queue {}", self.tracts_queue.unfinished_tasks, self.tracts_queue.qsize())
  674 + self.coord_queue.join()
  675 + # self.visualization_queue.join()
  676 + if self.trigger_state:
  677 + self.trigger_queue.join()
  678 + if self.view_tracts:
  679 + self.coord_tracts_queue.join()
  680 + self.tracts_queue.join()
  681 +
  682 + # print("coord join unfinished: {}, queue {}", self.coord_queue.unfinished_tasks, self.coord_queue.qsize())
  683 + # print("vis join unfinished: {}, queue {}", self.visualization_queue.unfinished_tasks, self.visualization_queue.qsize())
  684 +
  685 + tooltip = wx.ToolTip(_("Start neuronavigation"))
  686 + btn_nav.SetToolTip(tooltip)
  687 +
  688 + # Enable all navigation buttons
  689 + choice_ref.Enable(True)
  690 + choice_trck.Enable(True)
  691 + for btn_c in self.btns_coord:
  692 + btn_c.Enable(True)
  693 +
  694 + # if self.trigger_state:
  695 + # self.trigger.stop()
  696 +
  697 + Publisher.sendMessage("Navigation status", nav_status=False, vis_status=vis_components)
  698 +
  699 + else:
  700 +
575 701 if np.isnan(self.fiducials).any():
576   - dlg.InvalidFiducials()
  702 + wx.MessageBox(_("Invalid fiducials, select all coordinates."), _("InVesalius 3"))
577 703 btn_nav.SetValue(False)
578 704  
579   - elif not self.trk_init[0]:
580   - dlg.NavigationTrackerWarning(0, 'choose')
  705 + elif not self.trk_init[0] or not self.tracker_id:
  706 + dlg.ShowNavigationTrackerWarning(0, 'choose')
  707 + errors = True
581 708  
582 709 else:
  710 + if self.event.is_set():
  711 + self.event.clear()
  712 +
  713 + # prepare GUI for navigation
  714 + Publisher.sendMessage("Navigation status", nav_status=True, vis_status=vis_components)
  715 + Publisher.sendMessage("Toggle Cross", id=const.SLICE_STATE_CROSS)
  716 + Publisher.sendMessage("Hide current mask")
583 717 tooltip = wx.ToolTip(_("Stop neuronavigation"))
584 718 btn_nav.SetToolTip(tooltip)
585 719  
586   - # Disable all navigation buttons
  720 + # disable all navigation buttons
587 721 choice_ref.Enable(False)
588 722 choice_trck.Enable(False)
589 723 for btn_c in self.btns_coord:
590 724 btn_c.Enable(False)
591 725  
592   - # fids_head_img = np.zeros([3, 3])
593   - # for ic in range(0, 3):
594   - # fids_head_img[ic, :] = np.asarray(db.flip_x_m(self.fiducials[ic, :]))
595   - #
596   - # m_head_aux, q_head_aux, m_inv_head_aux = db.base_creation(fids_head_img)
597   - # m_head = np.asmatrix(np.identity(4))
598   - # m_head[:3, :3] = m_head_aux[:3, :3]
599   -
600   - m, q1, minv = db.base_creation_old(self.fiducials[:3, :])
601   - n, q2, ninv = db.base_creation_old(self.fiducials[3:, :])
602   -
  726 + # fiducials matrix
603 727 m_change = tr.affine_matrix_from_points(self.fiducials[3:, :].T, self.fiducials[:3, :].T,
604 728 shear=False, scale=False)
605   -
606   - # coreg_data = [m_change, m_head]
607   -
  729 + # initialize spatial tracker parameters
608 730 tracker_mode = self.trk_init, self.tracker_id, self.ref_mode_id
609   - # FIXME: FRE is taking long to calculate so it updates on GUI delayed to navigation - I think its fixed
610   - # TODO: Exhibit FRE in a warning dialog and only starts navigation after user clicks ok
611   - fre = db.calculate_fre(self.fiducials, minv, n, q1, q2)
612   -
613   - self.txtctrl_fre.SetValue(str(round(fre, 2)))
614   - if fre <= 3:
615   - self.txtctrl_fre.SetBackgroundColour('GREEN')
616   - else:
617   - self.txtctrl_fre.SetBackgroundColour('RED')
618   -
619   - if self.trigger_state:
620   - self.trigger = trig.Trigger(nav_id)
621 731  
622   - Publisher.sendMessage("Navigation status", status=True)
623   - Publisher.sendMessage("Toggle Cross", id=const.SLICE_STATE_CROSS)
624   - Publisher.sendMessage("Hide current mask")
  732 + # compute fiducial registration error (FRE)
  733 + # this is the old way to compute the fre, left here to recheck if new works fine.
  734 + # fre = db.calculate_fre(self.fiducials, minv, n, q1, q2)
  735 + fre = db.calculate_fre_m(self.fiducials)
  736 + self.UpdateFRE(fre)
625 737  
626 738 if self.track_obj:
627   - if self.obj_reg_status:
  739 + # if object tracking is selected
  740 + if not self.obj_reg_status:
  741 + # check if object registration was performed
  742 + wx.MessageBox(_("Perform coil registration before navigation."), _("InVesalius 3"))
  743 + errors = True
  744 + else:
  745 + # if object registration was correctly performed continue with navigation
628 746 # obj_reg[0] is object 3x3 fiducial matrix and obj_reg[1] is 3x3 orientation matrix
629 747 obj_fiducials, obj_orients, obj_ref_mode, obj_name = self.obj_reg
630 748  
631   - if self.trk_init and self.tracker_id:
632   -
633   - coreg_data = [m_change, obj_ref_mode]
634   -
635   - if self.ref_mode_id:
636   - coord_raw = dco.GetCoordinates(self.trk_init, self.tracker_id, self.ref_mode_id)
637   - obj_data = db.object_registration(obj_fiducials, obj_orients, coord_raw, m_change)
638   - coreg_data.extend(obj_data)
639   -
640   - self.correg = dcr.CoregistrationObjectDynamic(coreg_data, nav_id, tracker_mode)
641   - else:
642   - coord_raw = np.array([None])
643   - obj_data = db.object_registration(obj_fiducials, obj_orients, coord_raw, m_change)
644   - coreg_data.extend(obj_data)
645   -
646   - self.correg = dcr.CoregistrationObjectStatic(coreg_data, nav_id, tracker_mode)
  749 + coreg_data = [m_change, obj_ref_mode]
647 750  
  751 + if self.ref_mode_id:
  752 + coord_raw = dco.GetCoordinates(self.trk_init, self.tracker_id, self.ref_mode_id)
648 753 else:
649   - dlg.NavigationTrackerWarning(0, 'choose')
650   -
651   - else:
652   - dlg.InvalidObjectRegistration()
653   -
654   - else:
655   - coreg_data = [m_change, 0]
656   - if self.ref_mode_id:
657   - # self.correg = dcr.CoregistrationDynamic_old(bases_coreg, nav_id, tracker_mode)
658   - self.correg = dcr.CoregistrationDynamic(coreg_data, nav_id, tracker_mode)
659   - else:
660   - self.correg = dcr.CoregistrationStatic(coreg_data, nav_id, tracker_mode)
  754 + coord_raw = np.array([None])
661 755  
662   - else:
663   - tooltip = wx.ToolTip(_("Start neuronavigation"))
664   - btn_nav.SetToolTip(tooltip)
665   -
666   - # Enable all navigation buttons
667   - choice_ref.Enable(True)
668   - choice_trck.Enable(True)
669   - for btn_c in self.btns_coord:
670   - btn_c.Enable(True)
671   -
672   - if self.trigger_state:
673   - self.trigger.stop()
  756 + obj_data = db.object_registration(obj_fiducials, obj_orients, coord_raw, m_change)
  757 + coreg_data.extend(obj_data)
674 758  
675   - self.correg.stop()
  759 + jobs_list.append(dcr.CoordinateCorregistrate(self.ref_mode_id, tracker_mode, coreg_data, self.coord_queue,
  760 + self.view_tracts, self.coord_tracts_queue,
  761 + self.event, self.sleep_nav))
676 762  
677   - Publisher.sendMessage("Navigation status", status=False)
  763 + else:
  764 + coreg_data = (m_change, 0)
  765 + jobs_list.append(dcr.CoordinateCorregistrateNoObject(self.ref_mode_id, tracker_mode, coreg_data,
  766 + self.coord_queue,
  767 + self.view_tracts, self.coord_tracts_queue,
  768 + self.event, self.sleep_nav))
  769 +
  770 + if not errors:
  771 + #TODO: Test the trigger thread
  772 + if self.trigger_state:
  773 + # self.trigger = trig.Trigger(nav_id)
  774 + jobs_list.append(trig.TriggerNew(self.trigger_queue, self.event, self.sleep_nav))
  775 +
  776 + if self.view_tracts:
  777 + # initialize Trekker parameters
  778 + slic = sl.Slice()
  779 + affine = slic.affine
  780 + affine_vtk = vtk_utils.numpy_to_vtkMatrix4x4(affine)
  781 + Publisher.sendMessage("Update marker offset state", create=True)
  782 + self.trk_inp = self.trekker, affine, self.seed_offset, self.n_tracts, self.seed_radius, self.n_threads
  783 + # print("Appending the tract computation thread!")
  784 + jobs_list.append(dti.ComputeTractsThread(self.trk_inp, affine_vtk,
  785 + self.coord_tracts_queue, self.tracts_queue,
  786 + self.event, self.sleep_nav))
  787 +
  788 + jobs_list.append(UpdateNavigationScene(vis_queues, vis_components,
  789 + self.event, self.sleep_nav))
  790 +
  791 + for jobs in jobs_list:
  792 + # jobs.daemon = True
  793 + jobs.start()
  794 + # del jobs
678 795  
679 796 def ResetImageFiducials(self):
680 797 for m in range(0, 3):
... ... @@ -699,6 +816,8 @@ class NeuronavigationPanel(wx.Panel):
699 816 Publisher.sendMessage('Update object registration')
700 817 Publisher.sendMessage('Update track object state', flag=False, obj_name=False)
701 818 Publisher.sendMessage('Delete all markers')
  819 + Publisher.sendMessage("Update marker offset state", create=False)
  820 + Publisher.sendMessage("Remove tracts")
702 821 # TODO: Reset camera initial focus
703 822 Publisher.sendMessage('Reset cam clipping range')
704 823  
... ... @@ -830,8 +949,7 @@ class ObjectRegistrationPanel(wx.Panel):
830 949 def UpdateTrackerInit(self, nav_prop):
831 950 self.nav_prop = nav_prop
832 951  
833   - def UpdateNavigationStatus(self, status):
834   - nav_status = status
  952 + def UpdateNavigationStatus(self, nav_status, vis_status):
835 953 if nav_status:
836 954 self.checkrecordcoords.Enable(1)
837 955 self.checktrack.Enable(0)
... ... @@ -898,37 +1016,50 @@ class ObjectRegistrationPanel(wx.Panel):
898 1016 pass
899 1017  
900 1018 else:
901   - dlg.NavigationTrackerWarning(0, 'choose')
  1019 + dlg.ShowNavigationTrackerWarning(0, 'choose')
902 1020  
903 1021 def OnLinkLoad(self, event=None):
904   - filename = dlg.ShowLoadRegistrationDialog()
  1022 + filename = dlg.ShowLoadSaveDialog(message=_(u"Load object registration"),
  1023 + wildcard=_("Registration files (*.obr)|*.obr"))
  1024 + # data_dir = os.environ.get('OneDrive') + r'\data\dti_navigation\baran\pilot_20200131'
  1025 + # coil_path = 'magstim_coil_dell_laptop.obr'
  1026 + # filename = os.path.join(data_dir, coil_path)
905 1027  
906   - if filename:
907   - data = np.loadtxt(filename, delimiter='\t')
908   - self.obj_fiducials = data[:, :3]
909   - self.obj_orients = data[:, 3:]
910   -
911   - text_file = open(filename, "r")
912   - header = text_file.readline().split('\t')
913   - text_file.close()
914   -
915   - self.obj_name = header[1]
916   - self.obj_ref_mode = int(header[-1])
917   -
918   - self.checktrack.Enable(1)
919   - Publisher.sendMessage('Update object registration',
920   - data=(self.obj_fiducials, self.obj_orients, self.obj_ref_mode, self.obj_name))
921   - Publisher.sendMessage('Update status text in GUI', label=_("Ready"))
922   - self.checktrack.SetValue(True)
923   - Publisher.sendMessage('Update track object state', flag=True, obj_name=self.obj_name)
924   - Publisher.sendMessage('Change camera checkbox', status=False)
925   - wx.MessageBox(_("Object file successfully loaded"), _("Load"))
  1028 + try:
  1029 + if filename:
  1030 + #TODO: Improve method to read the file, using "with" similar to OnLoadParameters
  1031 + data = np.loadtxt(filename, delimiter='\t')
  1032 + self.obj_fiducials = data[:, :3]
  1033 + self.obj_orients = data[:, 3:]
  1034 +
  1035 + text_file = open(filename, "r")
  1036 + header = text_file.readline().split('\t')
  1037 + text_file.close()
  1038 +
  1039 + self.obj_name = header[1]
  1040 + self.obj_ref_mode = int(header[-1])
  1041 +
  1042 + self.checktrack.Enable(1)
  1043 + self.checktrack.SetValue(True)
  1044 + Publisher.sendMessage('Update object registration',
  1045 + data=(self.obj_fiducials, self.obj_orients, self.obj_ref_mode, self.obj_name))
  1046 + Publisher.sendMessage('Update status text in GUI',
  1047 + label=_("Object file successfully loaded"))
  1048 + Publisher.sendMessage('Update track object state', flag=True, obj_name=self.obj_name)
  1049 + Publisher.sendMessage('Change camera checkbox', status=False)
  1050 + # wx.MessageBox(_("Object file successfully loaded"), _("Load"))
  1051 + except:
  1052 + wx.MessageBox(_("Object registration file incompatible."), _("InVesalius 3"))
  1053 + Publisher.sendMessage('Update status text in GUI', label="")
926 1054  
927 1055 def ShowSaveObjectDialog(self, evt):
928 1056 if np.isnan(self.obj_fiducials).any() or np.isnan(self.obj_orients).any():
929 1057 wx.MessageBox(_("Digitize all object fiducials before saving"), _("Save error"))
930 1058 else:
931   - filename = dlg.ShowSaveRegistrationDialog("object_registration.obr")
  1059 + filename = dlg.ShowLoadSaveDialog(message=_(u"Save object registration as..."),
  1060 + wildcard=_("Registration files (*.obr)|*.obr"),
  1061 + style=wx.FD_SAVE | wx.FD_OVERWRITE_PROMPT,
  1062 + default_filename="object_registration.obr", save_ext="obr")
932 1063 if filename:
933 1064 hdr = 'Object' + "\t" + utils.decode(self.obj_name, const.FS_ENCODE) + "\t" + 'Reference' + "\t" + str('%d' % self.obj_ref_mode)
934 1065 data = np.hstack([self.obj_fiducials, self.obj_orients])
... ... @@ -969,8 +1100,8 @@ class MarkersPanel(wx.Panel):
969 1100 self.tgt_flag = self.tgt_index = None
970 1101 self.nav_status = False
971 1102  
972   - self.marker_colour = (1.0, 1.0, 0.)
973   - self.marker_size = 3
  1103 + self.marker_colour = const.MARKER_COLOUR
  1104 + self.marker_size = const.MARKER_SIZE
974 1105  
975 1106 # Change marker size
976 1107 spin_size = wx.SpinCtrl(self, -1, "", size=wx.Size(40, 23))
... ... @@ -1045,7 +1176,8 @@ class MarkersPanel(wx.Panel):
1045 1176 self.Update()
1046 1177  
1047 1178 def __bind_events(self):
1048   - Publisher.subscribe(self.UpdateCurrentCoord, 'Co-registered points')
  1179 + # Publisher.subscribe(self.UpdateCurrentCoord, 'Co-registered points')
  1180 + Publisher.subscribe(self.UpdateCurrentCoord, 'Update cross position')
1049 1181 Publisher.subscribe(self.OnDeleteSingleMarker, 'Delete fiducial marker')
1050 1182 Publisher.subscribe(self.OnDeleteAllMarkers, 'Delete all markers')
1051 1183 Publisher.subscribe(self.OnCreateMarker, 'Create marker')
... ... @@ -1055,8 +1187,8 @@ class MarkersPanel(wx.Panel):
1055 1187 self.current_coord = position[:]
1056 1188 #self.current_angle = pubsub_evt.data[1][3:]
1057 1189  
1058   - def UpdateNavigationStatus(self, status):
1059   - if not status:
  1190 + def UpdateNavigationStatus(self, nav_status, vis_status):
  1191 + if not nav_status:
1060 1192 sleep(0.5)
1061 1193 #self.current_coord[3:] = 0, 0, 0
1062 1194 self.nav_status = False
... ... @@ -1073,6 +1205,9 @@ class MarkersPanel(wx.Panel):
1073 1205 menu_id.AppendSeparator()
1074 1206 target_menu = menu_id.Append(1, _('Set as target'))
1075 1207 menu_id.Bind(wx.EVT_MENU, self.OnMenuSetTarget, target_menu)
  1208 + # TODO: Create the remove target option so the user can disable the target without removing the marker
  1209 + # target_menu_rem = menu_id.Append(3, _('Remove target'))
  1210 + # menu_id.Bind(wx.EVT_MENU, self.OnMenuRemoveTarget, target_menu_rem)
1076 1211  
1077 1212 target_menu.Enable(True)
1078 1213 self.PopupMenu(menu_id)
... ... @@ -1089,10 +1224,10 @@ class MarkersPanel(wx.Panel):
1089 1224 if evt == 'TARGET':
1090 1225 id_label = evt
1091 1226 else:
1092   - id_label = dlg.EnterMarkerID(self.lc.GetItemText(list_index, 4))
  1227 + id_label = dlg.ShowEnterMarkerID(self.lc.GetItemText(list_index, 4))
1093 1228 if id_label == 'TARGET':
1094 1229 id_label = ''
1095   - dlg.InvalidTargetID()
  1230 + wx.MessageBox(_("Invalid TARGET ID."), _("InVesalius 3"))
1096 1231 self.lc.SetItem(list_index, 4, id_label)
1097 1232 # Add the new ID to exported list
1098 1233 if len(self.list_coord[list_index]) > 8:
... ... @@ -1122,34 +1257,29 @@ class MarkersPanel(wx.Panel):
1122 1257 Publisher.sendMessage('Disable or enable coil tracker', status=True)
1123 1258 self.OnMenuEditMarkerId('TARGET')
1124 1259 self.tgt_flag = True
1125   - dlg.NewTarget()
  1260 + wx.MessageBox(_("New target selected."), _("InVesalius 3"))
1126 1261  
1127 1262 def OnMenuSetColor(self, evt):
1128 1263 index = self.lc.GetFocusedItem()
1129   - cdata = wx.ColourData()
1130   - cdata.SetColour(wx.Colour(self.list_coord[index][6]*255,self.list_coord[index][7]*255,self.list_coord[index][8]*255))
1131   - dlg = wx.ColourDialog(self, data=cdata)
1132   - dlg.GetColourData().SetChooseFull(True)
1133   - if dlg.ShowModal() == wx.ID_OK:
1134   - self.r, self.g, self.b = dlg.GetColourData().GetColour().Get(includeAlpha=False)
1135   - r = float(self.r) / 255.0
1136   - g = float(self.g) / 255.0
1137   - b = float(self.b) / 255.0
1138   - dlg.Destroy()
1139   - color = [r,g,b]
1140   -
1141   - Publisher.sendMessage('Set new color', index=index, color=color)
1142   -
1143   - self.list_coord[index][6] = r
1144   - self.list_coord[index][7] = g
1145   - self.list_coord[index][8] = b
  1264 +
  1265 + color_current = [self.list_coord[index][n] * 255 for n in range(6, 9)]
  1266 +
  1267 + color_new = dlg.ShowColorDialog(color_current=color_current)
  1268 +
  1269 + if color_new:
  1270 + assert len(color_new) == 3
  1271 + for n, col in enumerate(color_new):
  1272 + self.list_coord[index][n+6] = col/255.0
  1273 +
  1274 + Publisher.sendMessage('Set new color', index=index, color=color_new)
1146 1275  
1147 1276 def OnDeleteAllMarkers(self, evt=None):
1148 1277 if self.list_coord:
1149 1278 if evt is None:
1150 1279 result = wx.ID_OK
1151 1280 else:
1152   - result = dlg.DeleteAllMarkers()
  1281 + # result = dlg.DeleteAllMarkers()
  1282 + result = dlg.ShowConfirmationDialog(msg=_("Remove all markers? Cannot be undone."))
1153 1283  
1154 1284 if result == wx.ID_OK:
1155 1285 self.list_coord = []
... ... @@ -1162,7 +1292,7 @@ class MarkersPanel(wx.Panel):
1162 1292 self.tgt_flag = self.tgt_index = None
1163 1293 Publisher.sendMessage('Disable or enable coil tracker', status=False)
1164 1294 if not hasattr(evt, 'data'):
1165   - dlg.DeleteTarget()
  1295 + wx.MessageBox(_("Target deleted."), _("InVesalius 3"))
1166 1296  
1167 1297 def OnDeleteSingleMarker(self, evt=None, marker_id=None):
1168 1298 # OnDeleteSingleMarker is used for both pubsub and button click events
... ... @@ -1183,14 +1313,16 @@ class MarkersPanel(wx.Panel):
1183 1313 else:
1184 1314 index = None
1185 1315  
  1316 + #TODO: There are bugs when no marker is selected, test and improve
1186 1317 if index:
1187 1318 if self.tgt_flag and self.tgt_index == index[0]:
1188 1319 self.tgt_flag = self.tgt_index = None
1189 1320 Publisher.sendMessage('Disable or enable coil tracker', status=False)
1190   - dlg.DeleteTarget()
  1321 + wx.MessageBox(_("No data selected."), _("InVesalius 3"))
  1322 +
1191 1323 self.DeleteMarker(index)
1192 1324 else:
1193   - dlg.NoMarkerSelected()
  1325 + wx.MessageBox(_("Target deleted."), _("InVesalius 3"))
1194 1326  
1195 1327 def DeleteMarker(self, index):
1196 1328 for i in reversed(index):
... ... @@ -1213,12 +1345,16 @@ class MarkersPanel(wx.Panel):
1213 1345 self.CreateMarker(self.current_coord, self.marker_colour, self.marker_size)
1214 1346  
1215 1347 def OnLoadMarkers(self, evt):
1216   - filepath = dlg.ShowLoadMarkersDialog()
  1348 + filename = dlg.ShowLoadSaveDialog(message=_(u"Load markers"),
  1349 + wildcard=_("Markers files (*.mks)|*.mks"))
  1350 + # data_dir = os.environ.get('OneDrive') + r'\data\dti_navigation\baran\pilot_20200131'
  1351 + # marker_path = 'markers.mks'
  1352 + # filename = os.path.join(data_dir, marker_path)
1217 1353  
1218   - if filepath:
  1354 + if filename:
1219 1355 try:
1220 1356 count_line = self.lc.GetItemCount()
1221   - content = [s.rstrip() for s in open(filepath)]
  1357 + content = [s.rstrip() for s in open(filename)]
1222 1358 for data in content:
1223 1359 target = None
1224 1360 line = [s for s in data.split()]
... ... @@ -1254,7 +1390,7 @@ class MarkersPanel(wx.Panel):
1254 1390 self.CreateMarker(coord, colour, size, line[7])
1255 1391 count_line += 1
1256 1392 except:
1257   - dlg.InvalidMarkersFile()
  1393 + wx.MessageBox(_("Invalid markers file."), _("InVesalius 3"))
1258 1394  
1259 1395 def OnMarkersVisibility(self, evt, ctrl):
1260 1396  
... ... @@ -1266,7 +1402,11 @@ class MarkersPanel(wx.Panel):
1266 1402 ctrl.SetLabel('Hide')
1267 1403  
1268 1404 def OnSaveMarkers(self, evt):
1269   - filename = dlg.ShowSaveMarkersDialog("markers.mks")
  1405 + filename = dlg.ShowLoadSaveDialog(message=_(u"Save markers as..."),
  1406 + wildcard=_("Marker files (*.mks)|*.mks"),
  1407 + style=wx.FD_SAVE | wx.FD_OVERWRITE_PROMPT,
  1408 + default_filename="markers.mks", save_ext="mks")
  1409 +
1270 1410 if filename:
1271 1411 if self.list_coord:
1272 1412 text_file = open(filename, "w")
... ... @@ -1339,3 +1479,459 @@ class DbsPanel(wx.Panel):
1339 1479 except AttributeError:
1340 1480 default_colour = wx.SystemSettings_GetColour(wx.SYS_COLOUR_MENUBAR)
1341 1481  
  1482 +
  1483 +class TractographyPanel(wx.Panel):
  1484 +
  1485 + def __init__(self, parent):
  1486 + wx.Panel.__init__(self, parent)
  1487 + try:
  1488 + default_colour = wx.SystemSettings.GetColour(wx.SYS_COLOUR_MENUBAR)
  1489 + except AttributeError:
  1490 + default_colour = wx.SystemSettings_GetColour(wx.SYS_COLOUR_MENUBAR)
  1491 + self.SetBackgroundColour(default_colour)
  1492 +
  1493 + self.affine = None
  1494 + self.affine_vtk = None
  1495 + self.trekker = None
  1496 + self.n_tracts = const.N_TRACTS
  1497 + self.peel_depth = const.PEEL_DEPTH
  1498 + self.view_tracts = False
  1499 + self.seed_offset = const.SEED_OFFSET
  1500 + self.seed_radius = const.SEED_RADIUS
  1501 + self.sleep_nav = const.SLEEP_NAVIGATION
  1502 + self.brain_opacity = const.BRAIN_OPACITY
  1503 + self.brain_peel = None
  1504 + self.brain_actor = None
  1505 + self.n_peels = const.MAX_PEEL_DEPTH
  1506 + self.p_old = np.array([[0., 0., 0.]])
  1507 + self.tracts_run = None
  1508 + self.trekker_cfg = const.TREKKER_CONFIG
  1509 + self.nav_status = False
  1510 +
  1511 + self.SetAutoLayout(1)
  1512 + self.__bind_events()
  1513 +
  1514 + # Button for creating new coil
  1515 + tooltip = wx.ToolTip(_("Load brain visualization"))
  1516 + btn_mask = wx.Button(self, -1, _("Load brain"), size=wx.Size(65, 23))
  1517 + btn_mask.SetToolTip(tooltip)
  1518 + btn_mask.Enable(1)
  1519 + btn_mask.Bind(wx.EVT_BUTTON, self.OnLinkBrain)
  1520 + # self.btn_new = btn_new
  1521 +
  1522 + # Button for import config coil file
  1523 + tooltip = wx.ToolTip(_("Load FOD"))
  1524 + btn_load = wx.Button(self, -1, _("Load FOD"), size=wx.Size(65, 23))
  1525 + btn_load.SetToolTip(tooltip)
  1526 + btn_load.Enable(1)
  1527 + btn_load.Bind(wx.EVT_BUTTON, self.OnLinkFOD)
  1528 + # self.btn_load = btn_load
  1529 +
  1530 + # Save button for object registration
  1531 + tooltip = wx.ToolTip(_(u"Load Trekker configuration parameters"))
  1532 + btn_load_cfg = wx.Button(self, -1, _(u"Configure"), size=wx.Size(65, 23))
  1533 + btn_load_cfg.SetToolTip(tooltip)
  1534 + btn_load_cfg.Enable(1)
  1535 + btn_load_cfg.Bind(wx.EVT_BUTTON, self.OnLoadParameters)
  1536 + # self.btn_load_cfg = btn_load_cfg
  1537 +
  1538 + # Create a horizontal sizer to represent button save
  1539 + line_btns = wx.BoxSizer(wx.HORIZONTAL)
  1540 + line_btns.Add(btn_load, 1, wx.LEFT | wx.TOP | wx.RIGHT, 4)
  1541 + line_btns.Add(btn_load_cfg, 1, wx.LEFT | wx.TOP | wx.RIGHT, 4)
  1542 + line_btns.Add(btn_mask, 1, wx.LEFT | wx.TOP | wx.RIGHT, 4)
  1543 +
  1544 + # Change peeling depth
  1545 + text_peel_depth = wx.StaticText(self, -1, _("Peeling depth (mm):"))
  1546 + spin_peel_depth = wx.SpinCtrl(self, -1, "", size=wx.Size(50, 23))
  1547 + spin_peel_depth.Enable(1)
  1548 + spin_peel_depth.SetRange(0, const.MAX_PEEL_DEPTH)
  1549 + spin_peel_depth.SetValue(const.PEEL_DEPTH)
  1550 + spin_peel_depth.Bind(wx.EVT_TEXT, partial(self.OnSelectPeelingDepth, ctrl=spin_peel_depth))
  1551 + spin_peel_depth.Bind(wx.EVT_SPINCTRL, partial(self.OnSelectPeelingDepth, ctrl=spin_peel_depth))
  1552 +
  1553 + # Change number of tracts
  1554 + text_ntracts = wx.StaticText(self, -1, _("Number tracts:"))
  1555 + spin_ntracts = wx.SpinCtrl(self, -1, "", size=wx.Size(50, 23))
  1556 + spin_ntracts.Enable(1)
  1557 + spin_ntracts.SetRange(1, 2000)
  1558 + spin_ntracts.SetValue(const.N_TRACTS)
  1559 + spin_ntracts.Bind(wx.EVT_TEXT, partial(self.OnSelectNumTracts, ctrl=spin_ntracts))
  1560 + spin_ntracts.Bind(wx.EVT_SPINCTRL, partial(self.OnSelectNumTracts, ctrl=spin_ntracts))
  1561 +
  1562 + # Change seed offset for computing tracts
  1563 + text_offset = wx.StaticText(self, -1, _("Seed offset (mm):"))
  1564 + spin_offset = wx.SpinCtrlDouble(self, -1, "", size=wx.Size(50, 23), inc = 0.1)
  1565 + spin_offset.Enable(1)
  1566 + spin_offset.SetRange(0, 100.0)
  1567 + spin_offset.SetValue(self.seed_offset)
  1568 + spin_offset.Bind(wx.EVT_TEXT, partial(self.OnSelectOffset, ctrl=spin_offset))
  1569 + spin_offset.Bind(wx.EVT_SPINCTRL, partial(self.OnSelectOffset, ctrl=spin_offset))
  1570 + # self.spin_offset = spin_offset
  1571 +
  1572 + # Change seed radius for computing tracts
  1573 + text_radius = wx.StaticText(self, -1, _("Seed radius (mm):"))
  1574 + spin_radius = wx.SpinCtrlDouble(self, -1, "", size=wx.Size(50, 23), inc=0.1)
  1575 + spin_radius.Enable(1)
  1576 + spin_radius.SetRange(0, 100.0)
  1577 + spin_radius.SetValue(self.seed_radius)
  1578 + spin_radius.Bind(wx.EVT_TEXT, partial(self.OnSelectRadius, ctrl=spin_radius))
  1579 + spin_radius.Bind(wx.EVT_SPINCTRL, partial(self.OnSelectRadius, ctrl=spin_radius))
  1580 + # self.spin_radius = spin_radius
  1581 +
  1582 + # Change sleep pause between navigation loops
  1583 + text_sleep = wx.StaticText(self, -1, _("Sleep (s):"))
  1584 + spin_sleep = wx.SpinCtrlDouble(self, -1, "", size=wx.Size(50, 23), inc=0.01)
  1585 + spin_sleep.Enable(1)
  1586 + spin_sleep.SetRange(0.01, 10.0)
  1587 + spin_sleep.SetValue(self.sleep_nav)
  1588 + spin_sleep.Bind(wx.EVT_TEXT, partial(self.OnSelectSleep, ctrl=spin_sleep))
  1589 + spin_sleep.Bind(wx.EVT_SPINCTRL, partial(self.OnSelectSleep, ctrl=spin_sleep))
  1590 +
  1591 + # Change opacity of brain mask visualization
  1592 + text_opacity = wx.StaticText(self, -1, _("Brain opacity:"))
  1593 + spin_opacity = wx.SpinCtrlDouble(self, -1, "", size=wx.Size(50, 23), inc=0.1)
  1594 + spin_opacity.Enable(0)
  1595 + spin_opacity.SetRange(0, 1.0)
  1596 + spin_opacity.SetValue(self.brain_opacity)
  1597 + spin_opacity.Bind(wx.EVT_TEXT, partial(self.OnSelectOpacity, ctrl=spin_opacity))
  1598 + spin_opacity.Bind(wx.EVT_SPINCTRL, partial(self.OnSelectOpacity, ctrl=spin_opacity))
  1599 + self.spin_opacity = spin_opacity
  1600 +
  1601 + # Create a horizontal sizer to threshold configs
  1602 + border = 1
  1603 + line_peel_depth = wx.BoxSizer(wx.HORIZONTAL)
  1604 + line_peel_depth.AddMany([(text_peel_depth, 1, wx.EXPAND | wx.GROW | wx.TOP | wx.RIGHT | wx.LEFT, border),
  1605 + (spin_peel_depth, 0, wx.ALL | wx.EXPAND | wx.GROW, border)])
  1606 +
  1607 + line_ntracts = wx.BoxSizer(wx.HORIZONTAL)
  1608 + line_ntracts.AddMany([(text_ntracts, 1, wx.EXPAND | wx.GROW | wx.TOP | wx.RIGHT | wx.LEFT, border),
  1609 + (spin_ntracts, 0, wx.ALL | wx.EXPAND | wx.GROW, border)])
  1610 +
  1611 + line_offset = wx.BoxSizer(wx.HORIZONTAL)
  1612 + line_offset.AddMany([(text_offset, 1, wx.EXPAND | wx.GROW | wx.TOP | wx.RIGHT | wx.LEFT, border),
  1613 + (spin_offset, 0, wx.ALL | wx.EXPAND | wx.GROW, border)])
  1614 +
  1615 + line_radius = wx.BoxSizer(wx.HORIZONTAL)
  1616 + line_radius.AddMany([(text_radius, 1, wx.EXPAND | wx.GROW | wx.TOP | wx.RIGHT | wx.LEFT, border),
  1617 + (spin_radius, 0, wx.ALL | wx.EXPAND | wx.GROW, border)])
  1618 +
  1619 + line_sleep = wx.BoxSizer(wx.HORIZONTAL)
  1620 + line_sleep.AddMany([(text_sleep, 1, wx.EXPAND | wx.GROW | wx.TOP | wx.RIGHT | wx.LEFT, border),
  1621 + (spin_sleep, 0, wx.ALL | wx.EXPAND | wx.GROW, border)])
  1622 +
  1623 + line_opacity = wx.BoxSizer(wx.HORIZONTAL)
  1624 + line_opacity.AddMany([(text_opacity, 1, wx.EXPAND | wx.GROW | wx.TOP | wx.RIGHT | wx.LEFT, border),
  1625 + (spin_opacity, 0, wx.ALL | wx.EXPAND | wx.GROW, border)])
  1626 +
  1627 + # Check box to enable tract visualization
  1628 + checktracts = wx.CheckBox(self, -1, _('Enable tracts'))
  1629 + checktracts.SetValue(False)
  1630 + checktracts.Enable(0)
  1631 + checktracts.Bind(wx.EVT_CHECKBOX, partial(self.OnEnableTracts, ctrl=checktracts))
  1632 + self.checktracts = checktracts
  1633 +
  1634 + # Check box to enable surface peeling
  1635 + checkpeeling = wx.CheckBox(self, -1, _('Peel surface'))
  1636 + checkpeeling.SetValue(False)
  1637 + checkpeeling.Enable(0)
  1638 + checkpeeling.Bind(wx.EVT_CHECKBOX, partial(self.OnShowPeeling, ctrl=checkpeeling))
  1639 + self.checkpeeling = checkpeeling
  1640 +
  1641 + border_last = 1
  1642 + line_checks = wx.BoxSizer(wx.HORIZONTAL)
  1643 + line_checks.Add(checktracts, 0, wx.ALIGN_LEFT | wx.RIGHT | wx.LEFT, border_last)
  1644 + line_checks.Add(checkpeeling, 0, wx.ALIGN_RIGHT | wx.RIGHT | wx.LEFT, border_last)
  1645 +
  1646 + # Add line sizers into main sizer
  1647 + border = 1
  1648 + border_last = 10
  1649 + main_sizer = wx.BoxSizer(wx.VERTICAL)
  1650 + main_sizer.Add(line_btns, 0, wx.BOTTOM | wx.ALIGN_CENTER_HORIZONTAL, border_last)
  1651 + main_sizer.Add(line_peel_depth, 0, wx.GROW | wx.EXPAND | wx.LEFT | wx.RIGHT | wx.TOP, border)
  1652 + main_sizer.Add(line_ntracts, 0, wx.GROW | wx.EXPAND | wx.LEFT | wx.RIGHT | wx.TOP, border)
  1653 + main_sizer.Add(line_offset, 0, wx.GROW | wx.EXPAND | wx.LEFT | wx.RIGHT | wx.TOP, border)
  1654 + main_sizer.Add(line_radius, 0, wx.GROW | wx.EXPAND | wx.LEFT | wx.RIGHT | wx.TOP, border)
  1655 + main_sizer.Add(line_sleep, 0, wx.GROW | wx.EXPAND | wx.LEFT | wx.RIGHT | wx.TOP, border)
  1656 + main_sizer.Add(line_opacity, 0, wx.GROW | wx.EXPAND | wx.LEFT | wx.RIGHT | wx.TOP, border)
  1657 + main_sizer.Add(line_checks, 0, wx.GROW | wx.EXPAND | wx.LEFT | wx.RIGHT | wx.TOP | wx.BOTTOM, border_last)
  1658 + main_sizer.Fit(self)
  1659 +
  1660 + self.SetSizer(main_sizer)
  1661 + self.Update()
  1662 +
  1663 + def __bind_events(self):
  1664 + Publisher.subscribe(self.OnCloseProject, 'Close project data')
  1665 + Publisher.subscribe(self.OnUpdateTracts, 'Update cross position')
  1666 + Publisher.subscribe(self.UpdateNavigationStatus, 'Navigation status')
  1667 +
  1668 + def OnSelectPeelingDepth(self, evt, ctrl):
  1669 + self.peel_depth = ctrl.GetValue()
  1670 + if self.checkpeeling.GetValue():
  1671 + actor = self.brain_peel.get_actor(self.peel_depth)
  1672 + Publisher.sendMessage('Update peel', flag=True, actor=actor)
  1673 +
  1674 + def OnSelectNumTracts(self, evt, ctrl):
  1675 + self.n_tracts = ctrl.GetValue()
  1676 + # self.tract.n_tracts = ctrl.GetValue()
  1677 + Publisher.sendMessage('Update number of tracts', data=self.n_tracts)
  1678 +
  1679 + def OnSelectOffset(self, evt, ctrl):
  1680 + self.seed_offset = ctrl.GetValue()
  1681 + # self.tract.seed_offset = ctrl.GetValue()
  1682 + Publisher.sendMessage('Update seed offset', data=self.seed_offset)
  1683 +
  1684 + def OnSelectRadius(self, evt, ctrl):
  1685 + self.seed_radius = ctrl.GetValue()
  1686 + # self.tract.seed_offset = ctrl.GetValue()
  1687 + Publisher.sendMessage('Update seed radius', data=self.seed_radius)
  1688 +
  1689 + def OnSelectSleep(self, evt, ctrl):
  1690 + self.sleep_nav = ctrl.GetValue()
  1691 + # self.tract.seed_offset = ctrl.GetValue()
  1692 + Publisher.sendMessage('Update sleep', data=self.sleep_nav)
  1693 +
  1694 + def OnSelectOpacity(self, evt, ctrl):
  1695 + self.brain_actor.GetProperty().SetOpacity(ctrl.GetValue())
  1696 + Publisher.sendMessage('Update peel', flag=True, actor=self.brain_actor)
  1697 +
  1698 + def OnShowPeeling(self, evt, ctrl):
  1699 + # self.view_peeling = ctrl.GetValue()
  1700 + if ctrl.GetValue():
  1701 + actor = self.brain_peel.get_actor(self.peel_depth)
  1702 + else:
  1703 + actor = None
  1704 + Publisher.sendMessage('Update peel', flag=ctrl.GetValue(), actor=actor)
  1705 +
  1706 + def OnEnableTracts(self, evt, ctrl):
  1707 + self.view_tracts = ctrl.GetValue()
  1708 + Publisher.sendMessage('Update tracts visualization', data=self.view_tracts)
  1709 + if not self.view_tracts:
  1710 + Publisher.sendMessage('Remove tracts')
  1711 + Publisher.sendMessage("Update marker offset state", create=False)
  1712 +
  1713 + def UpdateNavigationStatus(self, nav_status, vis_status):
  1714 + self.nav_status = nav_status
  1715 +
  1716 + def OnLinkBrain(self, event=None):
  1717 + Publisher.sendMessage('Update status text in GUI', label=_("Busy"))
  1718 + Publisher.sendMessage('Begin busy cursor')
  1719 + mask_path = dlg.ShowImportOtherFilesDialog(const.ID_TREKKER_MASK)
  1720 + img_path = dlg.ShowImportOtherFilesDialog(const.ID_TREKKER_IMG)
  1721 + # data_dir = os.environ.get('OneDriveConsumer') + '\\data\\dti'
  1722 + # mask_file = 'sub-P0_dwi_mask.nii'
  1723 + # mask_path = os.path.join(data_dir, mask_file)
  1724 + # img_file = 'sub-P0_T1w_biascorrected.nii'
  1725 + # img_path = os.path.join(data_dir, img_file)
  1726 +
  1727 + if not self.affine_vtk:
  1728 + slic = sl.Slice()
  1729 + self.affine = slic.affine
  1730 + self.affine_vtk = vtk_utils.numpy_to_vtkMatrix4x4(self.affine)
  1731 +
  1732 + try:
  1733 + self.brain_peel = brain.Brain(img_path, mask_path, self.n_peels, self.affine_vtk)
  1734 + self.brain_actor = self.brain_peel.get_actor(self.peel_depth)
  1735 + self.brain_actor.GetProperty().SetOpacity(self.brain_opacity)
  1736 + Publisher.sendMessage('Update peel', flag=True, actor=self.brain_actor)
  1737 + self.checkpeeling.Enable(1)
  1738 + self.checkpeeling.SetValue(True)
  1739 + self.spin_opacity.Enable(1)
  1740 + Publisher.sendMessage('Update status text in GUI', label=_("Brain model loaded"))
  1741 + except:
  1742 + wx.MessageBox(_("Unable to load brain mask."), _("InVesalius 3"))
  1743 +
  1744 + Publisher.sendMessage('End busy cursor')
  1745 +
  1746 + def OnLinkFOD(self, event=None):
  1747 + Publisher.sendMessage('Update status text in GUI', label=_("Busy"))
  1748 + Publisher.sendMessage('Begin busy cursor')
  1749 + filename = dlg.ShowImportOtherFilesDialog(const.ID_TREKKER_FOD)
  1750 + # Juuso
  1751 + # data_dir = os.environ.get('OneDriveConsumer') + '\\data\\dti'
  1752 + # FOD_path = 'sub-P0_dwi_FOD.nii'
  1753 + # Baran
  1754 + # data_dir = os.environ.get('OneDrive') + r'\data\dti_navigation\baran\pilot_20200131'
  1755 + # FOD_path = 'Baran_FOD.nii'
  1756 + # filename = os.path.join(data_dir, FOD_path)
  1757 +
  1758 + if not self.affine_vtk:
  1759 + slic = sl.Slice()
  1760 + self.affine = slic.affine
  1761 + self.affine_vtk = vtk_utils.numpy_to_vtkMatrix4x4(self.affine)
  1762 +
  1763 + # try:
  1764 +
  1765 + self.trekker = Trekker.initialize(filename.encode('utf-8'))
  1766 + self.trekker, n_threads = dti.set_trekker_parameters(self.trekker, self.trekker_cfg)
  1767 +
  1768 + self.checktracts.Enable(1)
  1769 + self.checktracts.SetValue(True)
  1770 + self.view_tracts = True
  1771 + Publisher.sendMessage('Update Trekker object', data=self.trekker)
  1772 + Publisher.sendMessage('Update number of threads', data=n_threads)
  1773 + Publisher.sendMessage('Update tracts visualization', data=1)
  1774 + Publisher.sendMessage('Update status text in GUI', label=_("Trekker initialized"))
  1775 + # except:
  1776 + # wx.MessageBox(_("Unable to initialize Trekker, check FOD and config files."), _("InVesalius 3"))
  1777 +
  1778 + Publisher.sendMessage('End busy cursor')
  1779 +
  1780 + def OnLoadParameters(self, event=None):
  1781 + import json
  1782 + filename = dlg.ShowLoadSaveDialog(message=_(u"Load Trekker configuration"),
  1783 + wildcard=_("JSON file (*.json)|*.json"))
  1784 + try:
  1785 + # Check if filename exists, read the JSON file and check if all parameters match
  1786 + # with the required list defined in the constants module
  1787 + # if a parameter is missing, raise an error
  1788 + if filename:
  1789 + with open(filename) as json_file:
  1790 + self.trekker_cfg = json.load(json_file)
  1791 + assert all(name in self.trekker_cfg for name in const.TREKKER_CONFIG)
  1792 + if self.trekker:
  1793 + self.trekker, n_threads = dti.set_trekker_parameters(self.trekker, self.trekker_cfg)
  1794 + Publisher.sendMessage('Update Trekker object', data=self.trekker)
  1795 + Publisher.sendMessage('Update number of threads', data=n_threads)
  1796 +
  1797 + Publisher.sendMessage('Update status text in GUI', label=_("Trekker config loaded"))
  1798 +
  1799 + except (AssertionError, json.decoder.JSONDecodeError):
  1800 + # Inform user that file is not compatible
  1801 + self.trekker_cfg = const.TREKKER_CONFIG
  1802 + wx.MessageBox(_("File incompatible, using default configuration."), _("InVesalius 3"))
  1803 + Publisher.sendMessage('Update status text in GUI', label="")
  1804 +
  1805 + def OnUpdateTracts(self, arg, position):
  1806 + """
  1807 + Minimal working version of tract computation. Updates when cross sends Pubsub message to update.
  1808 + Position refers to the coordinates in InVesalius 2D space. To represent the same coordinates in the 3D space,
  1809 + flip_x the coordinates and multiply the z coordinate by -1. This is all done in the flix_x function.
  1810 +
  1811 + :param arg: event for pubsub
  1812 + :param position: list or array with the x, y, and z coordinates in InVesalius space
  1813 + """
  1814 + # Minimal working version of tract computation
  1815 + # It updates when cross updates
  1816 + # pass
  1817 + if self.view_tracts and not self.nav_status:
  1818 + # print("Running during navigation")
  1819 + coord_flip = db.flip_x_m(position[:3])[:3, 0]
  1820 + dti.compute_tracts(self.trekker, coord_flip, self.affine, self.affine_vtk,
  1821 + self.n_tracts)
  1822 +
  1823 + def OnCloseProject(self):
  1824 + self.checktracts.SetValue(False)
  1825 + self.checktracts.Enable(0)
  1826 + self.checkpeeling.SetValue(False)
  1827 + self.checkpeeling.Enable(0)
  1828 +
  1829 + self.spin_opacity.SetValue(const.BRAIN_OPACITY)
  1830 + self.spin_opacity.Enable(0)
  1831 + Publisher.sendMessage('Update peel', flag=False, actor=self.brain_actor)
  1832 +
  1833 + self.peel_depth = const.PEEL_DEPTH
  1834 + self.n_tracts = const.N_TRACTS
  1835 +
  1836 + Publisher.sendMessage('Remove tracts')
  1837 +
  1838 +
  1839 +class QueueCustom(queue.Queue):
  1840 + """
  1841 + A custom queue subclass that provides a :meth:`clear` method.
  1842 + https://stackoverflow.com/questions/6517953/clear-all-items-from-the-queue
  1843 + Modified to a LIFO Queue type (Last-in-first-out). Seems to make sense for the navigation
  1844 + threads, as the last added coordinate should be the first to be processed.
  1845 + In the first tests in a short run, seems to increase the coord queue size considerably,
  1846 + possibly limiting the queue size is good.
  1847 + """
  1848 +
  1849 + def clear(self):
  1850 + """
  1851 + Clears all items from the queue.
  1852 + """
  1853 +
  1854 + with self.mutex:
  1855 + unfinished = self.unfinished_tasks - len(self.queue)
  1856 + if unfinished <= 0:
  1857 + if unfinished < 0:
  1858 + raise ValueError('task_done() called too many times')
  1859 + self.all_tasks_done.notify_all()
  1860 + self.unfinished_tasks = unfinished
  1861 + self.queue.clear()
  1862 + self.not_full.notify_all()
  1863 +
  1864 +
  1865 +class UpdateNavigationScene(threading.Thread):
  1866 +
  1867 + def __init__(self, vis_queues, vis_components, event, sle):
  1868 + """Class (threading) to update the navigation scene with all graphical elements.
  1869 +
  1870 + Sleep function in run method is used to avoid blocking GUI and more fluent, real-time navigation
  1871 +
  1872 + :param affine_vtk: Affine matrix in vtkMatrix4x4 instance to update objects position in 3D scene
  1873 + :type affine_vtk: vtkMatrix4x4
  1874 + :param visualization_queue: Queue instance that manage coordinates to be visualized
  1875 + :type visualization_queue: queue.Queue
  1876 + :param event: Threading event to coordinate when tasks as done and allow UI release
  1877 + :type event: threading.Event
  1878 + :param sle: Sleep pause in seconds
  1879 + :type sle: float
  1880 + """
  1881 +
  1882 + threading.Thread.__init__(self, name='UpdateScene')
  1883 + self.trigger_state, self.view_tracts = vis_components
  1884 + self.coord_queue, self.trigger_queue, self.tracts_queue = vis_queues
  1885 + self.sle = sle
  1886 + self.event = event
  1887 +
  1888 + def run(self):
  1889 + # count = 0
  1890 + while not self.event.is_set():
  1891 + got_coords = False
  1892 + try:
  1893 + coord, m_img, view_obj = self.coord_queue.get_nowait()
  1894 + got_coords = True
  1895 +
  1896 + # print('UpdateScene: get {}'.format(count))
  1897 +
  1898 + # use of CallAfter is mandatory otherwise crashes the wx interface
  1899 + if self.view_tracts:
  1900 + bundle, affine_vtk, coord_offset = self.tracts_queue.get_nowait()
  1901 + wx.CallAfter(Publisher.sendMessage, 'Remove tracts')
  1902 + wx.CallAfter(Publisher.sendMessage, 'Update tracts', flag=True, root=bundle,
  1903 + affine_vtk=affine_vtk)
  1904 + wx.CallAfter(Publisher.sendMessage, 'Update marker offset', coord_offset=coord_offset)
  1905 + self.tracts_queue.task_done()
  1906 +
  1907 + if self.trigger_state:
  1908 + trigger_on = self.trigger_queue.get_nowait()
  1909 + if trigger_on:
  1910 + wx.CallAfter(Publisher.sendMessage, 'Create marker')
  1911 + self.trigger_queue.task_done()
  1912 +
  1913 + # TODO: If using the view_tracts substitute the raw coord from the offset coordinate, so the user
  1914 + # see the red cross in the position of the offset marker
  1915 + wx.CallAfter(Publisher.sendMessage, 'Update cross position', arg=m_img, position=coord)
  1916 +
  1917 + if view_obj:
  1918 + wx.CallAfter(Publisher.sendMessage, 'Update object matrix', m_img=m_img, coord=coord)
  1919 +
  1920 + self.coord_queue.task_done()
  1921 + # print('UpdateScene: done {}'.format(count))
  1922 + # count += 1
  1923 +
  1924 + sleep(self.sle)
  1925 + except queue.Empty:
  1926 + if got_coords:
  1927 + self.coord_queue.task_done()
  1928 +
  1929 +
  1930 +class InputAttributes(object):
  1931 + # taken from https://stackoverflow.com/questions/2466191/set-attributes-from-dictionary-in-python
  1932 + def __init__(self, *initial_data, **kwargs):
  1933 + for dictionary in initial_data:
  1934 + for key in dictionary:
  1935 + setattr(self, key, dictionary[key])
  1936 + for key in kwargs:
  1937 + setattr(self, key, kwargs[key])
... ...
invesalius/project.py
... ... @@ -313,8 +313,9 @@ class Project(metaclass=Singleton):
313 313 self.spacing = project["spacing"]
314 314 if project.get("affine", ""):
315 315 self.affine = project["affine"]
  316 + # self.affine = project.get("affine")
316 317 Publisher.sendMessage('Update affine matrix',
317   - affine=self.affine, status=True)
  318 + affine=np.asarray(self.affine).reshape(4, 4), status=True)
318 319  
319 320 self.compress = project.get("compress", True)
320 321  
... ...