samplefmrspace.m

toolbox of BVQX, This is the access between BV and matlab. It will help you to analysis data from BV

function [y, c] = samplefmrspace(x, c, fmr, trf, imeth, conly)
% samplefmrspace  - samples a data slab at a given set of coordinates
%
% FORMAT:       [y, c] = samplefmrspace(x, c, fmr [, trf [, meth [, conly]])
%
% Input fields:
%
%       x       data in FMR resolution (first matching dims, so STC works)
%       c       Cx3 list of coordinates to sample at, or a cell array
%               suitable to pass as args into ndgrid (in BVsystem
%               coordinates, axes order as TAL, so that BVsys = 128 - Tal)
%       fmr     FMR object
%       trf     optional 1xT cell array of transformation files,
%               e.g. {ia, fa, acpc, tal}
%       meth    optional, 'nearest', {'linear'}, 'cubic', 'spline'
%       conly   if given, only return coordinates of sampling
%
% Output fields:
%
%       y       Cx1 data values at the C coordinates
%       c       1x3 cell array with Cx1 ordinates for interpn sampling

% Version:  v0.7a
% Build:    7081012
% Date:     Aug-10 2007, 12:15 PM CEST
% Author:   Jochen Weber, Brain Innovation, B.V., Maastricht, NL
% URL/Info: http://wiki.brainvoyager.com/BVQXtools

% argument check
if nargin < 3 || ...
    isempty(x) || ...
   (~isnumeric(x) && ~istransio(x)) || ...
    numel(size(x)) < 3 || ...
   ~isBVQXfile(fmr, 'fmr') || ...
    numel(fmr) ~= 1 || ...
   ((~iscell(c) || ...
      numel(c) ~= 3 || ...
      ~isa(c{1}, 'double') || ...
      ~isa(c{2}, 'double') || ...
      ~isa(c{3}, 'double') || ...
      any(isinf([c{1}(:); c{2}(:); c{3}(:)]) | isinf([c{1}(:); c{2}(:); c{3}(:)]))) && ...
    (~isa(c, 'double') || ...
      numel(size(c)) ~= 2 || ...
      size(c, 2) ~= 3))
    error( ...
        'BVQXtools:BadArgument', ...
        'Bad or missing argument.' ...
    );
end
szknown = false;
if iscell(c)
    szknown = [numel(c{1}), numel(c{2}), numel(c{3})];
    [cx{1:3}] = ndgrid(c{1}, c{2}, c{3});
    c = [cx{1}(:), cx{2}(:), cx{3}(:)];
end
if nargin < 4
    trf = {};
end
if any(isBVQXfile(trf))
    trfc = cell(1, numel(trf));
    for ci = 1:numel(trf)
        trfc{ci} = trf(ci);
    end
    trf = trfc;
end
if ~iscell(trf)
    error( ...
        'BVQXtools:BadArgument', ...
        'Bad trf argument.' ...
    );
end
if nargin < 5 || ...
   ~ischar(imeth) || ...
    isempty(imeth) || ...
   ~any(strcmpi(imeth(:)', {'nearest', 'linear', 'cubic', 'spline'}))
    imeth = 'linear';
else
    imeth = lower(imeth(:)');
end
if nargin < 6 || ...
   ~islogical(conly) || ...
    numel(conly) ~= 1
    conly = false;
end

% get fmr CoordinateFrame
fmrc = fmr.CoordinateFrame();

% get dimensions to work with
sd = size(x);
td = [fmrc.DimX, fmrc.DimY, fmrc.DimZ];
si = 1;
% ti = [1, 2, 3];
while si < numel(sd) && ...
    (sd(si) ~= td(1))
    si = si + 1;
end
if si > numel(sd)
    error( ...
        'BVQXtools:BadArgument', ...
        'Dimensions of dataslab and FMR mismatch.' ...
    );
end
% ti(1) = si;
si = si + 1;
while si < numel(sd) && ...
    (sd(si) ~= td(2))
    si = si + 1;
end
if si > numel(sd)
    error( ...
        'BVQXtools:BadArgument', ...
        'Dimensions of dataslab and FMR mismatch.' ...
    );
end
% ti(2) = si;
si = si + 1;
while si < numel(sd) && ...
    (sd(si) ~= td(3))
    si = si + 1;
end
if si > numel(sd)
    error( ...
        'BVQXtools:BadArgument', ...
        'Dimensions of dataslab and FMR mismatch.' ...
    );
end
% ti(3) = si;

% transform target TAL to source ACPC if necessary
for ci = numel(trf):-1:1
    if isBVQXfile(trf{ci}, 'tal')
        c = acpc2tal(c, trf{ci}, true);
        trf(ci) = [];
    end
end

% subtract center coordinate
if fmr.CoordinateSystem ~= 1
    c = 255 - c;
end
c = c - 127.5;

% add forth column ...
c(:, 4) = 1;

% inspect list of transformations with, supposedly, no tal transform
for ci = numel(trf):-1:1
    if ~isBVQXfile(trf{ci}, 'trf')
        error( ...
            'BVQXtools:BadArgument', ...
            'Invalid trf argument (particle).' ...
        );
    end
    if trf{ci}.TransformationType > 1
        c(:, [2, 3, 1, 4]) = c(:, [2, 3, 1, 4]) * trf{ci}.TFMatrix';
    elseif trf{ci}.AlignmentStep == 2
        c(:, [2, 3, 1, 4]) = c(:, [2, 3, 1, 4]) * inv(trf{ci}.TFMatrix)';
    else
        c = c(:, [2, 3, 1, 4]) * inv(trf{ci}.TFMatrix)';
    end
end

% prepare sampling
c = {c(:, 1) ./ fmrc.ResX + (fmrc.DimX + 1) / 2, ...
     c(:, 2) ./ fmrc.ResY + (fmrc.DimY + 1) / 2, ...
    -c(:, 3) ./ fmrc.ResZ + (fmrc.DimZ + 1) / 2};

% only cooordinates
if conly
    y = zsz(c{1});
    return;
end

% sample x at c
if numel(sd) == 3
    y = interpn(double(x), c{1}, c{2}, c{3}, imeth);
else
end

% remove NaNs
y(isnan(y(:))) = 0;

% reshape ?
if numel(szknown) == 3
    y = reshape(y, szknown);
end