2

I'm trying to clip a set of polylines to a polygon with geopandas but it's returning an empty set. Here's my code:

>>> import geopandas as gpd
>>> import matplotlib.pyplot as plt

>>> gdf_lines = gpd.read_file('by_christina_nhd_0.geojson')
>>> gdf_basin = gpd.read_file('christina_proj_4269.geojson')
>>> ax = gdf_lines.plot()
>>> gdf_basin.plot(ax=ax, color='white', linewidth=2, edgecolor='black')
>>> plt.show()

polylines with polygon overlaid

>>> gdf_lines.crs == gdf_basin.crs
True
>>>clipped_lines = gpd.clip(gdf_lines, gdf_basin)
>>>clipped_lines.empty
True

They obviously intersect, so what's up? Why isn't clip returning those lines within the polygon? Any ideas?

Here is a link to the geojson files

FYI - I've also tried this on QGIS and I get an empty set there too. It also has a warning about "invalid geometry." Not sure if that has anything to do with it.

1 Answer 1

1

Your polygon is indeed invalid. You can check it for example with ogrinfo

ogrinfo -sql "select isvalidreason(geometry) from christina_proj_4269" jsonpoly.json -dialect sqlite
...
OGRFeature(SELECT):0
  isvalidreason(geometry) (String) = Self-intersection[-75.5297803164871 39.7081441293768]

The viewer in Gist actually shows the error

enter image description here

Fix the data and try again.

1
  • That was it! Thank you
    – j sad
    Apr 14, 2021 at 20:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.