1""" 2========================= 3Getting started with DIPY 4========================= 5 6In diffusion MRI (dMRI) usually we use three types of files, a Nifti file with 7the diffusion weighted data, and two text files one with b-values and 8one with the b-vectors. 9 10In DIPY_ we provide tools to load and process these files and we also provide 11access to publicly available datasets for those who haven't acquired yet 12their own datasets. 13 14With the following commands we can download a dMRI dataset 15""" 16 17from dipy.data import fetch_sherbrooke_3shell 18fetch_sherbrooke_3shell() 19 20""" 21By default these datasets will go in the ``.dipy`` folder inside your home 22directory. Here is how you can access them. 23""" 24 25from os.path import expanduser, join 26home = expanduser('~') 27 28""" 29``dname`` holds the directory name where the 3 files are in. 30""" 31 32dname = join(home, '.dipy', 'sherbrooke_3shell') 33 34""" 35Here, we show the complete filenames of the 3 files 36""" 37 38fdwi = join(dname, 'HARDI193.nii.gz') 39 40print(fdwi) 41 42fbval = join(dname, 'HARDI193.bval') 43 44print(fbval) 45 46fbvec = join(dname, 'HARDI193.bvec') 47 48print(fbvec) 49 50""" 51``/home/username/.dipy/sherbrooke_3shell/HARDI193.nii.gz`` 52 53``/home/username/.dipy/sherbrooke_3shell/HARDI193.bval`` 54 55``/home/username/.dipy/sherbrooke_3shell/HARDI193.bvec`` 56 57Now, that we have their filenames we can start checking what these look like. 58 59Let's start first by loading the dMRI datasets. For this purpose, we 60use a python library called nibabel_ which enables us to read and write 61neuroimaging-specific file formats. 62""" 63 64from dipy.io.image import load_nifti 65data, affine, img = load_nifti(fdwi, return_img=True) 66 67""" 68``data`` is a 4D array where the first 3 dimensions are the i, j, k voxel 69coordinates and the last dimension is the number of non-weighted (S0s) and 70diffusion-weighted volumes. 71 72We can very easily check the size of ``data`` in the following way: 73""" 74 75print(data.shape) 76 77""" 78``(128, 128, 60, 193)`` 79 80We can also check the dimensions of each voxel in the following way: 81""" 82 83print(img.header.get_zooms()[:3]) 84 85""" 86``(2.0, 2.0, 2.0)`` 87 88We can quickly visualize the results using matplotlib_. For example, 89let's show here the middle axial slices of volume 0 and volume 10. 90""" 91 92import matplotlib.pyplot as plt 93 94axial_middle = data.shape[2] // 2 95plt.figure('Showing the datasets') 96plt.subplot(1, 2, 1).set_axis_off() 97plt.imshow(data[:, :, axial_middle, 0].T, cmap='gray', origin='lower') 98plt.subplot(1, 2, 2).set_axis_off() 99plt.imshow(data[:, :, axial_middle, 10].T, cmap='gray', origin='lower') 100plt.show() 101plt.savefig('data.png', bbox_inches='tight') 102 103""" 104.. figure:: data.png 105 :align: center 106 107 Showing the middle axial slice without (left) and with (right) diffusion 108 weighting. 109 110The next step is to load the b-values and b-vectors from the disk using 111the function ``read_bvals_bvecs``. 112""" 113 114from dipy.io import read_bvals_bvecs 115bvals, bvecs = read_bvals_bvecs(fbval, fbvec) 116 117""" 118In DIPY, we use an object called ``GradientTable`` which holds all the 119acquisition specific parameters, e.g. b-values, b-vectors, timings and others. 120To create this object you can use the function ``gradient_table``. 121""" 122 123from dipy.core.gradients import gradient_table 124gtab = gradient_table(bvals, bvecs) 125 126""" 127Finally, you can use ``gtab`` (the GradientTable object) to show some 128information about the acquisition parameters 129""" 130 131print(gtab.info) 132 133""" 134B-values shape (193,) 135 min 0.000000 136 max 3500.000000 137B-vectors shape (193, 3) 138 min -0.964050 139 max 0.999992 140 141You can also see the b-values using: 142""" 143 144print(gtab.bvals) 145 146""" 147 148:: 149 150 [ 0. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 151 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 152 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 153 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 154 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 155 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 156 1000. 1000. 1000. 1000. 1000. 2000. 2000. 2000. 2000. 2000. 157 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 158 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 159 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 160 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 161 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 162 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 2000. 3500. 163 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 164 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 165 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 166 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 167 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 168 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 3500. 169 3500. 3500. 3500.] 170 171Or, for example the 10 first b-vectors using: 172""" 173 174print(gtab.bvecs[:10, :]) 175 176""" 177 178:: 179 180 array([[ 0. , 0. , 0. ], 181 [ 0.999979 , -0.00504001, -0.00402795], 182 [ 0. , 0.999992 , -0.00398794], 183 [-0.0257055 , 0.653861 , -0.756178 ], 184 [ 0.589518 , -0.769236 , -0.246462 ], 185 [-0.235785 , -0.529095 , -0.815147 ], 186 [-0.893578 , -0.263559 , -0.363394 ], 187 [ 0.79784 , 0.133726 , -0.587851 ], 188 [ 0.232937 , 0.931884 , -0.278087 ], 189 [ 0.93672 , 0.144139 , -0.31903 ]]) 190 191``gtab`` can be used to tell what part of the data is the S0 volumes 192(volumes which correspond to b-values of 0). 193""" 194 195S0s = data[:, :, :, gtab.b0s_mask] 196 197""" 198Here, we had only 1 S0 as we can verify by looking at the dimensions of S0s 199""" 200 201print(S0s.shape) 202 203""" 204``(128, 128, 60, 1)`` 205 206Just, for fun let's save this in a new Nifti file. 207""" 208 209from dipy.io.image import save_nifti 210save_nifti('HARDI193_S0.nii.gz', S0s, affine) 211 212""" 213Now, that we learned how to load dMRI datasets we can start the analysis. 214See example :ref:`example_reconst_dti` to learn how to create FA maps. 215 216.. include:: ../links_names.inc 217 218""" 219