avw_img_read - read Analyze format data image (*.img) [ avw, machine ] = avw_img_read(filename) filename - filename of Analyze file (Modified by M. Sato) Do not change orientation of image data written in the file Returned values: avw.hdr - a struct with image data parameters. avw.img - a 3D matrix of image data (double precision). The 3D matrix of image is written in .img file as for z=1:NZ for y=1:NY for x=1:NX read(img(x,y,z)) end end end Default ANALYZE coordinate system is Left-handed LAS X-Y plane is Transverse X-Z plane is Coronal Y-Z plane is Sagittal X axis runs from patient right (low X) to patient Left (high X) Y axis runs from posterior (low Y) to Anterior (high Y) Z axis runs from inferior (low Z) to Superior (high Z) See also: avw_hdr_read (called by this function), avw_view, avw_write, avw_img_write, avw_flip
0001 function [ avw, machine ] = avw_img_read(fileprefix,IMGorient,machine) 0002 % avw_img_read - read Analyze format data image (*.img) 0003 % 0004 % [ avw, machine ] = avw_img_read(filename) 0005 % 0006 % filename - filename of Analyze file 0007 % 0008 % (Modified by M. Sato) 0009 % Do not change orientation of image data written in the file 0010 % 0011 % Returned values: 0012 % 0013 % avw.hdr - a struct with image data parameters. 0014 % avw.img - a 3D matrix of image data (double precision). 0015 % 0016 % The 3D matrix of image is written in .img file as 0017 % for z=1:NZ 0018 % for y=1:NY 0019 % for x=1:NX 0020 % read(img(x,y,z)) 0021 % end 0022 % end 0023 % end 0024 % 0025 % Default ANALYZE coordinate system is Left-handed LAS 0026 % 0027 % X-Y plane is Transverse 0028 % X-Z plane is Coronal 0029 % Y-Z plane is Sagittal 0030 % 0031 % X axis runs from patient right (low X) to patient Left (high X) 0032 % Y axis runs from posterior (low Y) to Anterior (high Y) 0033 % Z axis runs from inferior (low Z) to Superior (high Z) 0034 % 0035 % See also: avw_hdr_read (called by this function), 0036 % avw_view, avw_write, avw_img_write, avw_flip 0037 % 0038 0039 0040 % $Revision: 1332 $ $Date:: 2011-02-24 15:57:20 +0900#$ 0041 0042 % Licence: GNU GPL, no express or implied warranties 0043 % History: 05/2002, Darren.Weber@flinders.edu.au 0044 % The Analyze format is copyright 0045 % (c) Copyright, 1986-1995 0046 % Biomedical Imaging Resource, Mayo Foundation 0047 % 01/2003, Darren.Weber@flinders.edu.au 0048 % - adapted for matlab v5 0049 % - revised all orientation information and handling 0050 % after seeking further advice from AnalyzeDirect.com 0051 % 03/2003, Darren.Weber@flinders.edu.au 0052 % - adapted for -ve pixdim values (non standard Analyze) 0053 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0054 0055 if ~exist('fileprefix','var'), 0056 msg = sprintf('...no input fileprefix - see help avw_img_read\n\n'); 0057 error(msg); 0058 end 0059 if ~exist('IMGorient','var'), IMGorient = ''; end 0060 if ~exist('machine','var'), machine = 'ieee-le'; end 0061 0062 if findstr('.hdr',fileprefix), 0063 fileprefix = strrep(fileprefix,'.hdr',''); 0064 end 0065 if findstr('.img',fileprefix), 0066 fileprefix = strrep(fileprefix,'.img',''); 0067 end 0068 0069 % MAIN 0070 0071 % Read the file header 0072 [ avw, machine ] = avw_hdr_read(fileprefix,machine); 0073 0074 avw = read_image(avw,IMGorient,machine); 0075 0076 return 0077 0078 0079 0080 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0081 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0082 function [ avw ] = read_image(avw,IMGorient,machine) 0083 0084 fid = fopen(sprintf('%s.img',avw.fileprefix),'r',machine); 0085 if fid < 0, 0086 msg = sprintf('...cannot open file %s.img\n\n',avw.fileprefix); 0087 error(msg); 0088 end 0089 0090 % Modified by M. Sato 0091 %ver = '[$Revision: 1332 $]'; 0092 %fprintf('\nAVW_IMG_READ [v%s]\n',ver(12:16)); tic; 0093 0094 % short int bitpix; /* Number of bits per pixel; 1, 8, 16, 32, or 64. */ 0095 % short int datatype /* Datatype for this image set */ 0096 % /*Acceptable values for datatype are*/ 0097 % _NONE 0 0098 % _UNKNOWN 0 /*Unknown data type*/ 0099 % _BINARY 1 /*Binary ( 1 bit per voxel)*/ 0100 % _UNSIGNED_CHAR 2 /*Unsigned character ( 8 bits per voxel)*/ 0101 % _SIGNED_SHORT 4 /*Signed short (16 bits per voxel)*/ 0102 % _SIGNED_INT 8 /*Signed integer (32 bits per voxel)*/ 0103 % _FLOAT 16 /*Floating point (32 bits per voxel)*/ 0104 % _COMPLEX 32 /*Complex (64 bits per voxel; 2 float)/* 0105 % _DOUBLE 64 /*Double precision (64 bits per voxel)*/ 0106 % _RGB 128 /*A Red-Green-Blue datatype*/ 0107 % _ALL 255 /*Undocumented*/ 0108 0109 switch double(avw.hdr.dime.bitpix), 0110 case 1, precision = 'bit1'; 0111 case 8, precision = 'uchar'; 0112 case 16, precision = 'int16'; 0113 case 32, 0114 if isequal(avw.hdr.dime.datatype, 8), precision = 'int32'; 0115 else precision = 'single'; 0116 end 0117 case 64, precision = 'double'; 0118 otherwise, 0119 precision = 'uchar'; 0120 fprintf('...precision undefined in header, using ''uchar''\n'); 0121 end 0122 0123 % read the whole .img file into matlab (faster) 0124 fprintf('...reading %s Analyze %s image format.\n',machine,precision); 0125 0126 fseek(fid,0,'bof'); 0127 0128 % adjust for matlab version 0129 ver = version; 0130 ver = str2num(ver(1)); 0131 if ver < 6, 0132 tmp = fread(fid,inf,sprintf('%s',precision)); 0133 else, 0134 tmp = fread(fid,inf,sprintf('%s=>double',precision)); 0135 end 0136 fclose(fid); 0137 0138 % Update the global min and max values 0139 avw.hdr.dime.glmax = max(double(tmp)); 0140 avw.hdr.dime.glmin = min(double(tmp)); 0141 0142 0143 %--------------------------------------------------------------- 0144 % Now partition the img data into xyz 0145 % dim[0] Number of dimensions in database; usually 4. 0146 % dim[1] Image X dimension; number of pixels in an image row. 0147 % dim[2] Image Y dimension; number of pixel rows in slice. 0148 % dim[3] Volume Z dimension; number of slices in a volume. 0149 % dim[4] Time points; number of volumes in database. 0150 0151 PixelDim = double(avw.hdr.dime.dim(2)); 0152 RowDim = double(avw.hdr.dime.dim(3)); 0153 SliceDim = double(avw.hdr.dime.dim(4)); 0154 0155 PixelSz = double(avw.hdr.dime.pixdim(2)); 0156 RowSz = double(avw.hdr.dime.pixdim(3)); 0157 SliceSz = double(avw.hdr.dime.pixdim(4)); 0158 0159 % ---- NON STANDARD ANALYZE... 0160 0161 % Some Analyze files have been found to set -ve pixdim values, eg 0162 % the MNI template avg152T1_brain in the FSL etc/standard folder, 0163 % perhaps to indicate flipped orientation? If so, this code below 0164 % will NOT handle the flip correctly! 0165 if PixelSz < 0, 0166 warning('X pixdim < 0 !!! resetting to abs(avw.hdr.dime.pixdim(2))'); 0167 PixelSz = abs(PixelSz); 0168 avw.hdr.dime.pixdim(2) = single(PixelSz); 0169 end 0170 if RowSz < 0, 0171 warning('Y pixdim < 0 !!! resetting to abs(avw.hdr.dime.pixdim(3))'); 0172 RowSz = abs(RowSz); 0173 avw.hdr.dime.pixdim(3) = single(RowSz); 0174 end 0175 if SliceSz < 0, 0176 warning('Z pixdim < 0 !!! resetting to abs(avw.hdr.dime.pixdim(4))'); 0177 SliceSz = abs(SliceSz); 0178 avw.hdr.dime.pixdim(4) = single(SliceSz); 0179 end 0180 0181 % ---- END OF NON STANDARD ANALYZE 0182 0183 % --- Do not change orientation of image data 0184 % 0185 0186 % --- Modified by M. Sato 0187 0188 avw.img = reshape(tmp, [PixelDim, RowDim, SliceDim]); 0189 0190 % n = 1; 0191 % x = 1:PixelDim; 0192 % for z = 1:SliceDim, 0193 % for y = 1:RowDim, 0194 % % load Y row of X values into Z slice avw.img 0195 % avw.img(x,y,z) = tmp(n:n+(PixelDim-1)); 0196 % n = n + PixelDim; 0197 % end 0198 % end 0199 0200 % no need to rearrange avw.hdr.dime.dim or avw.hdr.dime.pixdim 0201 0202 % 0203 %t=toc; fprintf('...done (%5.2f sec).\n\n',t); 0204 0205 return 0206 0207 % From ANALYZE documentation... 0208 % 0209 % The ANALYZE coordinate system has an origin in the lower left 0210 % corner. That is, with the subject lying supine, the coordinate 0211 % origin is on the right side of the body (x), at the back (y), 0212 % and at the feet (z). This means that: 0213 % 0214 % +X increases from right (R) to left (L) 0215 % +Y increases from the back (posterior,P) to the front (anterior, A) 0216 % +Z increases from the feet (inferior,I) to the head (superior, S) 0217 % 0218 % The LAS orientation is the radiological convention, where patient 0219 % left is on the image right. The alternative neurological 0220 % convention is RAS (also Talairach convention). 0221 % 0222 % A major advantage of the Analzye origin convention is that the 0223 % coordinate origin of each orthogonal orientation (transverse, 0224 % coronal, and sagittal) lies in the lower left corner of the 0225 % slice as it is displayed. 0226 % 0227 % Orthogonal slices are numbered from one to the number of slices 0228 % in that orientation. For example, a volume (x, y, z) dimensioned 0229 % 128, 256, 48 has: 0230 % 0231 % 128 sagittal slices numbered 1 through 128 (X) 0232 % 256 coronal slices numbered 1 through 256 (Y) 0233 % 48 transverse slices numbered 1 through 48 (Z) 0234 % 0235 % Pixel coordinates are made with reference to the slice numbers from 0236 % which the pixels come. Thus, the first pixel in the volume is 0237 % referenced p(1,1,1) and not at p(0,0,0). 0238 % 0239 % Transverse slices are in the XY plane (also known as axial slices). 0240 % Sagittal slices are in the ZY plane. 0241 % Coronal slices are in the ZX plane. 0242 % 0243 %----------------------------------------------------------------------------