Metrics for tracks, where tracks are arrays of points
Select an arbitrary point along distance on the track (curve)
Parameters: | xyz : array-like shape (N,3)
distance : float
|
---|---|
Returns: | ap : array shape (3,)
|
Examples
>>> import numpy as np
>>> from dipy.tracking.metrics import arbitrarypoint, length
>>> theta=np.pi*np.linspace(0,1,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> ap=arbitrarypoint(xyz,length(xyz)/3)
Size of track in bytes.
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | b : int
|
Center of mass of streamline
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | com : array shape (3,)
|
Examples
>>> from dipy.tracking.metrics import center_of_mass
>>> center_of_mass([])
Traceback (most recent call last):
...
ValueError: xyz array cannot be empty
>>> center_of_mass([[1,1,1]])
array([ 1., 1., 1.])
>>> xyz = np.array([[0,0,0],[1,1,1],[2,2,2]])
>>> center_of_mass(xyz)
array([ 1., 1., 1.])
downsample for a specific number of points along the curve/track
Uses the length of the curve. It works in a similar fashion to midpoint and arbitrarypoint but it also reduces the number of segments of a track.
Parameters: | xyz : array-like shape (N,3)
n_pol : int
|
---|---|
Returns: | xyz2 : array shape (M,3)
|
Examples
>>> import numpy as np
>>> # a semi-circle
>>> theta=np.pi*np.linspace(0,1,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> xyz2=downsample(xyz,3)
>>> # a cosine
>>> x=np.pi*np.linspace(0,1,100)
>>> y=np.cos(theta)
>>> z=0*y
>>> xyz=np.vstack((x,y,z)).T
>>> xyz2=downsample(xyz,3)
>>> len(xyz2)
3
>>> xyz3=downsample(xyz,10)
>>> len(xyz3)
10
Parameters: | xyz : array, shape(N,3)
|
---|---|
Returns: | ep : array, shape(3,)
|
Examples
>>> from dipy.tracking.metrics import endpoint
>>> import numpy as np
>>> theta=np.pi*np.linspace(0,1,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> ep=endpoint(xyz)
>>> ep.any()==xyz[-1].any()
True
Frenet-Serret Space Curve Invariants
Calculates the 3 vector and 2 scalar invariants of a space curve defined by vectors r = (x,y,z). If z is omitted (i.e. the array xyz has shape (N,2), then the curve is only 2D (planar), but the equations are still valid.
Similar to http://www.mathworks.com/matlabcentral/fileexchange/11169
In the following equations the prime (') indicates differentiation with respect to the parameter s of a parametrised curve \mathbf{r}(s).
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | T : array shape (N,3)
N : array shape (N,3)
B : array shape (N,3)
k : array shape (N,1)
t : array shape (N,1)
|
Examples
Create a helix and calculate its tangent, normal, binormal, curvature and torsion
>>> from dipy.tracking import metrics as tm
>>> import numpy as np
>>> theta = 2*np.pi*np.linspace(0,2,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=theta/(2*np.pi)
>>> xyz=np.vstack((x,y,z)).T
>>> T,N,B,k,t=tm.frenet_serret(xyz)
Combine sets of size n from items
Parameters: | items : sequence n : int |
---|---|
Returns: | ic : iterator |
Examples
>>> from dipy.tracking.metrics import generate_combinations
>>> ic=generate_combinations(range(3),2)
>>> for i in ic: print(i)
[0, 1]
[0, 2]
[1, 2]
If any point of the track is inside a sphere of a specified center and radius return True otherwise False. Mathematicaly this can be simply described by |x-c|\le r where x a point c the center of the sphere and r the radius of the sphere.
Parameters: | xyz : array, shape (N,3)
center : array, shape (3,)
radius : float
|
---|---|
Returns: | tf : {True,False}
|
Examples
>>> from dipy.tracking.metrics import inside_sphere
>>> line=np.array(([0,0,0],[1,1,1],[2,2,2]))
>>> sph_cent=np.array([1,1,1])
>>> sph_radius = 1
>>> inside_sphere(line,sph_cent,sph_radius)
True
If a track intersects with a sphere of a specified center and radius return the points that are inside the sphere otherwise False. Mathematicaly this can be simply described by |x-c| \le r where x a point c the center of the sphere and r the radius of the sphere.
Parameters: | xyz : array, shape (N,3)
center : array, shape (3,)
radius : float
|
---|---|
Returns: | xyzn : array, shape(M,3)
|
Examples
>>> from dipy.tracking.metrics import inside_sphere_points
>>> line=np.array(([0,0,0],[1,1,1],[2,2,2]))
>>> sph_cent=np.array([1,1,1])
>>> sph_radius = 1
>>> inside_sphere_points(line,sph_cent,sph_radius)
array([[1, 1, 1]])
If any segment of the track is intersecting with a sphere of specific center and radius return True otherwise False
Parameters: | xyz : array, shape (N,3)
center : array, shape (3,)
radius : float
|
---|---|
Returns: | tf : {True, False}
>>> from dipy.tracking.metrics import intersect_sphere : >>> line=np.array(([0,0,0],[1,1,1],[2,2,2])) : >>> sph_cent=np.array([1,1,1]) : >>> sph_radius = 1 : >>> intersect_sphere(line,sph_cent,sph_radius) : True : |
Notes
The ray to sphere intersection method used here is similar with http://local.wasp.uwa.edu.au/~pbourke/geometry/sphereline/ http://local.wasp.uwa.edu.au/~pbourke/geometry/sphereline/source.cpp we just applied it for every segment neglecting the intersections where the intersecting points are not inside the segment
Euclidean length of track line
This will give length in mm if tracks are expressed in world coordinates.
Parameters: | xyz : array-like shape (N,3)
along : bool, optional
|
---|---|
Returns: | L : scalar or array shape (N-1,)
|
Examples
>>> from dipy.tracking.metrics import length
>>> xyz = np.array([[1,1,1],[2,3,4],[0,0,0]])
>>> expected_lens = np.sqrt([1+2**2+3**2, 2**2+3**2+4**2])
>>> length(xyz) == expected_lens.sum()
True
>>> len_along = length(xyz, along=True)
>>> np.allclose(len_along, expected_lens.cumsum())
True
>>> length([])
0
>>> length([[1, 2, 3]])
0
>>> length([], along=True)
array([0])
Return longest track or length sorted track indices in bundle
If sort == True, return the indices of the sorted tracks in the bundle, otherwise return the longest track.
Parameters: | bundle : sequence
sort : bool, optional
|
---|---|
Returns: | longest_or_indices : array
|
Examples
>>> from dipy.tracking.metrics import longest_track_bundle
>>> import numpy as np
>>> bundle = [np.array([[0,0,0],[2,2,2]]),np.array([[0,0,0],[4,4,4]])]
>>> longest_track_bundle(bundle)
array([[0, 0, 0],
[4, 4, 4]])
>>> longest_track_bundle(bundle, True)
array([0, 1]...)
magnitude of vector
Calculates the mean curvature of a curve
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | m : float
|
Examples
Create a straight line and a semi-circle and print their mean curvatures
>>> from dipy.tracking import metrics as tm
>>> import numpy as np
>>> x=np.linspace(0,1,100)
>>> y=0*x
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> m=tm.mean_curvature(xyz) #mean curvature straight line
>>> theta=np.pi*np.linspace(0,1,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> m=tm.mean_curvature(xyz) #mean curvature for semi-circle
Calculates the mean orientation of a curve
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | m : float
|
Midpoint of track
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | mp : array shape (3,)
|
Examples
>>> from dipy.tracking.metrics import midpoint
>>> midpoint([])
Traceback (most recent call last):
...
ValueError: xyz array cannot be empty
>>> midpoint([[1, 2, 3]])
array([1, 2, 3])
>>> xyz = np.array([[1,1,1],[2,3,4]])
>>> midpoint(xyz)
array([ 1.5, 2. , 2.5])
>>> xyz = np.array([[0,0,0],[1,1,1],[2,2,2]])
>>> midpoint(xyz)
array([ 1., 1., 1.])
>>> xyz = np.array([[0,0,0],[1,0,0],[3,0,0]])
>>> midpoint(xyz)
array([ 1.5, 0. , 0. ])
>>> xyz = np.array([[0,9,7],[1,9,7],[3,9,7]])
>>> midpoint(xyz)
array([ 1.5, 9. , 7. ])
Calculate distance from midpoint of a curve to arbitrary point p
Parameters: | xyz : array-like shape (N,3)
p : array shape (3,)
|
---|---|
Returns: | d : float
|
Examples
>>> import numpy as np
>>> from dipy.tracking.metrics import midpoint2point, midpoint
>>> theta=np.pi*np.linspace(0,1,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> dist=midpoint2point(xyz,np.array([0,0,0]))
We use PCA to calculate the 3 principal directions for a track
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | va : array_like
ve : array_like
|
Examples
>>> import numpy as np
>>> from dipy.tracking.metrics import principal_components
>>> theta=np.pi*np.linspace(0,1,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> va, ve = principal_components(xyz)
>>> np.allclose(va, [0.51010101, 0.09883545, 0])
True
Generate B-splines as documented in http://www.scipy.org/Cookbook/Interpolation
The scipy.interpolate packages wraps the netlib FITPACK routines (Dierckx) for calculating smoothing splines for various kinds of data and geometries. Although the data is evenly spaced in this example, it need not be so to use this routine.
Parameters: | xyz : array, shape (N,3)
s : float, optional
k : int, optional
nest : None or int, optional
|
---|---|
Returns: | xyzn : array, shape (M,3)
|
See also
scipy.interpolate.splprep, scipy.interpolate.splev
Examples
>>> import numpy as np
>>> t=np.linspace(0,1.75*2*np.pi,100)# make ascending spiral in 3-space
>>> x = np.sin(t)
>>> y = np.cos(t)
>>> z = t
>>> x+= np.random.normal(scale=0.1, size=x.shape) # add noise
>>> y+= np.random.normal(scale=0.1, size=y.shape)
>>> z+= np.random.normal(scale=0.1, size=z.shape)
>>> xyz=np.vstack((x,y,z)).T
>>> xyzn=spline(xyz,3,2,-1)
>>> len(xyzn) > len(xyz)
True
First point of the track
Parameters: | xyz : array, shape(N,3)
|
---|---|
Returns: | sp : array, shape(3,)
|
Examples
>>> from dipy.tracking.metrics import startpoint
>>> import numpy as np
>>> theta=np.pi*np.linspace(0,1,100)
>>> x=np.cos(theta)
>>> y=np.sin(theta)
>>> z=0*x
>>> xyz=np.vstack((x,y,z)).T
>>> sp=startpoint(xyz)
>>> sp.any()==xyz[0].any()
True
Total turning angle projected.
Project space curve to best fitting plane. Calculate the cumulative signed angle between each line segment and the previous one.
Parameters: | xyz : array-like shape (N,3)
|
---|---|
Returns: | a : scalar
|