function test_suite = test_surficial_io_niml_dset
% tests for AFNI NIML input/output
%
% # For CoSMoMVPA's copyright information and license terms, #
% # see the COPYING file distributed with CoSMoMVPA. #
try % assignment of 'localfunctions' is necessary in Matlab >= 2016
test_functions = localfunctions();
catch % no problem; early Matlab versions can use initTestSuite fine
end
initTestSuite;
function test_niml_dset_dataset_io()
if cosmo_skip_test_if_no_external('afni')
return
end
ds = get_niml_dataset();
fn = cosmo_make_temp_filename('_tmp', '.niml.dset');
cleaner = onCleanup(@()delete(fn));
formats = {'ascii', 'binary', 'binary.lsbfirst', 'binary.msbfirst'};
n_formats = numel(formats);
for k = 1:n_formats
format = formats{k};
niml_str = get_niml_dset_string(format);
fid = fopen(fn, 'w');
fwrite(fid, niml_str);
fclose(fid);
ds2 = cosmo_surface_dataset(fn);
assert_dataset_equal(ds, ds2);
cosmo_map2surface(ds, fn, 'encoding', format);
ds3 = cosmo_surface_dataset(fn);
assert_dataset_equal(ds, ds3);
% try mapping using functionality in AFNI Matlab library
s = struct();
s.data = ds.samples';
s.node_indices = ds.a.fdim.values{1}(ds.fa.node_indices) - 1;
s.stats = ds.sa.stats;
s.labels = ds.sa.labels;
afni_niml_writesimple(s, fn, format);
ds4 = cosmo_surface_dataset(fn);
assert_dataset_equal(ds, ds4);
% ensure that mapping to struct works
s2 = cosmo_map2surface(ds, '-niml_dset', 'encoding', format);
assertEqual(s2, s);
end
function assert_dataset_equal(x, y)
assertElementsAlmostEqual(x.samples, y.samples, 'absolute', 1e-4);
assertEqual(sort(fieldnames(x)), sort(fieldnames(y)));
assertEqual(x.fa, y.fa);
assertEqual(x.sa, y.sa);
assertEqual(x.a, y.a);
function ds = get_niml_dataset()
ds = cosmo_synthetic_dataset('ntargets', 2, 'nchunks', 1, ...
'type', 'surface');
ds.a.fdim.values{1} = 1 + [1 4 7 8 5 3];
ds.sa = struct();
ds.sa.stats = {'Ftest(3,4)'; 'Zscore()'};
% TODO: use these complicated labels
% ds.sa.labels={'label1<>"&''';'label2'};
ds.sa.labels = {'label1'; 'label2'};
function string = get_niml_dset_string(format)
string = [ ...
get_niml_dset_header() ...
get_niml_dset_data(format) ...
get_niml_dset_footer()
];
function string = get_niml_dset_header()
string = bytes_printf([ ...
'<AFNI_dataset\n' ...
' dset_type="Node_Bucket"\n' ...
' self_idcode="XYZ_QJCWHYMMSUIUXCMBUZJASSJT"\n' ...
' filename="data.niml.dset"\n' ...
' label="data.niml.dset"\n' ...
' ni_form="ni_group" >\n' ...
]);
function string = get_niml_dset_footer()
string = bytes_printf([ ...
'<AFNI_atr\n' ...
' atr_name="COLMS_RANGE"\n' ...
' ni_type="String"\n' ...
' ni_dimen="1" >\n' ...
'"-3.6849477 2.0316862 1 0;-1.3264908 2.3386572 2 4"'...
'</AFNI_atr>\n' ...
'<AFNI_atr\n' ...
' atr_name="COLMS_LABS"\n' ...
' ni_type="String"\n' ...
' ni_dimen="1" >\n' ...
'"label1;label2"</AFNI_atr>\n' ...
'<AFNI_atr\n' ...
' atr_name="COLMS_TYPE"\n' ...
' ni_type="String"\n' ...
' ni_dimen="1" >\n' ...
'"Generic_Float;Generic_Float"</AFNI_atr>\n' ...
'<AFNI_atr\n' ...
' atr_name="COLMS_STATSYM"\n' ...
' ni_type="String"\n' ...
' ni_dimen="1" >\n' ...
'"Ftest(3,4);Zscore()"</AFNI_atr>\n' ...
'</AFNI_dataset>\n' ...
]);
% TODO: use these complicated labels:
% '"label1<>"&';label2"</AFNI_atr>\n' ...
function bytes = get_niml_dset_data(format)
ds = get_niml_dataset();
data = ds.samples;
nodes_base0 = ds.a.fdim.values{1}(ds.fa.node_indices) - 1;
[nrows, ncols] = size(data);
row_sep = sprintf('\n');
col_sep = ' ';
switch lower(format)
case 'ascii'
data_pat = [repmat(['%.6f' col_sep], 1, nrows) row_sep];
bucket_data = uint8(sprintf(data_pat, data));
node_pat = ['%d' row_sep];
node_data = uint8(sprintf(node_pat, nodes_base0));
ni_form = '';
case {'binary', 'binary.lsbfirst', 'binary.msbfirst'}
data_single = single(data(:)');
node_int = int32(nodes_base0(:)');
[unused, unused, endian] = computer();
computer_format = sprintf('binary.%ssbfirst', lower(endian));
do_swap = ~strcmp(format, 'binary') && ...
~strcmp(format, computer_format);
if do_swap
if cosmo_wtf('is_octave')
% Octave 3.8 does not support byteswap for single
data_single_int = typecast(data_single, 'int32');
data_single_int_swap = swapbytes(data_single_int);
data_single = typecast(data_single_int_swap, 'single');
else
data_single = swapbytes(data_single);
end
node_int = swapbytes(node_int);
end
if strcmp(format, 'binary')
ni_form_format = computer_format;
else
ni_form_format = format;
end
bucket_data = typecast(data_single, 'uint8');
node_data = typecast(node_int, 'uint8');
ni_form = sprintf(' ni_form="%s"\n', ni_form_format);
otherwise
error('unsupported format %s', format);
end
bytes = [ ...
bytes_printf([ ...
'<SPARSE_DATA\n' ...
' data_type="Node_Bucket_data"\n' ...
' ni_type="2*float"\n' ...
ni_form ...
' ni_dimen="6" >'])...
bucket_data ...
bytes_printf([ ...
'</SPARSE_DATA>\n' ...
'<INDEX_LIST\n' ...
' data_type="Node_Bucket_node_indices"\n' ...
' sorted_node_def="No"\n' ...
' COLMS_RANGE="1 8 0 3"\n' ...
' COLMS_LABS="Node Indices"\n' ...
' COLMS_TYPE="Node_Index"\n' ...
' ni_type="int"\n' ...
ni_form ...
' ni_dimen="6" >']) ...
node_data ...
bytes_printf([ ...
'</INDEX_LIST>\n']) ...
];
function b = bytes_printf(varargin)
b = uint8(sprintf(varargin{:}));
