Another vectorisation question

[John Reed]

Hello,

I am looking at a nearest neighbor problem.  I have a rectangular grid and 
random coordinates that fall within the 2d grid.  I want to create an array 
that stores all the random coordinates that are nearest to each respective grid 
point.

I've taken two different approaches to solving this problem, one uses the find 
function, and the other the histc function, but in either case I don't know how 
to get rid of all the for loops and i would like to speed up the code.

Is there a way to vectorise either of these approaches?  or is there some 
totally different approach that would give me the answer I want without for 
loops?

The two different approaches I've tried are below.

Any help would be greatly appreciated,

John
..........

incr = 70;   % maximum x & y integer value on square unit grid
pt_num = 40000;
coord=incr*rand(pt_num,2);

% 1st Approach - find grid point that coordinate is closest to and assign to 
array 
% holding all coordinates closest to each respective grid point
tic;
for i = 0:incr
   for j = 0:incr
      x = [ i - 0.5; i + 0.5];
      y = [ j - 0.5; j + 0.5]; 

      indices = find( coord(:,1) >= x(1) & ...
                   coord(:,1) <  x(2) & ...         
                   coord(:,2) >= y(1) & ...
                   coord(:,2) <  y(2) );

% build array storing random coordinates indices associated with each grid 
coordinate
      num_pt_per_grid(i+1,j+1)=length(indices);
      if ( length(indices) >= 1 )
         coords_to_grid(i+1,j+1,1:length(indices)) = indices;
      end
   end
end
toc

fflush(stdout);

% 2nd approach - using histc group all coordinates closest to grid point and 
assign
% to array holding all coordinates closest to each grid point
tic;
% create edges for histogram
histo_edges = [-0.5:1:incr+.5];

[x_count, x_index] = histc(coord(:,1),histo_edges);
[y_count, y_index] = histc(coord(:,2),histo_edges);

% build array storing random coordinates indices associated with each grid 
coordinate
num_pt_per_grid2 = zeros(incr+1,incr+1);
for i=1:length(x_index)
   num_pt_per_grid2(x_index(i),y_index(i))++;
   
coords_to_grid2(x_index(i),y_index(i),num_pt_per_grid2(x_index(i),y_index(i))) 
= i;
end
toc