/usr/share/doc/dipy/examples/reconst_gqi.py is in python-dipy 0.10.1-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 | """
===============================================
Reconstruct with Generalized Q-Sampling Imaging
===============================================
We show how to apply Generalized Q-Sampling Imaging [Yeh2010]_
to diffusion MRI datasets. You can think of GQI as an analytical version of
DSI orientation distribution function (ODF) (Garyfallidis, PhD thesis, 2012).
First import the necessary modules:
"""
import numpy as np
from dipy.data import fetch_taiwan_ntu_dsi, read_taiwan_ntu_dsi, get_sphere
from dipy.reconst.gqi import GeneralizedQSamplingModel
from dipy.direction import peaks_from_model
"""
Download and read the data for this tutorial.
"""
fetch_taiwan_ntu_dsi()
img, gtab = read_taiwan_ntu_dsi()
"""
img contains a nibabel Nifti1Image object (data) and gtab contains a GradientTable
object (gradient information e.g. b-values). For example to read the b-values
it is possible to write print(gtab.bvals).
Load the raw diffusion data and the affine.
"""
data = img.get_data()
print('data.shape (%d, %d, %d, %d)' % data.shape)
"""
data.shape ``(96, 96, 60, 203)``
This dataset has anisotropic voxel sizes, therefore reslicing is necessary.
"""
affine = img.get_affine()
"""
Read the voxel size from the image header.
"""
voxel_size = img.get_header().get_zooms()[:3]
"""
Instantiate the Model and apply it to the data.
"""
gqmodel = GeneralizedQSamplingModel(gtab, sampling_length=3)
"""
The parameter `sampling_length` is used here to
Lets just use one slice only from the data.
"""
dataslice = data[:, :, data.shape[2] / 2]
mask = dataslice[..., 0] > 50
gqfit = gqmodel.fit(dataslice, mask=mask)
"""
Load an odf reconstruction sphere
"""
sphere = get_sphere('symmetric724')
"""
Calculate the ODFs with this specific sphere
"""
ODF = gqfit.odf(sphere)
print('ODF.shape (%d, %d, %d)' % ODF.shape)
"""
ODF.shape ``(96, 96, 724)``
Using peaks_from_model we can find the main peaks of the ODFs and other
properties.
"""
gqpeaks = peaks_from_model(model=gqmodel,
data=dataslice,
sphere=sphere,
relative_peak_threshold=.5,
min_separation_angle=25,
mask=mask,
return_odf=False,
normalize_peaks=True)
gqpeak_values = gqpeaks.peak_values
"""
gqpeak_indices show which sphere points have the maximum values.
"""
gqpeak_indices = gqpeaks.peak_indices
"""
It is also possible to calculate GFA.
"""
GFA = gqpeaks.gfa
print('GFA.shape (%d, %d)' % GFA.shape)
"""
With parameter `return_odf=True` we can obtain the ODF using gqpeaks.ODF
"""
gqpeaks = peaks_from_model(model=gqmodel,
data=dataslice,
sphere=sphere,
relative_peak_threshold=.5,
min_separation_angle=25,
mask=mask,
return_odf=True,
normalize_peaks=True)
"""
This ODF will be of course identical to the ODF calculated above as long as the same
data and mask are used.
"""
np.sum(gqpeaks.odf != ODF) == 0
"""
True
The advantage of using peaks_from_models is that it calculates the ODF only once and
saves it or deletes if it is not necessary to keep.
.. [Yeh2010] Yeh, F-C et al., Generalized Q-sampling imaging, IEEE
Transactions on Medical Imaging, vol 29, no 9, 2010.
.. include:: ../links_names.inc
"""
|