utilities_LW/mha_read_slice.m

   1 function S = mha_read_slice(info, num_slice)
   2 % Function for reading the volume of a Insight Meta-Image (.mha, .mhd) file
   3 %
   4 % volume = tk_read_volume(file-header)
   5 %
   6 % examples:
   7 % 1: info = mha_read_header()
   8 %    V = mha_read_volume(info);
   9 %    imshow(squeeze(V(:,:,round(end/2))),[]);
  10 %
  11 % 2: V = mha_read_volume('test.mha');
  12
  13 if(~isstruct(info)), info=mha_read_header(info); end
  14
  15
  16 switch(lower(info.DataFile))
  17     case 'local'
  18     otherwise
  19     % Seperate file
  20     info.Filename=fullfile(fileparts(info.Filename),info.DataFile);
  21 end
  22
  23 % Open file
  24 switch(info.ByteOrder(1))
  25     case ('true')
  26         fid=fopen(info.Filename,'rb','ieee-be');
  27     otherwise
  28         fid=fopen(info.Filename,'rb','ieee-le');
  29 end
  30
  31 switch(lower(info.DataFile))
  32     case 'local'
  33         % Skip header
  34         fseek(fid,info.HeaderSize,'bof');
  35     otherwise
  36         fseek(fid,0,'bof');
  37 end
  38
  39
  40 datasize=prod(info.Dimensions)*info.BitDepth/8; % in bytes
  41 slicesize=info.Dimensions(1)*info.Dimensions(2)*info.BitDepth/8;
  42 OffsetToSelectedSlice=(num_slice-1)*slicesize;
  43 fseek(fid, OffsetToSelectedSlice, 'cof');
  44
  45
  46 switch(info.CompressedData(1))
  47     case 'f'
  48         % Read the Data
  49         switch(info.DataType)
  50             case 'char'
  51                                         S = int8(fread(fid,slicesize,'char'));
  52             case 'uchar'
  53                                         S = uint8(fread(fid,slicesize,'uchar'));
  54             case 'short'
  55                                         S = int16(fread(fid,slicesize,'short'));
  56             case 'ushort'
  57                                         S = uint16(fread(fid,slicesize,'ushort'));
  58             case 'int'
  59                                         S = int32(fread(fid,slicesize,'int'));
  60             case 'uint'
  61                                         S = uint32(fread(fid,slicesize,'uint'));
  62             case 'float'
  63                                         S = single(fread(fid,info.Dimensions(1)*info.Dimensions(2),'float'));
  64             case 'double'
  65                                         S = double(fread(fid,slicesize,'double'));
  66         end
  67     case 't'
  68         switch(info.DataType)
  69             case 'char', DataType='int8';
  70             case 'uchar', DataType='uint8';
  71             case 'short', DataType='int16';
  72             case 'ushort', DataType='uint16';
  73             case 'int', DataType='int32';
  74             case 'uint', DataType='uint32';
  75             case 'float', DataType='single';
  76             case 'double', DataType='double';
  77         end
  78         Z  = fread(fid,inf,'uchar=>uint8');
  79         S = zlib_decompress(Z,DataType);
  80
  81 end
  82
  83 fclose(fid);
  84 S = reshape(S,[info.Dimensions(1), info.Dimensions(2)]);
  85
  86 function M = zlib_decompress(Z,DataType)
  87 import com.mathworks.mlwidgets.io.InterruptibleStreamCopier
  88 a=java.io.ByteArrayInputStream(Z);
  89 b=java.util.zip.InflaterInputStream(a);
  90 isc = InterruptibleStreamCopier.getInterruptibleStreamCopier;
  91 c = java.io.ByteArrayOutputStream;
  92 isc.copyStream(b,c);
  93 M=typecast(c.toByteArray,DataType);