Gradients and Spheres#

This example shows how you can create gradient tables and sphere objects using DIPY.

Usually, as we saw in Getting started with DIPY, you load your b-values and b-vectors from disk and then you can create your own gradient table. But this time let’s say that you are an MR physicist and you want to design a new gradient scheme or you are a scientist who wants to simulate many different gradient schemes.

Now let’s assume that you are interested in creating a multi-shell acquisition with 2-shells, one at b=1000 \(s/mm^2\) and one at b=2500 \(s/mm^2\). For both shells let’s say that we want a specific number of gradients (64) and we want to have the points on the sphere evenly distributed.

This is possible using the disperse_charges which is an implementation of electrostatic repulsion [Jones1999].

Let’s start by importing the necessary modules.

import numpy as np
from dipy.core.gradients import gradient_table
from dipy.core.sphere import disperse_charges, Sphere, HemiSphere
from dipy.viz import window, actor

We can first create some random points on a HemiSphere using spherical polar coordinates.

rng = np.random.default_rng()
n_pts = 64
theta = np.pi * rng.random(n_pts)
phi = 2 * np.pi * rng.random(n_pts)
hsph_initial = HemiSphere(theta=theta, phi=phi)

Next, we call disperse_charges which will iteratively move the points so that the electrostatic potential energy is minimized.

In hsph_updated we have the updated HemiSphere with the points nicely distributed on the hemisphere. Let’s visualize them.

# Enables/disables interactive visualization
interactive = False

scene = window.Scene()
scene.SetBackground(1, 1, 1)

scene.add(actor.point(hsph_initial.vertices, window.colors.red,
                      point_radius=0.05))
scene.add(actor.point(hsph_updated.vertices, window.colors.green,
                      point_radius=0.05))

window.record(scene, out_path='initial_vs_updated.png', size=(300, 300))
if interactive:
    window.show(scene)
gradients spheres

Illustration of electrostatic repulsion of red points which become green points.

We can also create a sphere from the hemisphere and show it in the following way.

sph = Sphere(xyz=np.vstack((hsph_updated.vertices, -hsph_updated.vertices)))

scene.clear()
scene.add(actor.point(sph.vertices, window.colors.green, point_radius=0.05))

window.record(scene, out_path='full_sphere.png', size=(300, 300))
if interactive:
    window.show(scene)
gradients spheres

Full sphere.

It is time to create the Gradients. For this purpose we will use the function gradient_table and fill it with the hsph_updated vectors that we created above.

vertices = hsph_updated.vertices
values = np.ones(vertices.shape[0])

We need two stacks of vertices, one for every shell, and we need two sets of b-values, one at 1000 \(s/mm^2\), and one at 2500 \(s/mm^2\), as we discussed previously.

bvecs = np.vstack((vertices, vertices))
bvals = np.hstack((1000 * values, 2500 * values))

We can also add some b0s. Let’s add one at the beginning and one at the end.

bvecs = np.insert(bvecs, (0, bvecs.shape[0]), np.array([0, 0, 0]), axis=0)
bvals = np.insert(bvals, (0, bvals.shape[0]), 0)

print(bvals)
print(bvecs)
[   0. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000.
 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000.
 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000.
 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000.
 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000. 1000.
 1000. 1000. 1000. 1000. 1000. 2500. 2500. 2500. 2500. 2500. 2500. 2500.
 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500.
 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500.
 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500.
 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500.
 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500. 2500.    0.]
[[ 0.          0.          0.        ]
 [-0.63996136 -0.44439965  0.62686395]
 [-0.75440897  0.60871441  0.24562956]
 [ 0.98294554 -0.14395517  0.11443327]
 [-0.18514446 -0.01717564  0.98256121]
 [-0.17136791 -0.94579123  0.275884  ]
 [ 0.07448775  0.44526699  0.89229417]
 [ 0.44912135  0.7311377   0.51354424]
 [ 0.88021744 -0.30340889  0.36491137]
 [ 0.0172223   0.91751007  0.39733949]
 [ 0.41930481 -0.40200301  0.81398836]
 [ 0.74877952  0.58356541  0.31429389]
 [-0.94831563 -0.26376272  0.17642761]
 [ 0.41410954 -0.89837931  0.14638274]
 [ 0.30459549  0.90434746  0.29896029]
 [-0.24093597  0.80997809  0.53468248]
 [-0.66565665  0.22575477  0.71129178]
 [ 0.57358735  0.7945475   0.19922806]
 [-0.45769216 -0.29507511  0.83871841]
 [ 0.43966987  0.30361358  0.84528646]
 [ 0.06501272  0.16615552  0.98395411]
 [-0.48054011  0.8575286   0.18364614]
 [-0.20228912  0.96198486  0.18347819]
 [-0.8435801  -0.2955378   0.44836371]
 [-0.16049732  0.6128623   0.77371856]
 [ 0.6989936   0.1631694   0.6962641 ]
 [-0.91194621  0.40706569  0.05149406]
 [ 0.67053713 -0.46543199  0.57771362]
 [ 0.12064139 -0.72301164  0.68022042]
 [ 0.84962835  0.32092277  0.4184976 ]
 [-0.1730404  -0.59255945  0.78672124]
 [ 0.67932634 -0.73321564  0.03017523]
 [ 0.1180453  -0.98296197  0.14089383]
 [ 0.10691311  0.75643611  0.64527048]
 [-0.42459631 -0.64326575  0.63712412]
 [-0.4818426  -0.79257411  0.37370307]
 [ 0.13915753 -0.17824287  0.97409684]
 [-0.14708067 -0.82789599  0.54125365]
 [ 0.20950968 -0.88100869  0.42418084]
 [-0.87723032  0.30385801  0.37166824]
 [ 0.3499303   0.04065808  0.935893  ]
 [ 0.32479054  0.56888321  0.755568  ]
 [-0.17332499 -0.31859342  0.93191023]
 [-0.69023219  0.48037198  0.54113056]
 [-0.74284626 -0.57792813  0.33790312]
 [ 0.41623113 -0.66251656  0.62275473]
 [ 0.55127069 -0.13630237  0.82311742]
 [ 0.56848304 -0.73523947  0.36912052]
 [ 0.91224815  0.02624222  0.40879661]
 [-0.6942519  -0.10853653  0.71150131]
 [-0.44846801  0.50062334  0.74044089]
 [ 0.12500718  0.99037073  0.0594896 ]
 [ 0.80271906 -0.54954735  0.23160273]
 [ 0.1267564  -0.47100083  0.87297826]
 [-0.98201323  0.05729855  0.17990803]
 [ 0.8732449   0.48622748  0.03203404]
 [-0.26209575  0.30970766  0.91399507]
 [-0.51089319  0.72505232  0.46183036]
 [-0.69355897 -0.71593672  0.08006601]
 [-0.42084449 -0.9041882   0.07303165]
 [ 0.62875702  0.47753322  0.61369914]
 [-0.46735473  0.05715018  0.88222073]
 [-0.87290169  0.01891337  0.48752941]
 [ 0.7663917  -0.15248476  0.62401295]
 [ 0.97351627  0.18840973  0.12949072]
 [-0.63996136 -0.44439965  0.62686395]
 [-0.75440897  0.60871441  0.24562956]
 [ 0.98294554 -0.14395517  0.11443327]
 [-0.18514446 -0.01717564  0.98256121]
 [-0.17136791 -0.94579123  0.275884  ]
 [ 0.07448775  0.44526699  0.89229417]
 [ 0.44912135  0.7311377   0.51354424]
 [ 0.88021744 -0.30340889  0.36491137]
 [ 0.0172223   0.91751007  0.39733949]
 [ 0.41930481 -0.40200301  0.81398836]
 [ 0.74877952  0.58356541  0.31429389]
 [-0.94831563 -0.26376272  0.17642761]
 [ 0.41410954 -0.89837931  0.14638274]
 [ 0.30459549  0.90434746  0.29896029]
 [-0.24093597  0.80997809  0.53468248]
 [-0.66565665  0.22575477  0.71129178]
 [ 0.57358735  0.7945475   0.19922806]
 [-0.45769216 -0.29507511  0.83871841]
 [ 0.43966987  0.30361358  0.84528646]
 [ 0.06501272  0.16615552  0.98395411]
 [-0.48054011  0.8575286   0.18364614]
 [-0.20228912  0.96198486  0.18347819]
 [-0.8435801  -0.2955378   0.44836371]
 [-0.16049732  0.6128623   0.77371856]
 [ 0.6989936   0.1631694   0.6962641 ]
 [-0.91194621  0.40706569  0.05149406]
 [ 0.67053713 -0.46543199  0.57771362]
 [ 0.12064139 -0.72301164  0.68022042]
 [ 0.84962835  0.32092277  0.4184976 ]
 [-0.1730404  -0.59255945  0.78672124]
 [ 0.67932634 -0.73321564  0.03017523]
 [ 0.1180453  -0.98296197  0.14089383]
 [ 0.10691311  0.75643611  0.64527048]
 [-0.42459631 -0.64326575  0.63712412]
 [-0.4818426  -0.79257411  0.37370307]
 [ 0.13915753 -0.17824287  0.97409684]
 [-0.14708067 -0.82789599  0.54125365]
 [ 0.20950968 -0.88100869  0.42418084]
 [-0.87723032  0.30385801  0.37166824]
 [ 0.3499303   0.04065808  0.935893  ]
 [ 0.32479054  0.56888321  0.755568  ]
 [-0.17332499 -0.31859342  0.93191023]
 [-0.69023219  0.48037198  0.54113056]
 [-0.74284626 -0.57792813  0.33790312]
 [ 0.41623113 -0.66251656  0.62275473]
 [ 0.55127069 -0.13630237  0.82311742]
 [ 0.56848304 -0.73523947  0.36912052]
 [ 0.91224815  0.02624222  0.40879661]
 [-0.6942519  -0.10853653  0.71150131]
 [-0.44846801  0.50062334  0.74044089]
 [ 0.12500718  0.99037073  0.0594896 ]
 [ 0.80271906 -0.54954735  0.23160273]
 [ 0.1267564  -0.47100083  0.87297826]
 [-0.98201323  0.05729855  0.17990803]
 [ 0.8732449   0.48622748  0.03203404]
 [-0.26209575  0.30970766  0.91399507]
 [-0.51089319  0.72505232  0.46183036]
 [-0.69355897 -0.71593672  0.08006601]
 [-0.42084449 -0.9041882   0.07303165]
 [ 0.62875702  0.47753322  0.61369914]
 [-0.46735473  0.05715018  0.88222073]
 [-0.87290169  0.01891337  0.48752941]
 [ 0.7663917  -0.15248476  0.62401295]
 [ 0.97351627  0.18840973  0.12949072]
 [ 0.          0.          0.        ]]

Both b-values and b-vectors look correct. Let’s now create the GradientTable.

gtab = gradient_table(bvals, bvecs)

scene.clear()

We can also visualize the gradients. Let’s color the first shell blue and the second shell cyan.

colors_b1000 = window.colors.blue * np.ones(vertices.shape)
colors_b2500 = window.colors.cyan * np.ones(vertices.shape)
colors = np.vstack((colors_b1000, colors_b2500))
colors = np.insert(colors, (0, colors.shape[0]), np.array([0, 0, 0]), axis=0)
colors = np.ascontiguousarray(colors)

scene.add(actor.point(gtab.gradients, colors, point_radius=100))

window.record(scene, out_path='gradients.png', size=(300, 300))
if interactive:
    window.show(scene)
gradients spheres

Diffusion gradients.

References#

[Jones1999]

Jones, DK. et al. Optimal strategies for measuring diffusion in anisotropic systems by magnetic resonance imaging, Magnetic Resonance in Medicine, vol 42, no 3, 515-525, 1999.

Total running time of the script: (0 minutes 3.973 seconds)

Gallery generated by Sphinx-Gallery