Intersection with Shapely and Fiona gives Key Error

Question

I have written a script for standalone python based on this example for doing intersections of linestrings (representing streets) and polygons (representing gridcells) using shapely/fiona. The polygons consist of square gridcells and the lines run all over the same space as the grid. My goal is to do an intersect operation similar to an intersect in e.g. ArcGIS: Get all the linestrings that are in the grid and if a linestring goes over a gridcell boundary cut it into two new linestrings, and then add the data of the respective gridcells to all the linestrings and place them to a new geodataframe. Unfortunately this script gives me a Key Error in the prop= line for everything after 'lin_length': lin_length (tested it one by one). What am I doing wrong?

This is my script so far:

from shapely.geometry import shape, mapping
import fiona

# Set filepaths    
fp_out_int = r"filepath_of_output_geodataframe"

# below filepaths are also named 'out' because of previous operations not shown in theis script
fp_out_grd = r"filepath_of_input_grid_geodataframe"
fp_out_str = r"filepath_of_input_street_geodataframe"



# schema of the new shapefile
schema =  {'geometry': 'LineString','properties': {'lin_length': 'float','gridcode': 'int','emissions': 'float','oneway': 'int'}}

# creation of the new shapefile with the intersection
with fiona.open(fp_out_int, 'w',driver='ESRI Shapefile', schema=schema) as output:
    for street in fiona.open(fp_out_str):
        for grid in fiona.open(fp_out_grd):
            if shape(street['geometry']).intersects(shape(grid['geometry'])):     
                lin_length = shape(street['geometry']).intersection(shape(grid['geometry'])).length
                prop = {'lin_length': lin_length, 'gridcode': grid['GRIDCODE'],'emissions': grid['emissions'], 'oneway': street['oneway']} 
                output.write({'geometry':mapping(shape(street['geometry']).intersection(shape(grid['geometry']))),'properties': prop})

This is the header of the street file:

id  osm_id  line_type   name    tunnel  bridge  oneway  ref z_order access  service line_class  geometry

This is the header of the grid file:

FID_de_box  ID  GRIDCODE    FID_select  OBJECTID    NUTS_ID STAT_LEVL_  AREA    LEN Shape_Leng  Shape_Area  EUROPA_ID   emissions   geometry

This script produces two shapefiles like the ones I use, just a lot smaller. I tested it and with this it gives me the same error.

from shapely.geometry import Polygon, LineString, shape, mapping
import geopandas as gpd
import numpy as np

# data for gridcells
pol1 =Polygon(([23.56101188709934, 38.07561887452059], [23.54981601299715, 38.07722115566696], [23.55161511615942, 38.0861348148835], [23.56281239269627, 38.08453228872965], [23.56101188709934, 38.07561887452059]))
pol2 =Polygon(([23.57220730394376, 38.07401528299584], [23.56101188709934, 38.07561887452059], [23.56281239269627, 38.08453228872965], [23.57400921171318, 38.08292845201315], [23.57220730394376, 38.07401528299584]))
pol3 =Polygon(([23.58340226320573, 38.0724103811776], [23.57220730394376, 38.07401528299584], [23.57400921171318, 38.08292845201315], [23.58520557288525, 38.08132330481892], [23.58340226320573, 38.0724103811776]))
pol4 =Polygon(([23.59459676456072, 38.07080416915088], [23.58340226320573, 38.0724103811776], [23.58520557288525, 38.08132330481892], [23.59640147588782, 38.07971684723199], [23.59459676456072, 38.07080416915088]))
pol5 =Polygon(([23.60579080768423, 38.06919664700065], [23.59459676456072, 38.07080416915088], [23.59640147588782, 38.07971684723199], [23.6075969203962, 38.07810907933746], [23.60579080768423, 38.06919664700065]))
# data for linestrings
lin1 =LineString(([23.55246416432669, 38.08328984787229], [23.55255284486228, 38.08461385329758], [23.55258947377914, 38.08514006917881], [23.55263934610301, 38.08585747627143], [23.55266943713539, 38.08628998247519], [23.55269047571238, 38.08667647203052], [23.55269843852039, 38.0867825869247], [23.55276943324026, 38.08772747886945], [23.55281452787935, 38.08822368753732], [23.5528445350927, 38.08856348989195], [23.55295182345333, 38.08898568635482], [23.5530637218607, 38.08932205212918], [23.55309171741729, 38.08939463941064], [23.55349756916894, 38.09044933428697]))
lin2 =LineString(([23.55247112130633, 38.08275122677435], [23.55256038857513, 38.08275684264947], [23.55276499083158, 38.08276706857134], [23.55344258388413, 38.08280218874566], [23.55453432677248, 38.08284845685117], [23.55510739749246, 38.08283236359708], [23.55568038439338, 38.08274963421275], [23.55628044484158, 38.08261627813326], [23.55687044700596, 38.08240547326844], [23.55748945055532, 38.08206022267674], [23.55777083104488, 38.08186299649506]))
lin3 =LineString(([23.55816084099953, 38.08171446917083], [23.5580297480339, 38.08186215830474], [23.55792178912102, 38.08202116300794], [23.55783864064153, 38.08218896870949]))
lin4 =LineString(([23.55240976577508, 38.08312992115972], [23.5523812673043, 38.08304115680514], [23.55232519237208, 38.08296102581079], [23.55224447464451, 38.0828950602328], [23.55214472999674, 38.08284803775601], [23.55203283158937, 38.08282297586553]))
lin5= LineString(([23.55203283158937, 38.08282297586553], [23.55074092885323, 38.0826831657206], [23.54840497625781, 38.08236029481034], [23.54598419880227, 38.08202937727305], [23.54427546402133, 38.08179158268001], [23.54340517101483, 38.08167649914944], [23.5429032626528, 38.08159083609902], [23.54245818359428, 38.081500982097], [23.54202182171505, 38.08134432432669], [23.54174614091968, 38.08121608120814]))

# Creating grid gdf
poly_gdf = gpd.GeoDataFrame(geometry=[pol1, pol2, pol3, pol4, pol5])

poly_gdf['emissions'] = None
poly_gdf['emissions'] = np.random.randint(1,25, poly_gdf.shape[0])

poly_gdf['GRIDCODE'] = None
poly_gdf['GRIDCODE'] = np.random.randint(10000,99999, poly_gdf.shape[0])

# Creating LineString gdf
line_gdf = gpd.GeoDataFrame(geometry=[lin1, lin2, lin3, lin4, lin5])

line_gdf['oneway'] = None
line_gdf['oneway'] = np.random.randint(0,1, poly_gdf.shape[0])


# save the data
poly_gdf.to_file(fp_out_grd)
line_gdf.to_file(fp_out_str)

This is the Error message i get:

Answer 1

The keys that are not being found are within the 'properties' dict nested within grid. So try changing that line like this:

prop = {'lin_length': lin_length, 'gridcode': grid['properties']['GRIDCODE'],'emissions': grid['properties']['emissions'], 'oneway': street['properties']['oneway']}