Merge pull request #591 from arokem/reuse_sh_apm

Reuse the CSD sh coefficients if you already have them.

AFQ/api.py

@@ -15,6 +15,7 @@
 from dipy.io.stateful_tractogram import StatefulTractogram, Space
 from dipy.io.gradients import read_bvals_bvecs
 from dipy.stats.analysis import afq_profile, gaussian_weights
+from dipy.reconst import shm
 from dipy.reconst.dki_micro import axonal_water_fraction
 
 from bids.layout import BIDSLayout
@@ -25,7 +26,6 @@
 from AFQ.models.dti import noise_from_b0
 from AFQ.models.dki import _fit as dki_fit
 from AFQ.models.csd import _fit as csd_fit
-from AFQ.models.csd import fit_anisotropic_power_map
 import AFQ.tractography as aft
 import dipy.reconst.dti as dpy_dti
 import dipy.reconst.dki as dpy_dki
@@ -786,9 +786,8 @@ def _anisotropic_power_map(self, row):
             row, '_model-CSD_APM.nii.gz')
         if not op.exists(pmap_file):
             dwi_data, gtab, img = self._get_data_gtab(row)
-            mask = self._brain_mask(row)
-            pmap = fit_anisotropic_power_map(
-                dwi_data, gtab, mask)
+            sh_coeff = nib.load(self._csd(row)).get_fdata()
+            pmap = shm.anisotropic_power(sh_coeff)
             pmap = nib.Nifti1Image(pmap, img.affine)
             self.log_and_save_nii(pmap, pmap_file)
             meta_fname = self._get_fname(row, '_model-CSD_APM.json')

AFQ/models/csd.py

@@ -127,49 +127,3 @@ def fit_csd(data_files, bval_files, bvec_files, mask=None, response=None,
     fname = op.join(out_dir, 'csd_sh_coeff.nii.gz')
     nib.save(nib.Nifti1Image(csdfit.shm_coeff, aff), fname)
     return fname
-
-
-def fit_anisotropic_power_map(dwi, gtab, mask=None):
-    """
-    Fits an anisotropic power map.
-
-    Parameters
-    ----------
-    dwi : str, ndarray, or nifti1image
-        Data to greate map with.
-
-    gtab : GradientTable
-        A GradientTable with all the gradient information.
-
-    mask : str or nifti1image, optional
-        mask to mask the data with.
-        Default: None.
-
-    Returns
-    -------
-    ndarray containing an anisotropic power map.
-    """
-
-    if isinstance(dwi, str):
-        dwi = nib.load(dwi)
-    if isinstance(dwi, nib.Nifti1Image):
-        dwi_data = dwi.get_fdata()
-    else:
-        dwi_data = dwi
-
-    if isinstance(mask, str):
-        mask = nib.load(mask)
-    mask = mask.get_fdata()
-
-    model = _model(gtab, dwi_data)
-    sphere = dpd.get_sphere('symmetric724')
-    peaks = csd.peaks_from_model(
-        model=model,
-        data=dwi_data,
-        sphere=sphere,
-        relative_peak_threshold=.5,
-        min_separation_angle=25,
-        mask=mask)
-    ap = shm.anisotropic_power(peaks.shm_coeff)
-
-    return ap