Polygon inside polygon inside polygon

Question

I have number of simple polygons which do not intersect each other but some polygons may be embedded into others.

For example:

+--------------------------------------------+
|                                            |
|   +----------------+          +--------+   |
|   |                |         /         |   |
|   |   +--------+   |        /          |   |
|   |   |        |   |       /  +----(2)-+   |
|   |   |        |   |      /   |            |
|   |   +----(3)-+   |     /    |   +---+    |
|   |                |    +-----+   |   |    |
|   +------------(2)-+              +(2)+    |
|                                            |
+----------------------------------------(1)-+

How to find out the "depth" of all the polygons? In other words, how to find out by how many polygons a polygon is encompassed by? The "depth" are the numbers in parentheses.

I could count how many times a point of a polygon is inside of all the other polygons but this has quadratic complexity. How to compute those depths faster?

Answer 1

Put your polygons into some kind of spatial lookup structure, for example an R-tree based on the minimum bounding rectangles for the polygons. Then you should be able to compute the containment relation that you're after in O(n log n). (Unless you're in a pathological case where many of the minimum bounding rectangles for your polygons overlap, which seems unlikely based on your description.)

Edited to add: of course, you don't rely on the R-tree to tell you if one polygon is inside another (it's based on minimum bounding rectangles so it can only give you an approximation). You use the R-tree to cheaply identify candidate inclusions which you then verify in the expensive way (checking that a point in one polygon is inside the other).

Answer 2

In your example you depict a non-convex polygon so the rtree approach that others suggested will not be effective by itself. But you could couple that with an additional check to see if a vertex of the polygon is inside the other whenever the you find a match with the rtree approach. This would reduce the number of times you have to do the "point in polygon" check.

Another approach is to take your first idea - i.e. to check a vertex of each polygon with everyother polygon - and modify it to cache the results that you already computed and reuse them to reduce time complexity. The approach desribed below essentially builds a tree structure out of the polygons themselves, rather than use an rtree. Its essentially a different tree structure that honors non-convex polygons.

1. Initially each polygon is a root node and in a tree of its own.
2. You have to loop through your set of root nodes and check if its inside another root node using the "point in polygon" test.
3. When you find a match, remove the smaller polygon from the set of root nodes and add it as a child node to the bigger polygon. You started building a meaningful tree structure and you are also reducing the size of the container you are iterating over (root nodes) so the future iterations will be faster.
4. Continue this until the number of root nodes stops changing.
5. Remember in the main nested loop you are only comparing if one root node is inside another. When you find a match (small node inside a big node) you can skip the subtree of the smaller node. But you do have to traverse the subtree of the bigger node and do the checks to find the right place for the smaller node in the hierarchy. But still you avoid a lot of checks compared to a simple nested loop with quadratic complexity.

Answer 3

(This approach follows a similar idea to @GarethRees 's: first, cheaply discover which pairs of polygons you don't need to check for inclusion. Once the number of pairs still needed to check is acceptable, do the exact (expensive) geometric check.)

It's easy to calculate for each polygon p the bounding rectangle, i.e. the left, right, top and bottom, so let's do that first. E.g. for left: p.L = min { x : (x,y) is a vertex of p } Time is linear in the number of points.

To avoid having to compare each polygon to each other one, you could divide up the area into a 2x2 grid. Suppose the width and height of the area are respectively given by Dx and Dy, then you can create nine sets top,bottom,left,right,topright,topleft,bottomright,bottomleft,rest and do the following:

for p in polygons:
    in_top    = p.B > Dy/2
    in_bottom = p.T < Dy/2
    in_left   = p.R < Dx/2
    in_right  = p.L > Dx/2 
    if in_top:
        if in_left:
            add p to topleft
        elsif in_right:
            add p to topright
        else:
            add p to top
    elsif in_bottom:
        if in_left:
            add p to bottomleft
        elsif in_right:
            add p to bottomright
        else:
            add p to bottom

    if in_right and not (in_top or in_bottom):
        add p to right
    elsif in_left and not (in_top or in_bottom):
        add p to left

    if not (in_top or in_bottom or in_left or in_right):
        add p to rest

This is again linear time. Each polygon has been binned into its most "tightly" containing sector. What have you gained by this? Well, you know for example that for any polygon p in left there can't possibly be any inclusion relationship with set right, so you don't need to compare them. Likewise between bottomleft and right, bottomleft and topleft, and so on. Here is what it would look like on your example:

                      Dx/2
+----------------------|---------------------+
|                      |                     |
|   +----------------+ |        +--------+   |
|   |                | |       /         |   |
|   |   +--------+   | |      /          |   |
|___|___|________|___|_|____ /__+===d(2)=+___|_ Dy/2
|   |   |        |   | |    /   |            |
|   |   +---b(3)-+   | |   /    |   +---+    |
|   |                | |  +-----+   |   |    |
|   +-----------c(2)-+ |            e(2)+    |
|                      |                     |
+----------------------|----------------a(1)-+

So rest = {a}, top = {}, bottom = {}, left = {b,c}, right = {d}

topleft = {}, topright = {}, bottomleft = {}, bottomright = {e}

So basically now you need to compare (with the expensive exact check) at most b to c, d to e, and a to all others -- in fact, if you order the checks in a clever way, you won't need to compare a to all others, because inclusion is transitive, so if you notice that c includes b, and a includes c, then you do not need to check if a includes b.

Another point is that you can apply the above reasoning recursively. Say for example the set topright is too big for your taste; you can then apply the same technique by further dividing up that subregion (just need to get the book-keeping right).

Answer 4

Seems to me you get get away with sorting the polygons, using a test of whether one is inside another as the comparison operator.

Step 1

Suppose we define the relation '<' between polygons as follows: A < B iff A is inside B. It so happens that if A < B and B < C, then A < C (i.e. if polygon A is inside B and B is inside C, then A must be inside C). Now, we have a strict weak ordering between arbitrary polygons.

[Edit: You'd need to use some sort of point-inside-non-convex-polygon test, but presumably you are doing that already.]

Step 2

Sort the polygons according to this comparison using your favorite comparison-based sorting algorithm. For instance, merge sort has a worst-case time complexity of O(nlogn) comparisons where n is the number of polygons.

[Edit: This is the important step, because it gets rid of the quadratic complexity.]

Step 3

Ensure that the "largest" (i.e. outermost) element is first on your sorted array of polygons. (Reverse the list if necessary to achieve this - it is linear on the number of polygons).

Step 4

Now the "largest" (i.e. outermost) polygon should be the first element.

[Edit: In fact, the polygons have been sorted according to their depth. However, two polygons that have the same depth may appear in different orders depending on whether the sort was stable. This doesn't matter to us; what we are interested in is the change in depth.]

We will now assign a depth to each polygon as follows. Firstly, initialize the depth of each one to 0 ([Edit: initialize to 1, according to the example]). Next, iterate through your sorted list, but this time compare each element p only to the next element p+1. If (p+1 < p) is true, then increment the depth of p+1. Else, set the depth of p+1 to be the same as the depth of p.

Answer 5

Step 1: Orient your polygons in the same direction, say counter-clockwise.

Step 2: For any point (x, y) for which you need to compute the "depth" for compute the total winding number. This can be done in a number of ways; the fastest one in practice is to compute the SIGNED number of intersections between the horizontal (or vertical) ray originating in (x, y).

In particular, the depth for each polygon would be the depth of any of its vertices.