_oriented_envelope.py 1.9 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556
  1. import math
  2. from itertools import islice
  3. import numpy as np
  4. import shapely
  5. from shapely.affinity import affine_transform
  6. def _oriented_envelope_min_area(geometry, **kwargs):
  7. """Compute the oriented envelope (minimum rotated rectangle).
  8. This is a fallback implementation for GEOS < 3.12 to have the correct
  9. minimum area behaviour.
  10. """
  11. if geometry is None:
  12. return None
  13. if geometry.is_empty:
  14. return shapely.from_wkt("POLYGON EMPTY")
  15. # first compute the convex hull
  16. hull = geometry.convex_hull
  17. try:
  18. coords = hull.exterior.coords
  19. except AttributeError: # may be a Point or a LineString
  20. return hull
  21. # generate the edge vectors between the convex hull's coords
  22. edges = (
  23. (pt2[0] - pt1[0], pt2[1] - pt1[1])
  24. for pt1, pt2 in zip(coords, islice(coords, 1, None))
  25. )
  26. def _transformed_rects():
  27. for dx, dy in edges:
  28. # compute the normalized direction vector of the edge
  29. # vector.
  30. length = math.sqrt(dx**2 + dy**2)
  31. ux, uy = dx / length, dy / length
  32. # compute the normalized perpendicular vector
  33. vx, vy = -uy, ux
  34. # transform hull from the original coordinate system to
  35. # the coordinate system defined by the edge and compute
  36. # the axes-parallel bounding rectangle.
  37. transf_rect = affine_transform(hull, (ux, uy, vx, vy, 0, 0)).envelope
  38. # yield the transformed rectangle and a matrix to
  39. # transform it back to the original coordinate system.
  40. yield (transf_rect, (ux, vx, uy, vy, 0, 0))
  41. # check for the minimum area rectangle and return it
  42. transf_rect, inv_matrix = min(_transformed_rects(), key=lambda r: r[0].area)
  43. return affine_transform(transf_rect, inv_matrix)
  44. _oriented_envelope_min_area_vectorized = np.frompyfunc(
  45. _oriented_envelope_min_area, 1, 1
  46. )