1. pna2015
  2. day10
In [1]:
# - 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
In [2]:
# - show figures inside the notebook
%matplotlib inline
In [3]:
# - import common modules
import numpy as np  # the Python array package
import matplotlib.pyplot as plt  # the Python plotting package
In [4]:
# - set gray colormap and nearest neighbor interpolation by default
plt.rcParams['image.cmap'] = 'gray'
plt.rcParams['image.interpolation'] = 'nearest'
In [5]:
# - import numpy.linalg with a shorter name
import numpy.linalg as npl
In [6]:
# - make numpy print out values in less verbose form
np.set_printoptions(suppress=True, precision=4)

What transforms did I add?

I took this image

  • ds114_sub009_highres.nii

and I made a new copy, in img, hdr format:

  • ds114_sub009_highres_moved.img

I then modified the affine in ds114_sub009_highres_moved.hdr to add an extra rotation and translation.

What extra translation and rotation did I add?

Hint: remember the rotation formulae:

$$ \begin{alignat}{1} R_x(\theta) &= \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \theta & -\sin \theta \\[3pt] 0 & \sin \theta & \cos \theta \\[3pt] \end{bmatrix} \\[6pt] R_y(\theta) &= \begin{bmatrix} \cos \theta & 0 & \sin \theta \\[3pt] 0 & 1 & 0 \\[3pt] -\sin \theta & 0 & \cos \theta \\ \end{bmatrix} \\[6pt] R_z(\theta) &= \begin{bmatrix} \cos \theta & -\sin \theta & 0 \\[3pt] \sin \theta & \cos \theta & 0\\[3pt] 0 & 0 & 1\\ \end{bmatrix} \end{alignat} $$

You might also want to compare against results from our rotations code:

In [7]:
# - get our functions to make rotation matrices
# If this gives an error - maybe you didn't set your PYTHONPATH? Let us know!
from rotations import x_rotmat, y_rotmat, z_rotmat

Final hint:

  • inverse sine, cosine are np.arcsin, np.arccos.
In [8]:
# Get the original image affine
import nibabel as nib
orig_img = nib.load('ds114_sub009_highres.nii')
print(orig_img.affine)

[[   0.9989   -0.0605    0.0109 -129.8257]
 [   0.0427    1.263     0.2336 -119.0906]
 [  -0.0215   -0.3028    0.9723 -143.4178]
 [   0.        0.        0.        1.    ]]
In [9]:
# Get the new image affine
moved_img = nib.load('ds114_sub009_highres_moved.img')
moved_img.affine
Out[9]:
array([[   0.9416,   -0.4311,   -0.0586,  -98.8336],
       [   0.336 ,    1.1887,    0.2264, -164.1377],
       [  -0.0215,   -0.3028,    0.9723, -158.4178],
       [   0.    ,    0.    ,    0.    ,    1.    ]])
In [10]:
# Work out what transform has been added in the new affine
extra = moved_img.affine.dot(npl.inv(orig_img.affine))
extra
Out[10]:
array([[  0.9553,  -0.2955,   0.    , -10.    ],
       [  0.2955,   0.9553,   0.    , -12.    ],
       [  0.    ,  -0.    ,   1.    , -15.    ],
       [  0.    ,   0.    ,   0.    ,   1.    ]])
In [11]:
# What rotation and translation has been applied?
z_angle = np.arccos(extra[0, 0])
translation = extra[:3, 3]
z_angle, translation
Out[11]:
(0.30000001985761038, array([-10., -12., -15.]))