I have recently started using Geopandas and am having a problem with the performance of intersects. The intersect tests whether land parcels have a road within (or across the boundary).

Parcel data (parcels) is a polygon geodataframe and roads (roads) are as a line geodataframe. I attempt to find out if the parcel interacts with a road with the follwoing code:

road = roads.unary_union
parcels['road_intersection'] = parcels.intersects(road)

Here I generate get a column 'road_intersection' with True/False values.

with 985 parcel polygons and 20854 road lines, this takes over a minute to generate a result.

I have previously tried with

parcels['road_intersection'] = parcels.intersects(roads)

but this did not seem to test all roads with all parcels.

Would anyone be able to provide any tips for correct use or optimisation?


For intersections in GeoPandas, it is better to use a spatial-join (see More Efficient Spatial join in Python without QGIS, ArcGIS, PostGIS, etc or How to efficiently determine which of thousands of polygons intersect with a linestring [duplicate]

sample image

import geopandas as gpd
parcels = gpd.read_file('parcels.shp')
roads = gpd.read_file('roads.shp')
intersections= gpd.sjoin(parcels, roads, how="left", op='intersects')
    parcel                  geometry                         index_right  road
0  Parcel 1  POLYGON ((-0.6824583866837387 0.78233034571062...  0.0       Road 1  
1  Parcel 2  POLYGON ((-0.09859154929577452 0.3239436619718...  NaN       Nan
2  Parcel 3  POLYGON ((-0.9103713188220229 -0.1062740076824...  1.0       Road 2  
3  Parcel 4  POLYGON ((0.2266325224071704 0.620998719590268...  0.0       Road 1

With your solution

road = roads.unary_union
parcels['road_intersection'] = parcels.intersects(road)
    parcel                  geometry                         index_left  road_intersection  
0  Parcel 1  POLYGON ((-0.6824583866837387 0.78233034571062...   0         True 
1  Parcel 2  POLYGON ((-0.09859154929577452 0.3239436619718...   1         False
2  Parcel 3  POLYGON ((-0.9103713188220229 -0.1062740076824...   2         True
3  Parcel 4  POLYGON ((0.2266325224071704 0.620998719590268...   3         True
  • parcels['road_intersection'] = parcels.intersects(road)
    – Clay
    Sep 12 '18 at 20:05
  • 1
    This is not an intersect. It is a spatial join where the matches are made at intersections. OP was asking about an intersect.
    – user32882
    Mar 27 '19 at 7:07

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.