8

Is it possible to create buffers (in my case, but polygons in general also) that would create a boundary at their point of contact ? (so that they don't overlap on each other).

Like this (from left to right)

enter image description here

I am trying to solve this question where I need to join points representing addresses to their nearby parcel without including points related to adjacent parcels or parcels on the other side of the street.

It seems to me that there are two main options: either using voronoi/thiessen diagram but based on polygons, not points, or creating many buffers (making sure they overlap with the address points outside on the street), and then making these new polygons crop each other at the medium of their intersection line (which would be the median line of the street for street facing vertices, and the property boundary line for the vertice between two adjacent parcels/plots of land).

Spatially join attributes from point to polygon side in QGIS3 with numerous contiguous polygons and numerous adjacent points

Philippe Morgan
  • 633
  • 4
  • 11
  • Basically a tool to "grow" the polygons until they find each other's limits. Not sure if this exists. – Philippe Morgan Sep 06 '20 at 19:53
  • I think I could use something called Thiessen polygons if this applies also to polygons and not just to points. – Philippe Morgan Sep 06 '20 at 20:13
  • 1
    May I ask why do you comment under your question instead of editing your question with nesting those comments or setting up a chat? – Taras Sep 07 '20 at 08:57
  • Good idea, I will do this, not very familiar with the forum yet. – Philippe Morgan Sep 07 '20 at 09:08
  • 3
    Good question!! – menes Sep 07 '20 at 10:11
  • 2
    Try to focus your efforts on 3 main details: 1) Work only with the outer edges, creating a combined one-way buffer; 2) Work only with the inner central horizontal edges and create a central line from them, which will separate the inner polygon; 3) Work with the inner vertical lines, increasing them to the buffer edges and cross them; 4) Collect the necessary buffers, combining them if necessary ... somehow.... – Cyril Mikhalchenko Sep 07 '20 at 18:30
  • I have solved the issue that prompted this question in another way (see linked question). That being said I agree that it could still be interesting to learn how to do that, I will think about it more! – Philippe Morgan Sep 09 '20 at 18:06
  • Not an answer yet, but one might want to isolate vertices closest to each other (like parcels facing each other), then use something like "concave hull", then "v.voronoi.skeleton" to find a mean line between the two lines, as described here: https://gis.stackexchange.com/questions/123731/create-mean-line-from-multiple-lines-using-qgis/271109 However I have no idea how to automate the fact that a mean line is only created between lines close to each other, not between all the lines of all features. – Philippe Morgan Sep 09 '20 at 18:51

1 Answers1

2

This is an approach to this problem using that strategy of the 'inner central horizontal edge'. In resume:

  1. Create buffers of the polygons.
  2. Intersection of buffers to get overlaping zones.
  3. Extract intersection coords to create the new medium edge.
  4. Clip buffers with lines and assign the parcel ID.

Here are the layers to try out.

Here is the code:

library(sf)
library(mapview)

read data


pol <- st_read('./data/Parcels.gpkg')
mapview(pol)

enter image description here

# make the buffer
polbuf <- st_buffer(pol, 10, 0)
mapview(pol)+mapview(st_buffer(pol, 10, 0), alpha.regions=0.1)

enter image description here

# get intersection points
p_intersect <- gIntersection(as(as_Spatial(polbuf)[1,], "SpatialLines"), as_Spatial(polbuf))
mapview(p_intersect)

extract intersection coordinates


adapted from: https://stackoverflow.com/questions/55106950/r-sf-intersect-lines-with-the-borders-of-multipolygons-extract-coordinates-of


line_1 <- st_as_sf(as(as_Spatial(polbuf), "SpatialLines"))
poly_1 <- st_as_sf(as(as_Spatial(polbuf), "SpatialLines"))
pnts   <- st_intersection(line_1, st_cast(poly_1, "MULTILINESTRING", group_or_split = FALSE))


catch points and get unique (as they will be duplicated...)


ps <- unique(pnts[st_geometry_type(pnts)=='MULTIPOINT',])
mapview(ps)+mapview(pol)+mapview(st_buffer(pol, 10, 0), alpha.regions=0.1)

enter image description here

# Create clipping lines with this points
# catched here: https://gis.stackexchange.com/questions/270725/r-sf-package-points-to-multiple-lines-with-st-cast
ls <- st_cast(ps, "LINESTRING")
mapview(ls)+mapview(ps)+mapview(pol)+mapview(st_buffer(pol, 10, 0), alpha.regions=0.1)

enter image description here

# split buffers with lines and consolidate
# cathed here: https://stackoverflow.com/questions/5726297/cut-polygons-using-contour-line-beneath-the-polygon-layers
library(rgeos)
lpi <- gIntersection(as_Spatial(polbuf), as_Spatial(ls))
blpi <- gBuffer(lpi, width = 0.000001)
dpi <- gDifference(as_Spatial(polbuf), blpi)
mapview(dpi)

enter image description here

# Asign parcel to each polygon
polbuf2 <- st_join(st_cast(st_as_sf(dpi), 'POLYGON'), pol)
mapview(pol)+mapview(polbuf2, alpha.regions = 0.2)

enter image description here

César Arquero Cabral
  • 3,439
  • 1
  • 16
  • 53