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