uniquely reporting intersections in a list of axis-aligned rectangles

Question

I have a list of coordinates that form axis-aligned 2D boxes (axis-oriented/iso-oriented rectangles, rects).
I want to see how many boxes intersect each other without repeating it. For example, I have boxes A, B, and C. If A intersects with B, B intersecting with A will still count as one intersection instead of two separate ones.

The way I'm thinking is to grab one box, compare it to all the rest and then throw it out since the comparison is done and I want no repeats. For example going back to A, B and C. If I'm on A and it intersects B and not C I have made its round and hence will not need to keep it in the loop. Hence once I'm on B it will not recheck A.

I can't think of a faster method. This is akin to finding duplicates using linear search and I think the worst-case complexity will be O(n^2). I don't see how to use binary search as it is possible that there are multiple intersections. I've been also looking at the sorting algorithm but I need one that won't go back and check the older values again.

Answer 1

You can solve this in O(n log n). Since you're trying to solve a static intersection problem in two dimensions, a common trick is to transform the problem into a dynamic intersection problem in one dimension (i.e., one space dimension becomes the time dimension).

Suppose you have a closed rectangle with lower left corner (x1, y1) and upper right corner (x2, y2). This rectangle becomes two "interval events":

Interval [y1, y2] inserted at time x1
Interval [y1, y2] removed after time x2

Transform all your rectangles into events in this way, and then sort by time (breaking ties with insertions coming before removals).

Now, you need a data structure that lets you add and remove intervals [A, B], and also query the number of intervals in the data structure intersecting [A, B]. Then you process the "interval events" in the sorted order, but keep a running sum of how many current intervals intersect [A, B] before inserting each interval [A, B].

One data structure to do this in O(log n) time per operation is with two balanced binary search trees: one holding beginning-points of intervals, and the other holding ending-points of intervals. You'll also need to augment each BST to be an Order Statistic Tree, to quickly query the number of points less than or equal to a certain value that are in the tree.

Then, finding how many intervals currently in your data structure intersect an arbitrary interval [A, B] is simple; that count is:

#(Intervals intersecting [A, B]) =

  #(values in the beginning-points tree that are <= B) 
- #(values in the ending-points tree that are < A)

which you can check is correct from the definition of two intervals intersecting: neither interval starts after the other one has ended.

You could also replace the order statistic tree with a data structure for prefix-sums, like a Fenwick tree, which requires much less code for the same complexity.

Answer 2

Sample implementation of kcsquared's algorithm. You didn't specify a language, so I chose Python since I'm familiar with it and it's short readable code. Rectangles have lower left coordinate (x,y) and upper right coordinate (X,Y).

def intersections(rects):
    events = []
    for A in rects:
        events.append((A.x, 'insert', A.y, A.Y))
        events.append((A.X, 'remove', A.y, A.Y))
    intersections = 0
    ys = SortedList()
    Ys = SortedList()
    for x, op, y, Y in sorted(events):
        if op == 'insert':
            intersections += ys.bisect_right(Y) - Ys.bisect_left(y)
            ys.add(y)
            Ys.add(Y)
        else:
            ys.remove(y)
            Ys.remove(Y)
    return intersections

Tests on lists of 5,000 random rectangles against an O(n²) reference implementation, showing the numbers of intersections and runtimes:

124257  5465 ms  reference
124257   127 ms  intersections

121166  5444 ms  reference
121166   124 ms  intersections

118980  5435 ms  reference
118980   124 ms  intersections

With 10,000 rectangles:

489617  22342 ms  reference
489617    292 ms  intersections

489346  22491 ms  reference
489346    296 ms  intersections

489990  22302 ms  reference
489990    290 ms  intersections

Full code:

def reference(rects):
    intersections = 0
    for A, B in combinations(rects, 2):
        if A.X >= B.x and A.x <= B.X and A.Y >= B.y and A.y <= B.Y:
            intersections += 1
    return intersections

def intersections(rects):
    events = []
    for A in rects:
        events.append((A.x, 'insert', A.y, A.Y))
        events.append((A.X, 'remove', A.y, A.Y))
    intersections = 0
    ys = SortedList()
    Ys = SortedList()
    for x, op, y, Y in sorted(events):
        if op == 'insert':
            intersections += ys.bisect_right(Y) - Ys.bisect_left(y)
            ys.add(y)
            Ys.add(Y)
        else:
            ys.remove(y)
            Ys.remove(Y)
    return intersections

from random import randint, randrange
from itertools import combinations
from timeit import default_timer as timer
from sortedcontainers import SortedList
from collections import namedtuple

Rect = namedtuple('Rect', 'x X y Y')

for _ in range(3):
    rects = [Rect(x, x + randint(1, 100), y, y + randint(1, 100))
             for _ in range(5000)
             for x, y in [(randrange(1000), randrange(1000))]]
    for func in reference, intersections:
        t0 = timer()
        result = func(rects)
        t1 = timer()
        print(result, '%4d ms' % ((t1 - t0) * 1e3), func.__name__)
    print()

Answer 3

You can somewhat reduce the average complexity, more precisely the computation time, but you'll always get a worst case in n²... Because, inherently, that's a O(n²) algorithm, since you have to test two by two every N boxes.

One way to reduce that computation time is to "attach", to each box, its circumscribed sphere. It's quite simple to compute: its center is the barycenter of the 8 vertexes of the box / 4 vertexes of the rectangle, and its radius is any distance from this barycenter to one vertex.

The trick is: given two boxes, if their circumscribed spheres aren't intersecting (distance between the two centers is greater than the sum of the two radiuses), then the boxes can not intersect. If the two spheres intersect, then the boxes may intersect.

This trick doesn't change at all the complexity, but it reduce A LOT the computation time, since it's way faster to test for spheres than for boxes (especially if they aren't parallel to axes...). And it will reduce the list size, up to empty it for any potential intersections. Also, you also have only a shortened (even empty) list of potential intersections for each box.

You can gain some computation time by using a list of boxes pairs to finely test, using squares of distances for comparisons, etc.: it's not so much, but in the end, it still saves some time. Then, you can compute the real boxes intersections with way less boxes pairs.

This would be a significative boost with a lot of 3D random boxes unaligned with axes, and it won't be with few 2D rectangles aligned with axes. In both cases, you'll always be with a O(n²) complexity (i.e. the worst case).