function ds_splits=cosmo_split(ds, split_by, dim, check)
% splits a dataset by unique values in (a) sample or feature attribute(s).
%
% cosmo_split(ds, split_by[, dim])
%
% Inputs:
% ds dataset struct
% split_by fieldname for the sample (if dim==1) or feature (if dim==2)
% attribute by which the dataset is split.
% It can also be a cell with fieldnames, in which case the
% dataset is split by the set of combinations from these
% fieldnames.
% dim 1 (split by samples; default) or 2 (split by features).
% check boolean indicating whether the input ds is checked as being
% a proper dataset (default: true)
%
% Returns:
% ds_splits 1xP cell, if there are P unique values for the (set of)
% attribute(s) indicated by split_by and dim.
%
% Examples:
% ds=cosmo_synthetic_dataset();
% %
% % split by targets
% splits=cosmo_split(ds,'targets');
% cosmo_disp(splits{2}.sa);
% %|| .targets
% %|| [ 2
% %|| 2
% %|| 2 ]
% %|| .chunks
% %|| [ 1
% %|| 2
% %|| 3 ]
% %
% % split by chunks
% splits=cosmo_split(ds,'chunks');
% cosmo_disp(splits{3}.sa);
% %|| .targets
% %|| [ 1
% %|| 2 ]
% %|| .chunks
% %|| [ 3
% %|| 3 ]
% %
% % split by chunks and targets
% splits=cosmo_split(ds,{'chunks','targets'});
% cosmo_disp(splits{5}.sa);
% %|| .targets
% %|| [ 1 ]
% %|| .chunks
% %|| [ 3 ]
%
% % take an MEEG time-freq dataset, and split by time and channel
% ds=cosmo_synthetic_dataset('type','timefreq','size','big');
% %
% % dataset has 11 channels, 7 frequencies and 5 time points
% cosmo_disp(ds.fa)
% %|| .chan
% %|| [ 1 2 3 ... 304 305 306 ]@1x10710
% %|| .freq
% %|| [ 1 1 1 ... 7 7 7 ]@1x10710
% %|| .time
% %|| [ 1 1 1 ... 5 5 5 ]@1x10710
% %
% % split by time and frequency. Since splitting is done on the feature
% % dimension, the third argument (with value 2) is mandatory
% splits=cosmo_split(ds,{'time','freq'},2);
% % there are 7 * 5 = 35 splits, each with 11 features
% numel(splits)
% %|| 35
% cosmo_disp(cellfun(@(x) size(x.samples,2),splits))
% %|| [ 306 306 306 ... 306 306 306 ]@1x35
% cosmo_disp(splits{18}.fa)
% %|| .chan
% %|| [ 1 2 3 ... 304 305 306 ]@1x306
% %|| .freq
% %|| [ 4 4 4 ... 4 4 4 ]@1x306
% %|| .time
% %|| [ 3 3 3 ... 3 3 3 ]@1x306
% %
% % using cosmo_stack brings the split elements together again
% humpty_dumpty=cosmo_stack(splits,2);
% cosmo_disp(humpty_dumpty.fa)
% %|| .chan
% %|| [ 1 2 3 ... 304 305 306 ]@1x10710
% %|| .freq
% %|| [ 1 1 1 ... 7 7 7 ]@1x10710
% %|| .time
% %|| [ 1 1 1 ... 5 5 5 ]@1x10710
%
% Note:
% - This function is like the inverse of cosmo_stack; if
%
% >> ds_splits=cosmo_split(ds, split_by, dim),
%
% produces output (i.e., does not throw an error), then using
%
% >> ds_humpty_dumpty=cosmo_stack(ds_splits,dim)
%
% means that ds and ds_humpty_dumpty contain the same data, except that
% the order of the data (in the rows [columns] of .samples, or
% .sa [.fa]) may be different if dim==1 [dim==2].
%
%
% See also: cosmo_stack, cosmo_slice
%
% # For CoSMoMVPA's copyright information and license terms, #
% # see the COPYING file distributed with CoSMoMVPA. #