From 9a245c67b96e654664fec9c730a119192317d4d8 Mon Sep 17 00:00:00 2001 From: frheault Date: Thu, 12 Dec 2024 10:25:51 -0500 Subject: [PATCH 1/2] Switch to map coordinates --- scilpy/connectivity/connectivity.py | 38 ++++++++++++++--------------- 1 file changed, 19 insertions(+), 19 deletions(-) diff --git a/scilpy/connectivity/connectivity.py b/scilpy/connectivity/connectivity.py index dd6b75901..15f259345 100644 --- a/scilpy/connectivity/connectivity.py +++ b/scilpy/connectivity/connectivity.py @@ -10,14 +10,18 @@ import h5py import nibabel as nib import numpy as np +from scipy.ndimage import map_coordinates from scilpy.image.labels import get_data_as_labels from scilpy.io.hdf5 import reconstruct_streamlines_from_hdf5 from scilpy.tractanalysis.reproducibility_measures import \ compute_bundle_adjacency_voxel from scilpy.tractanalysis.streamlines_metrics import compute_tract_counts_map +from scilpy.tractograms.streamline_operations import \ + resample_streamlines_num_points from scilpy.utils.metrics_tools import compute_lesion_stats + d = threading.local() @@ -31,7 +35,7 @@ def compute_triu_connectivity_from_labels(tractogram, data_labels, ---------- tractogram: StatefulTractogram, or list[np.ndarray] Streamlines. A StatefulTractogram input is recommanded. - When using directly with a list of streamlines, streamlinee must be in + When using directly with a list of streamlines, streamlines must be in vox space, corner origin. data_labels: np.ndarray The loaded nifti image. @@ -58,9 +62,11 @@ def compute_triu_connectivity_from_labels(tractogram, data_labels, # Vox space, corner origin # = we can get the nearest neighbor easily. # Coord 0 = voxel 0. Coord 0.9 = voxel 0. Coord 1 = voxel 1. - tractogram.to_vox() - tractogram.to_corner() - streamlines = tractogram.streamlines + sfs_2_pts = resample_streamlines_num_points(tractogram, 2) + sfs_2_pts.to_vox() + sfs_2_pts.to_center() + streamlines = sfs_2_pts.streamlines + else: streamlines = tractogram @@ -71,23 +77,17 @@ def compute_triu_connectivity_from_labels(tractogram, data_labels, .format(nb_labels)) matrix = np.zeros((nb_labels, nb_labels), dtype=int) - start_labels = [] - end_labels = [] - - for s in streamlines: - start = ordered_labels.index( - data_labels[tuple(np.floor(s[0, :]).astype(int))]) - end = ordered_labels.index( - data_labels[tuple(np.floor(s[-1, :]).astype(int))]) - start_labels.append(start) - end_labels.append(end) + labels = map_coordinates(data_labels, streamlines._data.T, + order=0, mode='nearest') + start_labels = labels[0::2] + end_labels = labels[1::2] - matrix[start, end] += 1 - if start != end: - matrix[end, start] += 1 + # sort each pair of labels for start to be smaller than end + start_labels, end_labels = zip(*[sorted(pair) for pair in + zip(start_labels, end_labels)]) - matrix = np.triu(matrix) + np.add.at(matrix, (start_labels, end_labels), 1) assert matrix.sum() == len(streamlines) # Rejecting background @@ -249,7 +249,7 @@ def compute_connectivity_matrices_from_hdf5( if compute_volume: measures_to_return['volume_mm3'] = np.count_nonzero(density) * \ - np.prod(voxel_sizes) + np.prod(voxel_sizes) if compute_streamline_count: measures_to_return['streamline_count'] = len(streamlines) From 3918ffb2c2def5bcbcccaf69a4d18f4de1964597 Mon Sep 17 00:00:00 2001 From: frheault Date: Thu, 12 Dec 2024 10:29:34 -0500 Subject: [PATCH 2/2] Modify doc for to_center() --- scilpy/connectivity/connectivity.py | 9 +++------ 1 file changed, 3 insertions(+), 6 deletions(-) diff --git a/scilpy/connectivity/connectivity.py b/scilpy/connectivity/connectivity.py index 15f259345..11ed51a39 100644 --- a/scilpy/connectivity/connectivity.py +++ b/scilpy/connectivity/connectivity.py @@ -36,7 +36,7 @@ def compute_triu_connectivity_from_labels(tractogram, data_labels, tractogram: StatefulTractogram, or list[np.ndarray] Streamlines. A StatefulTractogram input is recommanded. When using directly with a list of streamlines, streamlines must be in - vox space, corner origin. + vox space, center origin. data_labels: np.ndarray The loaded nifti image. keep_background: Bool @@ -59,9 +59,7 @@ def compute_triu_connectivity_from_labels(tractogram, data_labels, For each streamline, the label at ending point. """ if isinstance(tractogram, StatefulTractogram): - # Vox space, corner origin - # = we can get the nearest neighbor easily. - # Coord 0 = voxel 0. Coord 0.9 = voxel 0. Coord 1 = voxel 1. + # vox space, center origin: compatible with map_coordinates sfs_2_pts = resample_streamlines_num_points(tractogram, 2) sfs_2_pts.to_vox() sfs_2_pts.to_center() @@ -78,8 +76,7 @@ def compute_triu_connectivity_from_labels(tractogram, data_labels, matrix = np.zeros((nb_labels, nb_labels), dtype=int) - labels = map_coordinates(data_labels, streamlines._data.T, - order=0, mode='nearest') + labels = map_coordinates(data_labels, streamlines._data.T, order=0) start_labels = labels[0::2] end_labels = labels[1::2]