Merge pull request #1079 from CHrlS98/grid-intersections

[DOC] Add documentation for grid_intersections

docs/source/fake_files/voxel_boundary_intersection.py

@@ -0,0 +1,5 @@
+# -*- coding: utf-8 -*-
+
+
+def subdivide_streamlines_at_voxel_faces():
+    pass

docs/source/modules/scilpy.tractanalysis.rst

@@ -34,10 +34,10 @@ scilpy.tractanalysis.features module
    :show-inheritance:
 
 
-scilpy.tractanalysis.grid\_intersections module
+scilpy.tractanalysis.voxel\_boundary\_intersection module
 ------------------------------------------------------
 
-.. automodule:: scilpy.tractanalysis.grid_intersections
+.. automodule:: scilpy.tractanalysis.voxel_boundary_intersection
    :members:
    :undoc-members:
    :show-inheritance:

scilpy/tractanalysis/afd_along_streamlines.py

@@ -6,7 +6,8 @@
 from scipy.special import lpn
 
 from scilpy.reconst.utils import find_order_from_nb_coeff
-from scilpy.tractanalysis.grid_intersections import grid_intersections
+from scilpy.tractanalysis.voxel_boundary_intersection import\
+    subdivide_streamlines_at_voxel_faces
 
 
 def afd_map_along_streamlines(sft, fodf, fodf_basis, length_weighting,
@@ -94,9 +95,10 @@ def afd_and_rd_sums_along_streamlines(sft, fodf, fodf_basis,
     weight_map = np.zeros(shape=fodf_data.shape[:-1])
 
     p_matrix = np.eye(fodf_data.shape[3]) * legendre0_at_n
-    all_crossed_indices = grid_intersections(sft.streamlines)
-    for crossed_indices in all_crossed_indices:
-        segments = crossed_indices[1:] - crossed_indices[:-1]
+    all_split_streamlines =\
+        subdivide_streamlines_at_voxel_faces(sft.streamlines)
+    for split_streamlines in all_split_streamlines:
+        segments = split_streamlines[1:] - split_streamlines[:-1]
         seg_lengths = np.linalg.norm(segments, axis=1)
 
         # Remove points where the segment is zero.
@@ -112,7 +114,7 @@ def afd_and_rd_sums_along_streamlines(sft, fodf, fodf_basis,
         closest_vertex_indices = sorted_angles[:, 0]
 
         # Those starting points are used for the segment vox_idx computations
-        strl_start = crossed_indices[non_zero_lengths]
+        strl_start = split_streamlines[non_zero_lengths]
         vox_indices = (strl_start + (0.5 * segments)).astype(int)
 
         normalization_weights = np.ones_like(seg_lengths)

scilpy/tractanalysis/bingham_metric_along_streamlines.py

@@ -2,7 +2,8 @@
 
 import numpy as np
 from scilpy.reconst.bingham import bingham_to_peak_direction
-from scilpy.tractanalysis.grid_intersections import grid_intersections
+from scilpy.tractanalysis.voxel_boundary_intersection import\
+    subdivide_streamlines_at_voxel_faces
 
 
 def bingham_metric_map_along_streamlines(sft, bingham_coeffs,
@@ -73,9 +74,10 @@ def bingham_metric_sum_along_streamlines(sft, bingham_coeffs, metric,
     weight_map = np.zeros(metric.shape[:-1])
     min_cos_theta = np.cos(np.radians(max_theta))
 
-    all_crossed_indices = grid_intersections(sft.streamlines)
-    for crossed_indices in all_crossed_indices:
-        segments = crossed_indices[1:] - crossed_indices[:-1]
+    all_split_streamlines =\
+        subdivide_streamlines_at_voxel_faces(sft.streamlines)
+    for split_streamlines in all_split_streamlines:
+        segments = split_streamlines[1:] - split_streamlines[:-1]
         seg_lengths = np.linalg.norm(segments, axis=1)
 
         # Remove points where the segment is zero.
@@ -85,7 +87,7 @@ def bingham_metric_sum_along_streamlines(sft, bingham_coeffs, metric,
         seg_lengths = seg_lengths[non_zero_lengths]
 
         # Those starting points are used for the segment vox_idx computations
-        seg_start = crossed_indices[non_zero_lengths]
+        seg_start = split_streamlines[non_zero_lengths]
         vox_indices = (seg_start + (0.5 * segments)).astype(int)
 
         normalization_weights = np.ones_like(seg_lengths)

scilpy/tractanalysis/fixel_density.py

@@ -3,7 +3,8 @@
 import numpy as np
 
 from dipy.io.streamline import load_tractogram
-from scilpy.tractanalysis.grid_intersections import grid_intersections
+from scilpy.tractanalysis.voxel_boundary_intersection import\
+    subdivide_streamlines_at_voxel_faces
 
 
 def _fixel_density_parallel(args):
@@ -21,9 +22,10 @@ def _fixel_density_single_bundle(bundle, peaks, max_theta, dps_key):
 
     min_cos_theta = np.cos(np.radians(max_theta))
 
-    all_crossed_indices = grid_intersections(sft.streamlines)
-    for i, crossed_indices in enumerate(all_crossed_indices):
-        segments = crossed_indices[1:] - crossed_indices[:-1]
+    all_split_streamlines =\
+        subdivide_streamlines_at_voxel_faces(sft.streamlines)
+    for i, split_streamlines in enumerate(all_split_streamlines):
+        segments = split_streamlines[1:] - split_streamlines[:-1]
         seg_lengths = np.linalg.norm(segments, axis=1)
 
         # Remove points where the segment is zero.
@@ -33,7 +35,7 @@ def _fixel_density_single_bundle(bundle, peaks, max_theta, dps_key):
         seg_lengths = seg_lengths[non_zero_lengths]
 
         # Those starting points are used for the segment vox_idx computations
-        seg_start = crossed_indices[non_zero_lengths]
+        seg_start = split_streamlines[non_zero_lengths]
         vox_indices = (seg_start + (0.5 * segments)).astype(int)
 
         normalized_seg = np.reshape(segments / seg_lengths[..., None], (-1, 3))

scilpy/tractanalysis/mrds_along_streamlines.py

@@ -2,7 +2,8 @@
 
 import numpy as np
 
-from scilpy.tractanalysis.grid_intersections import grid_intersections
+from scilpy.tractanalysis.voxel_boundary_intersection import\
+    subdivide_streamlines_at_voxel_faces
 
 
 def mrds_metrics_along_streamlines(sft, mrds_pdds,
@@ -79,9 +80,10 @@ def mrds_metric_sums_along_streamlines(sft, mrds_pdds, metrics,
     weight_map = np.zeros(metrics[0].shape[:-1])
     min_cos_theta = np.cos(np.radians(max_theta))
 
-    all_crossed_indices = grid_intersections(sft.streamlines)
-    for crossed_indices in all_crossed_indices:
-        segments = crossed_indices[1:] - crossed_indices[:-1]
+    all_split_streamlines =\
+        subdivide_streamlines_at_voxel_faces(sft.streamlines)
+    for split_streamlines in all_split_streamlines:
+        segments = split_streamlines[1:] - split_streamlines[:-1]
         seg_lengths = np.linalg.norm(segments, axis=1)
 
         # Remove points where the segment is zero.
@@ -91,7 +93,7 @@ def mrds_metric_sums_along_streamlines(sft, mrds_pdds, metrics,
         seg_lengths = seg_lengths[non_zero_lengths]
 
         # Those starting points are used for the segment vox_idx computations
-        seg_start = crossed_indices[non_zero_lengths]
+        seg_start = split_streamlines[non_zero_lengths]
         vox_indices = (seg_start + (0.5 * segments)).astype(int)
 
         normalization_weights = np.ones_like(seg_lengths)

scilpy/tractanalysis/grid_intersections.pyx → scilpy/tractanalysis/voxel_boundary_intersection.pyx

@@ -29,7 +29,23 @@ cdef struct Pointers:
 @cython.boundscheck(False)
 @cython.wraparound(False)
 @cython.cdivision(True)
-def grid_intersections(streamlines):
+def subdivide_streamlines_at_voxel_faces(streamlines):
+    """
+    Cut streamlines segments into smaller segments such that a segment covering
+    more than one voxel is split into smaller segments that either end or start
+    at voxel boundaries.
+
+    Parameters
+    ----------
+    streamlines: list of ndarray
+        Streamlines coordinates in voxel space, corner origin.
+
+    Returns
+    -------
+    split_coordinates: list of ndarray
+        Updated streamline coordinates with added coordinate points
+        at voxel boundaries.
+    """
     cdef:
         cnp.npy_intp nb_streamlines = len(streamlines._lengths)
         cnp.npy_intp at_point = 0

setup.py

@@ -20,11 +20,14 @@ def get_extensions():
                            ['scilpy/tractograms/uncompress.pyx'])
     quick_tools = Extension('scilpy.tractanalysis.quick_tools',
                             ['scilpy/tractanalysis/quick_tools.pyx'])
-    grid_intersections = Extension('scilpy.tractanalysis.grid_intersections',
-                                   ['scilpy/tractanalysis/grid_intersections.pyx'])
-    streamlines_metrics = Extension('scilpy.tractanalysis.streamlines_metrics',
-                                    ['scilpy/tractanalysis/streamlines_metrics.pyx'])
-    return [uncompress, quick_tools, grid_intersections, streamlines_metrics]
+    voxel_boundary_intersection =\
+        Extension('scilpy.tractanalysis.voxel_boundary_intersection',
+                  ['scilpy/tractanalysis/voxel_boundary_intersection.pyx'])
+    streamlines_metrics =\
+        Extension('scilpy.tractanalysis.streamlines_metrics',
+                  ['scilpy/tractanalysis/streamlines_metrics.pyx'])
+    return [uncompress, quick_tools,
+            voxel_boundary_intersection, streamlines_metrics]
 
 
 class CustomBuildExtCommand(build_ext):