avw_img_write

 avw_img_write - write Analyze image files (*.img)
 
 avw_img_write(avw,fileprefix,[IMGorient],[machine],[verbose])
 
 avw.img    - a 3D matrix of image data (double precision).
 avw.hdr    - a struct with image data parameters.  If
              not empty, this function calls avw_hdr_write.
 
 fileprefix - a string, the filename without the .img
              extension. If empty, may use avw.fileprefix
 
 IMGorient - optional int, force writing of specified 
             orientation, with values:
 
   [],  if empty, will use avw.hdr.hist.orient field
    0,  transverse/axial unflipped (default, radiological)
    1,  coronal unflipped
    2,  sagittal unflipped
    3,  transverse/axial flipped, left to right
    4,  coronal flipped, anterior to posterior
    5,  sagittal flipped, superior to inferior
 
 This function will set avw.hdr.hist.orient and write the 
 image data in a corresponding order.  This function is 
 in alpha development, so it has not been exhaustively 
 tested (07/2003). See avw_img_read for more information 
 and documentation on the orientation option.  
 Orientations 3-5 are NOT recommended!  They are part 
 of the Analyze format, but only used in Analyze
 for faster raster graphics during movies.
 
 machine - a string, see machineformat in fread for details.
           The default here is 'ieee-le'.
 
 verbose - the default is to output processing information to the command
           window.  If verbose = 0, this will not happen.

 Tip: to change the data type, set avw.hdr.dime.datatype to:
 
     1    Binary             (  1 bit  per voxel)
     2    Unsigned character (  8 bits per voxel)
     4    Signed short       ( 16 bits per voxel)
     8    Signed integer     ( 32 bits per voxel)
    16    Floating point     ( 32 bits per voxel)
    32    Complex, 2 floats  ( 64 bits per voxel), not supported
    64    Double precision   ( 64 bits per voxel)
   128    Red-Green-Blue     (128 bits per voxel), not supported
 
 See also: avw_write, avw_hdr_write,
           avw_read, avw_hdr_read, avw_img_read, avw_view

0001 function avw_img_write(avw, fileprefix, IMGorient, machine, verbose)
0002 
0003 % avw_img_write - write Analyze image files (*.img)
0004 %
0005 % avw_img_write(avw,fileprefix,[IMGorient],[machine],[verbose])
0006 %
0007 % avw.img    - a 3D matrix of image data (double precision).
0008 % avw.hdr    - a struct with image data parameters.  If
0009 %              not empty, this function calls avw_hdr_write.
0010 %
0011 % fileprefix - a string, the filename without the .img
0012 %              extension. If empty, may use avw.fileprefix
0013 %
0014 % IMGorient - optional int, force writing of specified
0015 %             orientation, with values:
0016 %
0017 %   [],  if empty, will use avw.hdr.hist.orient field
0018 %    0,  transverse/axial unflipped (default, radiological)
0019 %    1,  coronal unflipped
0020 %    2,  sagittal unflipped
0021 %    3,  transverse/axial flipped, left to right
0022 %    4,  coronal flipped, anterior to posterior
0023 %    5,  sagittal flipped, superior to inferior
0024 %
0025 % This function will set avw.hdr.hist.orient and write the
0026 % image data in a corresponding order.  This function is
0027 % in alpha development, so it has not been exhaustively
0028 % tested (07/2003). See avw_img_read for more information
0029 % and documentation on the orientation option.
0030 % Orientations 3-5 are NOT recommended!  They are part
0031 % of the Analyze format, but only used in Analyze
0032 % for faster raster graphics during movies.
0033 %
0034 % machine - a string, see machineformat in fread for details.
0035 %           The default here is 'ieee-le'.
0036 %
0037 % verbose - the default is to output processing information to the command
0038 %           window.  If verbose = 0, this will not happen.
0039 %
0040 % Tip: to change the data type, set avw.hdr.dime.datatype to:
0041 %
0042 %     1    Binary             (  1 bit  per voxel)
0043 %     2    Unsigned character (  8 bits per voxel)
0044 %     4    Signed short       ( 16 bits per voxel)
0045 %     8    Signed integer     ( 32 bits per voxel)
0046 %    16    Floating point     ( 32 bits per voxel)
0047 %    32    Complex, 2 floats  ( 64 bits per voxel), not supported
0048 %    64    Double precision   ( 64 bits per voxel)
0049 %   128    Red-Green-Blue     (128 bits per voxel), not supported
0050 %
0051 % See also: avw_write, avw_hdr_write,
0052 %           avw_read, avw_hdr_read, avw_img_read, avw_view
0053 %
0054 
0055 % $Revision: 1.1 $ $Date: 2004/11/12 01:30:25 $
0056 
0057 % Licence:  GNU GPL, no express or implied warranties
0058 % History:  05/2002, Darren.Weber@flinders.edu.au
0059 %                    The Analyze format is copyright
0060 %                    (c) Copyright, 1986-1995
0061 %                    Biomedical Imaging Resource, Mayo Foundation
0062 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0063 
0064 
0065 %------------------------------------------------------------------------
0066 % Check inputs
0067 
0068 if ~exist('avw','var'),
0069     doc avw_img_write;
0070     error('...no input avw.');
0071 elseif isempty(avw),
0072     error('...empty input avw.');
0073 elseif ~isfield(avw,'img'),
0074     error('...empty input avw.img');
0075 end
0076 
0077 if ~exist('fileprefix','var'),
0078     if isfield(avw,'fileprefix'),
0079         if ~isempty(avw.fileprefix),
0080             fileprefix = avw.fileprefix;
0081         else
0082             fileprefix = [];
0083         end
0084     else
0085         fileprefix = [];
0086     end
0087 end
0088 if isempty(fileprefix),
0089     [fileprefix, pathname, filterindex] = uiputfile('*.hdr','Specify an output Analyze .hdr file');
0090     if pathname, cd(pathname); end
0091     if ~fileprefix,
0092         doc avw_img_write;
0093         error('no output .hdr file specified');
0094     end
0095 end
0096 
0097 if findstr('.hdr',fileprefix),
0098 %    fprintf('AVW_IMG_WRITE: Removing .hdr extension from ''%s''\n',fileprefix);
0099     fileprefix = strrep(fileprefix,'.hdr','');
0100 end
0101 if findstr('.img',fileprefix),
0102 %    fprintf('AVW_IMG_WRITE: Removing .img extension from ''%s''\n',fileprefix);
0103     fileprefix = strrep(fileprefix,'.img','');
0104 end
0105 
0106 if ~exist('IMGorient','var'), IMGorient = ''; end
0107 if ~exist('machine','var'), machine = 'ieee-le'; end
0108 if ~exist('verbose','var'), verbose = 1; end
0109 
0110 if isempty(IMGorient), IMGorient = ''; end
0111 if isempty(machine), machine = 'ieee-le'; end
0112 if isempty(verbose), verbose = 1; end
0113 
0114 
0115 
0116 %------------------------------------------------------------------------
0117 % MAIN
0118 if verbose,
0119     version = '[$Revision: 1.1 $]';
0120     fprintf('\nAVW_IMG_WRITE [v%s]\n',version(12:16));  tic;
0121 end
0122 
0123 fid = fopen(sprintf('%s.img',fileprefix),'w',machine);
0124 if fid < 0,
0125     msg = sprintf('Cannot open file %s.img\n',fileprefix);
0126     error(msg);
0127 else
0128     avw = write_image(fid,avw,fileprefix,IMGorient,machine,verbose);
0129 end
0130 
0131 if verbose,
0132     t=toc; fprintf('...done (%5.2f sec).\n\n',t);
0133 end
0134 
0135 % MUST write header after the image, to ensure any
0136 % orientation changes during image write are saved
0137 % in the header
0138 avw_hdr_write(avw,fileprefix,machine,verbose);
0139 
0140 return
0141 
0142 
0143 
0144 
0145 %-----------------------------------------------------------------------------------
0146 
0147 function avw = write_image(fid,avw,fileprefix,IMGorient,machine,verbose)
0148 
0149 % short int bitpix;    /* Number of bits per pixel; 1, 8, 16, 32, or 64. */
0150 % short int datatype   /* Datatype for this image set */
0151 % /*Acceptable values for datatype are*/
0152 % #define DT_NONE             0
0153 % #define DT_UNKNOWN          0    /*Unknown data type*/
0154 % #define DT_BINARY           1    /*Binary             ( 1 bit per voxel)*/
0155 % #define DT_UNSIGNED_CHAR    2    /*Unsigned character ( 8 bits per voxel)*/
0156 % #define DT_SIGNED_SHORT     4    /*Signed short       (16 bits per voxel)*/
0157 % #define DT_SIGNED_INT       8    /*Signed integer     (32 bits per voxel)*/
0158 % #define DT_FLOAT           16    /*Floating point     (32 bits per voxel)*/
0159 % #define DT_COMPLEX         32    /*Complex,2 floats   (64 bits per voxel)/*
0160 % #define DT_DOUBLE          64    /*Double precision   (64 bits per voxel)*/
0161 % #define DT_RGB            128    /*A Red-Green-Blue datatype*/
0162 % #define DT_ALL            255    /*Undocumented*/
0163 
0164 switch double(avw.hdr.dime.datatype),
0165 case   1,
0166     avw.hdr.dime.bitpix = int16( 1); precision = 'bit1';
0167 case   2,
0168     avw.hdr.dime.bitpix = int16( 8); precision = 'uchar';
0169 case   4,
0170     avw.hdr.dime.bitpix = int16(16); precision = 'int16';
0171 case   8,
0172     avw.hdr.dime.bitpix = int16(32); precision = 'int32';
0173 case  16,
0174     avw.hdr.dime.bitpix = int16(32); precision = 'single';
0175 case  32,
0176     error('...complex datatype not yet supported.\n');
0177 case  64,
0178     avw.hdr.dime.bitpix = int16(64); precision = 'double';
0179 case 128,
0180     error('...RGB datatype not yet supported.\n');
0181 otherwise
0182     warning('...unknown datatype, using type 16 (32 bit floats).\n');
0183     avw.hdr.dime.datatype = int16(16);
0184     avw.hdr.dime.bitpix = int16(32); precision = 'single';
0185 end
0186 
0187 
0188 % write the .img file, depending on the .img orientation
0189 if verbose,
0190     fprintf('...writing %s precision Analyze image (%s).\n',precision,machine);
0191 end
0192 
0193 fseek(fid,0,'bof');
0194 
0195 % The standard image orientation is axial unflipped
0196 if isempty(avw.hdr.hist.orient),
0197     msg = [ '...avw.hdr.hist.orient ~= 0.\n',...
0198             '   This function assumes the input avw.img is\n',...
0199             '   in axial unflipped orientation in memory.  This is\n',...
0200             '   created by the avw_img_read function, which converts\n',...
0201             '   any input file image to axial unflipped in memory.\n'];
0202     warning(msg)
0203 end
0204 
0205 if isempty(IMGorient),
0206     if verbose,
0207         fprintf('...no IMGorient specified, using avw.hdr.hist.orient value.\n');
0208     end
0209     IMGorient = double(avw.hdr.hist.orient);
0210 end
0211 
0212 if ~isfinite(IMGorient),
0213     if verbose,
0214         fprintf('...IMGorient is not finite!\n');
0215     end
0216     IMGorient = 99;
0217 end
0218 
0219 switch IMGorient,
0220     
0221 case 0, % transverse/axial unflipped
0222     
0223     % For the 'transverse unflipped' type, the voxels are stored with
0224     % Pixels in 'x' axis (varies fastest) - from patient right to left
0225     % Rows in   'y' axis                  - from patient posterior to anterior
0226     % Slices in 'z' axis                  - from patient inferior to superior
0227     
0228     if verbose, fprintf('...writing axial unflipped\n'); end
0229     
0230     avw.hdr.hist.orient = uint8(0);
0231     
0232     SliceDim = double(avw.hdr.dime.dim(4)); % z
0233     RowDim   = double(avw.hdr.dime.dim(3)); % y
0234     PixelDim = double(avw.hdr.dime.dim(2)); % x
0235     SliceSz  = double(avw.hdr.dime.pixdim(4));
0236     RowSz    = double(avw.hdr.dime.pixdim(3));
0237     PixelSz  = double(avw.hdr.dime.pixdim(2));
0238     
0239     x = 1:PixelDim;
0240     for z = 1:SliceDim,
0241         for y = 1:RowDim,
0242             fwrite(fid,avw.img(x,y,z),precision);
0243         end
0244     end
0245     
0246 case 1, % coronal unflipped
0247     
0248     % For the 'coronal unflipped' type, the voxels are stored with
0249     % Pixels in 'x' axis (varies fastest) - from patient right to left
0250     % Rows in   'z' axis                  - from patient inferior to superior
0251     % Slices in 'y' axis                  - from patient posterior to anterior
0252     
0253     if verbose, fprintf('...writing coronal unflipped\n'); end
0254     
0255     avw.hdr.hist.orient = uint8(1);
0256     
0257     SliceDim = double(avw.hdr.dime.dim(3)); % y
0258     RowDim   = double(avw.hdr.dime.dim(4)); % z
0259     PixelDim = double(avw.hdr.dime.dim(2)); % x
0260     SliceSz  = double(avw.hdr.dime.pixdim(3));
0261     RowSz    = double(avw.hdr.dime.pixdim(4));
0262     PixelSz  = double(avw.hdr.dime.pixdim(2));
0263     
0264     x = 1:PixelDim;
0265     for y = 1:SliceDim,
0266         for z = 1:RowDim,
0267             fwrite(fid,avw.img(x,y,z),precision);
0268         end
0269     end
0270     
0271 case 2, % sagittal unflipped
0272     
0273     % For the 'sagittal unflipped' type, the voxels are stored with
0274     % Pixels in 'y' axis (varies fastest) - from patient posterior to anterior
0275     % Rows in   'z' axis                  - from patient inferior to superior
0276     % Slices in 'x' axis                  - from patient right to left
0277     
0278     if verbose, fprintf('...writing sagittal unflipped\n'); end
0279     
0280     avw.hdr.hist.orient = uint8(2);
0281     
0282     SliceDim = double(avw.hdr.dime.dim(2)); % x
0283     RowDim   = double(avw.hdr.dime.dim(4)); % z
0284     PixelDim = double(avw.hdr.dime.dim(3)); % y
0285     SliceSz  = double(avw.hdr.dime.pixdim(2));
0286     RowSz    = double(avw.hdr.dime.pixdim(4));
0287     PixelSz  = double(avw.hdr.dime.pixdim(3));
0288     
0289     y = 1:PixelDim;
0290     for x = 1:SliceDim,
0291         for z = 1:RowDim,
0292             fwrite(fid,avw.img(x,y,z),precision);
0293         end
0294     end
0295     
0296 case 3, % transverse/axial flipped
0297     
0298     % For the 'transverse flipped' type, the voxels are stored with
0299     % Pixels in 'x' axis (varies fastest) - from patient right to left
0300     % Rows in   'y' axis                  - from patient anterior to posterior*
0301     % Slices in 'z' axis                  - from patient inferior to superior
0302     
0303     if verbose,
0304         fprintf('...writing axial flipped (+Y from Anterior to Posterior)\n');
0305     end
0306     
0307     avw.hdr.hist.orient = uint8(3);
0308     
0309     SliceDim = double(avw.hdr.dime.dim(4)); % z
0310     RowDim   = double(avw.hdr.dime.dim(3)); % y
0311     PixelDim = double(avw.hdr.dime.dim(2)); % x
0312     SliceSz  = double(avw.hdr.dime.pixdim(4));
0313     RowSz    = double(avw.hdr.dime.pixdim(3));
0314     PixelSz  = double(avw.hdr.dime.pixdim(2));
0315     
0316     x = 1:PixelDim;
0317     for z = 1:SliceDim,
0318         for y = RowDim:-1:1, % flipped in Y
0319             fwrite(fid,avw.img(x,y,z),precision);
0320         end
0321     end
0322     
0323 case 4, % coronal flipped
0324     
0325     % For the 'coronal flipped' type, the voxels are stored with
0326     % Pixels in 'x' axis (varies fastest) - from patient right to left
0327     % Rows in   'z' axis                  - from patient inferior to superior
0328     % Slices in 'y' axis                  - from patient anterior to posterior
0329     
0330     if verbose,
0331         fprintf('...writing coronal flipped (+Z from Superior to Inferior)\n');
0332     end
0333     
0334     avw.hdr.hist.orient = uint8(4);
0335     
0336     SliceDim = double(avw.hdr.dime.dim(3)); % y
0337     RowDim   = double(avw.hdr.dime.dim(4)); % z
0338     PixelDim = double(avw.hdr.dime.dim(2)); % x
0339     SliceSz  = double(avw.hdr.dime.pixdim(3));
0340     RowSz    = double(avw.hdr.dime.pixdim(4));
0341     PixelSz  = double(avw.hdr.dime.pixdim(2));
0342     
0343     x = 1:PixelDim;
0344     for y = 1:SliceDim,
0345         for z = RowDim:-1:1,
0346             fwrite(fid,avw.img(x,y,z),precision);
0347         end
0348     end
0349     
0350 case 5, % sagittal flipped
0351     
0352     % For the 'sagittal flipped' type, the voxels are stored with
0353     % Pixels in 'y' axis (varies fastest) - from patient posterior to anterior
0354     % Rows in   'z' axis                  - from patient superior to inferior
0355     % Slices in 'x' axis                  - from patient right to left
0356     
0357     if verbose,
0358         fprintf('...writing sagittal flipped (+Z from Superior to Inferior)\n');
0359     end
0360     
0361     avw.hdr.hist.orient = uint8(5);
0362     
0363     SliceDim = double(avw.hdr.dime.dim(2)); % x
0364     RowDim   = double(avw.hdr.dime.dim(4)); % z
0365     PixelDim = double(avw.hdr.dime.dim(3)); % y
0366     SliceSz  = double(avw.hdr.dime.pixdim(2));
0367     RowSz    = double(avw.hdr.dime.pixdim(4));
0368     PixelSz  = double(avw.hdr.dime.pixdim(3));
0369     
0370     y = 1:PixelDim;
0371     for x = 1:SliceDim,
0372         for z = RowDim:-1:1, % superior to inferior
0373             fwrite(fid,avw.img(x,y,z),precision);
0374         end
0375     end
0376     
0377 otherwise, % transverse/axial unflipped
0378     
0379     % For the 'transverse unflipped' type, the voxels are stored with
0380     % Pixels in 'x' axis (varies fastest) - from patient right to left
0381     % Rows in   'y' axis                  - from patient posterior to anterior
0382     % Slices in 'z' axis                  - from patient inferior to superior
0383     
0384     if verbose,
0385         fprintf('...unknown orientation specified, assuming default axial unflipped\n');
0386     end
0387     
0388     avw.hdr.hist.orient = uint8(0);
0389     
0390     SliceDim = double(avw.hdr.dime.dim(4)); % z
0391     RowDim   = double(avw.hdr.dime.dim(3)); % y
0392     PixelDim = double(avw.hdr.dime.dim(2)); % x
0393     SliceSz  = double(avw.hdr.dime.pixdim(4));
0394     RowSz    = double(avw.hdr.dime.pixdim(3));
0395     PixelSz  = double(avw.hdr.dime.pixdim(2));
0396     
0397     x = 1:PixelDim;
0398     for z = 1:SliceDim,
0399         for y = 1:RowDim,
0400             fwrite(fid,avw.img(x,y,z),precision);
0401         end
0402     end
0403     
0404 end
0405 
0406 fclose(fid);
0407 
0408 % Update the header
0409 avw.hdr.dime.dim(2:4) = int16([PixelDim,RowDim,SliceDim]);
0410 avw.hdr.dime.pixdim(2:4) = single([PixelSz,RowSz,SliceSz]);
0411 
0412 return