% [xdispl, ydispl, tdel, errim] = blockMatch(im,buffer,blksz,dispbits); 
%
%
% Inputs:
% -------
%  im: current image
%  buffer: buffer of M past frames with which to block match
%    (of size 128x128xM)
%  blksz: size in pixels of blocks (blksz=16 will use 16x16 blocks)
%  dispbits: number of bits with which to encode displacement vectors
%    (sets maximum block displacement)
%
% Outputs:
% --------
%  xdispl: matrix of size numblocks x numblocks of x-displacements
%  ydispl: matrix of size numblocks x numblocks of y-displacements
%  tdel: matrix of size numblocks x numblocks of time buffer indices
%  errim: error image between original and predicted images of size(im)  
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [xdispl, ydispl, tdel, errim] = blockMatch(im,buffer,blksz,dispbits); 

errim = 0.*im;

imsize = size(im,1);
nblocks = imsize/blksz;

bufflen = size(buffer,3);

xdispl = zeros(nblocks);
ydispl = zeros(nblocks);

% find the maximum = +/- displacement in one direction in pels:
disp_max = 2^dispbits/2;

% create a vector of the range of displacements:
%disp_vec = -disp_max+1:disp_max; 
disp_vec = [0:disp_max -disp_max+1:-1];
num_disp = length(disp_vec);

for yblock = 1:nblocks,
  for xblock = 1:nblocks,

    % get vectors representing pixel positions of current block:
    ypels = (yblock-1)*blksz+1:yblock*blksz;
    xpels = (xblock-1)*blksz+1:xblock*blksz;

    % get block in current image:
    block = im(ypels,xpels); 

    for yshift_n = 1:length(disp_vec),
      for xshift_n = 1:length(disp_vec),

        % get vectors representing pixel positions of shifted block:
        sh_ypels = ypels+disp_vec(yshift_n);
        sh_xpels = xpels+disp_vec(xshift_n); 

        % Handle image edges: first find pixels outside of image:
        ypels_outside = find(sh_ypels<1 | sh_ypels>imsize);
        xpels_outside = find(sh_xpels<1 | sh_xpels>imsize);

        if (0)
        % only use the pixels within the image for the blockmatch:
        sh_ypels = sh_ypels(find(sh_ypels>0 & sh_ypels<=imsize));
        sh_xpels = sh_xpels(find(sh_xpels>0 & sh_xpels<=imsize));
        end
       
        for j = 1:bufflen, 
          if (isempty(ypels_outside) & isempty(xpels_outside))
            testblock = buffer(sh_ypels,sh_xpels,j); 
  	    err(yshift_n,xshift_n,j) = sum2(abs(block-testblock));	
          else
            % set these to something negative so that later we
            % can set them higher than the max error 
	    err(yshift_n,xshift_n,j) = -10;
          end 
        end
        
      end
    end     
    err(err<0) = err(err<0)*-(maxn(err));
 
    % find the indices of the minimum displacement
    [yind,xind,tind] = ind2sub(size(err),find(err == minn(err)));   
    % and if there are multiple displacements
    yind = yind(find(tind==min(tind)));
    xind = xind(find(tind==min(tind))); 
    tind = tind(find(tind==min(tind)));
    yind = yind(1); xind = xind(1);
    tind = tind(1);

    ydispl(yblock,xblock) = disp_vec(yind); 
    xdispl(yblock,xblock) = disp_vec(xind);  
    tdel(yblock,xblock) = tind;

    sh_ypels = ypels+disp_vec(yind);
    sh_xpels = xpels+disp_vec(xind);

    errblock = block-buffer(sh_ypels,sh_xpels,tind);
    errim(ypels,xpels) = errblock;

  end
end
 
