Tuesday 2 April 2019

python - How to find the intersection areas of overlapping buffer zones in single shapefile?



Find the intersections of the buffer zones


I found an example where they tried finding the intersection of two different shapefiles -Intersect Shapefiles using Shapely. I've tried to adapt my code to do the same intersection, but I am only working with one shapefile with buffer features inside. The program is able to write the shapefiles, but the intersection aren't occurring, thus no data written to my output shapefile Anyone willing to help me understand where the problem lies or has an alternate way of finding intersection in a single shapefile?


# creation of the new shapefile with the intersection
for i in threshold:
with fiona.open("threshold_buffer_shape_"+str(i)+".shp","r") as layer1:
with fiona.open("threshold_buffer_shape_"+str(i)+".shp","r") as layer2:
# We copy schema and add the new property for the new resulting shp
schema=layer2.schema.copy()
schema['properties']['uGID'] = 'int:10'
# create an empty spatial index object

with fiona.open("intersection_region_"+str(i)+".shp","w","ESRI Shapefile",schema) as layer3:
index = rtree.index.Index()
# populate the spatial index
for feat1 in layer1:
fid = int(feat1['properties']['GID'])
geom1 = shape(feat1['geometry'])
index.insert(fid,geom1.bounds)

# get list of fids where bounding boxes intersect
for feat2 in layer2:

geom2=shape(feat2['geometry'])

for fid in list(index.intersection(geom2.bounds)):
if fid!= int(feat2['properties']['GID']):
feat1 = layer1[fid]
geom1 = shape(feat1['geometry'])
if geom1.intersects(geom2):
props = feat2['properties']
props['uGID'] = feat1['properties']['uGID']
layer3.write({

'properties':props,
'geometry':mapping(geom1.intersection(geom2))
})

enter image description here



Overlapping buffer regions





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...