vmp_voxelstats.m

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

function stats = vmp_VoxelStats(hfile, mapno, coords, ctype)
% VMP::VoxelStats  - retrieve voxel statistics from map
%
% FORMAT:       stats = vmp.VoxelStats(mapno, coords [, ctype])
%
% Input fields:
%
%       mapno       map number (in range [1..M])
%       coords      Nx3 coordinate list
%       ctype       coordinate system type, one of
%                   {'BV'}, 'BVTal', 'Tal'
%                   where BV is the native BV system (X->Y, Y->Z, Z->X)
%                   BVTal is the same, but uses Tal axes order
%                   Tal is just (127 - BVTal)
%
% Output fields:
%
%       stats       extracted statistics from Map

% Version:  v0.7b
% Build:    7090215
% Date:     Sep-15 2007, 3:31 AM CEST
% Author:   Jochen Weber, Brain Innovation, B.V., Maastricht, NL
% URL/Info: http://wiki.brainvoyager.com/BVQXtools

% argument check
if nargin < 3 || ...
    numel(hfile) ~= 1 || ...
   ~isBVQXfile(hfile, 'vmp') || ...
   ~isa(mapno, 'double') || ...
    numel(mapno) ~= 1 || ...
    isinf(mapno) || ...
    isnan(mapno) || ...
    mapno < 1 || ...
   ~isa(coords, 'double') || ...
    isempty(coords) || ...
    size(coords, 2) ~= 3 || ...
    any(isinf(coords(:)) | isnan(coords(:)) | coords(:) < -128 | coords(:) > 256)
    error( ...
        'BVQXfile:BadArgument', ...
        'Invalid call to %s.', ...
        mfilename ...
    );
end
bc = bvqxfile_getcont(hfile.L);
if mapno > bc.NrOfMaps
    error( ...
        'BVQXfile:BadArgument', ...
        'Map number out of bounds.' ...
    );
end
if nargin < 4 || ...
   ~ischar(ctype) || ...
   ~any(strcmpi(ctype(:)', {'bv', 'bvtal', 'tal'}))
    ctype = 'bv';
else
    ctype = lower(ctype(:)');
end

% only 1x1x1 resolution supported anyway!
coords = round(coords);
ncoord = size(coords, 1);
scoord = cumprod(size(bc.Map(1).VMPData));

% convert coords to internal space
if strcmp(ctype, 'tal')
    coords = 127 - coords;
    ctype = 'bvtal';
end
if strcmp(ctype, 'bvtal')
    coords = [coords(:,[2, 3]), coords(:, 1)];
end

% calculate real offsets into matrix
res = bc.Resolution;
x1 = bc.XStart - 1;
y1 = bc.YStart - 1;
z1 = bc.ZStart - 1;
coords = round( ...
    (coords - repmat([x1, y1, z1], [ncoord, 1])) / res);
coords = 1 + coords(:,1) + scoord(1) * coords(:, 2) + scoord(2) * coords(:, 3);
stats = bc.Map(fix(mapno)).VMPData(coords);