Fixing non-noded intersection problem using PostGIS
In my experience this problem is nearly always caused by:
- High precision in your coordinates (43.231499999999996), combined with
- Lines that are almost coincident but not identical
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 for ST_Difference
. You'll need to decide if the buffering and potential return of an empty polygon are acceptable in your situation.
CREATE OR REPLACE FUNCTION safe_isect(geom_a geometry, geom_b geometry)
RETURNS geometry AS
$$
BEGIN
RETURN ST_Intersection(geom_a, geom_b);
EXCEPTION
WHEN OTHERS THEN
BEGIN
RETURN ST_Intersection(ST_Buffer(geom_a, 0.0000001), ST_Buffer(geom_b, 0.0000001));
EXCEPTION
WHEN OTHERS THEN
RETURN ST_GeomFromText('POLYGON EMPTY');
END;
END
$$
LANGUAGE 'plpgsql' STABLE STRICT;
Another benefit with this approach is you can pinpoint the geometries that are actually causing your problems; just add some RAISE NOTICE
statements in the EXCEPTION
block to output WKT or something else that will help you track down the problem.
Through a lot of trial and error I eventually realized that the non-noded intersection
resulted from a self-intersection problem. I found a thread that suggested using ST_buffer(geom, 0)
can be used to fix the problem (though it does make it a lot slower overall). I then tried to use ST_MakeValid()
and when applied directly to the geometry before any other function. This seems to fix the problem robustly.
ipoint <- pg.spi.exec(
sprintf(
"SELECT
%3$s AS id,
st_intersection(ST_MakeValid(''SRID=''||st_srid(%2$s)||'';%4$s''::text), ST_MakeValid(''%5$s'', 0)) AS polygon
FROM %1$s
WHERE
ST_Intersects(ST_MakeValid(%2$s::text),ST_MakeValid(''SRID=''||st_srid(%2$s)||'';%4$s''));",
arg1,
arg2,
arg3,
curpoly,
buffer_set$ewkb[1]
)
)
I've marked this as the answer as it seems to be the only approach that fixes my problem.
I ran into this same problem (Postgres 9.1.4, PostGIS 2.1.1), and the only thing that worked for me was to wrap the geometry with a very small buffer.
SELECT ST_Intersection(
(SELECT geom FROM table1), ST_Union(ST_Buffer(geom, 0.0000001))
) FROM table2
ST_MakeValid
didn't work for me, nor did the combination of ST_Node
and ST_Dump
. The buffer didn't seem to result in any degradation in performance, but if I made it any smaller I still received a non-noded intersection error.
Ugly, but it works.
Update:
The ST_Buffer strategy seems to work well, but I ran into an issue where it produced errors when casting the geometry to geography. For example, if a point is originally at -90.0, and it is buffered by 0.0000001, it is now at -90.0000001, which is an invalid geography.
This meant that even though ST_IsValid(geom)
was t
, ST_Area(geom::geography)
returned NaN
for many features.
To avoid the non-noded intersection problem, while maintaining valid geography, I ended up using ST_SnapToGrid
like so
SELECT ST_Union(ST_MakeValid(ST_SnapToGrid(geom, 0.0001))) AS geom, common_id
FROM table
GROUP BY common_id;