利用matlab编写实现显示fmri切片slice图像 混合显示 不同侧面显示 可叠加t检验图显示 by DR. Rajeev Raizada

sagittal_slice_number = 20;

sagittal_slice_vals = subj_3danat(sagittal_slice_number,:,:);

%%%% Let's look at the variables in our workspace again

whos

%%%% Note that even though sagittal_slice_vals is just one slice,
% sagittal_slice_vals 1x63x46 23184 double array
%
% It turns out that the Matlab image command won't accept 3D matrices,
% so we need to get rid of this redundant one-voxel-wide dimension.
% Luckily, Matlab has a function that does precisely this, called squeeze.

sagittal_slice_2D = squeeze(sagittal_slice_vals);

%%% Get rid of the redundant 3rd dimension

%%%% Let's look at all the variables in our workspace
%%%% that start with "sagit"

whos sagit*

%%%% Squeeze worked! It made our sagittal slice 2D instead of 3D
%
% sagittal_slice_2D 63x46 23184 double array

%%%% Now we can plot it

figure(2); 

clf; 

imagesc(sagittal_slice_2D);

colormap gray; %%% Make it a grayscale plot

colorbar; %%% Show how the image intensities map onto the colormap
axis('image'); %%% Make the proportions correct
title('An initial attempt at showing the sagittal slice');

caxis([ 0 250 ]); % Keep the black value at 0, but set the white value to 250

%%%% That makes the image look a lot better!
%%%% Now let's redraw the colorbar to see what the new colour-scaling is.

colorbar;

rotated_sagittal_slice = rot90(sagittal_slice_2D);

figure(3); %%% Make a new figure window, Fig.3
clf; %%% Clear the figure

imagesc(rotated_sagittal_slice); 
%%% Show the image, with the brightness scaled 
%%% so that the darkest is black and lightest is white

colormap gray; %%% Make it a grayscale plot
caxis([ 0 250 ]); % Keep the black value at 0, but set the white value to 250 
colorbar; %%% Show how the image intensities map onto the colormap
axis('image'); %%% Make the proportions correct

title('Sagittal slice, rotated so that it is the right way up');

%%% For showing images, we often don't want to see the numbers
%%% on the x-axis and y-axis. We can turn these off using this command:

axis('off'); %%% Turn off the numbers running along the axes

figure(4); %%% Make a new figure, Fig.4
clf; %%% Clear the figure

for loop_counter = 1:9, %%% Go around 9 times, adding one to the 
%%% value of loop_counter each time.

coronal_slice_number = 7*loop_counter;
%%% Loop counter goes 1,2,3,...,9
%%% so this gives slices 7,14,21,...,63

coronal_slice_vals = subj_3danat(:,coronal_slice_number,:);
%%% The coronal slice is the 2nd dimension.
%%% So, read out that slice from the 3D subj_3danat matrix
%%% The ":" signs in the 1st and 3rd coordinate positions
%%% mean "give me all the x-vals and all the z-vals in that slice"

coronal_slice_2D = squeeze(coronal_slice_vals);
%%% We have to use squeeze to take out the redundant
%%% one-voxel-wide 2nd dimension, since imagesc will 
%%% not let us display a 3D matrix

rotated_coronal_slice = rot90(coronal_slice_2D);
%%% Roate the slice by 90 degrees, to make it the right way up

subplot(3,3,loop_counter);
%%% Make a 3x3 array of subplots, 
%%% and draw into the "loop_counter"-th one

imagesc(rotated_coronal_slice); %%% Make the image
colormap gray; %%% Make it grayscale

caxis([ 0 250 ]); 
%%% So that very intense voxels don't make everything
%%% else look dark, set the black-value to 0, 
%%% and set the white value to 250, like we did above.

axis('image'); %%% Make the proportions of the image correct

axis('off'); %%% Turn off the numbers on the x- and y-axes

end; %%% The end of this for-loop
%%% If loop_counter is less than 9, then add 1 to it,
%%% and then go back to the beginning. 
%%% If loop_counter is 9, then stop --- the looping is finished.

axial_slice_number = 18;

%%%% The matrix containing a brain-full of voxels 
%%%% with t-statistic numbers in them is called speech_Tmap

axial_slice_Tmap_vals = speech_Tmap(:,:,axial_slice_number);

%%%% Squeeze it, to get rid of the nuisance one-voxel wide 3rd dimension
%%%% in case there is one. (There might not always be, depending on
%%%% what type of slice you are taking. But it's a good idea to squeeze
%%%% all the time anyway, because it will make your programs more robust).

axial_slice_Tmap_2D = squeeze(axial_slice_Tmap_vals);

%%%% Now let's show this matrix as an image, in the way we did above.

figure(5);
clf;
imagesc(axial_slice_Tmap_2D);
axis('image'); %%% Make the proportions of the image correct
title('Raw unthresholded t-map');

%%%% For functional maps, a nice colormap to use is called "hot"
%%%% This makes the lowest numbers black, and then the colours
%%%% get "warmer" as the numbers increase: red, orange, yellow, white
%%%% To see a number of useful colormaps, type 
%%%% help graph3d 
%%%% into the Matlab command window.

colormap hot;

%%% Let's make a colorbar to show us which t-map values get shown
%%% as which colours.
colorbar;

Tmap_threshold = 3.5;

%%%% In Matlab, the way to figure out which voxels have values 
%%%% above this threshold is to make a matrix that has a 1 in every
%%%% above-threshold voxel, and a zero everywhere else.
%%%%
%%%% In other words, we want a 1 in every voxel where the statement
%%%% "this voxel's intensity is greater than the threshold" is true.
%%%%
%%%% Here's how we do that in Matlab:

True_or_false_map = ( axial_slice_Tmap_2D > Tmap_threshold );

%%%% The matrix will have a 1 in it in every voxel where the
%%%% bracketed expression ( axial_slice_Tmap_2D > Tmap_threshold )
%%%% is true, i.e. in all the above-threshold voxels, 
%%%% and a zero in all the voxels where it's false.
%%%% Making true-or-false matrices like this turns out
%%%% to be useful in all kinds of circumstances.

%%%% Let's have a look at our True_or_false_map:

figure(6);
clf;
imagesc(True_or_false_map);
axis('image');
colormap hot;
colorbar;
title('1 if t-map is above-threshold, 0 if it is below-threshold');

Above_threshold_Tmap = True_or_false_map .* axial_slice_Tmap_2D;

%%%% Let's look at our above-treshold T-map

figure(7);
clf;
imagesc(Above_threshold_Tmap);
axis('image');
colormap hot;
colorbar;
title('The thresholded t-map');

figure(8);
clf;

subplot(3,1,1); %%% Three rows of subplots, one column, draw in the first one
imagesc(axial_slice_Tmap_2D);
axis('image');
axis('off');
colormap hot;
colorbar;
title('Raw unthresholded t-map');

subplot(3,1,2); %%% The last argument is 2, meaning
%%% "Draw in the 2nd of the three subplots"
imagesc(True_or_false_map);
axis('image');
axis('off');
colormap hot;
colorbar;
title('1 if above-threshold, 0 if below');

subplot(3,1,3); %%% The last argument is 3, meaning
%%% "Draw in the 3rd of the three subplots"
imagesc(Above_threshold_Tmap);
axis('image');
axis('off');
colormap hot;
colorbar;
title('The thresholded t-map');

max_Tmap_value = max(max(axial_slice_Tmap_2D)); % Single biggest t-map score

%%%% Recall that the argument that we give to caxis is a two-element
%%%% vector, containing the maximum and minimum values we want 
%%%% for the colorbar.

colorbar_min_and_max_vals = [ 0 max_Tmap_value ];

%%%% Now we can plot our thresholded and unthresholded t-maps 
%%%% side-by-side, using exactly the same colour scaling in both.

figure(9);
clf;

subplot(2,1,1); %%% Two rows of subplots, one column, draw in the first one
imagesc(axial_slice_Tmap_2D);
axis('image');
axis('off');
colormap hot;

%%% Now set the colour scaling to what we want. 
%%% Do this before making the colorbar
caxis(colorbar_min_and_max_vals);

colorbar;
title('t-map showing all above-zero values');

subplot(2,1,2); %%% The last argument is 2, meaning
%%% "Draw in the 2nd of the two subplots"
imagesc(Above_threshold_Tmap);
axis('image');
axis('off');
colormap hot;
caxis(colorbar_min_and_max_vals); 
%%% This subplot will now have the same colour scaling as the first one
colorbar;
title('The thresholded t-map, with same colour-scaling as the plot above');

figure(10); 
clf; 

set(gcf,'doublebuffer','on'); 

for sagittal_slice_number = 1:53, 

%%% There are 53 sagittal slices, and this loop
%%% goes around 53 times, adding one to the 
%%% value of sagittal_slice_number each time.

sagittal_slice_vals = subj_3danat(sagittal_slice_number,:,:);

sagittal_slice_2D = squeeze(sagittal_slice_vals);

rotated_sagittal_slice = rot90(sagittal_slice_2D);


imagesc(rotated_sagittal_slice); 
colormap gray; 

caxis([ 0 250 ]); 

axis('image');

xlabel(['Sagittal slice number ' num2str(sagittal_slice_number) ]);

title('Press any key to show the next slice');

pause; 
%%% After the key press, the program continues,

end; %%% The end of this for-loop