Plotting large shapefiles with matplotlib
This is not a problem of Matplotlib but your script and the module you use for reading shapefiles
1) You know that there are points in the geometries of the Polygon shapefile thus eliminate try... except
2) you load and read the shapefile twice for x and y (memory)
for shape in shp.shapeRecords():
xy = [i for i in shape.shape.points[:]]
x = [i[0] for i in xy]
y = [i[1] for i in xy]
or directly
for shape in shp.shapeRecords():
xy = [i for i in shape.shape.points[:]]
x,y = zip(*[(j[0],j[1]) for j in xy])
3) You can also use the Geo_interface (look at Plot shapefile with matplotlib)
for shape in shp.shapeRecords():
poly = shape.shape.__geo_interface__
print(poly)
{'type': 'Polygon', 'coordinates': (((203602.55736766502, 89867.47994546698), (204061.86095852466, 89822.92064187612), (203983.02526755622, 89322.48538616339), (203684.82069737124, 89031.13609345393), (203280.35932631575, 89260.78788888374), (203184.3854416585, 89624.11759508614), (203602.55736766502, 89867.47994546698)),)}
And you have the GeoJSON representation of the geometry (Polygon). You can plot the Polygon as in the reference
The LinearRing of the Polygon
x,y = zip(*[(i[0],i[1]) for i in poly['coordinates'][0]])
And the nodes of the Polygon
4) The problem of Pyshp (shapefile) is that it loads the complete shapefile into memory and if the shapefile is too big...
You can use a generator (read the layer one feature by one feature)
def records(filename):
# generator
reader = shapefile.Reader(filename)
for sr in reader.shapeRecords():
geom = sr.shape.__geo_interface__
yield geom
features = records("a_polygon.shp")
features.next()
{'type': 'Polygon', 'coordinates': (((203602.55736766502, 89867.47994546698), (204061.86095852466, 89822.92064187612), (203983.02526755622, 89322.48538616339), (203684.82069737124, 89031.13609345393), (203280.35932631575, 89260.78788888374), (203184.3854416585, 89624.11759508614), (203602.55736766502, 89867.47994546698)),)}
Or directly
shapes = shapefile.Reader('a_polygon.shp')
shapes.iterShapes().next().__geo_interface__
{'type': 'Polygon', 'coordinates': (((203602.55736766502, 89867.47994546698), (204061.86095852466, 89822.92064187612), (203983.02526755622, 89322.48538616339), (203684.82069737124, 89031.13609345393), (203280.35932631575, 89260.78788888374), (203184.3854416585, 89624.11759508614), (203602.55736766502, 89867.47994546698)),)}
5) Or use a Python module that directly uses generators/iterators :Fiona
import fiona
shapes = fiona.open("a_polygon.shp")
first = shapes.next() # for for feature in shapes
print(first)
{'geometry': {'type': 'Polygon', 'coordinates': [[(203602.55736766502, 89867.47994546698), (204061.86095852466, 89822.92064187612), (203983.02526755622, 89322.48538616339), (203684.82069737124, 89031.13609345393), (203280.35932631575, 89260.78788888374), (203184.3854416585, 89624.11759508614), (203602.55736766502, 89867.47994546698)]]}, 'type': 'Feature', 'id': '0', 'properties': OrderedDict([(u'id', None)])}
print(first['geometry']['coordinates']
[[(203602.55736766502, 89867.47994546698), (204061.86095852466, 89822.92064187612), (203983.02526755622, 89322.48538616339), (203684.82069737124, 89031.13609345393), (203280.35932631575, 89260.78788888374), (203184.3854416585, 89624.11759508614), (203602.55736766502, 89867.47994546698)]]
You can use geopandas for plotting as discussed in this answer.
You can also use pyshp as in following code
from descartes import PolygonPatch
import shapefile
sf=shapefile.Reader('shapefile')
fig = plt.figure()
ax = fig.gca()
for poly in sf.shapes():
poly_geo=poly.__geo_interface__
ax.add_patch(PolygonPatch(poly_geo, fc='#ffffff', ec='#000000', alpha=0.5, zorder=2 ))
ax.axis('scaled')
plt.show()
Using geopandas
, the code would look like:
import geopandas
df = geopandas.read_file(shpFilePath)
df.plot()
and this should also be faster (at least starting from geopandas 0.3.0)