cosmo fmri reorient skl

function ds = cosmo_fmri_reorient(ds, new_orient)
    % Change the orientation of an fmri dataset
    %
    % ds_reorient=cosmo_fmri_reorient(ds, new_orient)
    %
    % Inputs
    %     ds                fmri-dataset
    %     new_orient        new orientation for the dataset (see below for a
    %                       full list)
    %
    % Example:
    %     ds=cosmo_synthetic_dataset();
    %     cosmo_disp(ds.a);
    %     %|| .fdim
    %     %||   .labels
    %     %||     { 'i'  'j'  'k' }
    %     %||   .values
    %     %||     { [ 1         2         3 ]  [ 1         2 ]  [ 1 ] }
    %     %|| .vol
    %     %||   .mat
    %     %||       [ 2         0         0        -3
    %     %||         0         2         0        -3
    %     %||         0         0         2        -3
    %     %||         0         0         0         1 ]
    %     %||   .dim
    %     %||     [ 3         2         1 ]
    %     %||     .xform
    %     %||       'scanner_anat'
    %     cosmo_disp(ds.fa);
    %     %|| .i
    %     %||   [ 1         2         3         1         2         3 ]
    %     %|| .j
    %     %||   [ 1         1         1         2         2         2 ]
    %     %|| .k
    %     %||   [ 1         1         1         1         1         1 ]
    %     ds_reorient=cosmo_fmri_reorient(ds,'AIR');
    %     cosmo_disp(ds_reorient.a);
    %     %|| .fdim
    %     %||   .labels
    %     %||     { 'i'  'j'  'k' }
    %     %||   .values
    %     %||     { [ 1         2 ]  [ 1 ]  [ 1         2         3 ] }
    %     %|| .vol
    %     %||   .mat
    %     %||     [   0         0        -2         5
    %     %||        -2         0         0         3
    %     %||         0         2         0        -3
    %     %||         0         0         0         1 ]
    %     %||   .dim
    %     %||     [ 2         1         3 ]
    %     %||     .xform
    %     %||       'scanner_anat'
    %     cosmo_disp(ds_reorient.fa);
    %     %|| .i
    %     %||   [ 1         2         1         2         1         2 ]
    %     %|| .j
    %     %||   [ 1         1         1         1         1         1 ]
    %     %|| .k
    %     %||   [ 1         1         2         2         3         3 ]
    %
    %     % Many orientations are invalid, for example
    %     ds=cosmo_synthetic_dataset();
    %     cosmo_reorient(ds,'ALR');
    %     error('illegal orientation');
    %
    % Notes:
    %   - there are 3!*3^2 valid orientations, these are:
    %         'SAR'  'SAL'  'SPR'  'SPL'  'IAR'  'IAL'  'IPR'  'IPL'
    %         'SRA'  'SLA'  'SRP'  'SLP'  'IRA'  'ILA'  'IRP'  'ILP'
    %         'ASR'  'ASL'  'PSR'  'PSL'  'AIR'  'AIL'  'PIR'  'PIL'
    %         'ARS'  'ALS'  'PRS'  'PLS'  'ARI'  'ALI'  'PRI'  'PLI'
    %         'RAS'  'LAS'  'RPS'  'LPS'  'RAI'  'LAI'  'RPI'  'LPI'
    %         'RSA'  'LSA'  'RSP'  'LSP'  'RIA'  'LIA'  'RIP'  'LIP'
    %     For example, 'LPI' (used in Talairach/MNI) means that
    %       * the first dimension goes from left to right
    %       * the second dimension goes from posterior to anterior
    %       * the third dimension goes from inferior to superior
    %   - this function chances the orientation information by adjusting
    %     information in .fa and .a.fdim; contents of .samples remains
    %     unchanged.
    %
    % #   For CoSMoMVPA's copyright information and license terms,   #
    % #   see the COPYING file distributed with CoSMoMVPA.           #

    cosmo_check_dataset(ds, 'fmri');

    ds_orient = cosmo_fmri_orientation(ds);
    upper_new_orient = upper(new_orient);
    [perm, flip] = get_transform(ds_orient, upper_new_orient);

    fmri_dim_labels = {'i', 'j', 'k'};

    ds_fa = ds.fa;

    % initialize output
    fa = ds_fa;
    mat = zeros(4);
    mat(4, 4) = 1;
    dim_size = zeros(1, 3);
    dim_values = cell(3, 1);

    for dim = 1:3
        dim_label = fmri_dim_labels{perm(dim)};
        v = ds_fa.(dim_label);

        [unused, idx, unused, dim_name, values] = cosmo_dim_find(ds, dim_label);
        nvalues = max(values);
        if flip(dim)
            v = nvalues + 1 - v;
            mat(perm(dim), 4) = nvalues + 1;
            mat(perm(dim), dim) = -1;
        else
            mat(perm(dim), dim) = 1;
        end

        dim_size(dim) = nvalues;

        new_dim_label = fmri_dim_labels{dim};
        fa.(new_dim_label) = v;
        dim_values{dim} = 1:nvalues;
    end

    ds.a.fdim.values = dim_values;
    ds.a.vol.mat = ds.a.vol.mat * (mat);
    ds.a.vol.dim = dim_size;
    ds.fa = fa;

    ind = sub2ind(dim_size, ds.fa.i, ds.fa.j, ds.fa.k);
    [foo, i] = sort(ind);
    ds = cosmo_slice(ds, i, 2, false);

    assert(isequal(cosmo_fmri_orientation(ds), upper_new_orient));

function [perm, flip] = get_transform(src, trg)
    is_valid = true;

    labs = ['LR'; 'PA'; 'IS'];
    perm = zeros(1, 3);
    flip = false(1, 3);
    for dim = 1:3
        [src_i, src_j] = find(src(dim) == labs);
        [trg_i, trg_j] = find(bsxfun(@eq, trg', labs(src_i, :)));
        if numel(trg_i) ~= 1
            is_valid = false;
            break

        end
        perm(trg_i) = dim;
        flip(trg_i) = src_j ~= trg_j;
    end

    if ~isequal(sort(perm), 1:3)
        is_valid = false;
    end

    if ~is_valid
        error(['illegal target orientation; run\n'...
               '   help %s\nto see a list of valid orientations'], ...
              mfilename());
    end