1"""
2==============================================================
3Crossing invariant fiber response function with FORECAST model
4==============================================================
5
6We show how to obtain a voxel specific response function in the form of
7axially symmetric tensor and the fODF using the FORECAST model from
8[Anderson2005]_ , [Kaden2016]_ and [Zucchelli2017]_.
9
10First import the necessary modules:
11"""
12
13import numpy as np
14import matplotlib.pyplot as plt
15from dipy.reconst.forecast import ForecastModel
16from dipy.viz import actor, window
17from dipy.data import fetch_cenir_multib, read_cenir_multib, get_sphere
18
19"""
20Download and read the data for this tutorial.
21Our implementation of FORECAST requires multi-shell data.
22fetch_cenir_multib() provides data acquired using the shells at b-values 1000,
232000, and 3000 (see MAPMRI example for more information on this dataset).
24"""
25
26fetch_cenir_multib(with_raw=False)
27
28bvals = [1000, 2000, 3000]
29img, gtab = read_cenir_multib(bvals)
30data = np.asarray(img.dataobj)
31
32"""
33Let us consider only a single slice for the FORECAST fitting
34"""
35
36data_small = data[18:87, 51:52, 10:70]
37mask_small = data_small[..., 0] > 1000
38
39"""
40Instantiate the FORECAST Model.
41
42"sh_order" is the spherical harmonics order used for the fODF.
43
44dec_alg is the spherical deconvolution algorithm used for the FORECAST basis
45fitting, in this case we used the Constrained Spherical Deconvolution (CSD)
46algorithm.
47"""
48
49fm = ForecastModel(gtab, sh_order=6, dec_alg='CSD')
50
51"""
52Fit the FORECAST to the data
53"""
54
55f_fit = fm.fit(data_small, mask_small)
56
57"""
58Calculate the crossing invariant tensor indices [Kaden2016]_ : the parallel
59diffusivity, the perpendicular diffusivity, the fractional anisotropy and the
60mean diffusivity.
61"""
62
63d_par = f_fit.dpar
64d_perp = f_fit.dperp
65fa = f_fit.fractional_anisotropy()
66md = f_fit.mean_diffusivity()
67
68"""
69Show the indices and save them in FORECAST_indices.png.
70"""
71
72fig = plt.figure(figsize=(6, 6))
73ax1 = fig.add_subplot(2, 2, 1, title='parallel diffusivity')
74ax1.set_axis_off()
75ind = ax1.imshow(d_par[:, 0, :].T, interpolation='nearest',
76                 origin='lower', cmap=plt.cm.gray)
77plt.colorbar(ind, shrink=0.6)
78ax2 = fig.add_subplot(2, 2, 2, title='perpendicular diffusivity')
79ax2.set_axis_off()
80ind = ax2.imshow(d_perp[:, 0, :].T, interpolation='nearest',
81                 origin='lower', cmap=plt.cm.gray)
82plt.colorbar(ind, shrink=0.6)
83ax3 = fig.add_subplot(2, 2, 3, title='fractional anisotropy')
84ax3.set_axis_off()
85ind = ax3.imshow(fa[:, 0, :].T, interpolation='nearest',
86                 origin='lower', cmap=plt.cm.gray)
87plt.colorbar(ind, shrink=0.6)
88ax4 = fig.add_subplot(2, 2, 4, title='mean diffusivity')
89ax4.set_axis_off()
90ind = ax4.imshow(md[:, 0, :].T, interpolation='nearest',
91                 origin='lower', cmap=plt.cm.gray)
92plt.colorbar(ind, shrink=0.6)
93plt.savefig('FORECAST_indices.png', dpi=300, bbox_inches='tight')
94
95"""
96.. figure:: FORECAST_indices.png
97   :align: center
98
99   **FORECAST scalar indices**.
100
101"""
102
103"""
104Load an ODF reconstruction sphere
105"""
106
107sphere = get_sphere('repulsion724')
108
109"""
110Compute the fODFs.
111"""
112
113odf = f_fit.odf(sphere)
114print('fODF.shape (%d, %d, %d, %d)' % odf.shape)
115
116"""
117Display a part of the fODFs
118"""
119
120odf_actor = actor.odf_slicer(odf[16:36, :, 30:45], sphere=sphere,
121                             colormap='plasma', scale=0.6)
122odf_actor.display(y=0)
123odf_actor.RotateX(-90)
124scene = window.Scene()
125scene.add(odf_actor)
126window.record(scene, out_path='fODFs.png', size=(600, 600), magnification=4)
127
128"""
129.. figure:: fODFs.png
130   :align: center
131
132   **Fiber Orientation Distribution Functions, in a small ROI of the brain**.
133
134References
135----------
136
137.. [Anderson2005] Anderson A. W., "Measurement of Fiber Orientation
138       Distributions Using High Angular Resolution Diffusion Imaging", Magnetic
139       Resonance in Medicine, 2005.
140
141.. [Kaden2016] Kaden E. et al., "Quantitative Mapping of the Per-Axon Diffusion
142       Coefficients in Brain White Matter", Magnetic Resonance in
143       Medicine, 2016.
144
145.. [Zucchelli2017] Zucchelli E. et al., "A generalized SMT-based framework for
146       Diffusion MRI microstructural model estimation", MICCAI Workshop
147       on Computational DIFFUSION MRI (CDMRI), 2017.
148
149.. include:: ../links_names.inc
150
151"""
152