3

In trying to answer Writing selected features to shapefile using OGR and Python? I ran into problems with geopandas.GeoSeries intersects methods.

I have the following line and polygon:

import pandas as pd
import geopandas
import pickle

line = pickle.loads(
    b'\x80\x03cgeopandas.geodataframe\nGeoDataFrame\nq\x00)\x81q\x01}q\x02(X\x05\x00\x00\x00_dataq\x03cpandas.core.internals.managers\nBlockManager\nq\x04)\x81q\x05(]q\x06(cpandas.core.indexes.base\n_new_Index\nq\x07cpandas.core.indexes.base\nIndex\nq\x08}q\t(X\x04\x00\x00\x00dataq\ncnumpy.core.multiarray\n_reconstruct\nq\x0bcnumpy\nndarray\nq\x0cK\x00\x85q\rC\x01bq\x0e\x87q\x0fRq\x10(K\x01K\x02\x85q\x11cnumpy\ndtype\nq\x12X\x02\x00\x00\x00O8q\x13K\x00K\x01\x87q\x14Rq\x15(K\x03X\x01\x00\x00\x00|q\x16NNNJ\xff\xff\xff\xffJ\xff\xff\xff\xffK?tq\x17b\x89]q\x18(X\x03\x00\x00\x00FIDq\x19X\x08\x00\x00\x00geometryq\x1aetq\x1bbX\x04\x00\x00\x00nameq\x1cNu\x86q\x1dRq\x1eh\x07cpandas.core.indexes.range\nRangeIndex\nq\x1f}q (h\x1cNX\x05\x00\x00\x00startq!K\x00X\x04\x00\x00\x00stopq"K\x01X\x04\x00\x00\x00stepq#K\x01u\x86q$Rq%e]q&(h\x0bh\x0cK\x00\x85q\'h\x0e\x87q(Rq)(K\x01K\x01K\x01\x86q*h\x12X\x02\x00\x00\x00i8q+K\x00K\x01\x87q,Rq-(K\x03X\x01\x00\x00\x00<q.NNNJ\xff\xff\xff\xffJ\xff\xff\xff\xffK\x00tq/b\x89C\x08\x00\x00\x00\x00\x00\x00\x00\x00q0tq1bh\x0bh\x0cK\x00\x85q2h\x0e\x87q3Rq4(K\x01K\x01K\x01\x86q5h\x15\x89]q6cshapely.geometry.linestring\nLineString\nq7)Rq8CI\x01\x02\x00\x00\x00\x04\x00\x00\x00\xd0\x15\xb3\x1c\xf9\xb9\xf2\xbf\xd2\xd7<\x95B\xfa\xda?`\xf3T\n\xe9k\xee\xbf\xc8\xb0ZoA\x19\xc6\xbf\x12/\t\x08\'\xe2\xe5\xbf\x88\x97\x04\x84\x13\xf1\xc2?\x96#?W\xccr\xe4\xbf\xe8\x8aY\x8e\xfc\\\xe1?q9batq:be]q;(h\x07h\x08}q<(h\nh\x0bh\x0cK\x00\x85q=h\x0e\x87q>Rq?(K\x01K\x01\x85q@h\x15\x89]qAh\x19atqBbh\x1cNu\x86qCRqDh\x07h\x08}qE(h\nh\x0bh\x0cK\x00\x85qFh\x0e\x87qGRqH(K\x01K\x01\x85qIh\x15\x89]qJh\x1aatqKbh\x1cNu\x86qLRqMe}qNX\x06\x00\x00\x000.14.1qO}qP(X\x04\x00\x00\x00axesqQh\x06X\x06\x00\x00\x00blocksqR]qS(}qT(X\x06\x00\x00\x00valuesqUh)X\x08\x00\x00\x00mgr_locsqVcbuiltins\nslice\nqWK\x00K\x01K\x01\x87qXRqYu}qZ(hUh4hVhWK\x01K\x02K\x01\x87q[Rq\\ueustq]bX\x04\x00\x00\x00_typq^X\t\x00\x00\x00dataframeq_X\t\x00\x00\x00_metadataq`]qa(X\x03\x00\x00\x00crsqbX\x15\x00\x00\x00_geometry_column_nameqcehb}qdX\x04\x00\x00\x00initqeX\t\x00\x00\x00epsg:4326qfshch\x1aub.'
)

poly = pickle.loads(
    b'\x80\x03cgeopandas.geodataframe\nGeoDataFrame\nq\x00)\x81q\x01}q\x02(X\x05\x00\x00\x00_dataq\x03cpandas.core.internals.managers\nBlockManager\nq\x04)\x81q\x05(]q\x06(cpandas.core.indexes.base\n_new_Index\nq\x07cpandas.core.indexes.base\nIndex\nq\x08}q\t(X\x04\x00\x00\x00dataq\ncnumpy.core.multiarray\n_reconstruct\nq\x0bcnumpy\nndarray\nq\x0cK\x00\x85q\rC\x01bq\x0e\x87q\x0fRq\x10(K\x01K\x02\x85q\x11cnumpy\ndtype\nq\x12X\x02\x00\x00\x00O8q\x13K\x00K\x01\x87q\x14Rq\x15(K\x03X\x01\x00\x00\x00|q\x16NNNJ\xff\xff\xff\xffJ\xff\xff\xff\xffK?tq\x17b\x89]q\x18(X\x03\x00\x00\x00FIDq\x19X\x08\x00\x00\x00geometryq\x1aetq\x1bbX\x04\x00\x00\x00nameq\x1cNu\x86q\x1dRq\x1eh\x07cpandas.core.indexes.range\nRangeIndex\nq\x1f}q (h\x1cNX\x05\x00\x00\x00startq!K\x00X\x04\x00\x00\x00stopq"K\x06X\x04\x00\x00\x00stepq#K\x01u\x86q$Rq%e]q&(h\x0bh\x0cK\x00\x85q\'h\x0e\x87q(Rq)(K\x01K\x01K\x06\x86q*h\x12X\x02\x00\x00\x00i8q+K\x00K\x01\x87q,Rq-(K\x03X\x01\x00\x00\x00<q.NNNJ\xff\xff\xff\xffJ\xff\xff\xff\xffK\x00tq/b\x89C0\x00\x00\x00\x00\x00\x00\x00\x00\x01\x00\x00\x00\x00\x00\x00\x00\x02\x00\x00\x00\x00\x00\x00\x00\x03\x00\x00\x00\x00\x00\x00\x00\x04\x00\x00\x00\x00\x00\x00\x00\x05\x00\x00\x00\x00\x00\x00\x00q0tq1bh\x0bh\x0cK\x00\x85q2h\x0e\x87q3Rq4(K\x01K\x01K\x06\x86q5h\x15\x89]q6(cshapely.geometry.polygon\nPolygon\nq7)Rq8CM\x01\x03\x00\x00\x00\x01\x00\x00\x00\x04\x00\x00\x00^\x12\x10N\xc4K\xf2\xbf\xe8\x8aY\x8e\xfc\\\xd1?\xb7\xbf\xc7\xee\xd0\xf6\xe7\xbf\x0e\x8c:\xd3\xb8\x81\xd1?\x84\x13\xf1\x92\x80p\xf2\xbf`\xf3T\n\xe9k\xce\xbf^\x12\x10N\xc4K\xf2\xbf\xe8\x8aY\x8e\xfc\\\xd1?q9bh7)Rq:CM\x01\x03\x00\x00\x00\x01\x00\x00\x00\x04\x00\x00\x00\xf1\x92\x80p"^\xf2\xbf`\xf3T\n\xe9k\xce\xbfJ@8\x11/\t\xe8\xbf\x9c\x88\x97\x04\x84\x13\xd1?\xfe=v\x87\xb6\xbf\xe7\xbf\xa8\xf5\x16\x94a\xb5\xce\xbf\xf1\x92\x80p"^\xf2\xbf`\xf3T\n\xe9k\xce\xbfq;bh7)Rq<CM\x01\x03\x00\x00\x00\x01\x00\x00\x00\x04\x00\x00\x00\xb7\xbf\xc7\xee\xd0\xf6\xe7\xbf\x9c\x88\x97\x04\x84\x13\xd1? \xbbC\xdb\xdfc\xd7\xbf\x9c\x88\x97\x04\x84\x13\xd1?\xd8<\x95B\xfa\x9a\xe7\xbf\xf8\xf7\xd8\x1d\xda\xfe\xce\xbf\xb7\xbf\xc7\xee\xd0\xf6\xe7\xbf\x9c\x88\x97\x04\x84\x13\xd1?q=bh7)Rq>CM\x01\x03\x00\x00\x00\x01\x00\x00\x00\x04\x00\x00\x00\xfa\xb9b\x96#?\xe7\xbf`\xf3T\n\xe9k\xce\xbf\xd4\xb8\x81Qg\x1a\xd7\xbf\xde\x822\xac\xd6[\xd0?\xd4\xb8\x81Qg\x1a\xd7\xbf@\xfa\x9a\xa7RH\xcf\xbf\xfa\xb9b\x96#?\xe7\xbf`\xf3T\n\xe9k\xce\xbfq?bh7)Rq@CM\x01\x03\x00\x00\x00\x01\x00\x00\x00\x04\x00\x00\x00\x88\xb6\xbf\xc7\xee\xd0\xd6\xbf\xf8\xf7\xd8\x1d\xda\xfe\xce\xbf\x88\xb6\xbf\xc7\xee\xd0\xd6\xbf*\x85\xf45O\xa5\xd0?\xa0\xe5\xc8\xcf\x15\xb3\xac?4\x8d\x1b\x18u\xa6\xd1?\x88\xb6\xbf\xc7\xee\xd0\xd6\xbf\xf8\xf7\xd8\x1d\xda\xfe\xce\xbfqAbh7)RqBCM\x01\x03\x00\x00\x00\x01\x00\x00\x00\x04\x00\x00\x00<\xb4\xfd=v\x87\xd6\xbf\xa8\xf5\x16\x94a\xb5\xce\xbf\x00\xca\xb0ZoA\xa9?\x04\x84\x13\xf1\x92\x80\xd0?@\xa5\x90\xbe\xe6\xa9\xa4?@\xfa\x9a\xa7RH\xcf\xbf<\xb4\xfd=v\x87\xd6\xbf\xa8\xf5\x16\x94a\xb5\xce\xbfqCbetqDbe]qE(h\x07h\x08}qF(h\nh\x0bh\x0cK\x00\x85qGh\x0e\x87qHRqI(K\x01K\x01\x85qJh\x15\x89]qKh\x19atqLbh\x1cNu\x86qMRqNh\x07h\x08}qO(h\nh\x0bh\x0cK\x00\x85qPh\x0e\x87qQRqR(K\x01K\x01\x85qSh\x15\x89]qTh\x1aatqUbh\x1cNu\x86qVRqWe}qXX\x06\x00\x00\x000.14.1qY}qZ(X\x04\x00\x00\x00axesq[h\x06X\x06\x00\x00\x00blocksq\\]q](}q^(X\x06\x00\x00\x00valuesq_h)X\x08\x00\x00\x00mgr_locsq`cbuiltins\nslice\nqaK\x00K\x01K\x01\x87qbRqcu}qd(h_h4h`haK\x01K\x02K\x01\x87qeRqfueustqgbX\x04\x00\x00\x00_typqhX\t\x00\x00\x00dataframeqiX\t\x00\x00\x00_metadataqj]qk(X\x03\x00\x00\x00crsqlX\x15\x00\x00\x00_geometry_column_nameqmehl}qnX\x04\x00\x00\x00initqoX\t\x00\x00\x00epsg:4326qpshmh\x1aub.'
)

Which when plotted together look like this:

poly.plot(ax=line.plot(), color='red')

Line and polygons intersecting

I want to calculate and return which polygons in the triangular strip intersects with the line. I tried both of the following:

poly.intersects(line)
line.intersects(poly)

Both return True for the first item but False for the remaining 5 polygons:

0     True
1    False
2    False
3    False
4    False
5    False
dtype: bool

I had expected that it would return True for polygons 0, 1, and 2 and False for the rest but clearly this is not the case. From the GeoSeries intersects documentation I thought it would return a result based on each item in the Series.

I tried to instead iterate over each row:

for row in poly.iterrows():
    row = geopandas.GeoSeries(row)
    line.intersects(row)

But this returns the following error:

AttributeError: 'int' object has no attribute 'is_empty'

How should I calculate which polygons the line intersects with, so that my returned result is something like the following?

0     True
1     True
2     True
3    False
4    False
5    False
3

Not sure but I think you are running out of lines after the first comparision of polygon1-line1, second comparison will be polygon2-nothing.

If you for each polygon compare all lines (in your case one) you will get the results you want:

poly.geometry.map(lambda x: x.intersects(line.geometry.any()))

0     True
1     True
2     True
3    False
4    False
5    False

It is probably better to use Spatial Join.

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.