77
Convert the DataFrame's content (e.g. Lat and Lon columns) into appropriate Shapely geometries first and then use them together with the original DataFrame to create a GeoDataFrame.
from geopandas import GeoDataFrame
from shapely.geometry import Point
geometry = [Point(xy) for xy in zip(df.Lon, df.Lat)]
df = df.drop(['Lon', 'Lat'], axis=1)
crs = {'init': '...
25
Using Panda's to_sql method and SQLAlchemy you can store a dataframe in Postgres. And since you're storing a Geodataframe, GeoAlchemy will handle the geom column for you.
Here's a code sample:
# Imports
from geoalchemy2 import Geometry, WKTElement
from sqlalchemy import *
import pandas as pd
import geopandas as gpd
# Creating SQLAlchemy's engine to use
...
22
The question is about Fiona and Shapely and the other answer using GeoPandas requires to also know Pandas. Moreover GeoPandas uses Fiona to read/write shapefiles.
I do not question here the utility of GeoPandas, but you can do it directly with Fiona using the standard module itertools, specially with the command groupby ("In a nutshell, groupby takes an ...
20
it seems that this is the right way to do that right now:
rdf = gpd.GeoDataFrame( pd.concat( dataframesList, ignore_index=True) )
19
The examples provided are for executing the scripts in a Jupyter/IPython notebooks environment.
In a normal Python environment, you need to import matplotlib to show the image
import geopandas as gpd
world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))
world.plot()
import matplotlib.pyplot as plt
plt.show()
17
You can use the cx method on a geodataframe to select rows within a bounding box. For your example frames:
xmin, ymin, xmax, ymax = df_sussex.total_bounds
sac_sussex = df_sac.cx[xmin:xmax, ymin:ymax]
From http://geopandas.org/indexing.html:
In addition to the standard pandas methods, GeoPandas also provides coordinate based indexing with the cx indexer, ...
16
adding the argument op='within' in the sjoin function speeds up the point-in-polygon operation dramatically.
Default value is op='intersects', which I guess would also lead to correct result, but is 100 to 1000 times slower.
15
You are on the right track and the geopandas GeoDataFrame is a good choice for rasterization over Fiona. Fiona is a great toolset, but I think that the DataFrame is better suited to shapefiles and geometries than nested dictionaries.
import geopandas as gpd
import rasterio
from rasterio import features
Set up your filenames
shp_fn = '...
15
I think I found an interim solution, which I'm posting in case it's useful for anyone:
import pandas as pd
import numpy as np
from geopandas import GeoDataFrame
from shapely.geometry import Point, LineString
# Zip the coordinates into a point object and convert to a GeoDataFrame
geometry = [Point(xy) for xy in zip(df.lon, df.lat)]
df = GeoDataFrame(df, ...
15
You can directly use the Shapely function Nearest points (the geometries of the GeoSeries are Shapely geometries):
from shapely.ops import nearest_points
# unary union of the gpd2 geomtries
pts3 = gpd2.geometry.unary_union
def near(point, pts=pts3):
# find the nearest point and return the corresponding Place value
nearest = gpd2.geometry == ...
15
p1 = Point((1,2))
p2 = Point((5,6))
df = pd.DataFrame({'a': [11,22]})
gdf = gpd.GeoDataFrame(df, geometry = [p1,p2])
gdf
#out:
# a geometry
#0 11 POINT (1 2)
#1 22 POINT (5 6)
You can directly assign the buffer as a new geometry column to your GeoDataFrame:
gdf['geometry'] = gdf.geometry.buffer(2)
#out:
# a geometry
#0 11 POLYGON ((3 2, 2....
14
Shapely uses a cartesian plane system for computing geometries (distance = euclidean distance)
Shapely does not support coordinate system transformations. All operations on two or more features presume that the features exist in the same Cartesian plane.
GeoPandas uses Fiona to read shapefiles (and others) and Pyproj for cartographic transformations.
...
14
One-liners! Plus some performance pointers for big-data people.
Given a pandas.DataFrame that has x Longitude and y Latitude like so:
df.head()
x y
0 229.617902 -73.133816
1 229.611157 -73.141299
2 229.609825 -73.142795
3 229.607159 -73.145782
4 229.605825 -73.147274
Let's convert the pandas.DataFrame into a geopandas.GeoDataFrame as ...
14
If the crs of the GeoDataFrame is known (EPSG:4326 unit=degree, here), you don't need Shapely, nor pyproj in your script because GeoPandas uses them).
import geopandas as gpd
test = gpd.read_file("test_wgs84.shp")
print test.crs
test.head(2)
Now copy your GeoDataFrame and change the projection to a Cartesian system (EPSG:3857, unit= m as in the answer of ...
13
Ran into this problem myself. If you want the x and y as separate GeoDataFrame columns, then this works nicely:
gdf["x"] = gdf.centroid.map(lambda p: p.x)
gdf["y"] = gdf.centroid.map(lambda p: p.y)
Starting with GeoPandas 0.3.0, you can use the provided x and y properties instead:
gdf["x"] = gdf.centroid.x
gdf["y"] = gdf.centroid.y
12
in case of csv, it probably would be easier to read it with pandas and then convert it to geopandas Dataframe
import pandas as pd
import geopandas as gp
from shapely.geometry import Point
stations = pd.read_csv('../data/stations.csv')
stations['geometry'] = stations.apply(lambda z: Point(z.X, z.Y), axis=1)
stations = gp.GeoDataFrame(stations)
12
You can pass the json directly to the GeoDataFrame constructor:
import geopandas as gpd
import requests
data = requests.get("https://data.cityofnewyork.us/api/geospatial/arq3-7z49?method=export&format=GeoJSON")
gdf = gpd.GeoDataFrame(data.json())
gdf.head()
Outputs:
features type
0 {'type': '...
11
I highly recommend GeoPandas for dealing with large assortments of features and performing bulk operations.
It extends Pandas dataframes, and uses shapely under the hood.
from geopandas import GeoSeries, GeoDataFrame
# define your directories and file names
dir_input = '/path/to/your/file/'
name_in = 'cb_2013_us_county_20m.shp'
dir_output = '/path/to/...
10
If you have large dataframes, I've found that scipy's cKDTree spatial index .query method returns very fast results for nearest neighbor searches. As it uses a spatial index it's orders of magnitude faster than looping though the dataframe and then finding the minimum of all distances. It is also faster than using shapely's nearest_points with RTree (the ...
10
To write to GeoJSON:
dataframe.to_file("output.json", driver="GeoJSON")
To write to GeoPackage:
dataframe.to_file("output.gpkg", driver="GPKG")
Documentation is here, though somewhat sparse.
9
Leaving the rest below, but the main thing was accessing the geometry properly. If iterating over rows, e.g. for index, row in zones.iterrows(): you can simply use row.geometry.centroid.x and row.geometry.centroid.y. Geometry is a special column included in a GeoDataFrame, so every row has a geometry attribute.
You are accessing that attribute, which ...
9
I just experimented with this - maybe in GeoPandas 0.2.1 and Pandas 0.20.3 it is a bit more concise:
gdf = pd.concat([gdf1, gdf2])
gdf is automatically created as a GeoDataFrame. Of course if there is a chance of conflicting indices you'll want to keep the 'ignore_index=True' parameter.
9
I added that recipe to the rasterio documentation. Since it was such a simple shape, in this case I just unzipped the coords in the single record contained by the shapefile. That is, x, y = zip(*features[0]['coordinates'][0]), and then just plot.
More generally, I use PolygonPatch from descartes, and matplotlib.collections.
import fiona
import rasterio
...
9
The GeoDataFrame
import geopandas as gpd
g1 = gpd.GeoDataFrame.from_file("poly_intersect.shp")
g1.shape
(4, 3)
1) You can use the itertools module
a) If you want to merge the intersections of the overlapping polygons
import itertools
geoms = g1['geometry'].tolist()
intersection_iter = gpd.GeoDataFrame(gpd.GeoSeries([poly[0].intersection(poly[1]) for poly ...
9
Simple dissolve by GeoDataFrame field (Aggregation with dissolve)
import geopandas as gpd
world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))
continents = world.dissolve('continent')
continents.head(3)
But if we use as_index=False from Pandas DataFrame.groupby
continents2 = world.dissolve('continent', as_index=False)
continents2.head(3)
...
8
Geopandas plot takes a colormap parameter
The parameter is called cmap. There's an alias for it called colormap but it should not be used anymore.
GeoDataFrame.plot(column=None, cmap=None, alpha=0.5, categorical=False, legend=False, axes=None)
The pandas documentation mentions colormaps such as
cmap='cubehelix'
cmap='Greens'
cmap='winter'
8
I have come across this behavior before.
You need to explicitly pass the well known text (crs_wkt) string to the to_file() method. The string will then get passed to fiona.open(), which writes out the .prj file.
Using your sample code, doing something like this:
ws = r"D:\temp_se"
prj_file = gpd.datasets.get_path('naturalearth_lowres').replace(".shp","....
8
As I don't know your data, I give you a solution with my data (with x,y, z and a colum to test < 30). If I use your solution
import geopandas as gpd
import numpy as np
import pandas as pd
numpy_point_array = np.array([[202104.271187,90516.656257,170.520004272, 45],[202139.659561,90516.656257,170.740005493, 15],[202175.047935,90516.656257,170.809997559, ...
8
The fiona.listlayers() function returns a list of names of layers in a dataset.
import fiona
fiona.listlayers('NYCFutureHighTideWithSLR.gdb')
Any of the elements of the list can be used as a value of the layer keyword argument for gpd.read_file().
The integer index of a list element may also be used. If the layers of your dataset are ['layer_a', 'layer_b']...
8
The type of env here is a Shapely Polygon. In this line
envgdf['geometry'] = env
You're trying to assign a Polygon to a Geometry column. You can instead create a Geoseries from the Polygon and create a Geodataframe based on that. Here's the updated code:
import sys
import geopandas as gpd
shp = (sys.argv[1])
gdf = gpd.read_file(shp)
union = gdf....
Only top voted, non community-wiki answers of a minimum length are eligible
