function cosmo_plot_slices(data, dim, slice_step, slice_start, slice_stop)
% Plots a set of slices from a dataset, nifti image, or 3D data array
%
% cosmo_plot_slices(data[, dim][, slice_step][, slice_start][, slice_stop])
%
% Inputs:
% data an fmri dataset (e.g., from cosmo_fmri_dataset), or a 3D
% array with data. data should contain samples from a single
% volume (sample) only.
% dim dimension according to which slices are plotted
% (default: 3). Values between 1 and 3 allows using a
% saggital, axial or coronal view; but the mapping between
% the 3 numbers and 3 views depends on the particular
% orientation of the data.
% slice_step step between slices (default: 1). If negative then
% -slice_step indicates the number of slices
% slice_start the index of the first slice to plot (default: 1).
% slice_stop the index of the last slice to plot (default: the number of
% slices in the dim-th dimension).
%
% Examples:
% % plot an fMRI dataset struct with default options
% cosmo_plot_slices(ds)
%
% % plot an fMRI dataset struct along the second spatial dimension
% cosmo_plot_slices(ds, 2)
%
% % plot a random gaussian 3D array along the first dimension
% cosmo_plot_slices(randn([40,50,20]),1)
%
% % plot an fMRI dataset struct along the default spatial dimension
% % every 5-th slice
% cosmo_plot_slices(ds, [], 5)
%
% % plot an fMRI dataset struct along the third spatial dimension
% % with about 12 slices
% cosmo_plot_slices(ds, 3, -12)
%
% % plot an fMRI dataset struct along the third spatial dimension
% % with about 12 slices, starting at slice 10 and stopping at slice 25
% cosmo_plot_slices(ds, 3, -12, 10, 25)
%
% Notes:
% - Using this function only really makes sense for fMRI-like data.
% - This function does not provide a consistent orientation for slices,
% as this depends on the voxel-to-world transformation matrix, which is
% completely ignored in this function. Thus left-right and top-down
% swaps can occur. Different datasets may provide different views, for
% example dim=1 may give a saggital view if the dataset comes from one
% program and an axial view if it comes from another program.
%
% # For CoSMoMVPA's copyright information and license terms, #
% # see the COPYING file distributed with CoSMoMVPA. #
if nargin<2 || isempty(dim), dim=3; end
if nargin<3 || isempty(slice_step), slice_step=-20; end
if nargin<4 || isempty(slice_start), slice_start=1; end
if nargin<5 || isempty(slice_stop), slice_stop=[]; end % set later
data3D=get_data3D(data);
% get min and max values across the entire volume
data_lin=data3D(:);
mn=min(data_lin);
mx=max(data_lin);
% shift it so that we can walk over the first dimension
data_sh=shiftdim(data3D, dim-1);
if isempty(slice_stop)
slice_stop=size(data_sh,1);
end
if slice_step<0
nslices=-slice_step;
slice_step=ceil((slice_stop-slice_start+1)/nslices);
end
% determine which slices to show
slice_idxs=slice_start:slice_step:slice_stop;
nslices=numel(slice_idxs);
plot_ratio=.8; % ratio between number of rows and colums
nrows=ceil(sqrt(nslices)*plot_ratio);
ncols=ceil(nslices/nrows);
% use header depending on dim
header_labels={'i','j','k'};
% order of slices and whether the slice should be transposed
xorder=[-1 -1 1];
yorder=[-1 1 -1];
do_transpose=[true false true];
for k=1:nslices
slice_idx=slice_idxs(k);
slice=squeeze(data_sh(slice_idx,:,:));
if xorder(dim)<0
slice=slice(end:-1:1,:);
end
if yorder(dim)<0
slice=slice(:,end:-1:1);
end
if do_transpose(dim)
slice=slice';
end
subplot(nrows, ncols, k);
imagesc(slice, [mn, mx]);
title(sprintf('%s = %d', header_labels{dim}, slice_idx));
end
function data3D=get_data3D(data)
if isstruct(data)
cosmo_check_dataset(data,'fmri');
data4D=cosmo_unflatten(data);
sz=size(data4D);
if sz(1)>1
error(['expected single volume data, but found %d '...
'volumes. To select a single volume, use '...
'cosmo_slice. For example, to show the %d-th '...
'volume from a dataset struct ds, use:\n\n '...
'cosmo_plot_slices(cosmo_slice(ds,%d))\n'],...
sz(1),sz(1),sz(1));
end
if numel(sz)>4
error('data must be 4D');
end
data3D=reshape(data4D, sz(2:end));
elseif isnumeric(data) || islogical(data)
ndim=numel(size(data));
if ndim>3
error('expected 3D image - did you select a single volume?');
end
data3D=data;
else
error('illegal input: expected dataset or 3D array');
end