Fast Buffer for Shapely Point

fast_buffer.md

It's common to use the shapely.stretree.STRtree.query to find what's near a point. Often the point is buffered to get whatever is roughly within a radius of the point. For example,

dist = 123.456
geoms_near_point = tree.query(pt.buffer(dist))

What's really being returned here are any geometries in the tree whose bounding boxes intersect with the point's bounding box. By default, pt.buffer(dist) returns a buffer with 66 points that make a circle. That is 16 for each quandrent. Really only 4 are needed though because only the bounding box is being used.

So for a slight efficiency gain, set the resolution parameter to 1 so that a polygon with just 4 points is created.

dist = 123.456
geoms_near_point = tree.query(pt.buffer(dist, resolution=1))

This will have the exact same behavior as the previous example, while being slightly faster. How much faster? Probably not enough to really matter. But the following script gives some numbers.

from shapely.geometry import Point
from timeit import default_timer

iterations = 100000
zero_pt = Point(0, 0)
	
start = default_timer()
for i in range(iterations):
    zero_pt.buffer(1)
stop = default_timer()
print(f"buffer(1)              : {stop - start}")

start = default_timer()
for i in range(iterations):
    zero_pt.buffer(1, resolution=1)
stop = default_timer()
print(f"buffer(1, resolution=1): {stop - start}")

buffer(1)              : 2.8301283220062032
buffer(1, resolution=1): 1.9343692009570077

However, if the Point object is just a means to a box, using shapely.geometry.box is about 2x faster.

"""
Create a bunch of shapely boxes.
"""
import random

from shapely.geometry import Point, Polygon, box

fjds = Polygon, box, Point

# def point_buffer(pt: Point, buffer: float):
#     x, y = pt.x, pt.y
#
#     return Polygon(
#         [
#             (x - buffer, y - buffer),
#             (x + buffer, y - buffer),
#             (x + buffer, y + buffer),
#             (x - buffer, y + buffer),
#             (x - buffer, y - buffer),
#         ]
#     )


# def point_buffer(x: float, y: float, buffer: float):
#     return Polygon(
#         [
#             (x - buffer, y - buffer),
#             (x + buffer, y - buffer),
#             (x + buffer, y + buffer),
#             (x - buffer, y + buffer),
#             (x - buffer, y - buffer),
#         ]
#     )


# def point_buffer(xy: tuple[float, float], buffer: float):
#     return box(xy[0] - buffer, xy[1] - buffer, xy[0] + buffer, xy[1] + buffer)

def point_buffer(xy: tuple[float, float], buffer: float):
    return Point(xy).buffer(buffer, resolution=1)


def create_lots(n):
    random.seed(0)
    for i in range(n):
        x = random.random()
        y = random.random()
        p = point_buffer((x, y), 0.01)


if __name__ == "__main__":
    import timeit

    print(timeit.timeit("create_lots(100_000)", setup="from __main__ import create_lots", number=5))