Как наложить график колчана на график геоданных в Python

Я хотел бы наложить колчанный график направления ветра на базовую карту в моем блокноте jupyter. У меня есть фреймворк pandas, который включает столбцы: | Широта | Долгота | Предполагаемый истинный ветер |

Я уже использовал геопанды для создания фрейма геоданных и построения данных трека gps на базовой карте osm, используя контекст (код ниже). Я также смог объединить широту и долготу, чтобы получить среднее значение True Wind Inferred (направление ветра) для «прямоугольника» на карте. Тем не менее, я не нашел никаких примеров того, как построить график колчана для разбитого по ячейкам True Wind Inferred в ящиках. Я пока построил только график рассеяния, но цветовая карта плохо визуализирует данные о направлении.

Импорт:

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()

диаграмма рассеяния

Я хотел бы изменить тип графика с разброса с цветами на график колчана со стрелками, представляющими направление ветра по компасу.


person BenCaldwell    schedule 25.03.2019    source источник
comment
Вы смотрели колчан matplotlib? matplotlib.org/api/_as_gen/matplotlib.pyplot .quiver.html   -  person Julia    schedule 25.03.2019
comment
Я пытался воспроизвести это, но установка geopandas / fiona / gdal и т. Д. В моем env не стоит целого дня ...   -  person Josh Friedlander    schedule 25.03.2019
comment
@Julia В итоге я использовал колчан matplotlib. Сложнее всего было нанести его на те же оси, что и карта. Я решил это сейчас и опубликовал в ответе.   -  person BenCaldwell    schedule 26.03.2019


Ответы (1)


Я разобрался. 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()

колчан на карте

person BenCaldwell    schedule 26.03.2019