Reproducing the FA map
FouierL opened this issue · 0 comments
Hello, thank you for this amazing research and open source codes.
Sorry to take up your time, I have a question not related to this open source code but related to your DDM2 work.
When I was reproducing the FA map in the article, I realized that using only the noisy data didn't come out with the matching figure in the paper.
my code is:
import numpy as np
import dipy.reconst.dki as dki
import dipy.reconst.dti as dti
from dipy.core.gradients import gradient_table
from dipy.data import get_fnames
from dipy.io.gradients import read_bvals_bvecs
from dipy.io.image import load_nifti
from dipy.segment.mask import median_otsu
from dipy.viz.plotting import compare_maps
from scipy.ndimage import gaussian_filter
import matplotlib.pyplot as plt
from warnings import warn
import cv2
import numpy as np
hardi_fname, hardi_bval_fname, hardi_bvec_fname = get_fnames(
'stanford_hardi')
b0_size = 10
data, affine = load_nifti(hardi_fname)
bvals, bvecs = read_bvals_bvecs(hardi_bval_fname, hardi_bvec_fname)
gtab = gradient_table(bvals, bvecs)
maskdata, mask = median_otsu(data, vol_idx=[0, 1], median_radius=4, numpass=2,
autocrop=False, dilate=1)
slice=40
data = data[:, :, slice:slice+1]
mask = mask[:, :, slice:slice+1]
tenmodel = dti.TensorModel(gtab)
tenfit = tenmodel.fit(data, mask=mask)
fits = [tenfit]
maps = ['fa', 'md', 'ad', 'rd']
fit_labels = ['DTI']
for i in range(len(maps)):
attr = getattr(fits[0], maps[i])
attr=(attr - np.min(attr)) / (np.max(attr) - np.min(attr))
plt.imshow(attr,cmap='jet')
plt.axis('off')
plt.savefig('{}.png'.format('origin'+maps[i]))
plt.close()
the map in paper:
the map my code give:
Thank you for taking the time. Good luck with your research.