20

My question is an extension of Vertical lines in a polygon shapefile. Kindly refer to that question first.

What you will see is a method of generating vertical lines with respect to the bounding box, at user-defined spacing. I understand that OGR, Fiona, Shapely etc. can be used to do the next step of clipping, but I do not understand their utilization.

How do I read one line of a polygon shapefile? Every application that uses Shapely shows how to generate the LineString, Point or Polygon but never to read an existing shapefile

Kindly assist me with at-least a skeleton structure so I can build on it.

38

1) read your shapefile with Fiona, PyShp, ogr or ...using the geo_interface protocol (GeoJSON):

with Fiona

import fiona
shape = fiona.open("my_shapefile.shp")
print shape.schema
{'geometry': 'LineString', 'properties': OrderedDict([(u'FID', 'float:11')])}
#first feature of the shapefile
first = shape.next()
print first # (GeoJSON format)
{'geometry': {'type': 'LineString', 'coordinates': [(0.0, 0.0), (25.0, 10.0), (50.0, 50.0)]}, 'type': 'Feature', 'id': '0', 'properties': OrderedDict([(u'FID', 0.0)])}

with PyShp

import shapefile
shape = shapefile.Reader("my_shapefile.shp")
#first feature of the shapefile
feature = shape.shapeRecords()[0]
first = feature.shape.__geo_interface__  
print first # (GeoJSON format)
{'type': 'LineString', 'coordinates': ((0.0, 0.0), (25.0, 10.0), (50.0, 50.0))}

with ogr:

from osgeo import ogr
file = ogr.Open("my_shapefile.shp")
shape = file.GetLayer(0)
#first feature of the shapefile
feature = shape.GetFeature(0)
first = feature.ExportToJson()
print first # (GeoJSON format)
{"geometry": {"type": "LineString", "coordinates": [[0.0, 0.0], [25.0, 10.0], [50.0, 50.0]]}, "type": "Feature", "properties": {"FID": 0.0}, "id": 0}

2) conversion to Shapely geometry (with the shape function)

# now use the shape function of Shapely
from shapely.geometry import shape
shp_geom = shape(first['geometry']) # or shp_geom = shape(first) with PyShp)
print shp_geom
LINESTRING (0 0, 25 10, 50 50)
print type(shp_geom)
<class 'shapely.geometry.linestring.LineString'>

3) computations

4) save the resulting shapefile

  • 5
    I would add geopandas to the list: geopandas.read_file("my_shapefile.shp") – joris Aug 27 '17 at 7:30
  • As of GDAL 2.0, instead of osgeo.ogr.Open, use osgeo.gdal.OpenEx (details). – Kevin Nov 10 '17 at 11:28
  • with ogr, I first had to define the json as a json to be able to process it further with shapely: 'first = json.loads(first)' – Leo Sep 12 at 7:41
9

I find geopandas as the best performer here. Code:

import geopandas as gpd
shapefile = gpd.read_file("shapefile.shp")
print(shapefile)

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.