Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions docs/history.rst
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@ Latest
- ENH: Added :class:`pyproj.enums.GridAvailabilityUse` enum for TransformerGroup ``grid_check`` kwarg
- ENH: Added always_xy kwarg to :meth:`pyproj.transformer.Transformer.from_pipeline`
- DOC: Update fiona CRS compatibility for fiona 1.9+ (issue #1360)
- DOC: Describe the size limit and negative results of :meth:`pyproj.Geod.geometry_area_perimeter` concretely (issue #1592)
- BUG: Clear CONTEXT_THREAD_KEY when destroying PJ_CONTEXT (pull #1541)
- BUG: Default skew angle for CF grid mapping oblique mercator to 90 (issue #1506)
- BLD: update build-time dependencies
Expand Down
31 changes: 27 additions & 4 deletions pyproj/geod.py
Original file line number Diff line number Diff line change
Expand Up @@ -1031,15 +1031,24 @@ def geometry_area_perimeter(

.. note:: lats should be in the range [-90 deg, 90 deg].

.. note:: | There are a few limitations :
| - only works with areas up to half the size of the globe ;
| - certain large polygons may return negative values.
.. note:: The area is reduced modulo the total area of the ellipsoid,
into the range ``(-A/2, A/2]``, where ``A`` is that total area
(about ``5.10e14`` m\\ :sup:`2` on WGS84, so ``A/2`` is about
``2.55e14`` m\\ :sup:`2`). A region covering more than half the
ellipsoid is therefore reported as the negated area of its
complement: a polygon enclosing 99% of the globe returns about
-1% of the total area, not +99%. Size is the only reason a
correctly oriented polygon comes back negative, and "large"
means precisely "more than half the ellipsoid" -- shape,
proximity to a pole and crossing the antimeridian do not
matter on their own.

.. warning:: The area returned is signed with counter-clockwise (CCW) traversal
being treated as positive. For polygons, holes should use the
opposite traversal to the exterior (if the exterior is CCW, the
holes/interiors should be CW). You can use `shapely.ops.orient` to
modify the orientation.
modify the orientation. This is the other reason a result can be
negative, and unlike the size limit above it applies at any size.

If it is a Polygon, it will return the area and exterior perimeter.
It will subtract the area of the interior holes.
Expand All @@ -1063,6 +1072,20 @@ def geometry_area_perimeter(
>>> f"{poly_area:.0f} {poly_perimeter:.0f}"
'944373881400 3979008'

A polygon enclosing more than half the ellipsoid wraps around to a
negative value, even though its vertices are ordered counter-clockwise:

>>> nearly_global = Polygon(
... LineString([
... Point(-179, -89), Point(179, -89), Point(179, 89), Point(-179, 89)
... ])
... )
>>> geod.geometry_area_perimeter(nearly_global)[0] < 0
True

Its magnitude is the area of the thin strip left outside the polygon, so
adding the total area of the ellipsoid recovers the enclosed area.


Parameters
----------
Expand Down
5 changes: 3 additions & 2 deletions test/test_doctest_wrapper.py
Original file line number Diff line number Diff line change
Expand Up @@ -29,8 +29,9 @@ def test_doctests():
try:
import shapely.geometry # noqa: F401 pylint: disable=unused-import
except (ImportError, OSError):
# missing shapely
expected_failure_count = 6
# missing shapely; these are the examples in Geod.geometry_length and
# Geod.geometry_area_perimeter that need a shapely geometry
expected_failure_count = 8

# if the below line fails, doctests have failed
assert failure_count == expected_failure_count, (
Expand Down