6

Let's say there is a table 'polygon_a' and a table 'polygon_b'. Each table contains three polygon features. The polygon features are spatially equal.

CREATE TABLE polygon_a (
gid serial NOT NULL,
geom geometry(polygon, your_SRID),
CONSTRAINT polygon_a_pkey PRIMARY KEY (gid)
);

CREATE TABLE polygon_b (
gid serial NOT NULL,
geom geometry(polygon, your_SRID),
CONSTRAINT polygon_b_pkey PRIMARY KEY (gid)
);

enter image description here

Now I've created a view 'equals_false' to select all features from table 'polygon_a' where ST_Equals(polygon_a.geom, polygon_b.geom) is false. After using the QGIS 'Split features' tool to split a feature from table 'polygon_a' the view 'equals_false' returns two results.

CREATE VIEW equals_false AS SELECT
polygon_a.gid,
polygon_a.geom
FROM
polygon_a LEFT JOIN polygon_b
ON ST_Equals(polygon_a.geom, polygon_b.geom)
WHERE polygon_b.gid IS NULL;

enter image description here enter image description here

But after using the 'Merge selected features' tool to merge the splitted features 'equals false' still returns a result. What's the reason of this behaviour?

enter image description here

  • Query ST_AsText(geom) directly from PostGIS and you will see Well Known Text presentation of the geometries. It looks like the merge tool creates a geometry that is not exactly the same as the one you splitted. – user30184 Feb 13 '16 at 12:26
6

When you split a polygon with a line, nodes are added to the boundary at the point where the line intersects the polygon shell. Because the coordinates of these intersection points can't be exactly represented with floating point numbers (except in special cases), the boundary actually deflects slightly at these points. When these two halves are dissolved back together, the added nodes remain, and the reconstituted polygon is ever-so-slightly different from the original.

  • Is there a way to remove the added nodes using PostGIS? – Lunar Sea Feb 13 '16 at 15:45
1

As mentioned by 'dbaston' ST_Equals() doesn't return the expected results because of the added nodes. But using the hausdorff distance to compare the geometries is working.

CREATE VIEW equals_false AS SELECT DISTINCT
polygon_a.gid, polygon_a.geom
FROM
polygon_a, polygon_b
WHERE ST_Intersects(polygon_a.geom, polygon_b.geom)
AND NOT ST_Touches(polygon_a.geom, polygon_b.geom)
AND ST_HausdorffDistance(polygon_a.geom, polygon_b.geom) > 0.1;

I'm just not sure about the proper threshold.

1

As answered here, use ST_Snap to get exact noding of input datasets required to get expected overlay results, with something like:

SELECT ST_Equals(ST_Snap(A, B, tol), ST_Snap(B, A, tol))
FROM ...

using a tol tolerance to suit your dataset's precision.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.