Unfortunately, coincident geometries are the norm in our case. The two in the example I gave are the results of previous differences with a third overlapping polygon(s).
I tried Ben's suggestion first and got a compile error about a missing include file (this is all on a x86_64 Linux install using a very recent svn version). Thought you might want to know. geos/include/geos/geom/BinaryOp.h:54:44: error: geos/precision/GeometrySnapper.h: No such file or directory GeometrySnapper.h is listed in noinst_HEADERS of source/headers/geos/precision/Makefile Just to see what would happen I copied the file to include/geos/precision and it compiled, linked, and ran successfully. Stuart Ben Jubb <[EMAIL PROTECTED]> Sent by: [EMAIL PROTECTED] 09/19/2007 01:54 PM Please respond to GEOS Development List <geos-devel@geos.refractions.net> To GEOS Development List <geos-devel@geos.refractions.net> cc Subject Re: [geos-devel] no outgoing dirEdge found error There is a way to make this particular operation succeed in GEOS, which is to use the BinaryOp method instead of difference(). For example, the following code works (in MSVC8): #include <string> #include <iostream> #include <geos/geom/PrecisionModel.h> #include <geos/geom/GeometryFactory.h> #include <geos/geom/Geometry.h> #include <geos/geom/MultiPolygon.h> #include <geos/io/WKTReader.h> #include <geos/geom/BinaryOp.h> #include <geos/operation/overlay/OverlayOp.h> using namespace geos; int main (int argc, char *argv[]) { using namespace operation::overlay; typedef std::auto_ptr< geom::Geometry > GeomAutoPtr; std::string mp1s = "MULTIPOLYGON (((88.0541300495810049 -40.3914336120775417," " 88.0013586636740683 -40.3952580652168507," " 87.9856522556763849 -40.3125682141688557," " 88.0410955427405639 -40.3084853410231361," " 88.0447770149226727 -40.3320733430391556," " 88.0541300495810049 -40.3914336120775417))," " ((87.9569549228361041 -40.1614847724049824," " 87.9856522556763849 -40.3125682141688557," " 87.9594803530699494 -40.3144955269957563," " 87.9484023890450430 -40.2597855714132962," " 87.9236385364479673 -40.1638760511710515," " 87.9569549228361041 -40.1614847724049824)))"; std::string mp2s = "MULTIPOLYGON (((87.7184147523609425 -40.3314878875577705," " 87.7106409229361930 -40.2916861830969211," " 87.6363928113999862 -39.9074815079801084," " 87.8588765594887349 -39.8933226300382628," " 88.3245648150878964 -39.8616404104367064," " 88.3252910982550503 -39.8653834273614933," " 88.4058355104796334 -40.2774282017416283," " 88.4066631207306273 -40.2815647501607259," " 87.9574816223991007 -40.3146427145853394," " 87.7184147523609425 -40.3314878875577705)))"; try { geos::geom::PrecisionModel *model = new geos::geom::PrecisionModel(geos::geom::PrecisionModel::FLOATING); geos::geom::GeometryFactory *factory = new geos::geom::GeometryFactory(model); geos::io::WKTReader *wkt = new geos::io::WKTReader(); geos::geom::MultiPolygon *mp1 = (geos::geom::MultiPolygon *)wkt->read(mp1s); geos::geom::MultiPolygon *mp2 = (geos::geom::MultiPolygon *)wkt->read(mp2s); std::cout << "Is valid of one = " << mp1->isValid() << std::endl; std::cout << "Is valid of two = " << mp2->isValid() << std::endl; //geos::geom::Geometry *diff = mp1->difference(mp2); GeomAutoPtr gRealRes = BinaryOp(mp1, mp2, overlayOp(OverlayOp::opDIFFERENCE)); std::cout << mp1->toString() << std::endl; std::cout << mp2->toString() << std::endl; //std::cout << diff->toString() << std::endl; std::cout << gRealRes.get()->toString() << std::endl; } catch (std::exception const &se) { std::cout << "ERROR - " << se.what() << std::endl; } catch (...) { std::cout << "General error" << std::endl; } } BInaryOp does the same difference op, but when that fails it tries again after shifting the geometry close to the origin. Presumably this gains a few bits of precision in the intermediate results, that allows the operation to succeed. This method is perhaps not as robust tho. b Martin Davis wrote: > This computation works fine in the current version of JTS. This may > indicate a porting bug in GEOS, or a slightly different choice of > tolerances. > > This is a classic case which causes robustness problems - overlaying > two geometries with linework which is almost coincident. > > The only suggestion I can make is to reduce the precision of your data > slightly in cases which fail. Presumably your source data isn't > actually accurate to 16 decimal places? It could be a while before > we can look at the GEOS code for this. > > Stuart C Sides wrote: >> >> Hi GEOS devels, >> >> We are writing a C++ application which is looking for areas of >> overlap between 1000's of polygons. The app is being >> developed on Linux with g++ 4.1.0. >> >> We first tried the stable release of GEOS 2.2.3 and received this >> error after processing several polygons: >> >> TopologyException: found non-noded intersection between 87.9595 >> -40.3145, 87.9484 -40.2598 and 87.9857 -40.3126, 87.9575 -40.3146 >> 87.9595 -40.3145 >> >> While reviewing the list archives we found references to major code >> changes that might fix some of these errors, so we >> upgraded to RC3 and eventually to an SVN version from 2007-09-18. We >> got through a lot more polygons without errors, >> but finally received this error: >> >> TopologyException: no outgoing dirEdge found 87.957 -40.1615 >> >> Note: this coordinate is from the first segment of the second polygon >> of the first multi-polygon below. >> >> In the original code, before I converted the multipolygons to strings >> to create the example, we received this error: >> >> TopologyException: no outgoing dirEdge found 87.7184 -40.3315 >> >> Note: this coordinate is from the first segment of the second >> multi-polygon below: >> >> I've boiled the error down to just a few lines of code, and would >> appreciate any suggestions. The error is thrown at the >> difference. >> >> Thanks >> Stuart >> >> >> >> #include <string> >> #include <iostream> >> >> #include <geos/geom/PrecisionModel.h> >> #include <geos/geom/GeometryFactory.h> >> #include <geos/geom/Geometry.h> >> #include <geos/geom/MultiPolygon.h> >> #include <geos/io/WKTReader.h> >> >> int main (int argc, char *argv[]) >> { >> std::string mp1s = "MULTIPOLYGON (((88.0541300495810049 >> -40.3914336120775417," >> " 88.0013586636740683 -40.3952580652168507," >> " 87.9856522556763849 -40.3125682141688557," >> " 88.0410955427405639 -40.3084853410231361," >> " 88.0447770149226727 -40.3320733430391556," >> " 88.0541300495810049 -40.3914336120775417))," >> " ((87.9569549228361041 -40.1614847724049824," >> " 87.9856522556763849 -40.3125682141688557," >> " 87.9594803530699494 -40.3144955269957563," >> " 87.9484023890450430 -40.2597855714132962," >> " 87.9236385364479673 -40.1638760511710515," >> " 87.9569549228361041 -40.1614847724049824)))"; >> >> std::string mp2s = "MULTIPOLYGON (((87.7184147523609425 >> -40.3314878875577705," >> " 87.7106409229361930 -40.2916861830969211," >> " 87.6363928113999862 -39.9074815079801084," >> " 87.8588765594887349 -39.8933226300382628," >> " 88.3245648150878964 -39.8616404104367064," >> " 88.3252910982550503 -39.8653834273614933," >> " 88.4058355104796334 -40.2774282017416283," >> " 88.4066631207306273 -40.2815647501607259," >> " 87.9574816223991007 -40.3146427145853394," >> " 87.7184147523609425 -40.3314878875577705)))"; >> >> try { >> geos::geom::PrecisionModel *model = >> new >> geos::geom::PrecisionModel(geos::geom::PrecisionModel::FLOATING); >> geos::geom::GeometryFactory *factory = new >> geos::geom::GeometryFactory(model); >> >> geos::io::WKTReader *wkt = new geos::io::WKTReader(); >> geos::geom::MultiPolygon *mp1 = (geos::geom::MultiPolygon >> *)wkt->read(mp1s); >> geos::geom::MultiPolygon *mp2 = (geos::geom::MultiPolygon >> *)wkt->read(mp2s); >> >> std::cout << "Is valid of one = " << mp1->isValid() << std::endl; >> std::cout << "Is valid of two = " << mp2->isValid() << std::endl; >> >> geos::geom::Geometry *diff = mp1->difference(mp2); >> >> std::cout << mp1->toString() << std::endl; >> std::cout << mp2->toString() << std::endl; >> std::cout << diff->toString() << std::endl; >> >> } >> catch (std::exception const &se) { >> std::cout << "ERROR - " << se.what() << std::endl; >> } >> catch (...) { >> std::cout << "General error" << std::endl; >> } >> >> } >> ------------------------------------------------------------------------ >> >> _______________________________________________ >> geos-devel mailing list >> geos-devel@geos.refractions.net >> http://geos.refractions.net/mailman/listinfo/geos-devel >> > _______________________________________________ geos-devel mailing list geos-devel@geos.refractions.net http://geos.refractions.net/mailman/listinfo/geos-devel
_______________________________________________ geos-devel mailing list geos-devel@geos.refractions.net http://geos.refractions.net/mailman/listinfo/geos-devel