function [img] = rri_load_img(fname)
%
% syntax: [img] = rri_load_img(fname)
%
%    This script is used to read the img file specified by fname.
%    The corresponding header file (with name of image.hdr) is
%    read to determine the data type and dimensions.
%
%
%  W. Chau (May 27, 99)

%  determine the dimensions and data type, as well as the data byte ordering 

[fpath,iname] = fileparts(fname);
hfile = fullfile(fpath,[iname,'.hdr']);

mf = 'ieee-le';
fid = fopen(hfile,'r',mf);    % try little-endian byte ordering

  [fn,perm,mf]=fopen(fid);

  fseek(fid,42,'bof');
  [dims]=fread(fid,3,'int16');		

  fseek(fid,70,'bof');
  [datatype]=fread(fid,1,'int16');  % 2=8 bits   4=16 bits

  % Test whether the machine and data byte ordering are the same 
  % It is a hack, datatype cannot be bigger than 500

  if datatype >= 512	% machine and data byte ordering is different

     fclose(fid);

     %  reread the data again

     mf = 'ieee-be';  % use little-endian byte ordering 
     fid = fopen(hfile,'r',mf);

     fseek(fid,42,'bof');
     [dims]=fread(fid,3,'int16');		

     fseek(fid,70,'bof');
     [datatype]=fread(fid,1,'int16');  % 2=8 bits   4=16 bits

  end

fclose(fid);


% read the image data 

num_voxels = prod(dims);
ifile = fullfile(fpath,[iname,'.img']);
fid = fopen(ifile,'r',mf);
  switch datatype 
     case 2, img = fread(fid,num_voxels,'uchar');
     case 4, img = fread(fid,num_voxels,'int16');
     case 8, img = fread(fid,num_voxels,'int32');
     case 16, img = fread(fid,num_voxels,'float');
     otherwise,  disp('Image data type is not supported!')
  end
fclose(fid);

img=reshape(img,[dims(1) dims(2) 1 dims(3)]);