Well, the part I'm referring to is the 1. and 2 (Compute self-intersection nodes for A and B) as far as I understand. Traced through from https://github.com/libgeos/geos/blob/8bc1ea0f8019bcb83c60126e7883ae428ea7c677/src/geomgraph/GeometryGraph.cpp#L393 to https://github.com/libgeos/geos/blob/8bc1ea0f8019bcb83c60126e7883ae428ea7c677/src/geomgraph/index/SimpleMCSweepLineIntersector.cpp#L55 to https://github.com/libgeos/geos/blob/8bc1ea0f8019bcb83c60126e7883ae428ea7c677/src/geomgraph/index/SimpleMCSweepLineIntersector.cpp#L126
Increased a counter here https://github.com/libgeos/geos/blob/8bc1ea0f8019bcb83c60126e7883ae428ea7c677/src/geomgraph/index/ SimpleMCSweepLineIntersector.cpp#L161 and put a breakpoint here https://github.com/libgeos/geos/blob/8bc1ea0f8019bcb83c60126e7883ae428ea7c677/src/geomgraph/index/SimpleMCSweepLineIntersector.cpp#L166 counter's value at the end of the computeintersections was n*n, and breakpoint never hit, so a lot of checks were done without any useful work. Not really sure how "representative" my polygons are, but they are squares with zig-zag boundaries /\/\/\/\/\/\/\/\ \ / / \ \ / / \ \ / \/\/\/\/\/\/\/ don't know if that will come out right, but here's the snippet I used to check this: Geometry* getlargeGeo(double x, double y) { int i; const int n = 1000; GeometryFactory::Ptr factory = GeometryFactory::create(); CoordinateSequence* coords = factory->getCoordinateSequenceFactory()->create((size_t)0); for (i = 0; i < n; i++) coords->add(Coordinate(x + i * 10, y + (i % 2) * 3)); Coordinate last = coords->getAt(coords->size() - 1); for (i = 0; i < n; i++) coords->add(Coordinate(last.x + (i % 2) * 3, last.y + i * 10)); last = coords->getAt(coords->size() - 1); for (i = 0; i < n; i++) coords->add(Coordinate(last.x - i * 10, last.y + (i % 2) * 3)); last = coords->getAt(coords->size() - 1); for (i = 0; i < n; i++) coords->add(Coordinate(last.x + (i % 2) * 3, last.y - i * 10)); coords->add(coords->getAt(0)); LinearRing *shell = factory->createLinearRing(coords); return factory->createPolygon(shell, NULL); }; ... Geometry* lg1 = getlargeGeo(0, 0); Geometry* lg2 = getlargeGeo(100, 50); Geometry* res3 = lg1->intersection(lg2);
_______________________________________________ geos-devel mailing list geos-devel@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/geos-devel