# - compatibility with Python 3
from __future__ import print_function # print('me') instead of print 'me'
from __future__ import division # 1/2 == 0.5, not 0
# - show figures inside the notebook
%matplotlib inline
# - import common modules
import numpy as np # the Python array package
import matplotlib.pyplot as plt # the Python plotting package
# - set gray colormap and nearest neighbor interpolation by default
plt.rcParams['image.cmap'] = 'gray'
plt.rcParams['image.interpolation'] = 'nearest'
First we need the number of 'slices' (length of the third image dimension):
# - load image to get parameters
import nibabel as nib
img = nib.load('fds114_sub009_t2r1.nii')
num_slices = img.shape[2]
num_slices
Next we need the TR from the image header
# - get the TR from this image
TR = img.header.get_zooms()[-1]
TR
We need the time to acquire all but the last slice (don't ask!). SPM calls this TA (time of acqusition?):
# - calculate TA
time_for_one_slice = TR / num_slices
TA = TR - time_for_one_slice
TA
Next we need the acquisition order. This is a list acq_order of length num_slices where acq_order[i] gives the slice index of the slice acquired at position i in time.
This is for MATLAB, so the array indices start at 1 rather than 0 as for Python.
So, if acq_order is [1, 3, 5, 7, 9, 2, 4, 6, 8, 10] that means that the first slice acquired was (1-based) slice index 1, followed by slice index 3, followed by slice index 5, and so on.
For our image, the scanner collected the odd index slices first, starting at the bottom, and then came back and collected the even index slices:
# - generate acq_order list
odd = range(1, num_slices+1, 2)
even = range(2, num_slices+1, 2)
acq_order = list(odd) + list(even)
print(acq_order)
Lastly, SPM slice timing also asks for a reference slice. This is the (1-based) index of the slice to which the others will be matched in terms of time. For example, if we want to use the interpolation to match the times of all the slices to the acquisition time of the first slice, the reference slice would be 1 (bottom slice, with 1-based indexing).
Now we can use the SPM interface to run the slice timing.
Now we are going to use the nipype interface to run the same thing.
First, let's delete the image that was the result of slice timing via the SPM GUI:
# - import the routines for working with the operating system
import os
# Delete file if it exists
if os.path.exists('afds114_sub009_t2r1.nii'):
os.unlink('afds114_sub009_t2r1.nii') # delete file
# - import our initilization of nipype / MATLAB interface
import nipype_settings
# - import SliceTiming from nipype SPM interfaces
from nipype.interfaces.spm import SliceTiming
Your job now is to create the SliceTiming SPM batch job from the parameters you need.
Try looking for help on the SliceTiming object for clues how to do this.
# Initialize the SliceTiming batch object, and fill parameters
# Run the slice timing on the fixed image `fds114_sub009_t2r1.nii`
# Run the batch job
Let's check we did write a new slice-time corrected file:
# - check for file created by nipype slice timing
assert os.path.exists('afds114_sub009_t2r1.nii')