1from copy import deepcopy 2from warnings import warn 3import types 4 5from scipy.spatial.distance import cdist 6 7import numpy as np 8from nibabel.affines import apply_affine 9from nibabel.streamlines import ArraySequence as Streamlines 10from dipy.tracking.streamlinespeed import (compress_streamlines, length, 11 set_number_of_points) 12from dipy.tracking.distances import bundles_distances_mdf 13import dipy.tracking.utils as ut 14from dipy.core.geometry import dist_to_corner 15from dipy.core.interpolation import (interpolate_vector_3d, 16 interpolate_scalar_3d) 17 18 19def unlist_streamlines(streamlines): 20 """ Return the streamlines not as a list but as an array and an offset 21 22 Parameters 23 ---------- 24 streamlines: sequence 25 26 Returns 27 ------- 28 points : array 29 offsets : array 30 31 """ 32 33 points = np.concatenate(streamlines, axis=0) 34 offsets = np.zeros(len(streamlines), dtype='i8') 35 36 curr_pos = 0 37 for (i, s) in enumerate(streamlines): 38 39 prev_pos = curr_pos 40 curr_pos += s.shape[0] 41 points[prev_pos:curr_pos] = s 42 offsets[i] = curr_pos 43 44 return points, offsets 45 46 47def relist_streamlines(points, offsets): 48 """ Given a representation of a set of streamlines as a large array and 49 an offsets array return the streamlines as a list of shorter arrays. 50 51 Parameters 52 ----------- 53 points : array 54 offsets : array 55 56 Returns 57 ------- 58 streamlines: sequence 59 """ 60 61 streamlines = [] 62 63 streamlines.append(points[0: offsets[0]]) 64 65 for i in range(len(offsets) - 1): 66 streamlines.append(points[offsets[i]: offsets[i + 1]]) 67 68 return streamlines 69 70 71def center_streamlines(streamlines): 72 """ Move streamlines to the origin 73 74 Parameters 75 ---------- 76 streamlines : list 77 List of 2D ndarrays of shape[-1]==3 78 79 Returns 80 ------- 81 new_streamlines : list 82 List of 2D ndarrays of shape[-1]==3 83 inv_shift : ndarray 84 Translation in x,y,z to go back in the initial position 85 86 """ 87 center = np.mean(np.concatenate(streamlines, axis=0), axis=0) 88 return [s - center for s in streamlines], center 89 90 91def deform_streamlines(streamlines, 92 deform_field, 93 stream_to_current_grid, 94 current_grid_to_world, 95 stream_to_ref_grid, 96 ref_grid_to_world): 97 """ Apply deformation field to streamlines 98 99 Parameters 100 ---------- 101 streamlines : list 102 List of 2D ndarrays of shape[-1]==3 103 deform_field : 4D numpy array 104 x,y,z displacements stored in volume, shape[-1]==3 105 stream_to_current_grid : array, (4, 4) 106 transform matrix voxmm space to original grid space 107 current_grid_to_world : array (4, 4) 108 transform matrix original grid space to world coordinates 109 stream_to_ref_grid : array (4, 4) 110 transform matrix voxmm space to new grid space 111 ref_grid_to_world : array(4, 4) 112 transform matrix new grid space to world coordinates 113 114 Returns 115 ------- 116 new_streamlines : list 117 List of the transformed 2D ndarrays of shape[-1]==3 118 """ 119 120 if deform_field.shape[-1] != 3: 121 raise ValueError("Last dimension of deform_field needs shape==3") 122 123 stream_in_curr_grid = transform_streamlines(streamlines, 124 stream_to_current_grid) 125 displacements = values_from_volume(deform_field, stream_in_curr_grid, 126 np.eye(4)) 127 stream_in_world = transform_streamlines(stream_in_curr_grid, 128 current_grid_to_world) 129 new_streams_in_world = [sum(d, s) for d, s in zip(displacements, 130 stream_in_world)] 131 new_streams_grid = transform_streamlines(new_streams_in_world, 132 np.linalg.inv(ref_grid_to_world)) 133 new_streamlines = transform_streamlines(new_streams_grid, 134 np.linalg.inv(stream_to_ref_grid)) 135 return new_streamlines 136 137 138def transform_streamlines(streamlines, mat, in_place=False): 139 """ Apply affine transformation to streamlines 140 141 Parameters 142 ---------- 143 streamlines : Streamlines 144 Streamlines object 145 mat : array, (4, 4) 146 transformation matrix 147 in_place : bool 148 If True then change data in place. 149 Be careful changes input streamlines. 150 151 Returns 152 ------- 153 new_streamlines : Streamlines 154 Sequence transformed 2D ndarrays of shape[-1]==3 155 """ 156 # using new Streamlines API 157 if isinstance(streamlines, Streamlines): 158 if in_place: 159 streamlines._data = apply_affine(mat, streamlines._data) 160 return streamlines 161 new_streamlines = streamlines.copy() 162 new_streamlines._data = apply_affine(mat, new_streamlines._data) 163 return new_streamlines 164 # supporting old data structure of streamlines 165 return [apply_affine(mat, s) for s in streamlines] 166 167 168def select_random_set_of_streamlines(streamlines, select, rng=None): 169 """ Select a random set of streamlines 170 171 Parameters 172 ---------- 173 streamlines : Steamlines 174 Object of 2D ndarrays of shape[-1]==3 175 176 select : int 177 Number of streamlines to select. If there are less streamlines 178 than ``select`` then ``select=len(streamlines)``. 179 180 rng : RandomState 181 Default None. 182 183 Returns 184 ------- 185 selected_streamlines : list 186 187 Notes 188 ----- 189 The same streamline will not be selected twice. 190 """ 191 len_s = len(streamlines) 192 if rng is None: 193 rng = np.random.RandomState() 194 index = rng.choice(len_s, min(select, len_s), replace=False) 195 if isinstance(streamlines, Streamlines): 196 return streamlines[index] 197 return [streamlines[i] for i in index] 198 199 200def select_by_rois(streamlines, affine, rois, include, mode=None, 201 tol=None): 202 """Select streamlines based on logical relations with several regions of 203 interest (ROIs). For example, select streamlines that pass near ROI1, 204 but only if they do not pass near ROI2. 205 206 Parameters 207 ---------- 208 streamlines : list 209 A list of candidate streamlines for selection 210 affine : array_like (4, 4) 211 The mapping from voxel coordinates to streamline points. 212 The voxel_to_rasmm matrix, typically from a NIFTI file. 213 rois : list or ndarray 214 A list of 3D arrays, each with shape (x, y, z) corresponding to the 215 shape of the brain volume, or a 4D array with shape (n_rois, x, y, 216 z). Non-zeros in each volume are considered to be within the region 217 include : array or list 218 A list or 1D array of boolean values marking inclusion or exclusion 219 criteria. If a streamline is near any of the inclusion ROIs, it 220 should evaluate to True, unless it is also near any of the exclusion 221 ROIs. 222 mode : string, optional 223 One of {"any", "all", "either_end", "both_end"}, where a streamline is 224 associated with an ROI if: 225 226 "any" : any point is within tol from ROI. Default. 227 228 "all" : all points are within tol from ROI. 229 230 "either_end" : either of the end-points is within tol from ROI 231 232 "both_end" : both end points are within tol from ROI. 233 tol : float 234 Distance (in the units of the streamlines, usually mm). If any 235 coordinate in the streamline is within this distance from the center 236 of any voxel in the ROI, the filtering criterion is set to True for 237 this streamline, otherwise False. Defaults to the distance between 238 the center of each voxel and the corner of the voxel. 239 240 Notes 241 ----- 242 The only operation currently possible is "(A or B or ...) and not (X or Y 243 or ...)", where A, B are inclusion regions and X, Y are exclusion regions. 244 245 Returns 246 ------- 247 generator 248 Generates the streamlines to be included based on these criteria. 249 250 See also 251 -------- 252 :func:`dipy.tracking.utils.near_roi` 253 :func:`dipy.tracking.utils.reduce_rois` 254 255 Examples 256 -------- 257 >>> streamlines = [np.array([[0, 0., 0.9], 258 ... [1.9, 0., 0.]]), 259 ... np.array([[0., 0., 0], 260 ... [0, 1., 1.], 261 ... [0, 2., 2.]]), 262 ... np.array([[2, 2, 2], 263 ... [3, 3, 3]])] 264 >>> mask1 = np.zeros((4, 4, 4), dtype=bool) 265 >>> mask2 = np.zeros_like(mask1) 266 >>> mask1[0, 0, 0] = True 267 >>> mask2[1, 0, 0] = True 268 >>> selection = select_by_rois(streamlines, np.eye(4), [mask1, mask2], 269 ... [True, True], 270 ... tol=1) 271 >>> list(selection) # The result is a generator 272 [array([[ 0. , 0. , 0.9], 273 [ 1.9, 0. , 0. ]]), array([[ 0., 0., 0.], 274 [ 0., 1., 1.], 275 [ 0., 2., 2.]])] 276 >>> selection = select_by_rois(streamlines, np.eye(4), [mask1, mask2], 277 ... [True, False], 278 ... tol=0.87) 279 >>> list(selection) 280 [array([[ 0., 0., 0.], 281 [ 0., 1., 1.], 282 [ 0., 2., 2.]])] 283 >>> selection = select_by_rois(streamlines, np.eye(4), [mask1, mask2], 284 ... [True, True], 285 ... mode="both_end", 286 ... tol=1.0) 287 >>> list(selection) 288 [array([[ 0. , 0. , 0.9], 289 [ 1.9, 0. , 0. ]])] 290 >>> mask2[0, 2, 2] = True 291 >>> selection = select_by_rois(streamlines, np.eye(4), [mask1, mask2], 292 ... [True, True], 293 ... mode="both_end", 294 ... tol=1.0) 295 >>> list(selection) 296 [array([[ 0. , 0. , 0.9], 297 [ 1.9, 0. , 0. ]]), array([[ 0., 0., 0.], 298 [ 0., 1., 1.], 299 [ 0., 2., 2.]])] 300 """ 301 # This calculates the maximal distance to a corner of the voxel: 302 dtc = dist_to_corner(affine) 303 if tol is None: 304 tol = dtc 305 elif tol < dtc: 306 w_s = "Tolerance input provided would create gaps in your" 307 w_s += " inclusion ROI. Setting to: %s" % dist_to_corner 308 warn(w_s) 309 tol = dtc 310 include_roi, exclude_roi = ut.reduce_rois(rois, include) 311 include_roi_coords = np.array(np.where(include_roi)).T 312 x_include_roi_coords = apply_affine(affine, include_roi_coords) 313 exclude_roi_coords = np.array(np.where(exclude_roi)).T 314 x_exclude_roi_coords = apply_affine(affine, exclude_roi_coords) 315 316 if mode is None: 317 mode = "any" 318 for sl in streamlines: 319 include = ut.streamline_near_roi(sl, x_include_roi_coords, tol=tol, 320 mode=mode) 321 exclude = ut.streamline_near_roi(sl, x_exclude_roi_coords, tol=tol, 322 mode=mode) 323 if include & ~exclude: 324 yield sl 325 326 327def cluster_confidence(streamlines, max_mdf=5, subsample=12, power=1, 328 override=False): 329 """ Computes the cluster confidence index (cci), which is an 330 estimation of the support a set of streamlines gives to 331 a particular pathway. 332 333 Ex: A single streamline with no others in the dataset 334 following a similar pathway has a low cci. A streamline 335 in a bundle of 100 streamlines that follow similar 336 pathways has a high cci. 337 338 See: Jordan et al. 2017 339 (Based on streamline MDF distance from Garyfallidis et al. 2012) 340 341 Parameters 342 ---------- 343 streamlines : list of 2D (N, 3) arrays 344 A sequence of streamlines of length N (# streamlines) 345 max_mdf : int 346 The maximum MDF distance (mm) that will be considered a 347 "supporting" streamline and included in cci calculation 348 subsample: int 349 The number of points that are considered for each streamline 350 in the calculation. To save on calculation time, each 351 streamline is subsampled to subsampleN points. 352 power: int 353 The power to which the MDF distance for each streamline 354 will be raised to determine how much it contributes to 355 the cci. High values of power make the contribution value 356 degrade much faster. E.g., a streamline with 5mm MDF 357 similarity contributes 1/5 to the cci if power is 1, but 358 only contributes 1/5^2 = 1/25 if power is 2. 359 override: bool, False by default 360 override means that the cci calculation will still occur even 361 though there are short streamlines in the dataset that may alter 362 expected behaviour. 363 364 Returns 365 ------- 366 Returns an array of CCI scores 367 368 References 369 ---------- 370 [Jordan17] Jordan K. Et al., Cluster Confidence Index: A Streamline-Wise 371 Pathway Reproducibility Metric for Diffusion-Weighted MRI Tractography, 372 Journal of Neuroimaging, vol 28, no 1, 2017. 373 374 [Garyfallidis12] Garyfallidis E. et al., QuickBundles a method for 375 tractography simplification, Frontiers in Neuroscience, 376 vol 6, no 175, 2012. 377 378 """ 379 380 # error if any streamlines are shorter than 20mm 381 lengths = list(length(streamlines)) 382 if min(lengths) < 20 and not override: 383 raise ValueError('Short streamlines found. We recommend removing them.' 384 ' To continue without removing short streamlines set' 385 ' override=True') 386 387 # calculate the pairwise MDF distance between all streamlines in dataset 388 subsamp_sls = set_number_of_points(streamlines, subsample) 389 390 cci_score_mtrx = np.zeros([len(subsamp_sls)]) 391 392 for i, sl in enumerate(subsamp_sls): 393 mdf_mx = bundles_distances_mdf([subsamp_sls[i]], subsamp_sls) 394 if (mdf_mx == 0).sum() > 1: 395 raise ValueError('Identical streamlines. CCI calculation invalid') 396 mdf_mx_oi = (mdf_mx > 0) & (mdf_mx < max_mdf) & ~ np.isnan(mdf_mx) 397 mdf_mx_oi_only = mdf_mx[mdf_mx_oi] 398 cci_score = np.sum(np.divide(1, np.power(mdf_mx_oi_only, power))) 399 cci_score_mtrx[i] = cci_score 400 401 return cci_score_mtrx 402 403 404def _orient_by_roi_generator(out, roi1, roi2): 405 """ 406 Helper function to `orient_by_rois` 407 408 Performs the inner loop separately. This is needed, because functions with 409 `yield` always return a generator 410 """ 411 for idx, sl in enumerate(out): 412 dist1 = cdist(sl, roi1, 'euclidean') 413 dist2 = cdist(sl, roi2, 'euclidean') 414 min1 = np.argmin(dist1, 0) 415 min2 = np.argmin(dist2, 0) 416 if min1[0] > min2[0]: 417 yield sl[::-1] 418 else: 419 yield sl 420 421 422def _orient_by_roi_list(out, roi1, roi2): 423 """ 424 Helper function to `orient_by_rois` 425 426 Performs the inner loop separately. This is needed, because functions with 427 `yield` always return a generator. 428 429 Flips the streamlines in place (as needed) and returns a reference to the 430 updated list. 431 """ 432 for idx, sl in enumerate(out): 433 dist1 = cdist(sl, roi1, 'euclidean') 434 dist2 = cdist(sl, roi2, 'euclidean') 435 min1 = np.argmin(dist1, 0) 436 min2 = np.argmin(dist2, 0) 437 if min1[0] > min2[0]: 438 out[idx][:, 0] = sl[::-1][:, 0] 439 out[idx][:, 1] = sl[::-1][:, 1] 440 out[idx][:, 2] = sl[::-1][:, 2] 441 return out 442 443 444def orient_by_rois(streamlines, affine, roi1, roi2, in_place=False, 445 as_generator=False): 446 """Orient a set of streamlines according to a pair of ROIs 447 448 Parameters 449 ---------- 450 streamlines : list or generator 451 List or generator of 2d arrays of 3d coordinates. Each array contains 452 the xyz coordinates of a single streamline. 453 affine : array_like (4, 4) 454 The mapping from voxel coordinates to streamline points. 455 The voxel_to_rasmm matrix, typically from a NIFTI file. 456 roi1, roi2 : ndarray 457 Binary masks designating the location of the regions of interest, or 458 coordinate arrays (n-by-3 array with ROI coordinate in each row). 459 in_place : bool 460 Whether to make the change in-place in the original list 461 (and return a reference to the list), or to make a copy of the list 462 and return this copy, with the relevant streamlines reoriented. 463 Default: False. 464 as_generator : bool 465 Whether to return a generator as output. Default: False 466 467 Returns 468 ------- 469 streamlines : list or generator 470 The same 3D arrays as a list or generator, but reoriented with respect 471 to the ROIs 472 473 Examples 474 -------- 475 >>> streamlines = [np.array([[0, 0., 0], 476 ... [1, 0., 0.], 477 ... [2, 0., 0.]]), 478 ... np.array([[2, 0., 0.], 479 ... [1, 0., 0], 480 ... [0, 0, 0.]])] 481 >>> roi1 = np.zeros((4, 4, 4), dtype=bool) 482 >>> roi2 = np.zeros_like(roi1) 483 >>> roi1[0, 0, 0] = True 484 >>> roi2[1, 0, 0] = True 485 >>> orient_by_rois(streamlines, np.eye(4), roi1, roi2) 486 [array([[ 0., 0., 0.], 487 [ 1., 0., 0.], 488 [ 2., 0., 0.]]), array([[ 0., 0., 0.], 489 [ 1., 0., 0.], 490 [ 2., 0., 0.]])] 491 492 """ 493 # If we don't already have coordinates on our hands: 494 if len(roi1.shape) == 3: 495 roi1 = np.asarray(np.where(roi1.astype(bool))).T 496 if len(roi2.shape) == 3: 497 roi2 = np.asarray(np.where(roi2.astype(bool))).T 498 499 roi1 = apply_affine(affine, roi1) 500 roi2 = apply_affine(affine, roi2) 501 502 if as_generator: 503 if in_place: 504 w_s = "Cannot return a generator when in_place is set to True" 505 raise ValueError(w_s) 506 return _orient_by_roi_generator(streamlines, roi1, roi2) 507 508 # If it's a generator on input, we may as well generate it 509 # here and now: 510 if isinstance(streamlines, types.GeneratorType): 511 out = Streamlines(streamlines) 512 513 elif in_place: 514 out = streamlines 515 else: 516 # Make a copy, so you don't change the output in place: 517 out = deepcopy(streamlines) 518 519 return _orient_by_roi_list(out, roi1, roi2) 520 521 522def _orient_by_sl_generator(out, std_array, fgarray): 523 """Helper function that implements the generator version of this """ 524 for idx, sl in enumerate(fgarray): 525 dist_direct = np.sum(np.sqrt(np.sum((sl - std_array) ** 2, -1))) 526 dist_flipped = np.sum(np.sqrt(np.sum((sl[::-1] - std_array) ** 2, -1))) 527 if dist_direct > dist_flipped: 528 yield out[idx][::-1] 529 else: 530 yield out[idx] 531 532 533def _orient_by_sl_list(out, std_array, fgarray): 534 """Helper function that implements the sequence version of this""" 535 for idx, sl in enumerate(fgarray): 536 dist_direct = np.sum(np.sqrt(np.sum((sl - std_array) ** 2, -1))) 537 dist_flipped = np.sum(np.sqrt(np.sum((sl[::-1] - std_array) ** 2, -1))) 538 if dist_direct > dist_flipped: 539 out[idx][:, 0] = out[idx][::-1][:, 0] 540 out[idx][:, 1] = out[idx][::-1][:, 1] 541 out[idx][:, 2] = out[idx][::-1][:, 2] 542 return out 543 544 545def orient_by_streamline(streamlines, standard, n_points=12, in_place=False, 546 as_generator=False): 547 """ 548 Orient a bundle of streamlines to a standard streamline. 549 550 Parameters 551 ---------- 552 streamlines : Streamlines, list 553 The input streamlines to orient. 554 standard : Streamlines, list, or ndarrray 555 This provides the standard orientation according to which the 556 streamlines in the provided bundle should be reoriented. 557 n_points: int, optional 558 The number of samples to apply to each of the streamlines. 559 in_place : bool 560 Whether to make the change in-place in the original input 561 (and return a reference), or to make a copy of the list 562 and return this copy, with the relevant streamlines reoriented. 563 Default: False. 564 as_generator : bool 565 Whether to return a generator as output. Default: False 566 567 Returns 568 ------- 569 Streamlines : with each individual array oriented to be as similar as 570 possible to the standard. 571 572 """ 573 # Start by resampling, so that distance calculation is easy: 574 fgarray = set_number_of_points(streamlines, n_points) 575 std_array = set_number_of_points([standard], n_points) 576 577 if as_generator: 578 if in_place: 579 w_s = "Cannot return a generator when in_place is set to True" 580 raise ValueError(w_s) 581 return _orient_by_sl_generator(streamlines, std_array, fgarray) 582 583 # If it's a generator on input, we may as well generate it 584 # here and now: 585 if isinstance(streamlines, types.GeneratorType): 586 out = Streamlines(streamlines) 587 588 elif in_place: 589 out = streamlines 590 else: 591 # Make a copy, so you don't change the output in place: 592 out = deepcopy(streamlines) 593 594 return _orient_by_sl_list(out, std_array, fgarray) 595 596 597def _extract_vals(data, streamlines, affine, threedvec=False): 598 """ 599 Helper function for use with `values_from_volume`. 600 601 Parameters 602 ---------- 603 data : 3D or 4D array 604 Scalar (for 3D) and vector (for 4D) values to be extracted. For 4D 605 data, interpolation will be done on the 3 spatial dimensions in each 606 volume. 607 608 streamlines : ndarray or list 609 If array, of shape (n_streamlines, n_nodes, 3) 610 If list, len(n_streamlines) with (n_nodes, 3) array in 611 each element of the list. 612 613 affine : array_like (4, 4) 614 The mapping from voxel coordinates to streamline points. 615 The voxel_to_rasmm matrix, typically from a NIFTI file. 616 617 threedvec : bool 618 Whether the last dimension has length 3. This is a special case in 619 which we can use :func:`dipy.core.interpolate.interpolate_vector_3d` 620 for the interploation of 4D volumes without looping over the elements 621 of the last dimension. 622 623 Returns 624 --------- 625 array or list (depending on the input) : values interpolate to each 626 coordinate along the length of each streamline 627 """ 628 data = data.astype(float) 629 if (isinstance(streamlines, list) or 630 isinstance(streamlines, types.GeneratorType) or 631 isinstance(streamlines, Streamlines)): 632 streamlines = transform_streamlines(streamlines, 633 np.linalg.inv(affine)) 634 vals = [] 635 for sl in streamlines: 636 if threedvec: 637 vals.append(list(interpolate_vector_3d( 638 data, sl.astype(float))[0])) 639 else: 640 vals.append(list(interpolate_scalar_3d( 641 data, sl.astype(float))[0])) 642 643 elif isinstance(streamlines, np.ndarray): 644 sl_shape = streamlines.shape 645 sl_cat = streamlines.reshape(sl_shape[0] * 646 sl_shape[1], 3).astype(float) 647 648 inv_affine = np.linalg.inv(affine) 649 sl_cat = (np.dot(sl_cat, inv_affine[:3, :3]) + 650 inv_affine[:3, 3]) 651 652 # So that we can index in one operation: 653 if threedvec: 654 vals = np.array(interpolate_vector_3d(data, sl_cat)[0]) 655 else: 656 vals = np.array(interpolate_scalar_3d(data, sl_cat)[0]) 657 vals = np.reshape(vals, (sl_shape[0], sl_shape[1], -1)) 658 if vals.shape[-1] == 1: 659 vals = np.reshape(vals, vals.shape[:-1]) 660 else: 661 raise RuntimeError("Extracting values from a volume ", 662 "requires streamlines input as an array, ", 663 "a list of arrays, or a streamline generator.") 664 665 return vals 666 667 668def values_from_volume(data, streamlines, affine): 669 """Extract values of a scalar/vector along each streamline from a volume. 670 671 Parameters 672 ---------- 673 data : 3D or 4D array 674 Scalar (for 3D) and vector (for 4D) values to be extracted. For 4D 675 data, interpolation will be done on the 3 spatial dimensions in each 676 volume. 677 678 streamlines : ndarray or list 679 If array, of shape (n_streamlines, n_nodes, 3) 680 If list, len(n_streamlines) with (n_nodes, 3) array in 681 each element of the list. 682 683 affine : array_like (4, 4) 684 The mapping from voxel coordinates to streamline points. 685 The voxel_to_rasmm matrix, typically from a NIFTI file. 686 687 Returns 688 --------- 689 array or list (depending on the input) : values interpolate to each 690 coordinate along the length of each streamline. 691 692 Notes 693 ----- 694 Values are extracted from the image based on the 3D coordinates of the 695 nodes that comprise the points in the streamline, without any interpolation 696 into segments between the nodes. Using this function with streamlines that 697 have been resampled into a very small number of nodes will result in very 698 few values. 699 """ 700 data = np.asarray(data) 701 if len(data.shape) == 4: 702 if data.shape[-1] == 3: 703 return _extract_vals(data, streamlines, affine, threedvec=True) 704 if isinstance(streamlines, types.GeneratorType): 705 streamlines = Streamlines(streamlines) 706 vals = [] 707 for ii in range(data.shape[-1]): 708 vals.append(_extract_vals(data[..., ii], streamlines, affine)) 709 710 if isinstance(vals[-1], np.ndarray): 711 return np.swapaxes(np.array(vals), 2, 1).T 712 else: 713 new_vals = [] 714 for sl_idx in range(len(streamlines)): 715 sl_vals = [] 716 for ii in range(data.shape[-1]): 717 sl_vals.append(vals[ii][sl_idx]) 718 new_vals.append(np.array(sl_vals).T) 719 return new_vals 720 721 elif len(data.shape) == 3: 722 return _extract_vals(data, streamlines, affine) 723 else: 724 raise ValueError("Data needs to have 3 or 4 dimensions") 725 726 727def nbytes(streamlines): 728 return streamlines._data.nbytes / 1024. ** 2 729