I'm using a PL/R
function and PostGIS
to generate voronoi polygons around a set of points. The function that I am using is defined here. When I use this function on a particular dataset I get the following error message:
Error : ERROR: R interpreter expression evaluation error
DETAIL: Error in pg.spi.exec(sprintf("SELECT %3$s AS id,
st_intersection('SRID='||st_srid(%2$s)||';%4$s'::text,'%5$s')
AS polygon FROM %1$s WHERE st_intersects(%2$s::text,'SRID='||st_srid(%2$s)||';%4$s');",
:error in SQL statement : Error performing intersection: TopologyException: found non-noded
intersection between LINESTRING (571304 310990, 568465 264611) and LINESTRING (568465
264611, 594406 286813) at 568465.05533706467 264610.82749605528
CONTEXT: In R support function pg.spi.exec In PL/R function r_voronoi
From examining this part of the error message:
Error performing intersection: TopologyException: found non-noded intersection between
LINESTRING (571304 310990, 568465 264611) and LINESTRING (568465 264611, 594406 286813)
at 568465.05533706467 264610.82749605528
This is what the problem listed above looks like:
I initially thought that this message might be caused by the existence of identical points, and tried to solve this using the st_translate()
function, used in the following way:
ST_Translate(geom, random()*20, random()*20) as geom
This does fix the problem, but my concern is that I am now translating all points up to ~20m in the x/y direction. I also can't tell what an appropriate translation amount is need.
For example, in this dataset through trial and error a 20m * random number
is ok, but how can I tell if this needs to be bigger?
Based on the image above I think the problem is that the point is intersecting with the line while the algorithm is trying to intersect the point with a polygon. I'm not sure what I should be doing to ensure that the point is within a polygon, rather than intersecting with a line. The error is occurring on this line:
"SELECT
%3$s AS id,
st_intersection(''SRID=''||st_srid(%2$s)||'';%4$s''::text,''%5$s'') AS polygon
FROM
%1$s
WHERE
st_intersects(%2$s::text,''SRID=''||st_srid(%2$s)||'';%4$s'');"
I have read through What is a "non-noded intersection"? to try to better understand this problem.
Best Answer
In my experience this problem is nearly always caused by:
The "nudge" approach of the
ST_Buffer
solutions lets you get away with #2, but anything you can do to resolve these underlying causes, like snapping your geometry to a 1e-6 grid, will make your life easier. The buffered geometries are usually fine for intermediate calculations like overlap area, but you'll want to be careful about retaining them because they can make your close-but-not-quite problems worse in the long haul.PostgreSQL's exception-handling capability lets you write wrapper functions to handle these special cases, buffering only when needed. Here's an example for
ST_Intersection
; I use a similar function forST_Difference
. You'll need to decide if the buffering and potential return of an empty polygon are acceptable in your situation.Another benefit with this approach is you can pinpoint the geometries that are actually causing your problems; just add some
RAISE NOTICE
statements in theEXCEPTION
block to output WKT or something else that will help you track down the problem.