0

I have a set of {x, y} coordinates, for example:

xmin = 0;
xmax = 100;
ymin = 42;
ymax = 76;
numPoints = 10^5;
pointList = Table[{RandomReal[{xmin, xmax}], RandomReal[{ymin, ymax}]}, {numPoints}];

I would like to break up this set of points into N equal-area squares in the most efficient manner possible, and without binning a point into more than one squares, such that we return a list like the following:

decomposedPointList = 
  {{...Points In Box 1...},{...Points In Box 2...},...,{...Points In Box N...};

Where the intersection of any two sublists of points is always zero.

Is there a nice way to do this?

m_goldberg
  • 107,779
  • 16
  • 103
  • 257
T.T.
  • 353
  • 3
  • 6
  • Did you really mean to write pointList=Table[{RandomReal[{xmin,xmax}],RandomReal[ymin,ymax]},{i,1,numPoints}]; ? Should this not be pointList = Table[{RandomReal[{xmin, xmax}], RandomReal[{ymin, ymax}]}, {numPoints}]; so that the second entry is a point and not a list? also, your question is not clear, at least to me. – Nasser Aug 14 '13 at 03:36
  • @Nasser Thanks for that catch, I rewrote the line of code. – T.T. Aug 14 '13 at 03:38
  • @Nasser Basically I want to take my point like, and partition points into sublists based on these points falling into non-overlapping boxes (where we specify either the size of the boxes or the number of boxes). – T.T. Aug 14 '13 at 03:39
  • 1
    If I understand correctly this is a duplicate of Partitioning a 2D... – Kuba Aug 14 '13 at 04:50
  • @gpap I agree, I've posted my earlier than you have provided this link and I wasn't paying attention, why haven't you flagged this question as duplicate then? – Kuba Aug 14 '13 at 12:09
  • I thought I had. Done so now - not criticising you answering it, it's just if a user doesn't know what binning is, the question gets a variety of re-phrasings and it gets more difficult to find an answer – gpap Aug 14 '13 at 12:18
  • @gpap I understand, I also do not mind when people answer questions which are duplicates, it may be more helpful for OP, it's just those questions should be closed fast so later there is more handy database left. – Kuba Aug 14 '13 at 13:05

1 Answers1

0

This may not scale well but perhaps can be a start.

f[x1_, x2_, y1_, y2_, size_, pt_] := Module[
  {nx, ny, tab, np, dist},
  nx = (x2 - x1)/size;
  ny = (y2 - y1)/size;
  tab = Flatten[
    Table[{x1, y1} + 0.5 {i size, j size}, {i, 1, 2 nx, 2}, {j, 1, 
      2 ny, 2}], 1];
  np = First@Nearest[tab, pt];
  dist = EuclideanDistance[np, pt];
  If[dist <= size/2, {np, size},
   If[
    IntervalMemberQ[
      Interval[{np[[1]] - 0.5  size, np[[1]] + 0.5 size}], pt[[1]]] &&


     IntervalMemberQ[
      Interval[{np[[2]] - 0.5  size, np[[2]] + 0.5 size}], 
      pt[[2]]], {np, size}, {"not in region of interest",""}]
   ]]

Essentially partitions rectangular interest into sqaures of side 'size' and tabulates the centre. The final argument is the test point. There will be a unique nearest center. If the distance is size/2 then it is clearly in this square. The corner regions are dealt with intervals. Points outside region of interest will be coded.

I tested only on small test set 100 points of varying square sizes and seems to work but may not scale well. For example:

test = Thread[{RandomReal[{0, 100}, 100], RandomReal[{42, 76}, 100]}];
{StringForm[
     "Square centred at `1` of size `2`:", #[[1, 1]], #[[1, 
       2]]], #[[2]]} & /@ 
  Tally[f[0, 100, 42, 76, 1, #] & /@ test] // Grid
ubpdqn
  • 60,617
  • 3
  • 59
  • 148