Repository navigation
Conversation
When projecting a line, each segment was checked separately: whether its start point is in the domain, whether the segment is straight enough (projecting its mid-point), and whether it is covered by or disjoint from the domain (creating a LineString). These Python-level calls dominated the time for large geometries, e.g. about 445k LineStrings for the 10m Natural Earth ocean. Do the same checks for all segments at once before the loop. Segments that would be accepted as a whole in the first step of the bisection are handled directly, and only the others are still bisected. The straightness calculation is shared by both paths, so the result is unchanged. For the example from SciTools#2102 with a Mercator projection, this reduces the time from ~25 s to ~9 s. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
_attach_lines_to_boundary created a new LineString from all coordinates so far for every boundary point or line that was added to a ring. For rings with many pieces, e.g. the 10m Natural Earth ocean in Robinson where every boundary vertex is added separately, this took quadratic time (~16 s). Collect the coordinates of the pieces instead and create the LineString once the ring is closed. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
|
Thanks @headtr1ck. Points 1 and 2 seem to be independent of each other. Did you benchmark them individually? I'd be curious to know if one or other of them is providing a significant part of the speedup. Note that I'm not qualified to review the cython changes, so this one will have to wait for someone else. I would definitely be glad to see a speed up to contour plotting though, as I've noticed that is a bottleneck in some of my team's workflows. |
Yes, see the tables in the PR description. This PR speeds up projection altogether while the other one mainly benefits small cutouts. Edit: sorry I misunderstood. Let me have a look.
If you decide to merge it before I am happy to rebase this PR. |
Rationale
Projecting large geometries with
Projection.project_geometryis slow, which is a big part of #2102. For example, projecting the single 10m Natural Earth ocean polygon (444k vertices) into Mercator takes ~12 s, while transforming the same vertices takes 0.05 s. Profiling showed two causes, one per commit:Per-segment Python calls in
trace.pyx.project_lineartransforms all vertices at once, but then checks each segment separately. For every segment it:shapely.intersects_xy);transform, to check the segment is straight enough;LineStringto test whether the segment is covered by, or disjoint from, the domain.For the 10m ocean that's ~445k
LineStrings, ~440k point tests and ~890k single-point transforms.This PR does the same checks for all segments at once before the loop. Segments that
_project_segmentwould accept as a whole in the first step of its bisection are handled directly:Only the remaining segments are still bisected: ones that cross the boundary, contain points that can't be projected, or need resampling. The straightness calculation moves into a shared
is_straighthelper used by both paths, so the maths is unchanged.Quadratic time in
Projection._attach_lines_to_boundary. Each time a boundary point or another line was added to a ring, a newLineStringwas created from all coordinates so far. For the 10m ocean in Robinson, every boundary vertex is added separately to a ring with ~100k coordinates, which took ~16 s. The coordinates of the pieces are now collected and theLineStringis created once, when the ring is closed.The output is unchanged. I compared the WKB of the projected geometries before and after this PR for 24,624 cases:
PlateCarreeandGeodeticsources.All results are byte-identical, including the cases that raise errors. The images of the example from #2102 are pixel-identical as well.
Benchmarks
asv run --python=same, comparing main, #2752, this PR, and #2752 combined with this PR. The feature artist benchmarks come from #2752 and were copied in for the runs without it. Speed-ups are relative to main.feature_artist.DrawNaturalEarthFeature.time_draw(Natural Earth land, LambertConformal)feature_artist.DrawNaturalEarthFeaturesLazily.time_savefig: the example from #2102 (10m ocean, land and borders around Italy)No projection is needed for PlateCarree, so this PR doesn't change that case.
project_linear.Suite.time_project_linear(Natural Earth ocean)#2752 doesn't touch
project_geometry, so its column is just noise here.The two PRs complement each other:
Together, the regional 10m map goes from 11.4 s to 0.18 s, and the #2102 example in Mercator from 23.4 s to 3.6 s.
Implications
add_geometries, contours and anything else that goes throughproject_geometry.geom_fully_insidefast path, which only exists for (warped) rectangular projections, now matters much less, since the per-segment domain checks are cheap for every projection.project_geometryand the image tests cover both code paths.🤖 Generated with Claude Code