1""" 2===================== 3Gradients and Spheres 4===================== 5 6This example shows how you can create gradient tables and sphere objects using 7DIPY_. 8 9Usually, as we saw in :ref:`example_quick_start`, you load your b-values and 10b-vectors from disk and then you can create your own gradient table. But 11this time let's say that you are an MR physicist and you want to design a new 12gradient scheme or you are a scientist who wants to simulate many different 13gradient schemes. 14 15Now let's assume that you are interested in creating a multi-shell 16acquisition with 2-shells, one at b=1000 $s/mm^2$ and one at b=2500 $s/mm^2$. 17For both shells let's say that we want a specific number of gradients (64) and 18we want to have the points on the sphere evenly distributed. 19 20This is possible using the ``disperse_charges`` which is an implementation of 21electrostatic repulsion [Jones1999]_. 22""" 23 24import numpy as np 25from dipy.core.sphere import disperse_charges, Sphere, HemiSphere 26 27""" 28We can first create some random points on a ``HemiSphere`` using spherical polar 29coordinates. 30""" 31 32n_pts = 64 33theta = np.pi * np.random.rand(n_pts) 34phi = 2 * np.pi * np.random.rand(n_pts) 35hsph_initial = HemiSphere(theta=theta, phi=phi) 36 37""" 38Next, we call ``disperse_charges`` which will iteratively move the points so that 39the electrostatic potential energy is minimized. 40""" 41 42hsph_updated, potential = disperse_charges(hsph_initial, 5000) 43 44""" 45In ``hsph_updated`` we have the updated ``HemiSphere`` with the points nicely 46distributed on the hemisphere. Let's visualize them. 47""" 48 49from dipy.viz import window, actor 50 51# Enables/disables interactive visualization 52interactive = False 53 54scene = window.Scene() 55scene.SetBackground(1, 1, 1) 56 57scene.add(actor.point(hsph_initial.vertices, window.colors.red, 58 point_radius=0.05)) 59scene.add(actor.point(hsph_updated.vertices, window.colors.green, 60 point_radius=0.05)) 61 62print('Saving illustration as initial_vs_updated.png') 63window.record(scene, out_path='initial_vs_updated.png', size=(300, 300)) 64if interactive: 65 window.show(scene) 66 67""" 68.. figure:: initial_vs_updated.png 69 :align: center 70 71 Illustration of electrostatic repulsion of red points which become 72 green points. 73 74We can also create a sphere from the hemisphere and show it in the 75following way. 76 77""" 78 79sph = Sphere(xyz=np.vstack((hsph_updated.vertices, -hsph_updated.vertices))) 80 81scene.clear() 82scene.add(actor.point(sph.vertices, window.colors.green, point_radius=0.05)) 83 84print('Saving illustration as full_sphere.png') 85window.record(scene, out_path='full_sphere.png', size=(300, 300)) 86if interactive: 87 window.show(scene) 88 89""" 90.. figure:: full_sphere.png 91 :align: center 92 93 Full sphere. 94 95It is time to create the Gradients. For this purpose we will use the 96function ``gradient_table`` and fill it with the ``hsph_updated`` vectors that 97we created above. 98""" 99 100from dipy.core.gradients import gradient_table 101 102vertices = hsph_updated.vertices 103values = np.ones(vertices.shape[0]) 104 105""" 106We need two stacks of ``vertices``, one for every shell, and we need two sets 107of b-values, one at 1000 $s/mm^2$, and one at 2500 $s/mm^2$, as we discussed 108previously. 109""" 110 111bvecs = np.vstack((vertices, vertices)) 112bvals = np.hstack((1000 * values, 2500 * values)) 113 114""" 115We can also add some b0s. Let's add one at the beginning and one at the end. 116""" 117 118bvecs = np.insert(bvecs, (0, bvecs.shape[0]), np.array([0, 0, 0]), axis=0) 119bvals = np.insert(bvals, (0, bvals.shape[0]), 0) 120 121print(bvals) 122 123""" 124 125:: 126 127 [ 0. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 128 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 129 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 130 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 131 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 132 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 133 1000. 1000. 1000. 1000. 1000. 2500. 2500. 2500. 2500. 2500. 134 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 135 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 136 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 137 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 138 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 139 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 0.] 140 141""" 142 143print(bvecs) 144 145""" 146 147:: 148 149 [[ 0. 0. 0. ] 150 [-0.80451777 -0.16877559 0.56944355] 151 [ 0.32822557 -0.94355999 0.04430036] 152 [-0.23584135 -0.96241331 0.13468285] 153 [-0.39207424 -0.73505312 0.55314981] 154 [-0.32539386 -0.16751384 0.93062235] 155 [-0.82043195 -0.39411534 0.41420347] 156 [ 0.65741493 0.74947875 0.07802061] 157 [ 0.88853765 0.45303621 0.07251925] 158 [ 0.39638642 -0.15185138 0.90543855] 159 ... 160 [ 0.10175269 0.08197111 0.99142681] 161 [ 0.50577702 -0.37862345 0.77513476] 162 [ 0.42845026 0.40155296 0.80943535] 163 [ 0.26939707 0.81103868 0.51927014] 164 [-0.48938584 -0.43780086 0.75420946] 165 [ 0. 0. 0. ]] 166 167Both b-values and b-vectors look correct. Let's now create the 168``GradientTable``. 169""" 170 171gtab = gradient_table(bvals, bvecs) 172 173scene.clear() 174 175""" 176We can also visualize the gradients. Let's color the first shell blue and 177the second shell cyan. 178""" 179 180colors_b1000 = window.colors.blue * np.ones(vertices.shape) 181colors_b2500 = window.colors.cyan * np.ones(vertices.shape) 182colors = np.vstack((colors_b1000, colors_b2500)) 183colors = np.insert(colors, (0, colors.shape[0]), np.array([0, 0, 0]), axis=0) 184colors = np.ascontiguousarray(colors) 185 186scene.add(actor.point(gtab.gradients, colors, point_radius=100)) 187 188print('Saving illustration as gradients.png') 189window.record(scene, out_path='gradients.png', size=(300, 300)) 190if interactive: 191 window.show(scene) 192 193""" 194.. figure:: gradients.png 195 :align: center 196 197 Diffusion gradients. 198 199References 200---------- 201 202.. [Jones1999] Jones, DK. et al. Optimal strategies for measuring diffusion in 203 anisotropic systems by magnetic resonance imaging, Magnetic Resonance in 204 Medicine, vol 42, no 3, 515-525, 1999. 205 206.. include:: ../links_names.inc 207 208""" 209