如何在 python 中的地理数据框图上叠加箭袋图
How to overlay a quiver plot on a geodataframe plot in python
我想在我的 jupyter 笔记本中的底图上叠加一个风向颤动图。我有一个包含列的 pandas 数据框: |纬度 |经度 |真风推断 |
我已经使用 geopandas 创建地理数据框并使用上下文(下面的代码)在 osm 底图上绘制 gps 跟踪数据。我还能够对纬度和经度进行分类,以获得地图上 "box" 的平均真风推断(风向)。
但是,我还没有找到任何示例来说明如何在方框中绘制推断的装箱真风的颤抖图。到目前为止,我只绘制了散点图,但彩色地图不能很好地显示方向数据。
进口:
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
# plot inline graphics
%pylab inline
import os.path
from shapely.geometry import Point
import geopandas as gpd
import contextily as ctx
示例数据框:
df[['Latitude', 'Longitude', 'True Wind Inferred', 'coords']].head()
Latitude Longitude True Wind Inferred coords
0 -31.991899 115.848825 173.835559 POINT (115.848825 -31.991899)
1 -31.992036 115.848873 182.620880 POINT (115.848873 -31.992036)
2 -31.992181 115.848895 192.140276 POINT (115.848895 -31.992181)
3 -31.992308 115.848832 206.655730 POINT (115.848832 -31.992308)
4 -31.992430 115.848784 218.656646 POINT (115.848784 -31.99243)
合并数据帧:
step = 0.005
to_bin = lambda x: np.floor(x / step) * step
dfLocBin['latbin'] = df['Latitude'].map(to_bin)
dfLocBin['lonbin'] = df['Latitude'].map(to_bin)
dfLocBin = df.groupby(['latbin', 'lonbin'])[['True Wind Inferred']].mean()
dfLocBin.reset_index(inplace=True)
dfLocBin['coords'] = list(zip(dfLocBin['lonbin'], dfLocBin['latbin']))
dfLocBin['coords'] = dfLocBin['coords'].apply(Point)
dfLocBin.head()
latbin lonbin True Wind Inferred coords
0 -32.015 115.790 223.149075 POINT (115.79 -32.015)
1 -32.015 115.795 222.242870 POINT (115.795 -32.015)
2 -32.015 115.800 223.710092 POINT (115.8 -32.015)
3 -32.015 115.805 225.887096 POINT (115.805 -32.015)
4 -32.015 115.810 225.298059 POINT (115.81 -32.015)
并策划:
def add_basemap(ax, zoom, url='http://tile.stamen.com/terrain/tileZ/tileX/tileY.png'):
xmin, xmax, ymin, ymax = ax.axis()
basemap, extent = ctx.bounds2img(xmin, ymin, xmax, ymax, zoom=zoom, url=url)
ax.imshow(basemap, extent=extent, interpolation='bilinear')
# restore original x/y limits
ax.axis((xmin, xmax, ymin, ymax))
geo_df = gpd.GeoDataFrame(
dfLocBin, crs ={'init': 'epsg:4326'},
geometry = dfLocBin['coords']
).to_crs(epsg=3857)
ax = geo_df.plot(
figsize= (20, 20),
alpha = 1,
c=dfLocBin['True Wind Inferred']
)
add_basemap(ax, zoom=15, url=ctx.tile_providers.ST_TONER)
ax.set_axis_off()
plt.title('Binned True Wind Direction')
plt.show()
scatter plot
我想将绘图类型从带有颜色的散点图更改为带有代表风罗盘方向的箭头的箭袋图。
我解决了。箭袋的 X、Y 需要来自地理数据框的 geometry
才能在同一轴上正确绘制。地理数据框列如下所示:
geo_df.head()
latbin lonbin True Wind Inferred coords geometry
0 -32.014 115.798 220.492453 POINT (115.798 -32.014) POINT (12890574.39487949 -3765148.48502445)
1 -32.014 115.800 225.718756 POINT (115.8 -32.014) POINT (12890797.03386108 -3765148.48502445)
工作代码:
# bin the coordinates and plot a vector field
step = 0.002
to_bin = lambda x: np.floor(x / step) * step
df['latbin'] = df['Latitude'].map(to_bin)
df['lonbin'] = df['Longitude'].map(to_bin)
dfLocBin = df.groupby(['latbin', 'lonbin'])[['True Wind Inferred']].mean()
dfLocBin.reset_index(inplace=True)
dfLocBin['coords'] = list(zip(dfLocBin['lonbin'], dfLocBin['latbin']))
dfLocBin['coords'] = dfLocBin['coords'].apply(Point)
# ... turn them into geodataframe, and convert our
# epsg into 3857, since web map tiles are typically
# provided as such.
geo_df = gpd.GeoDataFrame(
dfLocBin, crs ={'init': 'epsg:4326'},
geometry = dfLocBin['coords']
).to_crs(epsg=3857)
# ... and make the plot
ax = geo_df.plot(
figsize= (20, 20),
alpha = 1
)
geo_df['X'] = geo_df['geometry'].x
geo_df['Y'] = geo_df['geometry'].y
geo_df['U'] = np.cos(np.radians(geo_df['True Wind Inferred']))
geo_df['V'] = np.sin(np.radians(geo_df['True Wind Inferred']))
ax.quiver(geo_df['X'],
geo_df['Y'],
geo_df['U'],
geo_df['V'],
color='deepskyblue')
add_basemap(ax, zoom=15, url=ctx.tile_providers.ST_TONER)
ax.set_axis_off()
plt.title('Binned True Wind Direction')
plt.show()
quiver on map
我想在我的 jupyter 笔记本中的底图上叠加一个风向颤动图。我有一个包含列的 pandas 数据框: |纬度 |经度 |真风推断 |
我已经使用 geopandas 创建地理数据框并使用上下文(下面的代码)在 osm 底图上绘制 gps 跟踪数据。我还能够对纬度和经度进行分类,以获得地图上 "box" 的平均真风推断(风向)。 但是,我还没有找到任何示例来说明如何在方框中绘制推断的装箱真风的颤抖图。到目前为止,我只绘制了散点图,但彩色地图不能很好地显示方向数据。
进口:
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
# plot inline graphics
%pylab inline
import os.path
from shapely.geometry import Point
import geopandas as gpd
import contextily as ctx
示例数据框:
df[['Latitude', 'Longitude', 'True Wind Inferred', 'coords']].head()
Latitude Longitude True Wind Inferred coords
0 -31.991899 115.848825 173.835559 POINT (115.848825 -31.991899)
1 -31.992036 115.848873 182.620880 POINT (115.848873 -31.992036)
2 -31.992181 115.848895 192.140276 POINT (115.848895 -31.992181)
3 -31.992308 115.848832 206.655730 POINT (115.848832 -31.992308)
4 -31.992430 115.848784 218.656646 POINT (115.848784 -31.99243)
合并数据帧:
step = 0.005
to_bin = lambda x: np.floor(x / step) * step
dfLocBin['latbin'] = df['Latitude'].map(to_bin)
dfLocBin['lonbin'] = df['Latitude'].map(to_bin)
dfLocBin = df.groupby(['latbin', 'lonbin'])[['True Wind Inferred']].mean()
dfLocBin.reset_index(inplace=True)
dfLocBin['coords'] = list(zip(dfLocBin['lonbin'], dfLocBin['latbin']))
dfLocBin['coords'] = dfLocBin['coords'].apply(Point)
dfLocBin.head()
latbin lonbin True Wind Inferred coords
0 -32.015 115.790 223.149075 POINT (115.79 -32.015)
1 -32.015 115.795 222.242870 POINT (115.795 -32.015)
2 -32.015 115.800 223.710092 POINT (115.8 -32.015)
3 -32.015 115.805 225.887096 POINT (115.805 -32.015)
4 -32.015 115.810 225.298059 POINT (115.81 -32.015)
并策划:
def add_basemap(ax, zoom, url='http://tile.stamen.com/terrain/tileZ/tileX/tileY.png'):
xmin, xmax, ymin, ymax = ax.axis()
basemap, extent = ctx.bounds2img(xmin, ymin, xmax, ymax, zoom=zoom, url=url)
ax.imshow(basemap, extent=extent, interpolation='bilinear')
# restore original x/y limits
ax.axis((xmin, xmax, ymin, ymax))
geo_df = gpd.GeoDataFrame(
dfLocBin, crs ={'init': 'epsg:4326'},
geometry = dfLocBin['coords']
).to_crs(epsg=3857)
ax = geo_df.plot(
figsize= (20, 20),
alpha = 1,
c=dfLocBin['True Wind Inferred']
)
add_basemap(ax, zoom=15, url=ctx.tile_providers.ST_TONER)
ax.set_axis_off()
plt.title('Binned True Wind Direction')
plt.show()
scatter plot
我想将绘图类型从带有颜色的散点图更改为带有代表风罗盘方向的箭头的箭袋图。
我解决了。箭袋的 X、Y 需要来自地理数据框的 geometry
才能在同一轴上正确绘制。地理数据框列如下所示:
geo_df.head()
latbin lonbin True Wind Inferred coords geometry
0 -32.014 115.798 220.492453 POINT (115.798 -32.014) POINT (12890574.39487949 -3765148.48502445)
1 -32.014 115.800 225.718756 POINT (115.8 -32.014) POINT (12890797.03386108 -3765148.48502445)
工作代码:
# bin the coordinates and plot a vector field
step = 0.002
to_bin = lambda x: np.floor(x / step) * step
df['latbin'] = df['Latitude'].map(to_bin)
df['lonbin'] = df['Longitude'].map(to_bin)
dfLocBin = df.groupby(['latbin', 'lonbin'])[['True Wind Inferred']].mean()
dfLocBin.reset_index(inplace=True)
dfLocBin['coords'] = list(zip(dfLocBin['lonbin'], dfLocBin['latbin']))
dfLocBin['coords'] = dfLocBin['coords'].apply(Point)
# ... turn them into geodataframe, and convert our
# epsg into 3857, since web map tiles are typically
# provided as such.
geo_df = gpd.GeoDataFrame(
dfLocBin, crs ={'init': 'epsg:4326'},
geometry = dfLocBin['coords']
).to_crs(epsg=3857)
# ... and make the plot
ax = geo_df.plot(
figsize= (20, 20),
alpha = 1
)
geo_df['X'] = geo_df['geometry'].x
geo_df['Y'] = geo_df['geometry'].y
geo_df['U'] = np.cos(np.radians(geo_df['True Wind Inferred']))
geo_df['V'] = np.sin(np.radians(geo_df['True Wind Inferred']))
ax.quiver(geo_df['X'],
geo_df['Y'],
geo_df['U'],
geo_df['V'],
color='deepskyblue')
add_basemap(ax, zoom=15, url=ctx.tile_providers.ST_TONER)
ax.set_axis_off()
plt.title('Binned True Wind Direction')
plt.show()
quiver on map