From 77e77f9d579036d85d4af7f72e315ad0444d798a Mon Sep 17 00:00:00 2001 From: Elliott Sales de Andrade Date: Sat, 13 Oct 2018 04:07:07 -0400 Subject: [PATCH] Fix projected line strings getting truncated early. The limit to prevent an infinite loop attempting to split line segments "small enough" was applied to the entire line string, which meant the entire path could get truncated. Instead, make the limit apply to each segment individually. Since each segment is somewhat smaller, drop the limit a bit to compensate for having to go through more segments. This keeps things running in a reasonable time. --- lib/cartopy/_trace.cpp | 3 ++- lib/cartopy/tests/test_polygon.py | 45 +++++++++++++++++++++++++++++++ 2 files changed, 47 insertions(+), 1 deletion(-) diff --git a/lib/cartopy/_trace.cpp b/lib/cartopy/_trace.cpp index c8d2f4c95a..d7fb4728c9 100644 --- a/lib/cartopy/_trace.cpp +++ b/lib/cartopy/_trace.cpp @@ -563,7 +563,8 @@ void _project_segment(GEOSContextHandle_t handle, t_current = 0.0; state = get_state(p_current, gp_domain, handle); - while(t_current < 1.0 && lines.size() < 500) + size_t old_lines_size = lines.size(); + while(t_current < 1.0 && (lines.size() - old_lines_size) < 100) { //std::cerr << "Bisecting" << std::endl; #ifdef DEBUG diff --git a/lib/cartopy/tests/test_polygon.py b/lib/cartopy/tests/test_polygon.py index 7d4086d021..8d7900024f 100644 --- a/lib/cartopy/tests/test_polygon.py +++ b/lib/cartopy/tests/test_polygon.py @@ -18,6 +18,7 @@ from __future__ import (absolute_import, division, print_function) import numpy as np +import pytest import shapely.geometry as sgeom import shapely.wkt @@ -135,6 +136,50 @@ def test_project_previous_infinite_loop(self): src = ccrs.PlateCarree() src._attach_lines_to_boundary(multi_line_strings, True) + @pytest.mark.parametrize('proj', + [ccrs.InterruptedGoodeHomolosine, ccrs.Mollweide]) + def test_infinite_loop_bounds(self, proj): + # test a polygon which used to get stuck in an infinite loop but is now + # erroneously clipped. + # see https://github.com/SciTools/cartopy/issues/1131 + + # These names are for IGH; effectively the same for Mollweide. + bottom = [0., 70.] + right = [0., 90.] + top = [-180., 90.] + left = [-180., 70.] + verts = np.array([ + bottom, + right, + top, + left, + bottom, + ]) + bad_path = sgeom.Polygon(verts) + + target = proj() + source = ccrs.PlateCarree() + + projected = target.project_geometry(bad_path, source) + + # When transforming segments was broken, the resulting path did not + # close, and either filled most of the domain, or a smaller portion + # than it should. Check that the bounds match the individual points at + # the expected edges. + projected_left = target.transform_point(left[0], left[1], source) + assert projected.bounds[0] == pytest.approx(projected_left[0], + rel=target.threshold) + projected_bottom = target.transform_point(bottom[0], bottom[1], source) + assert projected.bounds[1] == pytest.approx(projected_bottom[1], + rel=target.threshold) + projected_right = target.transform_point(right[0], right[1], source) + assert projected.bounds[2] == pytest.approx(projected_right[0], + rel=target.threshold, + abs=1e-8) + projected_top = target.transform_point(top[0], top[1], source) + assert projected.bounds[3] == pytest.approx(projected_top[1], + rel=target.threshold) + def test_3pt_poly(self): projection = ccrs.OSGB() polygon = sgeom.Polygon([(-1000, -1000),