direction
¶
Module: direction.peaks
¶
|
Non Linear Direction Finder. |
|
Get the directions of odf peaks. |
|
Fit the model to data and computes peaks and metrics |
Reshape peaks for visualization. |
PeaksAndMetrics
¶
peak_directions_nl¶
- dipy.direction.peaks.peak_directions_nl(sphere_eval, relative_peak_threshold=0.25, min_separation_angle=25, sphere=<dipy.core.sphere.HemiSphere object>, xtol=1e-07)¶
Non Linear Direction Finder.
Parameters¶
- sphere_evalcallable
A function which can be evaluated on a sphere.
- relative_peak_thresholdfloat
Only return peaks greater than
relative_peak_threshold * m
where m is the largest peak.- min_separation_anglefloat in [0, 90]
The minimum distance between directions. If two peaks are too close only the larger of the two is returned.
- sphereSphere
A discrete Sphere. The points on the sphere will be used for initial estimate of maximums.
- xtolfloat
Relative tolerance for optimization.
Returns¶
- directionsarray (N, 3)
Points on the sphere corresponding to N local maxima on the sphere.
- valuesarray (N,)
Value of sphere_eval at each point on directions.
peak_directions¶
- dipy.direction.peaks.peak_directions(odf, sphere, relative_peak_threshold=0.5, min_separation_angle=25, is_symmetric=True)¶
Get the directions of odf peaks.
Peaks are defined as points on the odf that are greater than at least one neighbor and greater than or equal to all neighbors. Peaks are sorted in descending order by their values then filtered based on their relative size and spacing on the sphere. An odf may have 0 peaks, for example if the odf is perfectly isotropic.
Parameters¶
- odf1d ndarray
The odf function evaluated on the vertices of sphere
- sphereSphere
The Sphere providing discrete directions for evaluation.
- relative_peak_thresholdfloat in [0., 1.]
Only peaks greater than
min + relative_peak_threshold * scale
are kept, wheremin = max(0, odf.min())
andscale = odf.max() - min
.- min_separation_anglefloat in [0, 90]
The minimum distance between directions. If two peaks are too close only the larger of the two is returned.
- is_symmetricbool, optional
If True, v is considered equal to -v.
Returns¶
- directions(N, 3) ndarray
N vertices for sphere, one for each peak
- values(N,) ndarray
peak values
- indices(N,) ndarray
peak indices of the directions on the sphere
Notes¶
If the odf has any negative values, they will be clipped to zeros.
peaks_from_model¶
- dipy.direction.peaks.peaks_from_model(model, data, sphere, relative_peak_threshold, min_separation_angle, mask=None, return_odf=False, return_sh=True, gfa_thr=0, normalize_peaks=False, sh_order=8, sh_basis_type=None, npeaks=5, B=None, invB=None, parallel=False, num_processes=None)¶
Fit the model to data and computes peaks and metrics
Parameters¶
- modela model instance
model will be used to fit the data.
- datandarray
Diffusion data.
- sphereSphere
The Sphere providing discrete directions for evaluation.
- relative_peak_thresholdfloat
Only return peaks greater than
relative_peak_threshold * m
where m is the largest peak.- min_separation_anglefloat in [0, 90] The minimum distance between
directions. If two peaks are too close only the larger of the two is returned.
- maskarray, optional
If mask is provided, voxels that are False in mask are skipped and no peaks are returned.
- return_odfbool
If True, the odfs are returned.
- return_shbool
If True, the odf as spherical harmonics coefficients is returned
- gfa_thrfloat
Voxels with gfa less than gfa_thr are skipped, no peaks are returned.
- normalize_peaksbool
If true, all peak values are calculated relative to max(odf).
- sh_orderint, optional
Maximum SH order in the SH fit. For sh_order, there will be
(sh_order + 1) * (sh_order + 2) / 2
SH coefficients (default 8).- sh_basis_type{None, ‘tournier07’, ‘descoteaux07’}
None
for the default DIPY basis,tournier07
for the Tournier 2007 [2] basis, anddescoteaux07
for the Descoteaux 2007 [1] basis (None
defaults todescoteaux07
).- npeaksint
Maximum number of peaks found (default 5 peaks).
- Bndarray, optional
Matrix that transforms spherical harmonics to spherical function
sf = np.dot(sh, B)
.- invBndarray, optional
Inverse of B.
- parallel: bool
If True, use multiprocessing to compute peaks and metric (default False). Temporary files are saved in the default temporary directory of the system. It can be changed using
import tempfile
andtempfile.tempdir = '/path/to/tempdir'
.- num_processes: int, optional
If parallel is True, the number of subprocesses to use (default multiprocessing.cpu_count()). If < 0 the maximal number of cores minus
num_processes + 1
is used (enter -1 to use as many cores as possible). 0 raises an error.
Returns¶
- pamPeaksAndMetrics
An object with
gfa
,peak_directions
,peak_values
,peak_indices
,odf
,shm_coeffs
as attributes
References¶
reshape_peaks_for_visualization¶
- dipy.direction.peaks.reshape_peaks_for_visualization(peaks)¶
Reshape peaks for visualization.
Reshape and convert to float32 a set of peaks for visualisation with mrtrix or the fibernavigator.
Parameters¶
- peaks: nd array (…, N, 3) or PeaksAndMetrics object
The peaks to be reshaped and converted to float32.
Returns¶
peaks : nd array (…, 3*N)