Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

FIX: contourf nested contours #2225

Merged
merged 1 commit into from
Aug 4, 2023
Merged
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
15 changes: 10 additions & 5 deletions lib/cartopy/mpl/patch.py
Original file line number Diff line number Diff line change
Expand Up @@ -168,15 +168,20 @@ def path_to_geos(path, force_ccw=False):
geom = sgeom.LineString(path_verts)

# If geom is a Polygon and is contained within the last geom in
# collection, add it to its list of internal polygons, otherwise
# simply append it as a new external geom.
# collection, it usually needs to be an interior to that geom (e.g. a
# lake within a land mass). Sometimes there is a further geom within
# this interior (e.g. an island in a lake, or some instances of
# contours). This needs to be a new external geom in the collection.
if geom.is_empty:
pass
elif (len(collection) > 0 and
isinstance(collection[-1][0], sgeom.Polygon) and
isinstance(geom, sgeom.Polygon) and
collection[-1][0].contains(geom.exterior)):
collection[-1][1].append(geom.exterior)
if any(internal.contains(geom) for internal in collection[-1][1]):
collection.append((geom, []))
Copy link
Member Author

@rcomer rcomer Aug 4, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Although I’ve tested this across quite a few plots, I’m now questioning whether it is quite right. If we append the geom from the third level down to collection, could the next geom be at the second level down and so would need to check against something earlier than collection[-1]?

Copy link
Contributor

@greglucas greglucas Aug 4, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hm... that is a good question. I wonder if you could somehow keep track of the "last" geometry in the nested collections and sort of recursively drill down/out from there?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I was thinking of this as a tree problem with the path segments ordered in a "depth first" way. However looking at matplotlib/matplotlib#25247, I think the paths that used to be separate are now just concatenated together. So anything at the second level has to be directly preceded by either something else at the second level or something at the first level, and so my above hypothetical can't happen for contours at least.

else:
collection[-1][1].append(geom)
elif isinstance(geom, sgeom.Point):
other_result_geoms.append(geom)
else:
Expand All @@ -188,8 +193,8 @@ def path_to_geos(path, force_ccw=False):
geom_collection = []
for external_geom, internal_polys in collection:
if internal_polys:
# XXX worry about islands within lakes
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I assume that XXX is the same as TODO and I'm taking a punt that this is now covered.

geom = sgeom.Polygon(external_geom.exterior, internal_polys)
exteriors = [geom.exterior for geom in internal_polys]
geom = sgeom.Polygon(external_geom.exterior, exteriors)
else:
geom = external_geom

Expand Down
16 changes: 16 additions & 0 deletions lib/cartopy/tests/mpl/test_patch.py
Original file line number Diff line number Diff line change
Expand Up @@ -43,3 +43,19 @@ def test_polygon_with_interior_and_singularity(self):
geoms = cpatch.path_to_geos(p)
assert [type(geom) for geom in geoms] == [sgeom.Polygon, sgeom.Point]
assert len(geoms[0].interiors) == 1

def test_nested_polygons(self):
# A geometry with three nested squares.
vertices = [[0, 0], [0, 10], [10, 10], [10, 0], [0, 0],
[2, 2], [2, 8], [8, 8], [8, 2], [2, 2],
[4, 4], [4, 6], [6, 6], [6, 4], [4, 4]]
codes = [1, 2, 2, 2, 79, 1, 2, 2, 2, 79, 1, 2, 2, 2, 79]
p = Path(vertices, codes=codes)
geoms = cpatch.path_to_geos(p)

# The first square makes the first geometry with the second square as
# its interior. The third square is its own geometry with no interior.
assert len(geoms) == 2
assert all(type(geom) == sgeom.Polygon for geom in geoms)
assert len(geoms[0].interiors) == 1
assert len(geoms[1].interiors) == 0