Tuesday, 16 May 2017

postgis - Clip overlapping polygon in equal sharing size



Within single shapefile I am trying to clip overlapping polygons to their boundaries using centre line of overlapped area using Postgis or ogr2ogr.


Sample


I tried to use ogr2ogr command to remove slivers


ogr2ogr dissolved.shp input.shp -dialect sqlite \
-sql "select ST_union(ST_buffer(Geometry, 10)),common_attribute from input GROUP BY common_attribute"

But every polygon in my shapefile got only unique attribiutions so I cannot use ST_UNION.


I tried also to use ST_Difference postgreSQL database with the Postgis extension but also without success.


My main problem is that i can not find method to clip polygons on equal sharing size.



Answer




Here is a solution in postgis that fits the use-case (image) that you describe. There are some assumptions though:



  1. It assumes simple (non multi) polygons

  2. Every polygon is only intersected by 1 other polygon (you can fix that, but it requires more coding).

  3. There will be a straight line between the cutting points of the polygon (that may not be the center line that you are looking for)




WITH polygons AS (
SELECT 1 as id, ST_Buffer(ST_MakePoint(10,10),3) geom
UNION ALL

SELECT 2 as id, ST_Buffer(ST_MakePoint(5,10),3) geom
--Uncomment this to see what happens when multiple polygons overlap
--UNION ALL
--SELECT 3 as id, ST_Buffer(ST_MakePoint(7.5,7.5),3) geom
)
,cut_polygons AS (
SELECT a.id,
ST_LineMerge(
ST_Difference(ST_ExteriorRing(a.geom), b.geom)
) geom

FROM polygons a
LEFT JOIN polygons b ON (ST_Intersects(a.geom, b.geom) AND a.id != b.id)

)
SELECT id, ST_MakePolygon(ST_AddPoint(geom, ST_StartPoint(geom))) geom
FROM cut_polygons;

The trick is that the polygons are first converted to linestrings (ST_ExteriorRing) that are cookie-cut by an overlapping polygon (ST_Difference). The result is a linestring (ST_LineMerge) with a gap. When you fill that gap with a straight line (ST_AddPoint) you have a polygon (ST_MakePolygon) that is more or less split at the center of the overlap.


No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...