fmr_writeanalyzevol.m

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

function rvalue = fmr_WriteAnalyzeVol(hfile, volnum, filename)
% FMR::WriteAnalyzeVol  - write an Analyze image from one volume
%
% FORMAT:       [success] = fmr.WriteAnalyzeVol(volume, filename);
%
% Input fields:
%
%       volume      volume number to write
%       filename    analyze filename
%
% Output fields:
%
%       success     true if write was successful

% Version:  v0.7b
% Build:    7083010
% Date:     Aug-30 2007, 10:19 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, 'fmr')
    error( ...
        'BVQXfile:BadArgument', ...
        'Invalid call to %s.', ...
        mfilename ...
    );
end

% default rvalue
rvalue = false;

% get reference
sc = bvqxfile_getscont(hfile.L);
bc = sc.C;
ofile = sc.F;
[opath{1:2}] = fileparts(ofile);
ofile = opath{2};

% slice data loaded at all?
if isempty(bc.Slice) || ...
   ~isstruct(bc.Slice) || ...
   ~isfield(bc.Slice, 'STCData')
    error( ...
        'BVQXtools:InvalidCall', ...
        'Slice data not loaded, issue LoadSTC method first.' ...
    );
end

% check further arguments
if nargin < 2 || ...
   ~isa(volnum, 'double') || ...
    numel(volnum) ~= 1 || ...
    isnan(volnum) || ...
    isinf(volnum) || ...
    volnum ~= fix(volnum) || ...
    volnum < 1 || ...
    volnum > bc.NrOfVolumes 
    error( ...
        'BVQXfile:BadArgument', ...
        'Invalid or missing volnum argument for call to %s.', ...
        mfilename ...
    );
end
if nargin < 3 || ...
   ~ischar(filename) || ...
    isempty(filename) || ...
    numel(filename) ~= length(filename) || ...
    length(filename) < 5 || ...
  (~strcmpi(filename(end-3:end), '.hdr') && ...
   ~strcmpi(filename(end-3:end), '.img'))
    error( ...
        'BVQXfile:BadArgument', ...
        'Invalid or missing filename argument for call to %s.', ...
        mfilename ...
    );
end

% filename
filename = filename(:)';
filename = [filename(1:end-3) 'img'];
hdrfname = [filename(1:end-3) 'hdr'];
matfname = [filename(1:end-3) 'mat'];

% dimension and datatype
xpix = bc.ResolutionX;
ypix = bc.ResolutionY;
zpix = bc.NrOfSlices;
if (bc.DataStorageFormat == 1 && ...
    zpix ~= length(bc.Slice)) || ...
   (bc.DataStorageFormat == 2 && ...
    zpix ~= size(bc.Slice.STCData, 4))
    error( ...
        'BVQXfile:InvalidArraySize', ...
        'Invalid number of slices in header, cannot save data.' ...
    );
end
isiz = [xpix, ypix, zpix, 1];

% X/Y resolution
xres = bc.InplaneResolutionX;
yres = bc.InplaneResolutionY;
zres = bc.SliceThickness + bc.SliceGap;
tres = bc.TR / 1000;

% positional information
xvec = xres .* [bc.RowDirX; bc.RowDirY; bc.RowDirZ];
yvec = yres .* [bc.ColDirX; bc.ColDirY; bc.ColDirZ];
ovec = [ bc.Slice1CenterX; bc.Slice1CenterY; bc.Slice1CenterZ];
lvec = [ bc.SliceNCenterX; bc.SliceNCenterY; bc.SliceNCenterZ];
zspan = bc.NrOfSlices - 1;
zvec = (1 / zspan) .* (lvec - ovec);
svec = ovec - ...
    ((xpix + 1) / 2) * xvec - ...
    ((ypix + 1) / 2) * yvec;

% build mat
tmat = [[xvec, yvec, zvec], svec;  0, 0, 0, 1];

% radiological convention ?
if (isfield(bc, 'Convention') && ...
    bc.Convention(1)) || ...
   (isfield(bc, 'RadiologicalConvention') && ...
    strcmpi(bc.RadiologicalConvention(:)', 'yes'))
    tfmat = eye(4);
    tfmat(2, 2) = -1;
    tmat = tfmat * tmat;
end

% try volume creation
try
    hdr = BVQXfile('new:hdr');
    hdrc = bvqxfile_getcont(hdr.L);
catch
    error( ...
        'BVQXfile:InternalError', ...
        'Error creating Analyze header object: %s.', ...
        lasterr ...
    );
end

% set dims and data
hdrc.ImgDim.Dim(1:5) = [4, isiz];
hdrc.ImgDim.DataType = 4;
hdrc.ImgDim.BitsPerPixel = 16;
hdrc.ImgDim.PixSpacing(2:5) = [xres, yres, zres, tres];
hdrc.ImgDim.CalMaxDisplay = 32767;
hdrc.ImgDim.CalMinDisplay = 0;
hdrc.ImgDim.GLMax = 32767;
hdrc.ImgDim.GLMin = 0;
hdrc.DataHist.Description = sprintf('Volume %d of FMR %s', volnum, ofile);
hdrc.DataHist.ScanNumber = sprintf('%d', volnum);
hdrc.VoxelData = int16(zeros(isiz));

% which format
switch bc.DataStorageFormat
    
    % old format (BVQX <= 1.8.7)
    case {1}
        for sc = 1:zpix
            hdrc.VoxelData(:,:,sc,1) = bc.Slice(sc).STCData(:, :, volnum);
        end
        
    % new format (BVQX >= 1.9.9)
    case {2}
        hdrc.VoxelData(:, :, :) = squeeze(bc.Slice.STCData(:, :, volnum, :));
        
    % unsupported
    otherwise
        BVQXfile(0, 'clearobj', hdr.L);
        error( ...
            'BVQXfile:InvalidObject', ...
            'Unsupported DataStorageFormat of FMR.' ...
        );
end

% save hdr/img
try
    rvalue = false;
    aft_SaveAs(hdr, hdrfname);
    rvalue = hdr_SaveVoxelData(hdr);
catch
    warning( ...
        'BVQXfile:InternalError', ...
        'Error writing Analyze header/image %s.', ...
        hdrfname ...
    );
end
BVQXfile(0, 'clearobj', hdr.L);

% leave early ?
if ~rvalue
    return;
end

% save mat
eval('M=tmat;mat=tmat;save(matfname,''M'',''mat'',''-v6''', '');
if exist(matfname, 'file') ~= 2
    rvalue = false;
end