使用 Basemap 的 robin 投影和 contourf 数据问题
Issue with the robin projection and contourf data with Basemap
我正在使用底图库显示来自哥白尼程序的空间信息。
问题是我不知道如何在知更鸟投影上投影数据,但我用正交投影正确地做到了。
所以目前,我试过这个:
plt.ioff()
# adapt for location of datasources
filePath = '../data/grib/download.grib'
# load data
grbs = grb.open(filePath)
grbs.seek(0)
data, lats, lons = (None, None, None)
dataUnit = None
title = None
for g in grbs:
data, lats, lons = g.data()
name = g.name
level = g.level
pressureUnit = g.pressureUnits
date = g.validDate
dataUnit = g.units
title = name + ' at ' + str(level) + ' ' + str(pressureUnit) + ' [' + str(date) + ']'
print(title)
break
# mapPlot = Basemap(projection='ortho', lat_0=0, lon_0=0)
mapPlot = Basemap(projection='robin', lat_0=0, lon_0=0, resolution='l')
mapPlot.drawcoastlines(linewidth=0.25)
x, y = mapPlot(lons, lats)
mapPlot.contourf(x, y, data)
mapPlot.colorbar(location='bottom', format='%.1f', label=dataUnit)
plt.title(title)
plt.show()
正交投影工作正常。但是对于知更鸟投影,我有一个......有趣的模式。
我做错了什么?
所以我想办法了。我被误导了,但我看到的第一个例子。
这是我的代码:
import matplotlib
from mpl_toolkits.basemap import Basemap, shiftgrid
import matplotlib.pyplot as plt
import numpy as np
import pygrib as grb
# Get data
data = g['values']
lats = g['distinctLatitudes'] # 1D vector
lons = g['distinctLongitudes'] # 1D vector
# Useful information for late
name = g.name
level = str(g.level) + g.pressureUnits
date = g.validDate
dataUnit = g.units
# Parse the data
# Shit the data to start à -180. This is important to mark the data to start at -180°
data, lons = shiftgrid(180., data, lons, start=False) # shiftgrid
# Choose a representation (works with both)
# mapPlot = Basemap(projection='ortho', lat_0=0, lon_0=0)
mapPlot = Basemap(projection='robin', lat_0=0, lon_0=0)
mapPlot.drawcoastlines(linewidth=0.25)
# Convert the coordinates into the map projection
x, y = mapPlot(*np.meshgrid(lons, lats))
# Display data
map = mapPlot.contourf(x, y, data, levels=boundaries, cmap=plt.get_cmap('coolwarm'))
# Add what ever you want to your map.
mapPlot.nightshade(date, alpha=0.1)
# Legend
mapPlot.colorbar(map, label=dataUnit)
# Title
plt.title(name + ' at ' + str(level) + ' [' + str(date) + ']')
plt.show()
这就是 returns 我所期待的。
我正在使用底图库显示来自哥白尼程序的空间信息。 问题是我不知道如何在知更鸟投影上投影数据,但我用正交投影正确地做到了。
所以目前,我试过这个:
plt.ioff()
# adapt for location of datasources
filePath = '../data/grib/download.grib'
# load data
grbs = grb.open(filePath)
grbs.seek(0)
data, lats, lons = (None, None, None)
dataUnit = None
title = None
for g in grbs:
data, lats, lons = g.data()
name = g.name
level = g.level
pressureUnit = g.pressureUnits
date = g.validDate
dataUnit = g.units
title = name + ' at ' + str(level) + ' ' + str(pressureUnit) + ' [' + str(date) + ']'
print(title)
break
# mapPlot = Basemap(projection='ortho', lat_0=0, lon_0=0)
mapPlot = Basemap(projection='robin', lat_0=0, lon_0=0, resolution='l')
mapPlot.drawcoastlines(linewidth=0.25)
x, y = mapPlot(lons, lats)
mapPlot.contourf(x, y, data)
mapPlot.colorbar(location='bottom', format='%.1f', label=dataUnit)
plt.title(title)
plt.show()
正交投影工作正常。但是对于知更鸟投影,我有一个......有趣的模式。
我做错了什么?
所以我想办法了。我被误导了,但我看到的第一个例子。
这是我的代码:
import matplotlib
from mpl_toolkits.basemap import Basemap, shiftgrid
import matplotlib.pyplot as plt
import numpy as np
import pygrib as grb
# Get data
data = g['values']
lats = g['distinctLatitudes'] # 1D vector
lons = g['distinctLongitudes'] # 1D vector
# Useful information for late
name = g.name
level = str(g.level) + g.pressureUnits
date = g.validDate
dataUnit = g.units
# Parse the data
# Shit the data to start à -180. This is important to mark the data to start at -180°
data, lons = shiftgrid(180., data, lons, start=False) # shiftgrid
# Choose a representation (works with both)
# mapPlot = Basemap(projection='ortho', lat_0=0, lon_0=0)
mapPlot = Basemap(projection='robin', lat_0=0, lon_0=0)
mapPlot.drawcoastlines(linewidth=0.25)
# Convert the coordinates into the map projection
x, y = mapPlot(*np.meshgrid(lons, lats))
# Display data
map = mapPlot.contourf(x, y, data, levels=boundaries, cmap=plt.get_cmap('coolwarm'))
# Add what ever you want to your map.
mapPlot.nightshade(date, alpha=0.1)
# Legend
mapPlot.colorbar(map, label=dataUnit)
# Title
plt.title(name + ' at ' + str(level) + ' [' + str(date) + ']')
plt.show()
这就是 returns 我所期待的。