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
12 changes: 11 additions & 1 deletion lib/cartopy/crs.py
Original file line number Diff line number Diff line change
Expand Up @@ -1185,7 +1185,17 @@ def _rings_to_multi_polygon(self, rings, is_ccw):
exterior_rings = []
interior_rings = []
for ring in rings:
if ring.is_ccw != is_ccw:
# Use the shoelace signed-area to determine ring orientation.
# ring.is_ccw is documented as unreliable for self-intersecting
# rings (e.g. when a projected polygon spans exactly ±180°
# longitude and the boundary attachment creates a degenerate
# "tail"). The shoelace sum dot(x[:-1], y[1:]) - dot(x[1:], y[:-1])
# gives twice the signed area and any self-intersecting
# tail contributions cancel out.
coords = np.asarray(ring.coords)
x, y = coords[:, 0], coords[:, 1]
ring_is_ccw = np.dot(x[:-1], y[1:]) > np.dot(x[1:], y[:-1])
if ring_is_ccw != is_ccw:
interior_rings.append(ring)
else:
exterior_rings.append(ring)
Expand Down
21 changes: 21 additions & 0 deletions lib/cartopy/tests/test_polygon.py
Original file line number Diff line number Diff line change
Expand Up @@ -281,6 +281,27 @@ def test_project_degenerate_poly(self):
polygons = target.project_geometry(polygon, source)
assert isinstance(polygons, sgeom.MultiPolygon)

def test_full_width_band_not_inverted(self):
# A polygon spanning exactly ±180° longitude creates a projected ring
# with a degenerate self-intersecting "tail" at the ±180 seam.
# Shapely's is_ccw is unreliable for such rings; the shoelace signed
# area must be used instead so the ring is classified correctly as an
# interior (not an exterior) ring. See GH-2483.
north_tropic = sgeom.Polygon(
[(-180, 0), (180, 0), (180, 20), (-180, 20)])
proj = ccrs.Orthographic(central_longitude=0.0, central_latitude=90.0)
result = proj.project_geometry(north_tropic, ccrs.PlateCarree())

# Convert points we expect to be inside vs out and check the
# containment of the result.
src = ccrs.PlateCarree()
inside_band = sgeom.Point(proj.transform_point(0, 10, src))
outside_band = sgeom.Point(proj.transform_point(0, 45, src))
assert result.contains(inside_band), \
'(0°E, 10°N) should be inside the projected tropical band'
assert not result.contains(outside_band), \
'(0°E, 45°N) should be outside the projected tropical band'


class TestQuality:
def setup_class(self):
Expand Down
Loading