function values=extract(voxels,input_file_base,extension) 

%EXTRACT
%
% Extracts values of files at specified voxels.
%
% VALUES = EXTRACT ( VOXELS, INPUT_FILE_BASE, [EXTENSION] )
%
% VOXELS is a 3 column matrix of voxel locations in 'register' format
%
% INOUT_FILE_BASE is a matrix of file bases
%
% EXTENSION is a file extension to be added to INPUT_FILE_BASE.
% Default is [].
%
% VALUES is a #VOXELS x #FILEs matrix of extracted values.

if nargin<3
   extension=[]
end

numvoxels=size(voxels,1)
numfiles=size(input_file_base,1)

input_file=[deblank(input_file_base(1,:)) extension];
d=fmris_read_image(input_file,0,0);
numframes=d.dim(4)
values=zeros(numvoxels,numfiles,numframes);

for file=1:numfiles
   input_file=[deblank(input_file_base(file,:)) extension];
   d=fmris_read_image(input_file,0,0);
   if file==1
      i=voxels(:,1);
      if d.vox(1)<0
         i=d.dim(1)-1-i;
      end
      j=voxels(:,2);
      if d.vox(2)<0
         j=d.dim(2)-1-j;
      end
      k=voxels(:,3);
   end
   for m=1:numvoxels
      d=fmris_read_image(input_file,k(m)+1,1:numframes);
      values(m,file,:)=squeeze(d.data(i(m)+1,j(m)+1,1,:));
   end
end