Extracting Slices from a 3Dimensional MRI Data Set
The imtransform
and tformarray
functions can be used to interpolate and reslice a threedimensional MRI data set, providing a convenient way to view a volume of data..
Key concepts: 
Affine, custom, and composite transformations, 3D and 3D/2D transformations, input array permutations, resampling in selected dimensions 
Key functions 

Overview of Example
The example includes these steps:
Step 1: Load and View Horizontal MRI
This example uses the MRI data set that comes with MATLAB and that is used in the help examples for both montage and immovie. Loading MRI adds two variables to the workspace: D
(128by128by1by27, class uint8
) and a grayscale colormap, map
(89by3, class double).
D
comprises 27 128by128 horizontal slices from an MRI data scan of a human cranium. Values in D
range from 0 through 88, so the colormap is needed to generate a figure with a useful visual range. The dimensionality of D
makes it compatible with immovie
. The first two dimensions are spatial. The third dimension is the color dimension, with size 1 because it indexes into the color map. (size(D,3)
would be 3 for an RGB movie.) The fourth dimension is temporal (as with any movie), but in this particular case it is also spatial. So there are three spatial dimensions in D
and we can use imtransform
or tformarray
to convert the horizontal slices to sagittal slices (showing the view from the side of the head) or coronal (frontal) slices (showing the view from the front or back of the head).
The spatial dimensions of D
are ordered as follows:
An important factor is that the sampling intervals are not the same along the three dimensions: samples along the vertical dimension (4) are spaced 2.5 times more widely than along the horizontal dimensions.
Load the MRI data set and view the 27 horizontal slices as a movie, then as a montage. To avoid TRUESIZE
warnings when creating montages, the example uses iptgetpref
to turn off the warning.
truesizewarning = iptgetpref('TruesizeWarning'); iptsetpref('TruesizeWarning','off'); % Turn warning off load mri; figure; immovie(D,map); montage(D,map); title('Horizontal Slices');
Step 2: Extract Sagittal Slice from Horizontal Slices Using imtransform
We can construct a midsagittal slice from the MRI data by taking a subset of D
and transforming it to account for the different sampling intervals and the spatial orientation of the dimensions of D
.
The following statement extracts all the data needed for a midsagittal slice.
M1 = D(:,64,:,:); size(M1)
However we cannot view M1
as an image because it is 128by1by1by27. reshape
(or squeeze
) can convert M1
into a 128by27 image, viewable with imshow
.
M2 = reshape(M1,[128 27]); size(M2) figure, imshow(M2,map); title('Sagittal  Raw Data');
The dimensions in M2
are ordered as follows:
We can obtain a much more satisfying view by transforming M2
to change its orientation and increase the sampling along the vertical (inferiorsuperior) dimension by a factor of 2.5  making the sampling interval equal in all three spatial dimensions. We could do this in steps starting with a transpose, but the following affine transformation enables a singlestep transformation and more economical use of memory.
T0 = maketform('affine',[0 2.5; 1 0; 0 0]);
The upper 2by2 block of the matrix passed to maketform
[ 0 2.5 1 0 ]
combines the rotation and scaling. After transformation we have:
imtransform(M2,T0,'cubic')
would suffice to apply T
to M2
and provide good resolution while interpolating along the top to bottom direction. However, there is no need for cubic interpolation in the front to back direction, since no resampling will occur along (output) dimension 2. Therefore we specify nearestneighbor resampling in this dimension, with greater efficiency and identical results.
R2 = makeresampler({'cubic','nearest'},'fill'); M3 = imtransform(M2,T0,R2); figure, imshow(M3,map); title('Sagittal  IMTRANSFORM')
Step 3: Extract Sagittal Slice from the Horizontal Slices using tformarray
In this step we obtain the same result as step 2, but use tformarray
to go from three spatial dimensions to two in a single operation. Step 2 does start with an array having three spatial dimensions and end with an array having two spatial dimensions, but intermediate twodimensional images (M1
and M2
) pave the way for the call to imtransform
that creates M3
. These intermediate images are not necessary if we use tformarray
instead of imtransform
. imtransform
is very convenient for 2D to 2D transformations, but tformarray
supports ND to MD transformations, where M need not equal N.
Through its TDIMS_A
argument, tformarray
allows us to define a permutation for the input array. Since we want to create an image with:
and extract just a single sagittal plane via the original dimension 2, we specify tdims_a = [4 1 2]
. We create a tform
via composition starting with a 2D affine transformation that scales the (new) dimension 1 by a factor of 2.5 and adds a shift of 68.5 to keep the array coordinates positive. The second part of the composite is a custom transformation that extracts the 64th sagittal plane using a very simple INVERSE_FCN
:
Inverse Function: ipex003.m 
function U = ipex003( X, t ) 
Both T2
and Tc
below take a 3D input to a 2D output.
T1 = maketform('affine',[2.5 0; 0 1; 68.5 0]); T2 = maketform('custom',3,2,[],@ipex003,64); Tc = maketform('composite',T1,T2);
We use the same approach to resampling as before, but include a third dimension.
R3 = makeresampler({'cubic','nearest','nearest'},'fill');
tformarray
transforms the three spatial dimensions of D
to a 2D output in a single step. Our output image is 66by128, with the original 27 planes expanding to 66 in the vertical (inferiorsuperior) direction.
M4 = tformarray(D,Tc,R3,[4 1 2],[1 2],[66 128],[],0);
The result is identical to the previous output of imtransform
.
figure, imshow(M4,map); title('Sagittal  TFORMARRAY');
Step 4: Create Sagittal Slices and Show Them as a Movie
We create a 4D array (the third dimension is the color dimension) that can be used to generate a movie that goes from left to right, starts 30 planes in, skips every other plane, and has 35 frames in total. The transformed array has:
As in the previous step, we permute the input array using TDIMS_A = [4 1 2]
, again flipping and rescaling/resampling the vertical dimension. Our affine transformation is the same as the T1
above, except that we add a third dimension with a (3,3) element of 0.5 and (4,3) element of 14 chosen to map 30, 32, ... 98 to 1, 2, ..., 35. This centers our 35 frames on the midsagittal slice.
T3 = maketform('affine',[2.5 0 0; 0 1 0; 0 0 0.5; 68.5 0 14]);
In our call to tformarray
, TSIZE_B = [66 128 35]
now includes the 35 frames in the 4th, lefttoright dimension (which is the third transform dimension. The resampler remains the same.
S = tformarray(D,T3,R3,[4 1 2],[1 2 4],[66 128 35],[],0);
View the sagittal slices as a movie, then a montage (padding the array slightly to separate the elements of the montage).
figure; immovie(S,map); S2 = padarray(S,[6 0 0 0],0,'both'); montage(S2,map); title('Sagittal Slices');
Step 5: Create Coronal Slices and Show Them as a Movie
Constructing coronal slices is almost the same as constructing sagittal slices. We change TDIMS_A
from [4 1 2]
to [4 2 1]
. We create a series of 45 frames, starting 8 planes in and moving from back to front, skipping every other frame. The dimensions of the output array are ordered as follows:
T4 = maketform('affine',[2.5 0 0; 0 1 0; 0 0 0.5; 68.5 0 61]);
In our call to tformarray
, TSIZE_B = [66 128 48]
specifies the vertical, sidetoside, and fronttoback dimensions, respectively. The resampler remains the same.
C = tformarray(D,T4,R3,[4 2 1],[1 2 4],[66 128 45],[],0);
View the coronal slices as a movie, then a montage (padding the array slightly to separate the elements of the montage).
figure; immovie(C,map); C2 = padarray(C,[6 0 0 0],0,'both'); montage(C2,map); title('Coronal Slices');
Note that all array permutations and flips in steps 3, 4, and 5 were handled as part of the tformarray
operation.
% Restore preference for TruesizeWarning iptsetpref('TruesizeWarning',truesizewarning);