SurfStatReadVol1

Reads a single volumetric file in MINC, ANALYZE, NIFTI or AFNI format. 
 
 Usage: d = SurfStatReadVol1( file [, Z, T] ). 

 file = file name with extension .mnc, .img, .nii or .brik as above.
 Z    = vector of slices.
 T    = vector of times. 
 If Z and T are both 0, then just the header is read into d. 
 If Z and T are omitted, reads the entire volume.

0001 function d = SurfStatReadVol1( file, Z, T );
0002 
0003 %Reads a single volumetric file in MINC, ANALYZE, NIFTI or AFNI format.
0004 %
0005 % Usage: d = SurfStatReadVol1( file [, Z, T] ).
0006 %
0007 % file = file name with extension .mnc, .img, .nii or .brik as above.
0008 % Z    = vector of slices.
0009 % T    = vector of times.
0010 % If Z and T are both 0, then just the header is read into d.
0011 % If Z and T are omitted, reads the entire volume.
0012 
0013 if ~isstr(file)
0014     d = file;
0015     if isfield(d,'data') & nargin==3
0016         if length(size(d.data))==4 & Z~=0 & T~=0
0017             d.data=file.data(:,:,Z,T);
0018         end
0019         if length(size(d.data))==3 & length(Z)>1 & length(T)==1
0020             d.data=file.data(:,:,Z);
0021         end
0022         if length(size(d.data))==3 & length(Z)==1 & length(T)>1
0023             d.data=file.data(:,:,T);
0024         end
0025         if length(size(d.data))==3 & length(Z)==1 & length(T)==1 & Z>0 & T>0
0026             d.data=file.data(:,:,Z);
0027         end
0028     end
0029     return
0030 end
0031 
0032 file = deblank(file);
0033    
0034 [path,name,ext]=fileparts(file);
0035 if strcmp(ext,'.gz')
0036    if isunix
0037       unix(['gunzip ' file]);
0038    else
0039       ['Can''t gunzip on non-unix system']
0040       return
0041    end
0042    [path,name,ext]=fileparts(name)
0043 end 
0044 
0045 switch lower(ext)
0046 case '.mnc'
0047    if nargin == 3
0048       d = fmris_read_minc(file,Z,T);
0049    else
0050       d = fmris_read_minc(file);
0051    end
0052 case '.img'
0053    if nargin == 3
0054       d = fmris_read_analyze(file,Z,T);
0055    else
0056       d = fmris_read_analyze(file);
0057    end
0058    if isfield(d,'data') d.data(isnan(d.data))=0; end;
0059    d.origin=-d.origin.*d.vox(1:3);
0060 case '.brik'
0061    if nargin == 3
0062       d = fmris_read_afni(file,Z,T);
0063    else
0064       d = fmris_read_afni(file);
0065    end
0066 case '.nii'
0067    if nargin == 3
0068       d = fmris_read_nifti(file,Z,T);
0069    else
0070       d = fmris_read_nifti(file);
0071    end
0072 otherwise
0073    ['Unknown file extension']
0074 end
0075 
0076 d.parent_file=file;
0077 
0078 return
0079 end
0080 
0081 
0082 
0083 %%
0084 
0085 function [d]=fmris_read_afni(file,Z,T)
0086 
0087 [err, Info]=BrikInfo(file);
0088 
0089 d.dim     = [Info.DATASET_DIMENSIONS(1:3) Info.DATASET_RANK(2)];
0090 
0091 if nargin == 1 
0092    Z = 1:d.dim(3);
0093    T = 1:d.dim(4);
0094 end
0095 
0096 if (T(1)~=0)&(Z(1)~=0)
0097    Opt.Slices=Z;
0098    Opt.Frames=T;
0099    [err, d.data, Info, ErrMessage] = BrikLoad(file, Opt);
0100    d.calib        = [min(min(min(min(d.data)))) max(max(max(max(d.data))))];
0101 end
0102 
0103 d.vox                     = Info.DELTA;
0104 d.vox_units                = '';
0105 d.vox_offset            = 0;
0106 d.precision                = '';
0107 d.calib_units            = '';
0108 d.origin                 = Info.ORIGIN;
0109 d.descrip                = '';
0110 
0111 if isfield(Info,'WORSLEY_DF')
0112    df=Info.WORSLEY_DF;
0113 else
0114    df=[];
0115 end
0116 if ~isempty(df)
0117    d.df=df;
0118 end
0119 
0120 if isfield(Info,'WORSLEY_NCONJ')
0121    nconj=Info.WORSLEY_NCONJ;
0122 else
0123    nconj=[];
0124 end
0125 if ~isempty(nconj)
0126    d.nconj=nconj;
0127 end
0128 
0129 if isfield(Info,'WORSLEY_FWHM')
0130    fwhm=Info.WORSLEY_FWHM;
0131 else
0132    fwhm=[];
0133 end
0134 if ~isempty(fwhm)
0135    d.fwhm=fwhm;
0136 end
0137 
0138 return;
0139 end
0140 
0141 
0142 
0143 %%
0144 
0145 function    d=fmris_read_analyze(file,Z,T);
0146 
0147 % Read in 1 4D or >1 3D Analyze format Image Volumes
0148 %
0149 % Usage
0150 %
0151 % d=Read_Analyze(file);
0152 % Reads in all image data and attributes into the structure d.
0153 %
0154 % d=Read_Analyze(file,Z,T);
0155 % Reads in chosen planes and frames and attributes into the structure d.
0156 % Z and T are vectors e.g. Z=1:31;T=1:12; or Z=[24 26];T=[1 3 5];
0157 %
0158 % (c) Roger Gunn & John Aston
0159 
0160 % Get machine format on which the file was written:
0161 
0162 machineformat=getmachineformat(file);
0163 
0164 if ~isempty(machineformat)
0165    
0166    % Read Header Information
0167    d=Read_Analyze_Hdr(file,machineformat);
0168    d.file_name=file;
0169    
0170    % try to get file precision if it is unknown:
0171    if strcmp(d.precision,'Unknown')
0172       d.precision=getprecision(deblank(file),machineformat,d.dim,d.global);
0173    end
0174    
0175    % Read Image Data
0176    
0177    if ~isempty(d) & ~strcmp(d.precision,'Unknown') 
0178       
0179       if nargin==1 | strcmp(d.precision,'uint1') 
0180          
0181          % Read in Whole Analyze Volume
0182          fid = fopen(d.file_name,'r',machineformat);
0183          if fid > -1
0184             d.data=d.scale*reshape(fread(fid,prod(d.dim),d.precision),d.dim(1),d.dim(2),d.dim(3),d.dim(4));
0185             fclose(fid);
0186          else
0187             errordlg('Check Image File: Existence, Permissions ?','Read Error'); 
0188          end;
0189          
0190          if nargin==3
0191             if all(Z>0)&all(Z<=d.dim(3))&all(T>0)&all(T<=d.dim(4))
0192                d.data=d.data(:,:,Z,T);
0193                d.Z=Z;
0194                d.T=T;
0195             else
0196                errordlg('Incompatible Matrix Identifiers !','Read Error');  
0197             end
0198          end
0199          
0200       elseif nargin==3
0201          % Read in Chosen Planes and Frames
0202          if (T(1)~=0)|(Z(1)~=0)
0203             
0204             if all(Z>0)&all(Z<=d.dim(3))&all(T>0)&all(T<=d.dim(4))
0205                
0206                fid = fopen(d.file_name,'r',machineformat);
0207                if fid > -1
0208                   d.data=zeros(d.dim(1),d.dim(2),length(Z),length(T));
0209                   for t=1:length(T)
0210                      for z=1:length(Z)
0211                         status=fseek(fid,d.hdr.byte*((T(t)-1)*prod(d.dim(1:3))+(Z(z)-1)*prod(d.dim(1:2))),'bof');
0212                         d.data(:,:,z,t)=d.scale*fread(fid,[d.dim(1) d.dim(2)],d.precision);
0213                      end
0214                   end
0215                   d.Z=Z;
0216                   d.T=T;
0217                   fclose(fid);
0218                else
0219                   errordlg('Check Image File: Existence, Permissions ?','Read Error'); 
0220                end;
0221                
0222             else
0223                errordlg('Incompatible Matrix Identifiers !','Read Error'); 
0224             end;
0225             
0226          end;
0227       else
0228          errordlg('Unusual Number of Arguments','Read Error');
0229       end;
0230    else
0231       if strcmp(d.precision,'Unknown');
0232          errordlg('Unknown Data Type (Precision?)','Read Error');
0233       end
0234    end;
0235 else
0236    errordlg('Unknown Machine Format','Read Error'); 
0237 end
0238 
0239 % if there is no slice thickness, set it to 6mm:
0240 
0241 if d.vox(3)==0
0242    d.vox(3)=6;
0243 end
0244 
0245 return;
0246 end
0247 
0248 %%
0249 
0250 function machineformat=getmachineformat(file);
0251 
0252 % Get machine format by reading the d.hdr.dim(1) attribute and
0253 % making sure it is 1, 2, 3 or 4.
0254 
0255 machineformat=[];
0256 for mf='nlbdgcas'
0257    fid = fopen([file(1:(length(file)-3)) 'hdr'],'r',mf);
0258    if fid > -1
0259       fseek(fid,40,'bof');
0260       if any(fread(fid,1,'int16')==1:4)
0261          machineformat=mf;
0262          fclose(fid);
0263          break
0264       else
0265          fclose(fid);
0266       end
0267    else
0268       errordlg('Check Header File: Existence, Permissions ?','Read Error'); 
0269       break
0270    end
0271 end
0272 
0273 return
0274 end
0275 
0276 %%
0277 
0278 function precision=getprecision(file,machineformat,dim,range);
0279 
0280 % Get precision by reading a value from the middle of the .img file and
0281 % making sure it is within the global attribute
0282 
0283 precisions=['int8   '   
0284         'int16  '  
0285         'int32  '  
0286         'int64  '   
0287         'uint8  '   
0288         'uint16 '  
0289         'uint32 '  
0290         'uint64 ' 
0291         'single ' 
0292         'float32' 
0293         'double ' 
0294         'float64' ];
0295 nbytes=[1 2 4 8 1 2 4 8 4 4 8 8];
0296 middle_vol=dim(1)*dim(2)*floor(dim(3)/2)+dim(1)*round(dim(2)/2)+round(dim(1)/2);
0297 h=dir(file);
0298 n=ceil(h.bytes/prod(dim));
0299 
0300 fid = fopen(file,'r',machineformat);
0301 if fid > -1
0302    for i=1:size(precisions,1)
0303       if nbytes(i)==n
0304          status=fseek(fid,middle_vol*n,'bof');
0305          if status==0
0306             precision=deblank(precisions(i,:));
0307             x=fread(fid,10,precision);
0308             if all(range(1)<=x) & all(x<=range(2))
0309                return
0310             end
0311          end
0312       end
0313    end
0314 end
0315 errordlg('Check Header File: Existence, Permissions ?','Read Error'); 
0316 precision='Unknown';
0317 
0318 return
0319 end
0320 
0321 %%
0322 
0323 function d=Read_Analyze_Hdr(file,machineformat);
0324 
0325 % Read Analyze Header information into the structure d
0326 % Adapted from John Ashburners spm_hread.m
0327 
0328 fid  = fopen([file(1:(length(file)-3)) 'hdr'],'r',machineformat);
0329 
0330 if fid > -1
0331    
0332    % read (struct) header_key
0333    %---------------------------------------------------------------------------
0334    fseek(fid,0,'bof');
0335    
0336    d.hdr.sizeof_hdr         = fread(fid,1,'int32');
0337    d.hdr.data_type          = deblank(setstr(fread(fid,10,'char'))');
0338    d.hdr.db_name            = deblank(setstr(fread(fid,18,'char'))');
0339    d.hdr.extents            = fread(fid,1,'int32');
0340    d.hdr.session_error       = fread(fid,1,'int16');
0341    d.hdr.regular            = deblank(setstr(fread(fid,1,'char'))');
0342    d.hdr.hkey_un0            = deblank(setstr(fread(fid,1,'char'))');
0343    
0344    
0345    
0346    % read (struct) image_dimension
0347    %---------------------------------------------------------------------------
0348    fseek(fid,40,'bof');
0349    
0350    d.hdr.dim                = fread(fid,8,'int16');
0351 
0352    d.hdr.vox_units            = deblank(setstr(fread(fid,4,'char'))');
0353    d.hdr.cal_units            = deblank(setstr(fread(fid,8,'char'))');
0354    d.hdr.unused1            = fread(fid,1,'int16');
0355    d.hdr.datatype            = fread(fid,1,'int16');
0356    d.hdr.bitpix                = fread(fid,1,'int16');
0357    d.hdr.dim_un0            = fread(fid,1,'int16');
0358    d.hdr.pixdim                = fread(fid,8,'float');
0359    d.hdr.vox_offset            = fread(fid,1,'float');
0360    d.hdr.funused1            = fread(fid,1,'float');
0361    d.hdr.funused2            = fread(fid,1,'float');
0362    d.hdr.funused3            = fread(fid,1,'float');
0363    d.hdr.cal_max            = fread(fid,1,'float');
0364    d.hdr.cal_min            = fread(fid,1,'float');
0365    d.hdr.compressed            = fread(fid,1,'int32');
0366    d.hdr.verified            = fread(fid,1,'int32');
0367    d.hdr.glmax                = fread(fid,1,'int32');
0368    d.hdr.glmin                = fread(fid,1,'int32');
0369    
0370    % read (struct) data_history
0371    %---------------------------------------------------------------------------
0372    fseek(fid,148,'bof');
0373    
0374    d.hdr.descrip            = deblank(setstr(fread(fid,80,'char'))');
0375    d.hdr.aux_file            = deblank(setstr(fread(fid,24,'char'))');
0376    d.hdr.orient                = fread(fid,1,'char');
0377    d.hdr.origin                = fread(fid,5,'uint16');
0378    d.hdr.generated            = deblank(setstr(fread(fid,10,'char'))');
0379    d.hdr.scannum            = deblank(setstr(fread(fid,10,'char'))');
0380    d.hdr.patient_id            = deblank(setstr(fread(fid,10,'char'))');
0381    d.hdr.exp_date            = deblank(setstr(fread(fid,10,'char'))');
0382    d.hdr.exp_time            = deblank(setstr(fread(fid,10,'char'))');
0383    d.hdr.hist_un0            = deblank(setstr(fread(fid,3,'char'))');
0384    d.hdr.views                = fread(fid,1,'int32');
0385    d.hdr.vols_added            = fread(fid,1,'int32');
0386    d.hdr.start_field        = fread(fid,1,'int32');
0387    d.hdr.field_skip            = fread(fid,1,'int32');
0388    d.hdr.omax                = fread(fid,1,'int32');
0389    d.hdr.omin                = fread(fid,1,'int32');
0390    d.hdr.smax                = fread(fid,1,'int32');
0391    d.hdr.smin                = fread(fid,1,'int32');
0392    
0393    fclose(fid);
0394    
0395    % Put important information in main structure
0396    %---------------------------------------------------------------------------
0397    
0398    d.dim                      = d.hdr.dim(2:5)';
0399    vox                         = d.hdr.pixdim(2:5)';
0400    if     vox(4)==0 
0401       vox(4)=[];
0402    end
0403    d.vox                       = vox;
0404    d.vox_units               = d.hdr.vox_units;
0405    d.vox_offset                = d.hdr.vox_offset;
0406    scale                     = d.hdr.funused1;
0407    d.scale                       = ~scale + scale;
0408    d.global                    = [d.hdr.glmin d.hdr.glmax];
0409    d.calib                    = [d.hdr.cal_min d.hdr.cal_max];
0410    d.calib_units            = d.hdr.cal_units;
0411    d.origin                    = d.hdr.origin(1:3)';
0412    d.descrip                   = d.hdr.descrip(1:max(find(d.hdr.descrip)));
0413    
0414    if ~isnan(str2double(d.hdr.scannum));    d.df(1)  =str2double(d.hdr.scannum);    end
0415    if ~isnan(str2double(d.hdr.patient_id)); d.df(2)  =str2double(d.hdr.patient_id); end
0416    if ~isnan(str2double(d.hdr.exp_date));   d.fwhm(1)=str2double(d.hdr.exp_date);   end
0417    if ~isnan(str2double(d.hdr.exp_time));   d.fwhm(2)=str2double(d.hdr.exp_time);   end
0418    if ~isnan(str2double(d.hdr.hist_un0));   d.nconj  =str2double(d.hdr.hist_un0);   end
0419    
0420    switch d.hdr.datatype
0421    case 1
0422       d.precision     = 'uint1';
0423       d.hdr.byte     = 0;
0424    case 2
0425       d.precision     = 'uint8';
0426       d.hdr.byte     = 1;
0427    case 4
0428       d.precision     = 'int16';
0429       d.hdr.byte     = 2;
0430    case 8
0431       d.precision     = 'int32';
0432       d.hdr.byte     = 4;
0433    case 16
0434       d.precision     = 'float';
0435       d.hdr.byte     = 4;
0436    case 64
0437       d.precision     = 'double';
0438       d.hdr.byte     = 8;
0439    otherwise
0440       d.precision     = 'Unknown';
0441       d.hdr.byte     = 0;
0442    end
0443    
0444 else
0445    d=[];
0446    errordlg('Check Header File: Existence, Permissions ?','Read Error'); 
0447 end
0448 
0449 return
0450 end
0451 
0452 
0453 %%
0454 
0455 function [d]=fmris_read_minc(file,Z,T)
0456 
0457 fid=fopen(file);
0458 file=fopen(fid);
0459 fclose(fid);
0460 
0461 d.file_name = file;
0462 
0463 d.dim     = fliplr(miinquire(file,'imagesize')');
0464 d.dim     = d.dim + (d.dim == 0);
0465 
0466 if nargin == 1 
0467    Z = 1:d.dim(3);
0468    T = 1:d.dim(4);
0469 end
0470 
0471 if (T~=0)&(Z~=0)
0472    if d.dim(4)==1
0473       images     = mireadimages(file,Z-1);
0474    else
0475       n=length(T);
0476       images=[];
0477       for i=0:floor((n-1)/256)
0478          frames=T((i*256+1):min((i+1)*256,n));
0479          images     = [images mireadimages(file,Z-1,frames-1)];
0480       end
0481    end
0482    d.data         = reshape(images,[d.dim(1:2) length(Z) length(T)]);
0483    d.calib        = [min(min(min(min(d.data)))) max(max(max(max(d.data))))];
0484 end
0485 
0486 [x_vox, y_vox, z_vox] =  miinquire(file,'attvalue','xspace','step','attvalue','yspace','step','attvalue','zspace','step');
0487 d.vox                     = [x_vox y_vox z_vox];
0488 d.vox_units                = '';
0489 d.vox_offset            = 0;
0490 d.precision                = '';
0491 d.calib_units            = '';
0492 [x_origin, y_origin, z_origin] = miinquire(file,'attvalue','xspace','start','attvalue','yspace','start','attvalue','zspace','start');
0493 d.origin                 = [x_origin y_origin z_origin];
0494 d.descrip                = '';
0495 
0496 df=miinquire(file, 'attvalue', 'df', 'df');
0497 if ~isempty(df)
0498    d.df=df;
0499 end
0500 
0501 nconj=miinquire(file, 'attvalue', 'nconj', 'nconj');
0502 if ~isempty(nconj)
0503    d.nconj=nconj;
0504 end
0505 
0506 fwhm=miinquire(file, 'attvalue', 'fwhm', 'fwhm');
0507 if ~isempty(fwhm)
0508    d.fwhm=fwhm;
0509 end
0510 
0511 return;
0512 end
0513 
0514 
0515 %%
0516 
0517 function [d]=fmris_read_nifti(file,Z,T)
0518 
0519 d.file_name=file;
0520 machineformats='nlbdgcas';
0521 for i=1:length(machineformats)
0522   fid=fopen(file,'r',machineformats(i));
0523   if fid == -1;
0524      fprintf('Invalid file %s',file)
0525   end
0526   sizeof_hdr=fread(fid,1,'int');
0527   if sizeof_hdr==348;
0528      %Must have found a formt that works, so break and continue using this format
0529      break;
0530   else
0531     %Need to close if file is not native format
0532     %else if the file format is 's', then 7 stranded files are left orphaned and never closed.
0533     fclose(fid);
0534   end
0535 end
0536 data_type=fread(fid,10,'char');
0537 db_name=fread(fid,18,'char');
0538 extents=fread(fid,1,'int');
0539 session_error=fread(fid,1,'short');
0540 regular=char(fread(fid,1,'char')');
0541 dim_info=char(fread(fid,1,'char')');
0542 dim=fread(fid,8,'short');
0543 intent_p =fread(fid,2,'float');
0544 intent_q =fread(fid,2,'uint16');
0545 intent_code =fread(fid,1,'short');
0546 datatype=fread(fid,1,'short');
0547 bitpix=fread(fid,1,'short');
0548 slice_start=fread(fid,1,'short');
0549 pixdim=fread(fid,8,'float');
0550 vox_offset=fread(fid,1,'float');
0551 scl_slope =fread(fid,1,'float');
0552 scl_inter =fread(fid,1,'float');
0553 slice_end=fread(fid,1,'short');
0554 slice_code =char(fread(fid,1,'char')');
0555 xyzt_units =char(fread(fid,1,'char')');
0556 cal_max=fread(fid,1,'float');
0557 cal_min=fread(fid,1,'float');
0558 slice_duration=fread(fid,1,'float');
0559 toffset=fread(fid,1,'float');
0560 glmax=fread(fid,1,'int');
0561 glmin=fread(fid,1,'int');
0562 descrip=char(fread(fid,80,'char')');
0563 aux_file=char(fread(fid,24,'char')');
0564 qform_code =fread(fid,1,'short');
0565 sform_code =fread(fid,1,'short');
0566 quatern_b =fread(fid,1,'float');
0567 quatern_c =fread(fid,1,'float');
0568 quatern_d =fread(fid,1,'float');
0569 qoffset_x =fread(fid,1,'float');
0570 qoffset_y =fread(fid,1,'float');
0571 qoffset_z =fread(fid,1,'float');
0572 srow_x =fread(fid,4,'float');
0573 srow_y =fread(fid,4,'float');
0574 srow_z =fread(fid,4,'float');
0575 intent_name=char(fread(fid,16,'char')');
0576 magic =char(fread(fid,4,'char')');
0577 
0578 d.machineformat=machineformats(i);
0579 d.dim=ones(1,4);
0580 if dim(1)==5
0581    if dim(5)==1 
0582       dim=[4; dim(2:4); dim(6); zeros(3,1)]
0583    else
0584       dim
0585       fclose(fid);
0586       return
0587    end
0588 end
0589 d.dim(1:dim(1))=dim((1:dim(1))+1);
0590 d.vox=zeros(1,3);
0591 d.vox(1:dim(1))=pixdim((1:dim(1))+1);
0592 %Attempt to fill out more information for a complete nifti description
0593   d.vox_offset  = vox_offset;
0594   d.scale       = scl_slope;
0595   d.intercept   = scl_inter;
0596   d.global      = [glmin glmax];
0597   d.calib       = [cal_min cal_max];
0598   if qform_code>0;
0599     d.origin    = [qoffset_x qoffset_y qoffset_z];
0600   else
0601     d.origin    = [srow_x(4) srow_y(4) srow_z(4)];
0602   end
0603   d.descrip     = descrip;
0604 d.parent_file=file;
0605 
0606 if intent_p(1)>0; d.df(1)=intent_p(1); end;
0607 if intent_p(2)>0; d.df(2)=intent_p(2); end;
0608 
0609 intent_q=intent_q/2^16*100;
0610 if intent_q(1)>0; d.fwhm(1)=intent_q(1); end;
0611 if intent_q(2)>0; d.fwhm(2)=intent_q(2); end;
0612 
0613 if nargin<2
0614    Z=1:d.dim(3);
0615    T=1:d.dim(4);
0616 end
0617 
0618 if Z(1)==0 & T(1)==0
0619    fclose(fid);
0620    return
0621 end
0622 
0623 types=lower(['UINT8  ';'INT16  ';'INT32  ';'FLOAT32';'FLOAT64']);
0624 codes=[2; 4; 8; 16;  64];
0625 d.precision=[];
0626 for i=1:5
0627    if datatype==codes(i)
0628       d.precision=deblank(types(i,:));
0629    end
0630 end
0631 if isempty(d.precision)
0632    unknown_datatype=datatype
0633    fclose(fid);
0634    return
0635 end
0636 
0637 %d.byte=bitpix/8
0638 d.byte=4;
0639 if all(Z>0)&all(Z<=d.dim(3))&all(T>0)&all(T<=d.dim(4))
0640    d.data=zeros(d.dim(1),d.dim(2),length(Z),length(T));
0641    for t=1:length(T)
0642       for z=1:length(Z)
0643          position=d.byte*((T(t)-1)*prod(d.dim(1:3))+(Z(z)-1)*prod(d.dim(1:2)))+vox_offset;
0644          status=fseek(fid,position,'bof');
0645          d.data(:,:,z,t)=fread(fid,[d.dim(1) d.dim(2)],d.precision);
0646       end
0647    end
0648   fclose(fid);
0649 else
0650    Z
0651    T
0652    fclose(fid);
0653    return
0654 end
0655 
0656 if scl_slope~=0
0657    d.data=d.data*scl_slope+scl_inter;
0658 end
0659 
0660 
0661 return;
0662 end
0663