[dbaston]

This PR adds a GeometrySplitter class for splitting linear or polygonal features using points, lines, and/or polygons.

Curved geometries are not yet supported. A couple of pieces are missing for this:

• to split a curve with a point, we need to find the closest section of the curve to that point. This requires
  adding curve support to DistanceOp.
• to split a CurvePolygon with a line/curve, we need to figure out which of the resulting polygons are inside
  the original geometry. This is done with PIP tests, but the existing InteriorPointArea algorithm does not
  support curves. Either InteriorPointArea could be extended for this, or we could compute interior points
  by testing either (pseudo)random points or computing small offsets from segments/arcs. EDIT: added RandomPointsBuilder from JTS to handle this.

I've added applicable test cases from PostGIS and QGIS, most of which are passing.

Currently not passing:

• tests relying on missing curve support described above
• splitting a self-intersecting line with another line causes splits at the self-intersection points.
  PostGIS does the same but I understand why the QGIS behavior is desirable.
• a split of a nearly-collapsed polygon produces a different result from QGIS and PostGIS. I'm not sure
  the results from that software are desirable, though.

[pramsey]

Clearly we just need to move PostGIS and QGIS to the GEOS implementation ASAP and resolve the inconsistency problem 😄

[nyalldawson]

Thank you!! I would LOVE to drop QGIS' fragile and inefficient version of this! 🥳🥳

[m-kuhn]

This is very good news for everyone !

[jorisvandenbossche]

@dbaston I am testing this with Shapely, by using the new C API instead of the python homegrown splitting logic in shapely, and running the tests, and at the moment only noticing one issue / difference.

For the case of splitting a line with a line, when both lines overlap (so not just cross), Shapely at the moment detects this and raises an error disallowing such input.

With GEOS, the split function works:

$ geosop split -a 'LINESTRING (0 0, 1.5 1.5, 3 4)' -b 'LINESTRING (0 0, 15 15)'
GEOMETRYCOLLECTION (LINESTRING (0 0, 1.5 1.5), LINESTRING (1.5 1.5, 3 4))

I don't have a strong opinion about what is the best behaviour (the above output seems the logical result if you want to have it just work), but just raising this to be sure this was a conscious design decision (or might warrant a test).

[dbaston]

Thanks @jorisvandenbossche , I should have listed this difference this in the PR description but forgot to.

I don't know if either the GEOS or PostGIS behavior is a conscious design decision so much as a natural consequence of the implementation. PostGIS is taking advantage of a GEOSDifference behavior to perform the split, but that only produces the desired output if the intersection of the two inputs does not contain any lines. The GEOS implementation is not based on GEOSDifference and does not have this limitation, and I didn't see a reason to explicitly introduce it.

or might warrant a test

There are in fact a couple:

geos/tests/unit/operation/split/GeometrySplitterTest.cpp

Lines 371 to 393 in 7ce7194

	template<>
	template<>
	void object::test<30>()
	{
	set_test_name("Split line by overlapping line (1); PostGIS test 8.1");
	// Note: PostGIS produces an error here instead of a result
	testSplit("LINESTRING(0 0, 10 0)",
	"LINESTRING(5 0, 20 0)",
	"GEOMETRYCOLLECTION (LINESTRING (0 0, 5 0), LINESTRING (5 0, 10 0))");
	}
	template<>
	template<>
	void object::test<31>()
	{
	set_test_name("Split line by overlapping line (1); PostGIS test 8.2");
	// Note: PostGIS produces an error here instead of a result
	testSplit("LINESTRING(0 0, 10 0)",
	"LINESTRING(5 0, 8 0)",
	"GEOMETRYCOLLECTION (LINESTRING (0 0, 5 0), LINESTRING (5 0, 8 0), LINESTRING (8 0, 10 0))");
	}