sims
Module: sims.phantom
add_noise (vol[, snr, S0, noise_type])
|
Add noise of specified distribution to a 4D array. |
diff2eigenvectors (dx, dy, dz)
|
numerical derivatives 2 eigenvectors |
get_fnames ([name])
|
Provide full paths to example or test datasets. |
gradient_table (bvals[, bvecs, big_delta, ...])
|
A general function for creating diffusion MR gradients. |
orbital_phantom ([gtab, evals, func, t, ...])
|
Create a phantom based on a 3-D orbit f(t) -> (x,y,z) . |
single_tensor (gtab[, S0, evals, evecs, snr])
|
Simulate diffusion-weighted signals with a single tensor. |
vec2vec_rotmat (u, v)
|
rotation matrix from 2 unit vectors |
Module: sims.voxel
GradientTable (gradients[, big_delta, ...])
|
Diffusion gradient information |
add_noise (signal, snr, S0[, noise_type])
|
Add noise of specified distribution to the signal from a single voxel. |
all_tensor_evecs (e0)
|
Given the principle tensor axis, return the array of all eigenvectors column-wise (or, the rotation matrix that orientates the tensor). |
callaghan_perpendicular (q, radius)
|
Calculates the perpendicular diffusion signal E(q) in a cylinder of radius R using the Soderman model [1]. |
cylinders_and_ball_soderman (gtab, tau[, ...])
|
Calculates the three-dimensional signal attenuation E(q) originating from within a cylinder of radius R using the Soderman approximation [1]. |
dki_design_matrix (gtab)
|
Construct B design matrix for DKI. |
dki_signal (gtab, dt, kt[, S0, snr])
|
Simulated signal based on the diffusion and diffusion kurtosis tensors of a single voxel. |
dot (a, b[, out])
|
Dot product of two arrays. |
gaussian_parallel (q, tau[, D])
|
Calculates the parallel Gaussian diffusion signal. |
kurtosis_element (D_comps, frac, ind_i, ...)
|
Computes the diffusion kurtosis tensor element (with indexes i, j, k and l) based on the individual diffusion tensor components of a multicompartmental model. |
multi_tensor (gtab, mevals[, S0, angles, ...])
|
Simulate a Multi-Tensor signal. |
multi_tensor_dki (gtab, mevals[, S0, angles, ...])
|
Simulate the diffusion-weight signal, diffusion and kurtosis tensors based on the DKI model |
multi_tensor_msd (mf[, mevals, tau])
|
Simulate a Multi-Tensor rtop. |
multi_tensor_odf (odf_verts, mevals, angles, ...)
|
Simulate a Multi-Tensor ODF. |
multi_tensor_pdf (pdf_points, mevals, angles, ...)
|
Simulate a Multi-Tensor ODF. |
multi_tensor_rtop (mf[, mevals, tau])
|
Simulate a Multi-Tensor rtop. |
single_tensor (gtab[, S0, evals, evecs, snr])
|
Simulate diffusion-weighted signals with a single tensor. |
single_tensor_msd ([evals, tau])
|
Simulate a Multi-Tensor rtop. |
single_tensor_odf (r[, evals, evecs])
|
Simulated ODF with a single tensor. |
single_tensor_pdf (r[, evals, evecs, tau])
|
Simulated ODF with a single tensor. |
single_tensor_rtop ([evals, tau])
|
Simulate a Single-Tensor rtop. |
sphere2cart (r, theta, phi)
|
Spherical to Cartesian coordinates |
sticks_and_ball (gtab[, d, S0, angles, ...])
|
Simulate the signal for a Sticks & Ball model. |
vec2vec_rotmat (u, v)
|
rotation matrix from 2 unit vectors |
add_noise
-
dipy.sims.phantom.add_noise(vol, snr=1.0, S0=None, noise_type='rician')
Add noise of specified distribution to a 4D array.
- Parameters:
- volarray, shape (X,Y,Z,W)
Diffusion measurements in W directions at each (X, Y, Z)
voxel
position.
- snrfloat, optional
The desired signal-to-noise ratio. (See notes below.)
- S0float, optional
Reference signal for specifying snr (defaults to 1).
- noise_typestring, optional
The distribution of noise added. Can be either ‘gaussian’ for Gaussian
distributed noise, ‘rician’ for Rice-distributed noise (default) or
‘rayleigh’ for a Rayleigh distribution.
- Returns:
- volarray, same shape as vol
Volume with added noise.
Notes
SNR is defined here, following [1], as S0 / sigma
, where sigma
is
the standard deviation of the two Gaussian distributions forming the real
and imaginary components of the Rician noise distribution (see [2]).
References
[1]
Descoteaux, Angelino, Fitzgibbons and Deriche (2007) Regularized,
fast and robust q-ball imaging. MRM, 58: 497-510
[2]
Gudbjartson and Patz (2008). The Rician distribution of noisy MRI
data. MRM 34: 910-914.
Examples
>>> signal = np.arange(800).reshape(2, 2, 2, 100)
>>> signal_w_noise = add_noise(signal, snr=10, noise_type='rician')
diff2eigenvectors
-
dipy.sims.phantom.diff2eigenvectors(dx, dy, dz)
numerical derivatives 2 eigenvectors
get_fnames
-
dipy.sims.phantom.get_fnames(name='small_64D')
Provide full paths to example or test datasets.
- Parameters:
- namestr
the filename/s of which dataset to return, one of:
‘small_64D’ small region of interest nifti,bvecs,bvals 64 directions
‘small_101D’ small region of interest nifti, bvecs, bvals
101 directions
‘aniso_vox’ volume with anisotropic voxel size as Nifti
‘fornix’ 300 tracks in Trackvis format (from Pittsburgh
Brain Competition)
‘gqi_vectors’ the scanner wave vectors needed for a GQI acquisitions
of 101 directions tested on Siemens 3T Trio
‘small_25’ small ROI (10x8x2) DTI data (b value 2000, 25 directions)
‘test_piesno’ slice of N=8, K=14 diffusion data
‘reg_c’ small 2D image used for validating registration
‘reg_o’ small 2D image used for validation registration
‘cb_2’ two vectorized cingulum bundles
- Returns:
- fnamestuple
filenames for dataset
Examples
>>> import numpy as np
>>> from dipy.io.image import load_nifti
>>> from dipy.data import get_fnames
>>> fimg, fbvals, fbvecs = get_fnames('small_101D')
>>> bvals=np.loadtxt(fbvals)
>>> bvecs=np.loadtxt(fbvecs).T
>>> data, affine = load_nifti(fimg)
>>> data.shape == (6, 10, 10, 102)
True
>>> bvals.shape == (102,)
True
>>> bvecs.shape == (102, 3)
True
gradient_table
-
dipy.sims.phantom.gradient_table(bvals, bvecs=None, big_delta=None, small_delta=None, b0_threshold=50, atol=0.01, btens=None)
A general function for creating diffusion MR gradients.
It reads, loads and prepares scanner parameters like the b-values and
b-vectors so that they can be useful during the reconstruction process.
- Parameters:
- bvalscan be any of the four options
an array of shape (N,) or (1, N) or (N, 1) with the b-values.
a path for the file which contains an array like the above (1).
an array of shape (N, 4) or (4, N). Then this parameter is
considered to be a b-table which contains both bvals and bvecs. In
this case the next parameter is skipped.
a path for the file which contains an array like the one at (3).
- bvecscan be any of two options
an array of shape (N, 3) or (3, N) with the b-vectors.
a path for the file which contains an array like the previous.
- big_deltafloat
acquisition pulse separation time in seconds (default None)
- small_deltafloat
acquisition pulse duration time in seconds (default None)
- b0_thresholdfloat
All b-values with values less than or equal to bo_threshold are
considered as b0s i.e. without diffusion weighting.
- atolfloat
All b-vectors need to be unit vectors up to a tolerance.
- btenscan be any of three options
a string specifying the shape of the encoding tensor for all volumes
in data. Options: ‘LTE’, ‘PTE’, ‘STE’, ‘CTE’ corresponding to
linear, planar, spherical, and “cigar-shaped” tensor encoding.
Tensors are rotated so that linear and cigar tensors are aligned
with the corresponding gradient direction and the planar tensor’s
normal is aligned with the corresponding gradient direction.
Magnitude is scaled to match the b-value.
an array of strings of shape (N,), (N, 1), or (1, N) specifying
encoding tensor shape for each volume separately. N corresponds to
the number volumes in data. Options for elements in array: ‘LTE’,
‘PTE’, ‘STE’, ‘CTE’ corresponding to linear, planar, spherical, and
“cigar-shaped” tensor encoding. Tensors are rotated so that linear
and cigar tensors are aligned with the corresponding gradient
direction and the planar tensor’s normal is aligned with the
corresponding gradient direction. Magnitude is scaled to match the
b-value.
an array of shape (N,3,3) specifying the b-tensor of each volume
exactly. N corresponds to the number volumes in data. No rotation or
scaling is performed.
- Returns:
- gradientsGradientTable
A GradientTable with all the gradient information.
Notes
Often b0s (b-values which correspond to images without diffusion
weighting) have 0 values however in some cases the scanner cannot
provide b0s of an exact 0 value and it gives a bit higher values
e.g. 6 or 12. This is the purpose of the b0_threshold in the __init__.
We assume that the minimum number of b-values is 7.
B-vectors should be unit vectors.
Examples
>>> from dipy.core.gradients import gradient_table
>>> bvals = 1500 * np.ones(7)
>>> bvals[0] = 0
>>> sq2 = np.sqrt(2) / 2
>>> bvecs = np.array([[0, 0, 0],
... [1, 0, 0],
... [0, 1, 0],
... [0, 0, 1],
... [sq2, sq2, 0],
... [sq2, 0, sq2],
... [0, sq2, sq2]])
>>> gt = gradient_table(bvals, bvecs)
>>> gt.bvecs.shape == bvecs.shape
True
>>> gt = gradient_table(bvals, bvecs.T)
>>> gt.bvecs.shape == bvecs.T.shape
False
orbital_phantom
-
dipy.sims.phantom.orbital_phantom(gtab=None, evals=array([0.0015, 0.0004, 0.0004]), func=None, t=array([0., 0.00628947, 0.01257895, 0.01886842, 0.0251579, 0.03144737, 0.03773685, 0.04402632, 0.0503158, 0.05660527, 0.06289475, 0.06918422, 0.0754737, 0.08176317, 0.08805265, 0.09434212, 0.1006316, 0.10692107, 0.11321055, 0.11950002, 0.1257895, 0.13207897, 0.13836845, 0.14465792, 0.15094739, 0.15723687, 0.16352634, 0.16981582, 0.17610529, 0.18239477, 0.18868424, 0.19497372, 0.20126319, 0.20755267, 0.21384214, 0.22013162, 0.22642109, 0.23271057, 0.23900004, 0.24528952, 0.25157899, 0.25786847, 0.26415794, 0.27044742, 0.27673689, 0.28302637, 0.28931584, 0.29560531, 0.30189479, 0.30818426, 0.31447374, 0.32076321, 0.32705269, 0.33334216, 0.33963164, 0.34592111, 0.35221059, 0.35850006, 0.36478954, 0.37107901, 0.37736849, 0.38365796, 0.38994744, 0.39623691, 0.40252639, 0.40881586, 0.41510534, 0.42139481, 0.42768429, 0.43397376, 0.44026323, 0.44655271, 0.45284218, 0.45913166, 0.46542113, 0.47171061, 0.47800008, 0.48428956, 0.49057903, 0.49686851, 0.50315798, 0.50944746, 0.51573693, 0.52202641, 0.52831588, 0.53460536, 0.54089483, 0.54718431, 0.55347378, 0.55976326, 0.56605273, 0.57234221, 0.57863168, 0.58492115, 0.59121063, 0.5975001, 0.60378958, 0.61007905, 0.61636853, 0.622658, 0.62894748, 0.63523695, 0.64152643, 0.6478159, 0.65410538, 0.66039485, 0.66668433, 0.6729738, 0.67926328, 0.68555275, 0.69184223, 0.6981317, 0.70442118, 0.71071065, 0.71700013, 0.7232896, 0.72957907, 0.73586855, 0.74215802, 0.7484475, 0.75473697, 0.76102645, 0.76731592, 0.7736054, 0.77989487, 0.78618435, 0.79247382, 0.7987633, 0.80505277, 0.81134225, 0.81763172, 0.8239212, 0.83021067, 0.83650015, 0.84278962, 0.8490791, 0.85536857, 0.86165805, 0.86794752, 0.87423699, 0.88052647, 0.88681594, 0.89310542, 0.89939489, 0.90568437, 0.91197384, 0.91826332, 0.92455279, 0.93084227, 0.93713174, 0.94342122, 0.94971069, 0.95600017, 0.96228964, 0.96857912, 0.97486859, 0.98115807, 0.98744754, 0.99373702, 1.00002649, 1.00631597, 1.01260544, 1.01889491, 1.02518439, 1.03147386, 1.03776334, 1.04405281, 1.05034229, 1.05663176, 1.06292124, 1.06921071, 1.07550019, 1.08178966, 1.08807914, 1.09436861, 1.10065809, 1.10694756, 1.11323704, 1.11952651, 1.12581599, 1.13210546, 1.13839494, 1.14468441, 1.15097389, 1.15726336, 1.16355283, 1.16984231, 1.17613178, 1.18242126, 1.18871073, 1.19500021, 1.20128968, 1.20757916, 1.21386863, 1.22015811, 1.22644758, 1.23273706, 1.23902653, 1.24531601, 1.25160548, 1.25789496, 1.26418443, 1.27047391, 1.27676338, 1.28305286, 1.28934233, 1.29563181, 1.30192128, 1.30821075, 1.31450023, 1.3207897, 1.32707918, 1.33336865, 1.33965813, 1.3459476, 1.35223708, 1.35852655, 1.36481603, 1.3711055, 1.37739498, 1.38368445, 1.38997393, 1.3962634, 1.40255288, 1.40884235, 1.41513183, 1.4214213, 1.42771078, 1.43400025, 1.44028973, 1.4465792, 1.45286867, 1.45915815, 1.46544762, 1.4717371, 1.47802657, 1.48431605, 1.49060552, 1.496895, 1.50318447, 1.50947395, 1.51576342, 1.5220529, 1.52834237, 1.53463185, 1.54092132, 1.5472108, 1.55350027, 1.55978975, 1.56607922, 1.5723687, 1.57865817, 1.58494765, 1.59123712, 1.59752659, 1.60381607, 1.61010554, 1.61639502, 1.62268449, 1.62897397, 1.63526344, 1.64155292, 1.64784239, 1.65413187, 1.66042134, 1.66671082, 1.67300029, 1.67928977, 1.68557924, 1.69186872, 1.69815819, 1.70444767, 1.71073714, 1.71702662, 1.72331609, 1.72960557, 1.73589504, 1.74218451, 1.74847399, 1.75476346, 1.76105294, 1.76734241, 1.77363189, 1.77992136, 1.78621084, 1.79250031, 1.79878979, 1.80507926, 1.81136874, 1.81765821, 1.82394769, 1.83023716, 1.83652664, 1.84281611, 1.84910559, 1.85539506, 1.86168454, 1.86797401, 1.87426349, 1.88055296, 1.88684243, 1.89313191, 1.89942138, 1.90571086, 1.91200033, 1.91828981, 1.92457928, 1.93086876, 1.93715823, 1.94344771, 1.94973718, 1.95602666, 1.96231613, 1.96860561, 1.97489508, 1.98118456, 1.98747403, 1.99376351, 2.00005298, 2.00634246, 2.01263193, 2.01892141, 2.02521088, 2.03150035, 2.03778983, 2.0440793, 2.05036878, 2.05665825, 2.06294773, 2.0692372, 2.07552668, 2.08181615, 2.08810563, 2.0943951, 2.10068458, 2.10697405, 2.11326353, 2.119553, 2.12584248, 2.13213195, 2.13842143, 2.1447109, 2.15100038, 2.15728985, 2.16357932, 2.1698688, 2.17615827, 2.18244775, 2.18873722, 2.1950267, 2.20131617, 2.20760565, 2.21389512, 2.2201846, 2.22647407, 2.23276355, 2.23905302, 2.2453425, 2.25163197, 2.25792145, 2.26421092, 2.2705004, 2.27678987, 2.28307935, 2.28936882, 2.2956583, 2.30194777, 2.30823724, 2.31452672, 2.32081619, 2.32710567, 2.33339514, 2.33968462, 2.34597409, 2.35226357, 2.35855304, 2.36484252, 2.37113199, 2.37742147, 2.38371094, 2.39000042, 2.39628989, 2.40257937, 2.40886884, 2.41515832, 2.42144779, 2.42773727, 2.43402674, 2.44031622, 2.44660569, 2.45289516, 2.45918464, 2.46547411, 2.47176359, 2.47805306, 2.48434254, 2.49063201, 2.49692149, 2.50321096, 2.50950044, 2.51578991, 2.52207939, 2.52836886, 2.53465834, 2.54094781, 2.54723729, 2.55352676, 2.55981624, 2.56610571, 2.57239519, 2.57868466, 2.58497414, 2.59126361, 2.59755308, 2.60384256, 2.61013203, 2.61642151, 2.62271098, 2.62900046, 2.63528993, 2.64157941, 2.64786888, 2.65415836, 2.66044783, 2.66673731, 2.67302678, 2.67931626, 2.68560573, 2.69189521, 2.69818468, 2.70447416, 2.71076363, 2.71705311, 2.72334258, 2.72963206, 2.73592153, 2.742211, 2.74850048, 2.75478995, 2.76107943, 2.7673689, 2.77365838, 2.77994785, 2.78623733, 2.7925268, 2.79881628, 2.80510575, 2.81139523, 2.8176847, 2.82397418, 2.83026365, 2.83655313, 2.8428426, 2.84913208, 2.85542155, 2.86171103, 2.8680005, 2.87428998, 2.88057945, 2.88686892, 2.8931584, 2.89944787, 2.90573735, 2.91202682, 2.9183163, 2.92460577, 2.93089525, 2.93718472, 2.9434742, 2.94976367, 2.95605315, 2.96234262, 2.9686321, 2.97492157, 2.98121105, 2.98750052, 2.99379, 3.00007947, 3.00636895, 3.01265842, 3.0189479, 3.02523737, 3.03152684, 3.03781632, 3.04410579, 3.05039527, 3.05668474, 3.06297422, 3.06926369, 3.07555317, 3.08184264, 3.08813212, 3.09442159, 3.10071107, 3.10700054, 3.11329002, 3.11957949, 3.12586897, 3.13215844, 3.13844792, 3.14473739, 3.15102687, 3.15731634, 3.16360582, 3.16989529, 3.17618476, 3.18247424, 3.18876371, 3.19505319, 3.20134266, 3.20763214, 3.21392161, 3.22021109, 3.22650056, 3.23279004, 3.23907951, 3.24536899, 3.25165846, 3.25794794, 3.26423741, 3.27052689, 3.27681636, 3.28310584, 3.28939531, 3.29568479, 3.30197426, 3.30826374, 3.31455321, 3.32084268, 3.32713216, 3.33342163, 3.33971111, 3.34600058, 3.35229006, 3.35857953, 3.36486901, 3.37115848, 3.37744796, 3.38373743, 3.39002691, 3.39631638, 3.40260586, 3.40889533, 3.41518481, 3.42147428, 3.42776376, 3.43405323, 3.44034271, 3.44663218, 3.45292166, 3.45921113, 3.4655006, 3.47179008, 3.47807955, 3.48436903, 3.4906585, 3.49694798, 3.50323745, 3.50952693, 3.5158164, 3.52210588, 3.52839535, 3.53468483, 3.5409743, 3.54726378, 3.55355325, 3.55984273, 3.5661322, 3.57242168, 3.57871115, 3.58500063, 3.5912901, 3.59757958, 3.60386905, 3.61015852, 3.616448, 3.62273747, 3.62902695, 3.63531642, 3.6416059, 3.64789537, 3.65418485, 3.66047432, 3.6667638, 3.67305327, 3.67934275, 3.68563222, 3.6919217, 3.69821117, 3.70450065, 3.71079012, 3.7170796, 3.72336907, 3.72965855, 3.73594802, 3.7422375, 3.74852697, 3.75481644, 3.76110592, 3.76739539, 3.77368487, 3.77997434, 3.78626382, 3.79255329, 3.79884277, 3.80513224, 3.81142172, 3.81771119, 3.82400067, 3.83029014, 3.83657962, 3.84286909, 3.84915857, 3.85544804, 3.86173752, 3.86802699, 3.87431647, 3.88060594, 3.88689542, 3.89318489, 3.89947436, 3.90576384, 3.91205331, 3.91834279, 3.92463226, 3.93092174, 3.93721121, 3.94350069, 3.94979016, 3.95607964, 3.96236911, 3.96865859, 3.97494806, 3.98123754, 3.98752701, 3.99381649, 4.00010596, 4.00639544, 4.01268491, 4.01897439, 4.02526386, 4.03155334, 4.03784281, 4.04413228, 4.05042176, 4.05671123, 4.06300071, 4.06929018, 4.07557966, 4.08186913, 4.08815861, 4.09444808, 4.10073756, 4.10702703, 4.11331651, 4.11960598, 4.12589546, 4.13218493, 4.13847441, 4.14476388, 4.15105336, 4.15734283, 4.16363231, 4.16992178, 4.17621126, 4.18250073, 4.1887902, 4.19507968, 4.20136915, 4.20765863, 4.2139481, 4.22023758, 4.22652705, 4.23281653, 4.239106, 4.24539548, 4.25168495, 4.25797443, 4.2642639, 4.27055338, 4.27684285, 4.28313233, 4.2894218, 4.29571128, 4.30200075, 4.30829023, 4.3145797, 4.32086918, 4.32715865, 4.33344812, 4.3397376, 4.34602707, 4.35231655, 4.35860602, 4.3648955, 4.37118497, 4.37747445, 4.38376392, 4.3900534, 4.39634287, 4.40263235, 4.40892182, 4.4152113, 4.42150077, 4.42779025, 4.43407972, 4.4403692, 4.44665867, 4.45294815, 4.45923762, 4.4655271, 4.47181657, 4.47810604, 4.48439552, 4.49068499, 4.49697447, 4.50326394, 4.50955342, 4.51584289, 4.52213237, 4.52842184, 4.53471132, 4.54100079, 4.54729027, 4.55357974, 4.55986922, 4.56615869, 4.57244817, 4.57873764, 4.58502712, 4.59131659, 4.59760607, 4.60389554, 4.61018502, 4.61647449, 4.62276396, 4.62905344, 4.63534291, 4.64163239, 4.64792186, 4.65421134, 4.66050081, 4.66679029, 4.67307976, 4.67936924, 4.68565871, 4.69194819, 4.69823766, 4.70452714, 4.71081661, 4.71710609, 4.72339556, 4.72968504, 4.73597451, 4.74226399, 4.74855346, 4.75484294, 4.76113241, 4.76742188, 4.77371136, 4.78000083, 4.78629031, 4.79257978, 4.79886926, 4.80515873, 4.81144821, 4.81773768, 4.82402716, 4.83031663, 4.83660611, 4.84289558, 4.84918506, 4.85547453, 4.86176401, 4.86805348, 4.87434296, 4.88063243, 4.88692191, 4.89321138, 4.89950086, 4.90579033, 4.9120798, 4.91836928, 4.92465875, 4.93094823, 4.9372377, 4.94352718, 4.94981665, 4.95610613, 4.9623956, 4.96868508, 4.97497455, 4.98126403, 4.9875535, 4.99384298, 5.00013245, 5.00642193, 5.0127114, 5.01900088, 5.02529035, 5.03157983, 5.0378693, 5.04415878, 5.05044825, 5.05673772, 5.0630272, 5.06931667, 5.07560615, 5.08189562, 5.0881851, 5.09447457, 5.10076405, 5.10705352, 5.113343, 5.11963247, 5.12592195, 5.13221142, 5.1385009, 5.14479037, 5.15107985, 5.15736932, 5.1636588, 5.16994827, 5.17623775, 5.18252722, 5.1888167, 5.19510617, 5.20139564, 5.20768512, 5.21397459, 5.22026407, 5.22655354, 5.23284302, 5.23913249, 5.24542197, 5.25171144, 5.25800092, 5.26429039, 5.27057987, 5.27686934, 5.28315882, 5.28944829, 5.29573777, 5.30202724, 5.30831672, 5.31460619, 5.32089567, 5.32718514, 5.33347462, 5.33976409, 5.34605356, 5.35234304, 5.35863251, 5.36492199, 5.37121146, 5.37750094, 5.38379041, 5.39007989, 5.39636936, 5.40265884, 5.40894831, 5.41523779, 5.42152726, 5.42781674, 5.43410621, 5.44039569, 5.44668516, 5.45297464, 5.45926411, 5.46555359, 5.47184306, 5.47813254, 5.48442201, 5.49071148, 5.49700096, 5.50329043, 5.50957991, 5.51586938, 5.52215886, 5.52844833, 5.53473781, 5.54102728, 5.54731676, 5.55360623, 5.55989571, 5.56618518, 5.57247466, 5.57876413, 5.58505361, 5.59134308, 5.59763256, 5.60392203, 5.61021151, 5.61650098, 5.62279046, 5.62907993, 5.6353694, 5.64165888, 5.64794835, 5.65423783, 5.6605273, 5.66681678, 5.67310625, 5.67939573, 5.6856852, 5.69197468, 5.69826415, 5.70455363, 5.7108431, 5.71713258, 5.72342205, 5.72971153, 5.736001, 5.74229048, 5.74857995, 5.75486943, 5.7611589, 5.76744838, 5.77373785, 5.78002732, 5.7863168, 5.79260627, 5.79889575, 5.80518522, 5.8114747, 5.81776417, 5.82405365, 5.83034312, 5.8366326, 5.84292207, 5.84921155, 5.85550102, 5.8617905, 5.86807997, 5.87436945, 5.88065892, 5.8869484, 5.89323787, 5.89952735, 5.90581682, 5.9121063, 5.91839577, 5.92468524, 5.93097472, 5.93726419, 5.94355367, 5.94984314, 5.95613262, 5.96242209, 5.96871157, 5.97500104, 5.98129052, 5.98757999, 5.99386947, 6.00015894, 6.00644842, 6.01273789, 6.01902737, 6.02531684, 6.03160632, 6.03789579, 6.04418527, 6.05047474, 6.05676422, 6.06305369, 6.06934316, 6.07563264, 6.08192211, 6.08821159, 6.09450106, 6.10079054, 6.10708001, 6.11336949, 6.11965896, 6.12594844, 6.13223791, 6.13852739, 6.14481686, 6.15110634, 6.15739581, 6.16368529, 6.16997476, 6.17626424, 6.18255371, 6.18884319, 6.19513266, 6.20142214, 6.20771161, 6.21400108, 6.22029056, 6.22658003, 6.23286951, 6.23915898, 6.24544846, 6.25173793, 6.25802741, 6.26431688, 6.27060636, 6.27689583, 6.28318531]), datashape=(64, 64, 64, 65), origin=(32, 32, 32), scale=(25, 25, 25), angles=array([0., 0.2026834, 0.40536679, 0.60805019, 0.81073359, 1.01341699, 1.21610038, 1.41878378, 1.62146718, 1.82415057, 2.02683397, 2.22951737, 2.43220076, 2.63488416, 2.83756756, 3.04025096, 3.24293435, 3.44561775, 3.64830115, 3.85098454, 4.05366794, 4.25635134, 4.45903473, 4.66171813, 4.86440153, 5.06708493, 5.26976832, 5.47245172, 5.67513512, 5.87781851, 6.08050191, 6.28318531]), radii=array([0.2, 0.56, 0.92, 1.28, 1.64, 2.]), S0=100.0, snr=None)
Create a phantom based on a 3-D orbit f(t) -> (x,y,z)
.
- Parameters:
- gtabGradientTable
Gradient table of measurement directions.
- evalsarray, shape (3,)
Tensor eigenvalues.
- funcuser defined function f(t)->(x,y,z)
It could be desirable for -1=<x,y,z <=1
.
If None creates a circular orbit.
- tarray, shape (K,)
Represents time for the orbit. Default is
np.linspace(0, 2 * np.pi, 1000)
.
- datashapearray, shape (X,Y,Z,W)
Size of the output simulated data
- origintuple, shape (3,)
Define the center for the volume
- scaletuple, shape (3,)
Scale the function before applying to the grid
- anglesarray, shape (L,)
Density angle points, always perpendicular to the first eigen vector
Default np.linspace(0, 2 * np.pi, 32).
- radiiarray, shape (M,)
Thickness radii. Default np.linspace(0.2, 2, 6)
.
angles and radii define the total thickness options
- S0double, optional
Maximum simulated signal. Default 100.
- snrfloat, optional
The signal to noise ratio set to apply Rician noise to the data.
Default is to not add noise at all.
- Returns:
- dataarray, shape (datashape)
Examples
>>> def f(t):
... x = np.sin(t)
... y = np.cos(t)
... z = np.linspace(-1, 1, len(x))
... return x, y, z
>>> data = orbital_phantom(func=f)
single_tensor
-
dipy.sims.phantom.single_tensor(gtab, S0=1, evals=None, evecs=None, snr=None)
Simulate diffusion-weighted signals with a single tensor.
- Parameters:
- gtabGradientTable
Table with information of b-values and gradient directions g.
Note that if gtab has a btens attribute, simulations will be performed
according to the given b-tensor B information.
- S0double, optional
Strength of signal in the presence of no diffusion gradient (also
called the b=0
value).
- evals(3,) ndarray, optional
Eigenvalues of the diffusion tensor. By default, values typical for
prolate white matter are used.
- evecs(3, 3) ndarray, optional
Eigenvectors of the tensor. You can also think of this as a rotation
matrix that transforms the direction of the tensor. The eigenvectors
need to be column wise.
- snrfloat, optional
Signal to noise ratio, assuming Rician noise. None implies no noise.
- Returns:
- S(N,) ndarray
- Simulated signal:
S(b, g) = S_0 e^(-b g^T R D R.T g)
, if gtab.tens=None
S(B) = S_0 e^(-B:D)
, if gtab.tens information is given
References
[1]
M. Descoteaux, “High Angular Resolution Diffusion MRI: from Local
Estimation to Segmentation and Tractography”, PhD thesis,
University of Nice-Sophia Antipolis, p. 42, 2008.
[2]
E. Stejskal and J. Tanner, “Spin diffusion measurements: spin echos
in the presence of a time-dependent field gradient”, Journal of
Chemical Physics, nr. 42, pp. 288–292, 1965.
vec2vec_rotmat
-
dipy.sims.phantom.vec2vec_rotmat(u, v)
rotation matrix from 2 unit vectors
u, v being unit 3d vectors return a 3x3 rotation matrix R than aligns u to
v.
In general there are many rotations that will map u to v. If S is any
rotation using v as an axis then R.S will also map u to v since (S.R)u =
S(Ru) = Sv = v. The rotation R returned by vec2vec_rotmat leaves fixed the
perpendicular to the plane spanned by u and v.
The transpose of R will align v to u.
- Parameters:
- uarray, shape(3,)
- varray, shape(3,)
- Returns:
- Rarray, shape(3,3)
Examples
>>> import numpy as np
>>> from dipy.core.geometry import vec2vec_rotmat
>>> u=np.array([1,0,0])
>>> v=np.array([0,1,0])
>>> R=vec2vec_rotmat(u,v)
>>> np.dot(R,u)
array([ 0., 1., 0.])
>>> np.dot(R.T,v)
array([ 1., 0., 0.])
-
class dipy.sims.voxel.GradientTable(gradients, big_delta=None, small_delta=None, b0_threshold=50, btens=None)
Bases: object
Diffusion gradient information
- Parameters:
- gradientsarray_like (N, 3)
Diffusion gradients. The direction of each of these vectors corresponds
to the b-vector, and the length corresponds to the b-value.
- b0_thresholdfloat
Gradients with b-value less than or equal to b0_threshold are
considered as b0s i.e. without diffusion weighting.
Notes
The GradientTable object is immutable. Do NOT assign attributes.
If you have your gradient table in a bval & bvec format, we recommend
using the factory function gradient_table
- Attributes:
- gradients(N,3) ndarray
diffusion gradients
- bvals(N,) ndarray
The b-value, or magnitude, of each gradient direction.
- qvals: (N,) ndarray
The q-value for each gradient direction. Needs big and small
delta.
- bvecs(N,3) ndarray
The direction, represented as a unit vector, of each gradient.
- b0s_mask(N,) ndarray
Boolean array indicating which gradients have no diffusion
weighting, ie b-value is close to 0.
- b0_thresholdfloat
Gradients with b-value less than or equal to b0_threshold are
considered to not have diffusion weighting.
- btens(N,3,3) ndarray
The b-tensor of each gradient direction.
Methods
b0s_mask |
|
bvals |
|
bvecs |
|
gradient_strength |
|
qvals |
|
tau |
|
-
__init__(gradients, big_delta=None, small_delta=None, b0_threshold=50, btens=None)
Constructor for GradientTable class
-
b0s_mask()
-
bvals()
-
bvecs()
-
gradient_strength()
-
property info
-
qvals()
-
tau()
add_noise
-
dipy.sims.voxel.add_noise(signal, snr, S0, noise_type='rician')
Add noise of specified distribution to the signal from a single voxel.
- Parameters:
- signal1-d ndarray
The signal in the voxel.
- snrfloat
The desired signal-to-noise ratio. (See notes below.)
If snr is None, return the signal as-is.
- S0float
Reference signal for specifying snr.
- noise_typestring, optional
The distribution of noise added. Can be either ‘gaussian’ for Gaussian
distributed noise, ‘rician’ for Rice-distributed noise (default) or
‘rayleigh’ for a Rayleigh distribution.
- Returns:
- signalarray, same shape as the input
Signal with added noise.
Notes
SNR is defined here, following [1], as S0 / sigma
, where sigma
is
the standard deviation of the two Gaussian distributions forming the real
and imaginary components of the Rician noise distribution (see [2]).
References
[1]
Descoteaux, Angelino, Fitzgibbons and Deriche (2007) Regularized,
fast and robust q-ball imaging. MRM, 58: 497-510
[2]
Gudbjartson and Patz (2008). The Rician distribution of noisy MRI
data. MRM 34: 910-914.
Examples
>>> signal = np.arange(800).reshape(2, 2, 2, 100)
>>> signal_w_noise = add_noise(signal, 10., 100., noise_type='rician')
all_tensor_evecs
-
dipy.sims.voxel.all_tensor_evecs(e0)
Given the principle tensor axis, return the array of all
eigenvectors column-wise (or, the rotation matrix that orientates the
tensor).
- Parameters:
- e0(3,) ndarray
Principle tensor axis.
- Returns:
- evecs(3,3) ndarray
Tensor eigenvectors, arranged column-wise.
callaghan_perpendicular
-
dipy.sims.voxel.callaghan_perpendicular(q, radius)
Calculates the perpendicular diffusion signal E(q) in a cylinder of
radius R using the Soderman model [1]. Assumes that the pulse length
is infinitely short and the diffusion time is infinitely long.
- Parameters:
- qarray, shape (N,)
q-space value in 1/mm
- radiusfloat
cylinder radius in mm
- Returns:
- Earray, shape (N,)
signal attenuation
References
[1]
(1,2)
Söderman, Olle, and Bengt Jönsson. “Restricted diffusion in
cylindrical geometry.” Journal of Magnetic Resonance, Series A
117.1 (1995): 94-97.
cylinders_and_ball_soderman
-
dipy.sims.voxel.cylinders_and_ball_soderman(gtab, tau, radii=(0.005, 0.005), D=0.0007, S0=1.0, angles=((0, 0), (90, 0)), fractions=(35, 35), snr=20)
Calculates the three-dimensional signal attenuation E(q) originating
from within a cylinder of radius R using the Soderman approximation [1].
The diffusion signal is assumed to be separable perpendicular and parallel
to the cylinder axis [2].
This function is basically an extension of the ball and stick model.
Setting the radius to zero makes them equivalent.
- Parameters:
- gtabGradientTable
Signal measurement directions.
- taufloat
diffusion time in s
- radiiarray-like, optional
cylinder radius in mm
- Dfloat, optional
diffusion constant
- S0float, optional
Unweighted signal value.
- anglesarray (K, 2) or (K, 3), optional
List of K polar angles (in degrees) for the sticks or array of K
sticks as unit vectors.
- fractionsarray-like, optional
Percentage of each stick. Remainder to 100 specifies isotropic
component.
- snrfloat, optional
Signal to noise ratio, assuming Rician noise. If set to None, no
noise is added.
- Returns:
- Earray, shape (N,)
signal attenuation
References
[1]
(1,2)
Söderman, Olle, and Bengt Jönsson. “Restricted diffusion in
cylindrical geometry.” Journal of Magnetic Resonance, Series A
117.1 (1995): 94-97.
[2]
Assaf, Yaniv, et al. “New modeling and experimental framework to
characterize hindered and restricted water diffusion in brain white
matter.” Magnetic Resonance in Medicine 52.5 (2004): 965-978.
dki_design_matrix
-
dipy.sims.voxel.dki_design_matrix(gtab)
Construct B design matrix for DKI.
- Parameters:
- gtabGradientTable
Measurement directions.
- Returns:
- Barray (N, 22)
Design matrix or B matrix for the DKI model
B[j, :] = (Bxx, Bxy, Bzz, Bxz, Byz, Bzz,
Bxxxx, Byyyy, Bzzzz, Bxxxy, Bxxxz,
Bxyyy, Byyyz, Bxzzz, Byzzz, Bxxyy,
Bxxzz, Byyzz, Bxxyz, Bxyyz, Bxyzz,
BlogS0)
dki_signal
-
dipy.sims.voxel.dki_signal(gtab, dt, kt, S0=150, snr=None)
Simulated signal based on the diffusion and diffusion kurtosis
tensors of a single voxel. Simulations are preformed assuming the DKI
model.
- Parameters:
- gtabGradientTable
Measurement directions.
- dt(6,) ndarray
Elements of the diffusion tensor.
- kt(15, ) ndarray
Elements of the diffusion kurtosis tensor.
- S0float (optional)
Strength of signal in the presence of no diffusion gradient.
- snrfloat (optional)
Signal to noise ratio, assuming Rician noise. None implies no noise.
- Returns:
- S(N,) ndarray
Simulated signal based on the DKI model:
\[S=S_{0}e^{-bD+\frac{1}{6}b^{2}D^{2}K}\]
References
[1]
R. Neto Henriques et al., “Exploring the 3D geometry of the
diffusion kurtosis tensor - Impact on the development of robust
tractography procedures and novel biomarkers”, NeuroImage (2015)
111, 85-99.
dot
-
dipy.sims.voxel.dot(a, b, out=None)
Dot product of two arrays. Specifically,
If both a and b are 1-D arrays, it is inner product of vectors
(without complex conjugation).
If both a and b are 2-D arrays, it is matrix multiplication,
but using matmul()
or a @ b
is preferred.
If either a or b is 0-D (scalar), it is equivalent to
multiply()
and using numpy.multiply(a, b)
or a * b
is
preferred.
If a is an N-D array and b is a 1-D array, it is a sum product over
the last axis of a and b.
If a is an N-D array and b is an M-D array (where M>=2
), it is a
sum product over the last axis of a and the second-to-last axis of
b:
dot(a, b)[i,j,k,m] = sum(a[i,j,:] * b[k,:,m])
It uses an optimized BLAS library when possible (see numpy.linalg).
- Parameters:
- aarray_like
First argument.
- barray_like
Second argument.
- outndarray, optional
Output argument. This must have the exact kind that would be returned
if it was not used. In particular, it must have the right type, must be
C-contiguous, and its dtype must be the dtype that would be returned
for dot(a,b). This is a performance feature. Therefore, if these
conditions are not met, an exception is raised, instead of attempting
to be flexible.
- Returns:
- outputndarray
Returns the dot product of a and b. If a and b are both
scalars or both 1-D arrays then a scalar is returned; otherwise
an array is returned.
If out is given, then it is returned.
- Raises:
- ValueError
If the last dimension of a is not the same size as
the second-to-last dimension of b.
See also
vdot
Complex-conjugating dot product.
tensordot
Sum products over arbitrary axes.
einsum
Einstein summation convention.
matmul
‘@’ operator as method with out parameter.
linalg.multi_dot
Chained dot product.
Examples
Neither argument is complex-conjugated:
>>> np.dot([2j, 3j], [2j, 3j])
(-13+0j)
For 2-D arrays it is the matrix product:
>>> a = [[1, 0], [0, 1]]
>>> b = [[4, 1], [2, 2]]
>>> np.dot(a, b)
array([[4, 1],
[2, 2]])
>>> a = np.arange(3*4*5*6).reshape((3,4,5,6))
>>> b = np.arange(3*4*5*6)[::-1].reshape((5,4,6,3))
>>> np.dot(a, b)[2,3,2,1,2,2]
499128
>>> sum(a[2,3,2,:] * b[1,2,:,2])
499128
gaussian_parallel
-
dipy.sims.voxel.gaussian_parallel(q, tau, D=0.0007)
Calculates the parallel Gaussian diffusion signal.
- Parameters:
- qarray, shape (N,)
q-space value in 1/mm
- taufloat
diffusion time in s
- Dfloat, optional
diffusion constant
- Returns:
- Earray, shape (N,)
signal attenuation
kurtosis_element
-
dipy.sims.voxel.kurtosis_element(D_comps, frac, ind_i, ind_j, ind_k, ind_l, DT=None, MD=None)
Computes the diffusion kurtosis tensor element (with indexes i, j, k
and l) based on the individual diffusion tensor components of a
multicompartmental model.
- Parameters:
- D_comps(K,3,3) ndarray
Diffusion tensors for all K individual compartment of the
multicompartmental model.
- frac[float]
Percentage of the contribution of each tensor. The sum of fractions
should be equal to 100%.
- ind_iint
Element’s index i (0 for x, 1 for y, 2 for z)
- ind_jint
Element’s index j (0 for x, 1 for y, 2 for z)
- ind_kint
Element’s index k (0 for x, 1 for y, 2 for z)
- ind_l: int
Elements index l (0 for x, 1 for y, 2 for z)
- DT(3,3) ndarray (optional)
Voxel’s global diffusion tensor.
- MDfloat (optional)
Voxel’s global mean diffusivity.
- Returns:
- wijklfloat
kurtosis tensor element of index i, j, k, l
Notes
wijkl is calculated using equation 8 given in [1]
References
[1]
R. Neto Henriques et al., “Exploring the 3D geometry of the
diffusion kurtosis tensor - Impact on the development of robust
tractography procedures and novel biomarkers”, NeuroImage (2015)
111, 85-99.
multi_tensor
-
dipy.sims.voxel.multi_tensor(gtab, mevals, S0=1.0, angles=((0, 0), (90, 0)), fractions=(50, 50), snr=20)
Simulate a Multi-Tensor signal.
- Parameters:
- gtabGradientTable
Table with information of b-values and gradient directions.
Note that if gtab has a btens attribute, simulations will be performed
according to the given b-tensor information.
- mevalsarray (K, 3)
each tensor’s eigenvalues in each row
- S0float, optional
Unweighted signal value (b0 signal).
- anglesarray (K, 2) or (K, 3), optional
List of K tensor directions in polar angles (in degrees) or unit
vectors
- fractionsarray-like, optional
Percentage of the contribution of each tensor. The sum of fractions
should be equal to 100%.
- snrfloat, optional
Signal to noise ratio, assuming Rician noise. If set to None, no
noise is added.
- Returns:
- S(N,) ndarray
Simulated signal.
- sticks(M,3)
Sticks in cartesian coordinates.
Examples
>>> import numpy as np
>>> from dipy.sims.voxel import multi_tensor
>>> from dipy.data import get_fnames
>>> from dipy.core.gradients import gradient_table
>>> from dipy.io.gradients import read_bvals_bvecs
>>> fimg, fbvals, fbvecs = get_fnames('small_101D')
>>> bvals, bvecs = read_bvals_bvecs(fbvals, fbvecs)
>>> gtab = gradient_table(bvals, bvecs)
>>> mevals=np.array(([0.0015, 0.0003, 0.0003],[0.0015, 0.0003, 0.0003]))
>>> e0 = np.array([1, 0, 0.])
>>> e1 = np.array([0., 1, 0])
>>> S = multi_tensor(gtab, mevals)
multi_tensor_dki
-
dipy.sims.voxel.multi_tensor_dki(gtab, mevals, S0=1.0, angles=((90.0, 0.0), (90.0, 0.0)), fractions=(50, 50), snr=20)
Simulate the diffusion-weight signal, diffusion and kurtosis tensors
based on the DKI model
- Parameters:
- gtabGradientTable
- mevalsarray (K, 3)
eigenvalues of the diffusion tensor for each individual compartment
- S0float (optional)
Unweighted signal value (b0 signal).
- anglesarray (K,2) or (K,3) (optional)
List of K tensor directions of the diffusion tensor of each compartment
in polar angles (in degrees) or unit vectors
- fractionsfloat (K,) (optional)
Percentage of the contribution of each tensor. The sum of fractions
should be equal to 100%.
- snrfloat (optional)
Signal to noise ratio, assuming Rician noise. If set to None, no
noise is added.
- Returns:
- S(N,) ndarray
Simulated signal based on the DKI model.
- dt(6,)
elements of the diffusion tensor.
- kt(15,)
elements of the kurtosis tensor.
Notes
Simulations are based on multicompartmental models which assumes that
tissue is well described by impermeable diffusion compartments
characterized by their only diffusion tensor. Since simulations are based
on the DKI model, coefficients larger than the fourth order of the signal’s
taylor expansion approximation are neglected.
References
[1]
R. Neto Henriques et al., “Exploring the 3D geometry of the
diffusion kurtosis tensor - Impact on the development of robust
tractography procedures and novel biomarkers”, NeuroImage (2015)
111, 85-99.
Examples
>>> import numpy as np
>>> from dipy.sims.voxel import multi_tensor_dki
>>> from dipy.data import get_fnames
>>> from dipy.core.gradients import gradient_table
>>> from dipy.io.gradients import read_bvals_bvecs
>>> fimg, fbvals, fbvecs = get_fnames('small_64D')
>>> bvals, bvecs = read_bvals_bvecs(fbvals, fbvecs)
>>> bvals_2s = np.concatenate((bvals, bvals * 2), axis=0)
>>> bvecs_2s = np.concatenate((bvecs, bvecs), axis=0)
>>> gtab = gradient_table(bvals_2s, bvecs_2s)
>>> mevals = np.array([[0.00099, 0, 0],[0.00226, 0.00087, 0.00087]])
>>> S, dt, kt = multi_tensor_dki(gtab, mevals)
multi_tensor_msd
-
dipy.sims.voxel.multi_tensor_msd(mf, mevals=None, tau=0.025330295910584444)
Simulate a Multi-Tensor rtop.
- Parameters:
- mfsequence of floats, bounded [0,1]
Percentages of the fractions for each tensor.
- mevalssequence of 1D arrays,
Eigen-values for each tensor. By default, values typical for prolate
white matter are used.
- taufloat,
diffusion time. By default the value that makes q=sqrt(b).
- Returns:
- msdfloat,
Mean square displacement.
References
[1]
Cheng J., “Estimation and Processing of Ensemble Average Propagator
and Its Features in Diffusion MRI”, PhD Thesis, 2012.
multi_tensor_odf
-
dipy.sims.voxel.multi_tensor_odf(odf_verts, mevals, angles, fractions)
Simulate a Multi-Tensor ODF.
- Parameters:
- odf_verts(N,3) ndarray
Vertices of the reconstruction sphere.
- mevalssequence of 1D arrays,
Eigen-values for each tensor.
- anglessequence of 2d tuples,
Sequence of principal directions for each tensor in polar angles
or cartesian unit coordinates.
- fractionssequence of floats,
Percentages of the fractions for each tensor.
- Returns:
- ODF(N,) ndarray
Orientation distribution function.
Examples
Simulate a MultiTensor ODF with two peaks and calculate its exact ODF.
>>> import numpy as np
>>> from dipy.sims.voxel import multi_tensor_odf, all_tensor_evecs
>>> from dipy.data import default_sphere
>>> vertices, faces = default_sphere.vertices, default_sphere.faces
>>> mevals = np.array(([0.0015, 0.0003, 0.0003],[0.0015, 0.0003, 0.0003]))
>>> angles = [(0, 0), (90, 0)]
>>> odf = multi_tensor_odf(vertices, mevals, angles, [50, 50])
multi_tensor_pdf
-
dipy.sims.voxel.multi_tensor_pdf(pdf_points, mevals, angles, fractions, tau=0.025330295910584444)
Simulate a Multi-Tensor ODF.
- Parameters:
- pdf_points(N, 3) ndarray
Points to evaluate the PDF.
- mevalssequence of 1D arrays,
Eigen-values for each tensor. By default, values typical for prolate
white matter are used.
- anglessequence,
Sequence of principal directions for each tensor in polar angles
or cartesian unit coordinates.
- fractionssequence of floats,
Percentages of the fractions for each tensor.
- taufloat,
diffusion time. By default the value that makes q=sqrt(b).
- Returns:
- pdf(N,) ndarray,
Probability density function of the water displacement.
References
[1]
Cheng J., “Estimation and Processing of Ensemble Average Propagator
and its Features in Diffusion MRI”, PhD Thesis, 2012.
multi_tensor_rtop
-
dipy.sims.voxel.multi_tensor_rtop(mf, mevals=None, tau=0.025330295910584444)
Simulate a Multi-Tensor rtop.
- Parameters:
- mfsequence of floats, bounded [0,1]
Percentages of the fractions for each tensor.
- mevalssequence of 1D arrays,
Eigen-values for each tensor. By default, values typical for prolate
white matter are used.
- taufloat,
diffusion time. By default the value that makes q=sqrt(b).
- Returns:
- rtopfloat,
Return to origin probability.
References
[1]
Cheng J., “Estimation and Processing of Ensemble Average Propagator
and Its Features in Diffusion MRI”, PhD Thesis, 2012.
single_tensor
-
dipy.sims.voxel.single_tensor(gtab, S0=1, evals=None, evecs=None, snr=None)
Simulate diffusion-weighted signals with a single tensor.
- Parameters:
- gtabGradientTable
Table with information of b-values and gradient directions g.
Note that if gtab has a btens attribute, simulations will be performed
according to the given b-tensor B information.
- S0double, optional
Strength of signal in the presence of no diffusion gradient (also
called the b=0
value).
- evals(3,) ndarray, optional
Eigenvalues of the diffusion tensor. By default, values typical for
prolate white matter are used.
- evecs(3, 3) ndarray, optional
Eigenvectors of the tensor. You can also think of this as a rotation
matrix that transforms the direction of the tensor. The eigenvectors
need to be column wise.
- snrfloat, optional
Signal to noise ratio, assuming Rician noise. None implies no noise.
- Returns:
- S(N,) ndarray
- Simulated signal:
S(b, g) = S_0 e^(-b g^T R D R.T g)
, if gtab.tens=None
S(B) = S_0 e^(-B:D)
, if gtab.tens information is given
References
[1]
M. Descoteaux, “High Angular Resolution Diffusion MRI: from Local
Estimation to Segmentation and Tractography”, PhD thesis,
University of Nice-Sophia Antipolis, p. 42, 2008.
[2]
E. Stejskal and J. Tanner, “Spin diffusion measurements: spin echos
in the presence of a time-dependent field gradient”, Journal of
Chemical Physics, nr. 42, pp. 288–292, 1965.
single_tensor_msd
-
dipy.sims.voxel.single_tensor_msd(evals=None, tau=0.025330295910584444)
Simulate a Multi-Tensor rtop.
- Parameters:
- evals1D arrays,
Eigen-values for the tensor. By default, values typical for prolate
white matter are used.
- taufloat,
diffusion time. By default the value that makes q=sqrt(b).
- Returns:
- msdfloat,
Mean square displacement.
References
[1]
Cheng J., “Estimation and Processing of Ensemble Average Propagator
and Its Features in Diffusion MRI”, PhD Thesis, 2012.
single_tensor_odf
-
dipy.sims.voxel.single_tensor_odf(r, evals=None, evecs=None)
Simulated ODF with a single tensor.
- Parameters:
- r(N,3) or (M,N,3) ndarray
Measurement positions in (x, y, z), either as a list or on a grid.
- evals(3,)
Eigenvalues of diffusion tensor. By default, use values typical for
prolate white matter.
- evecs(3, 3) ndarray
Eigenvectors of the tensor, written column-wise. You can also think
of these as the rotation matrix that determines the orientation of
the diffusion tensor.
- Returns:
- ODF(N,) ndarray
The diffusion probability at r
after time tau
.
References
[1]
Aganj et al., “Reconstruction of the Orientation Distribution
Function in Single- and Multiple-Shell q-Ball Imaging Within
Constant Solid Angle”, Magnetic Resonance in Medicine, nr. 64,
pp. 554–566, 2010.
single_tensor_pdf
-
dipy.sims.voxel.single_tensor_pdf(r, evals=None, evecs=None, tau=0.025330295910584444)
Simulated ODF with a single tensor.
- Parameters:
- r(N,3) or (M,N,3) ndarray
Measurement positions in (x, y, z), either as a list or on a grid.
- evals(3,)
Eigenvalues of diffusion tensor. By default, use values typical for
prolate white matter.
- evecs(3, 3) ndarray
Eigenvectors of the tensor. You can also think of these as the
rotation matrix that determines the orientation of the diffusion
tensor.
- taufloat,
diffusion time. By default the value that makes q=sqrt(b).
- Returns:
- pdf(N,) ndarray
The diffusion probability at r
after time tau
.
References
[1]
Cheng J., “Estimation and Processing of Ensemble Average Propagator
and Its Features in Diffusion MRI”, PhD Thesis, 2012.
single_tensor_rtop
-
dipy.sims.voxel.single_tensor_rtop(evals=None, tau=0.025330295910584444)
Simulate a Single-Tensor rtop.
- Parameters:
- evals1D arrays,
Eigen-values for the tensor. By default, values typical for prolate
white matter are used.
- taufloat,
diffusion time. By default the value that makes q=sqrt(b).
- Returns:
- rtopfloat,
Return to origin probability.
References
[1]
Cheng J., “Estimation and Processing of Ensemble Average Propagator
and Its Features in Diffusion MRI”, PhD Thesis, 2012.
sphere2cart
-
dipy.sims.voxel.sphere2cart(r, theta, phi)
Spherical to Cartesian coordinates
This is the standard physics convention where theta is the
inclination (polar) angle, and phi is the azimuth angle.
Imagine a sphere with center (0,0,0). Orient it with the z axis
running south-north, the y axis running west-east and the x axis
from posterior to anterior. theta (the inclination angle) is the
angle to rotate from the z-axis (the zenith) around the y-axis,
towards the x axis. Thus the rotation is counter-clockwise from the
point of view of positive y. phi (azimuth) gives the angle of
rotation around the z-axis towards the y axis. The rotation is
counter-clockwise from the point of view of positive z.
Equivalently, given a point P on the sphere, with coordinates x, y,
z, theta is the angle between P and the z-axis, and phi is
the angle between the projection of P onto the XY plane, and the X
axis.
Geographical nomenclature designates theta as ‘co-latitude’, and phi
as ‘longitude’
- Parameters:
- rarray_like
radius
- thetaarray_like
inclination or polar angle
- phiarray_like
azimuth angle
- Returns:
- xarray
x coordinate(s) in Cartesion space
- yarray
y coordinate(s) in Cartesian space
- zarray
z coordinate
Notes
See these pages:
for excellent discussion of the many different conventions
possible. Here we use the physics conventions, used in the
wikipedia page.
Derivations of the formulae are simple. Consider a vector x, y, z of
length r (norm of x, y, z). The inclination angle (theta) can be
found from: cos(theta) == z / r -> z == r * cos(theta). This gives
the hypotenuse of the projection onto the XY plane, which we will
call Q. Q == r*sin(theta). Now x / Q == cos(phi) -> x == r *
sin(theta) * cos(phi) and so on.
We have deliberately named this function sphere2cart
rather than
sph2cart
to distinguish it from the Matlab function of that
name, because the Matlab function uses an unusual convention for the
angles that we did not want to replicate. The Matlab function is
trivial to implement with the formulae given in the Matlab help.
sticks_and_ball
-
dipy.sims.voxel.sticks_and_ball(gtab, d=0.0015, S0=1.0, angles=((0, 0), (90, 0)), fractions=(35, 35), snr=20)
Simulate the signal for a Sticks & Ball model.
- Parameters:
- gtabGradientTable
Signal measurement directions.
- dfloat, optional
Diffusivity value.
- S0float, optional
Unweighted signal value.
- anglesarray (K, 2) or (K, 3), optional
List of K polar angles (in degrees) for the sticks or array of K
sticks as unit vectors.
- fractionsarray-like, optional
Percentage of each stick. Remainder to 100 specifies isotropic
component.
- snrfloat, optional
Signal to noise ratio, assuming Rician noise. If set to None, no
noise is added.
- Returns:
- S(N,) ndarray
Simulated signal.
- sticks(M,3)
Sticks in cartesian coordinates.
References
[1]
Behrens et al., “Probabilistic diffusion
tractography with multiple fiber orientations: what can we gain?”,
Neuroimage, 2007.
vec2vec_rotmat
-
dipy.sims.voxel.vec2vec_rotmat(u, v)
rotation matrix from 2 unit vectors
u, v being unit 3d vectors return a 3x3 rotation matrix R than aligns u to
v.
In general there are many rotations that will map u to v. If S is any
rotation using v as an axis then R.S will also map u to v since (S.R)u =
S(Ru) = Sv = v. The rotation R returned by vec2vec_rotmat leaves fixed the
perpendicular to the plane spanned by u and v.
The transpose of R will align v to u.
- Parameters:
- uarray, shape(3,)
- varray, shape(3,)
- Returns:
- Rarray, shape(3,3)
Examples
>>> import numpy as np
>>> from dipy.core.geometry import vec2vec_rotmat
>>> u=np.array([1,0,0])
>>> v=np.array([0,1,0])
>>> R=vec2vec_rotmat(u,v)
>>> np.dot(R,u)
array([ 0., 1., 0.])
>>> np.dot(R.T,v)
array([ 1., 0., 0.])