www.gusucode.com > images 案例代码 matlab源码程序 > images/MRISliceExample.m
%% Exploring Slices from a 3-Dimensional MRI Data Set % This example shows how to explore a volume of data by extracting slices % through a three-dimensional MRI data set using |imtransform| and |tformarray| % functions. % Copyright 1993-2013 The MathWorks, Inc. %% Step 1: Load and View Horizontal MRI % This example uses the MRI data set that comes with MATLAB(R) and that is % used in the help examples for both |montage| and |immovie|. Loading % |mri.mat| adds two variables to the workspace: |D| (128-by-128-by-1-by-27, % class uint8) and a grayscale colormap, |map| (89-by-3, class double). % % |D| comprises 27 128-by-128 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 |montage|. 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 image sequence.) The fourth % dimension is temporal (as with any image sequence), 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: % % * Dimension 1: Front to back of head (rostral/anterior to caudal/posterior) % * Dimension 2: Left to right of head % * Dimension 4: Bottom to top of head (inferior to superior). % % 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 montage. load mri; montage(D,map) title('Horizontal Slices'); %% Step 2: Extract Sagittal Slice from Horizontal Slices Using IMTRANSFORM % We can construct a mid-sagittal 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 128-by-1-by-1-by-27. % |reshape| (or |squeeze|) can convert |M1| into a 128-by-27 image that % is 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: % % * Dimension 1: Front to back of head (rostral to caudal) % * Dimension 2: Bottom to top of head (inferior to superior). % % We can obtain a much more satisfying view by transforming |M2| to change % its orientation and increase the sampling along the vertical % (inferior-superior) 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 single-step transformation and more economical use of memory. T0 = maketform('affine',[0 -2.5; 1 0; 0 0]); %% % The upper 2-by-2 block of the matrix passed to maketform, |[0 -2.5;1 0]|, % combines the rotation and scaling. After transformation we have: % % * Dimension 1: Top to bottom of head (superior to inferior). % * Dimension 2: Front to back of head (rostral to caudal) %% % The call % % 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 % nearest-neighbor 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 two-dimensional 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 2-D to 2-D % transformations, but |tformarray| supports N-D to M-D 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: % % * Dimension 1: Superior to inferior (original dimension 4, reversed) % * Dimension 2: Caudal to rostral (original dimension 1) % % 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 2-D affine transformation |T1| 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 |T2| % that extracts the 64th sagittal plane using a very simple |INVERSE_FCN|. T1 = maketform('affine',[-2.5 0; 0 1; 68.5 0]); inverseFcn = @(X,t) [X repmat(t.tdata,[size(X,1) 1])]; T2 = maketform('custom',3,2,[],inverseFcn,64); Tc = maketform('composite',T1,T2); %% % Note that |T2| and |Tc| take a 3-D input to a 2-D input. %% % 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 2-D output % in a single step. Our output image is 66-by-128, with the original 27 planes % expanding to 66 in the vertical (inferior-superior) 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 and Display Sagittal Slices % We create a 4-D array (the third dimension is the color dimension) that can be % used to generate an image sequence that goes from left to right, starts 30 % planes in, skips every other plane, and has 35 frames in total. The % transformed array has: % % * Dimension 1: Top to bottom (superior to inferior) % * Dimension 2: Front to back (rostral to caudal) % * Dimension 4: Left to right. %% % 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 mid-sagittal 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, left-to-right 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 montage (padding the array slightly to separate % the elements of the montage). S2 = padarray(S,[6 0 0 0],0,'both'); figure, montage(S2,map) title('Sagittal Slices'); %% Step 5: Create and Display Coronal Slices % 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: % % * Dimension 1: Top to bottom (superior to inferior) % * Dimension 2: Left to right % * Dimension 4: Back to front (caudal to rostral). 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, % side-to-side, and front-to-back dimensions, respectively. The resampler % remains the same. C = tformarray(D,T4,R3,[4 2 1],[1 2 4],[66 128 45],[],0); %% % Note that all array permutations and flips in steps 3, 4, and 5 were handled % as part of the |tformarray| operation. %% % View the coronal slices as a montage (padding the array slightly to separate % the elements of the montage). C2 = padarray(C,[6 0 0 0],0,'both'); figure, montage(C2,map) title('Coronal Slices'); %%