CODEWITHSUNDEEP
pythonshapefilegeopandas
asheets
asheets

Reputation: 860

How to convert a shapefile to a complete list of latitude and longitude points

I'm trying to convert a shapefile to a list of latitude and longitude points that represent every point defined by the shapefile. Reading the file with geopandas and using the .plot() function displays these points as a graph, but I would like the raw points. I attempted to iterate through the polygons in geopandas .geometry and store all the points within the polygons. I plotted these points to test if they gave an accurate representation of the area, but they do not. I did all this with the following code:

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

def geoToList(geodataframe):
    points = []
    for s in geodataframe.geometry:iq
        s = str(s)
        s = re.sub('[^0-9., ]+', '', s).split(',')
        s = map(lambda x: x.strip(), s)
        s = map(lambda x: (float(x.split()[0]), float(x.split()[1])), s)
        points.extend(list(s))   
    return points

habitat = gpd.read_file('desktop/species_19377/species_19377.shp')
#borough = borough.to_crs(epsg=4326)

points = geoToList(habitat)
x = [point[0] for point in points]
y = [point[1] for point in points]

plt.scatter(x, y)
plt.show() #representation of the points in all polygons
habitat.plot() #representtation of the points I want

I would like some function that returns a list of points that could be plotted and look identical to the output of habitat.plot()

My next idea is to store the graphs as images and assign the pixel values latitude and longitude values based on the scale of the graph, but I'm sure this is more complicated than it needs to be.

Any help would be appreciated!

Upvotes: 1

Views: 5512

Answers (1)

joris
joris

Reputation: 139162

To extract all points from a set of Polygons / MultiPolygons, you could do something like this:

from shapely.geometry import MultiPolygon

def points_from_polygons(polygons):
    points = []
    for mpoly in polygons:
        if isinstance(mpoly, MultiPolygon):
            polys = list(mpoly)
        else:
            polys = [mpoly]
        for polygon in polys:
            for point in polygon.exterior.coords:
                points.append(point)
            for interior in polygon.interiors:
                for point in interior.coords:
                    points.append(point)
    return points

points = points_from_polygons(habitat.geometry)
x = [point.x for point in points]
y = [point.y for point in points]

Upvotes: 1

Related Questions