5
$\begingroup$

I need solve the task of classification of the points regarding their position. So I use such a function:

selector[mas_, reg_] := Select[mas, 
                      RegionMember[
                           Rectangle[{First@reg - st/2, Last@reg - st/2}, {First@reg + st/2, Last@reg + st/2}], 
                           #[[2]]] &];

It checks if second elements of mas entries are inside the rectangle defined by point reg. Finally it returns the list of right mas entries. There are ~4000 rectangles with ~250000 points and it evidently works too long.

I try to make the compiled version of such function:

  cs=Compile[{mas,_Real,4},{reg,_Real,2},
       Select[mas, 
             RegionMember[
                 Rectangle[{First@reg - st/2, Last@reg - st/2}, {First@reg + st/2, Last@reg + st/2}],
                 #[[2]]] &]

]

This function produce an error like "The first argument should be a 4 rank tensor". My array has a complicated form like mas={{1,{123.4,654.3},{12.1,21.2}},{...},...}. Of course it is not a good matrix but what is a problem with an arbitrary shape array? I was trying make the packed array from the initial one but it does not help.

The question is how to compile a function with arbitrary shape array of Real arguments?

$\endgroup$
3
  • $\begingroup$ What are reg and st? $\endgroup$ Jan 30, 2022 at 15:19
  • $\begingroup$ Do all the rectangles have the same size? $\endgroup$ Jan 30, 2022 at 15:26
  • $\begingroup$ @HenrikSchumacher, st is an arbitrary number - the rectangle size. reg - is an arbitrary point (pair of numbers) that defines the position of the rectangle. Yes, the sizes are equal. $\endgroup$
    – Rom38
    Jan 31, 2022 at 3:28

1 Answer 1

5
$\begingroup$

You should not try to compile the brute-force algorithm. Instead you may use hierarchical data structures that allow the lookups in nearly linear time.

First read off all points pts= mas[[All,2]]. Then use Nearest to do the lookup (it employs a k-d-Tree). Simulating some data:

m = 4000;
n = 4000;
st = 0.5;
pts = RandomReal[{-1, 1}, {m, 2}];
rect = RandomReal[{-1 + st, 1 - st}, {n, 2}];
result = Nearest[pts -> "Index", rect, {\[Infinity], st/2}, DistanceFunction -> ChessboardDistance];

Now result[[k]] contains all the indices of the points in the k-th rectangle (and hence the indixec for the rows of your original array).

Here a test for correctness:

k = 13;

Graphics[{
  Point[pts],
  EdgeForm[Red], FaceForm[],
  Rectangle[rect[[k]] - st/2, rect[[k]] + st/2],
  Green, Point[pts[[result[[k]]]]]
  }]

enter image description here

For m = 250000 this takes two seconds on my machine.

$\endgroup$
2
  • $\begingroup$ It is very interesting approach. Thanks, will try! $\endgroup$
    – Rom38
    Jan 31, 2022 at 3:33
  • $\begingroup$ You're welcome! Btw., if you want to apply this to general axis-aligned rectangles then you merely have to scale the coodinates (by multiplying pts and rect with an appropriate diagonal matrix from the right) before calling Nearest. $\endgroup$ Jan 31, 2022 at 8:02

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.