function b = resamplenn( a, x, y )
%
%   b = resamplenn( a, x, y ) - resample  a at points given in x and y
%
%   a   -   input image
%   x   -   image of x coordinates to sample, same size as y
%   y   -   image of y coordinates to sample, same size as x
%   b   -   resampled image - same size as x and y
%

% get size of a then reshape
[ma,na] = size(a);
a = [a(:); 0];      % zero element for out of bounds
obounds = ma * na + 1;

% get size of b from x
[mb,nb] = size(x);
b = zeros(mb,nb);

% find grid coords from source image and offset of resample points
j = round(y(:));
k = round(x(:));
    
% get masks for out of bounds
mask1 = (1 <=  j    &  j    <= ma) & (1 <=  k    &  k    <= na);

% indices for points to use for bilinear interpolation
indx1 = j   + ma*(k-1);
    
% replace indices for out of bounds with index for zero element in 
indx1 = mask1 .* indx1 + ~mask1 * obounds;

% resample the image
b = a(indx1);
b = reshape( b, mb, nb );
