pygrib 数据在 359.5 到 360 度之间的 matplotlib 图上的白色区域
White area on matplotlib plot with pygrib data between 359.5 and 360 degrees
我尝试的是用 matplotlib 绘制 gfs 天气模型的输出，使用 pygrib 保存数据，这些数据保存在 grib 文件中。几乎一切正常，输出如下所示：
程序似乎没有通过使用 0 度的数据来缩小 359.5 和 360 度之间的差距。如果数据在常规列表或其他内容中，我会使用 0° 的数据并通过附加列表将其保存为 360°。我见过有人对非 pygrib 数据有同样的问题。
如果您知道如何更改 pygrib 数据（不幸的是，常规操作不适用于 pygrib 数据）或如何使 matplotlib 缩小差距，您真的会帮助我解决这个问题。也许函数 "addcyclic" 可以提供帮助，但我不知道如何。
编辑：我解决了问题，请看我的回答。
这里是产生问题的代码：
#!/usr/bin/python3
import os, sys, datetime, string
from abc import ABCMeta, abstractmethod
import numpy as np
import numpy.ma as ma
from scipy.ndimage.filters import minimum_filter, maximum_filter
import pygrib
from netCDF4 import Dataset
from pylab import *
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap, addcyclic, shiftgrid
import laplaceFilter
import mpl_util
class Plot(Basemap):
def __init__(self, basemapParams):
super().__init__(**basemapParams)
self.layers = []
def addLayer(self, layer):
self.layers.append(layer)
def plot(self, data):
for layer in self.layers:
layer.plot(self, data)
plt.title('Plot')
plt.show()
class Layer(metaclass=ABCMeta):
def __init__(self):
pass
@abstractmethod
def plot(self, plot, data):
return NotImplemented
class BackgroundLayer(Layer):
def __init__(self, bgtype, coords):
#possible bgtype values: borders, topo, both
self.bgtype = bgtype
self.lonStart = coords[0]
self.lonEnd = coords[1]
self.latStart = coords[2]
self.latEnd = coords[3]
def plot(self, plot, data):
[...]
def findSubsetIndices(self,min_lat,max_lat,min_lon,max_lon,lats,lons):
[...]
class LegendLayer(Layer):
def __init__(self):
pass
class GribDataLayer(Layer, metaclass=ABCMeta):
def __init__(self, varname, level, clevs, cmap, factor):
self.varname = varname
self.level = level
self.clevs = clevs
self.cmap = cmap
self.factor = factor
def plot(self, plot, data):
#depending on the height we want to use, we have to change the index
indexes = {1000:0, 2000:1, 3000:2, 5000:3, 7000:4, 10000:5, 15000:6, 20000:7, 25000:8, 30000:9,
35000:10, 40000:11, 45000:12, 50000:13, 55000:14, 60000:15, 65000:16, 70000:17,
75000:18, 80000:19, 85000:20, 90000:21, 92500:22, 95000:23, 97500:24, 100000:25, 0:0}
selecteddata = data.select(name = self.varname)[indexes[self.level]]
lats, lons = selecteddata.latlons()
layerdata = selecteddata.values*self.factor
x, y = plot(lons, lats) # compute map proj coordinates.
self.fillLayer(plot, x, y, layerdata, self.clevs, self.cmap)
@abstractmethod
def fillLayer(self, plot, x, y, layerdata, clevs, cmap):
return NotImplemented
class ContourLayer(GribDataLayer):
def __init__(self, varname, level, clevs, cmap, factor, linewidth=1.5, fontsize=15,
fmt="%3.1f", inline=0,labelcolor = 'k'):
self.linewidth = linewidth
self.fontsize = fontsize
self.fmt = fmt
self.inline = inline
self.labelcolor = labelcolor
super().__init__(varname, level, clevs, cmap, factor)
def fillLayer(self, plot, x, y, layerdata, clevs, cmap):
# contour data over the map.
cs = plot.contour(x,y,layerdata,clevs,colors = cmap,linewidths = self.linewidth)
plt.clabel(cs, clevs, fontsize = self.fontsize, fmt = self.fmt,
inline = self.inline, colors = self.labelcolor)
if self.varname == "Pressure reduced to MSL":
self.plotHighsLows(plot,layerdata,x,y)
def plotHighsLows(self,plot,layerdata,x,y):
[...]
class ContourFilledLayer(GribDataLayer):
def __init__(self, varname, level, clevs, cmap, factor, extend="both"):
self.extend = extend
super().__init__(varname, level, clevs, cmap, factor)
def fillLayer(self, plot, x, y, layerdata, clevs, cmap):
# contourfilled data over the map.
cs = plot.contourf(x,y,layerdata,levels=clevs,cmap=cmap,extend=self.extend)
#cbar = plot.colorbar.ColorbarBase(cs)
[...]
ger_coords = [4.,17.,46.,56.]
eu_coords = [-25.,57.,22.,70.]
### Choose Data
data = pygrib.open('gfs.t12z.mastergrb2f03')
### 500hPa Europe
coords = eu_coords
plot1 = Plot({"projection":"lcc","resolution":"h","rsphere":(6378137.00,6356752.3142), "area_thresh": 1000.,
"llcrnrlon":coords[0],"llcrnrlat":coords[2],"urcrnrlon":coords[1],"urcrnrlat":coords[3],
"lon_0":(coords[0]+coords[1])/2.,"lat_0":(coords[2]+coords[3])/2.})
clevs = range(480,600,4)
cmap = plt.cm.nipy_spectral
factor = .1
extend = "both"
level = 50000
layer1 = ContourFilledLayer('Geopotential Height', level, clevs, cmap, factor, extend)
clevs = [480.,552.,600.]
linewidth = 2.
fontsize = 14
fmt = "%d"
inline = 0
labelcolor = 'k'
layer2 = ContourLayer('Geopotential Height', level, clevs, 'k', factor, linewidth, fontsize, fmt, inline, labelcolor)
level = 0
clevs = range(800,1100,5)
factor = .01
linewidth = 1.5
inline = 0
labelcolor = 'k'
layer3 = ContourLayer('Pressure reduced to MSL', level, clevs, 'w', factor, linewidth, fontsize, fmt, inline, labelcolor)
plot1.addLayer(BackgroundLayer('borders', coords))
plot1.addLayer(layer1)
plot1.addLayer(layer2)
plot1.addLayer(layer3)
plot1.plot(data)
2个月后我自己解决了：
如果您的经度范围是 0 到 359.75，Matplotlib 不会填充该区域，因为从 matplotlib 的角度来看，它在那里结束。我通过划分数据然后堆叠来解决它。
selecteddata_all = data.select(name = "Temperature")[0]
selecteddata1, lats1, lons1 = selecteddata_all.data(lat1=20,lat2=60,lon1=335,lon2=360)
selecteddata2, lats2, lons2 = selecteddata_all.data(lat1=20,lat2=60,lon1=0,lon2=30)
lons = np.hstack((lons1,lons2))
lats = np.hstack((lats1,lats2))
selecteddata = np.hstack((selecteddata1,selecteddata2))
0° 左侧不再有白色区域。
如果你想绘制整个半球（0 到 359.75 度），我不知道是否有解决方法。
我自己 运行 也参与过几次，底图模块的 addcyclic 函数实际上工作得很好。 The basemap docs 把语法排好，用得很好。
就代码中的变量而言，您可以在 GribDataLayer class:
中乘以 self.factor 之前或之后添加循环点
layerdata, lons = addcyclic(layerdata, lons)
您也可以使用 np.append 并编写您自己的函数来完成同样的任务。它看起来像这样：
layerdata = np.append(layerdata,layerdata[...,0,None],axis=-1)
如果您的输入数据是二维的，那么上面的语法相当于选择第一个经度带中的所有数据（即 layerdata[:,0]）
layerdata = np.append(layerdata,layerdata[:,0,None],axis=-1)
希望对您有所帮助！
我尝试的是用 matplotlib 绘制 gfs 天气模型的输出，使用 pygrib 保存数据，这些数据保存在 grib 文件中。几乎一切正常，输出如下所示：
程序似乎没有通过使用 0 度的数据来缩小 359.5 和 360 度之间的差距。如果数据在常规列表或其他内容中，我会使用 0° 的数据并通过附加列表将其保存为 360°。我见过有人对非 pygrib 数据有同样的问题。 如果您知道如何更改 pygrib 数据（不幸的是，常规操作不适用于 pygrib 数据）或如何使 matplotlib 缩小差距，您真的会帮助我解决这个问题。也许函数 "addcyclic" 可以提供帮助，但我不知道如何。
编辑：我解决了问题，请看我的回答。
这里是产生问题的代码：
#!/usr/bin/python3
import os, sys, datetime, string
from abc import ABCMeta, abstractmethod
import numpy as np
import numpy.ma as ma
from scipy.ndimage.filters import minimum_filter, maximum_filter
import pygrib
from netCDF4 import Dataset
from pylab import *
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap, addcyclic, shiftgrid
import laplaceFilter
import mpl_util
class Plot(Basemap):
def __init__(self, basemapParams):
super().__init__(**basemapParams)
self.layers = []
def addLayer(self, layer):
self.layers.append(layer)
def plot(self, data):
for layer in self.layers:
layer.plot(self, data)
plt.title('Plot')
plt.show()
class Layer(metaclass=ABCMeta):
def __init__(self):
pass
@abstractmethod
def plot(self, plot, data):
return NotImplemented
class BackgroundLayer(Layer):
def __init__(self, bgtype, coords):
#possible bgtype values: borders, topo, both
self.bgtype = bgtype
self.lonStart = coords[0]
self.lonEnd = coords[1]
self.latStart = coords[2]
self.latEnd = coords[3]
def plot(self, plot, data):
[...]
def findSubsetIndices(self,min_lat,max_lat,min_lon,max_lon,lats,lons):
[...]
class LegendLayer(Layer):
def __init__(self):
pass
class GribDataLayer(Layer, metaclass=ABCMeta):
def __init__(self, varname, level, clevs, cmap, factor):
self.varname = varname
self.level = level
self.clevs = clevs
self.cmap = cmap
self.factor = factor
def plot(self, plot, data):
#depending on the height we want to use, we have to change the index
indexes = {1000:0, 2000:1, 3000:2, 5000:3, 7000:4, 10000:5, 15000:6, 20000:7, 25000:8, 30000:9,
35000:10, 40000:11, 45000:12, 50000:13, 55000:14, 60000:15, 65000:16, 70000:17,
75000:18, 80000:19, 85000:20, 90000:21, 92500:22, 95000:23, 97500:24, 100000:25, 0:0}
selecteddata = data.select(name = self.varname)[indexes[self.level]]
lats, lons = selecteddata.latlons()
layerdata = selecteddata.values*self.factor
x, y = plot(lons, lats) # compute map proj coordinates.
self.fillLayer(plot, x, y, layerdata, self.clevs, self.cmap)
@abstractmethod
def fillLayer(self, plot, x, y, layerdata, clevs, cmap):
return NotImplemented
class ContourLayer(GribDataLayer):
def __init__(self, varname, level, clevs, cmap, factor, linewidth=1.5, fontsize=15,
fmt="%3.1f", inline=0,labelcolor = 'k'):
self.linewidth = linewidth
self.fontsize = fontsize
self.fmt = fmt
self.inline = inline
self.labelcolor = labelcolor
super().__init__(varname, level, clevs, cmap, factor)
def fillLayer(self, plot, x, y, layerdata, clevs, cmap):
# contour data over the map.
cs = plot.contour(x,y,layerdata,clevs,colors = cmap,linewidths = self.linewidth)
plt.clabel(cs, clevs, fontsize = self.fontsize, fmt = self.fmt,
inline = self.inline, colors = self.labelcolor)
if self.varname == "Pressure reduced to MSL":
self.plotHighsLows(plot,layerdata,x,y)
def plotHighsLows(self,plot,layerdata,x,y):
[...]
class ContourFilledLayer(GribDataLayer):
def __init__(self, varname, level, clevs, cmap, factor, extend="both"):
self.extend = extend
super().__init__(varname, level, clevs, cmap, factor)
def fillLayer(self, plot, x, y, layerdata, clevs, cmap):
# contourfilled data over the map.
cs = plot.contourf(x,y,layerdata,levels=clevs,cmap=cmap,extend=self.extend)
#cbar = plot.colorbar.ColorbarBase(cs)
[...]
ger_coords = [4.,17.,46.,56.]
eu_coords = [-25.,57.,22.,70.]
### Choose Data
data = pygrib.open('gfs.t12z.mastergrb2f03')
### 500hPa Europe
coords = eu_coords
plot1 = Plot({"projection":"lcc","resolution":"h","rsphere":(6378137.00,6356752.3142), "area_thresh": 1000.,
"llcrnrlon":coords[0],"llcrnrlat":coords[2],"urcrnrlon":coords[1],"urcrnrlat":coords[3],
"lon_0":(coords[0]+coords[1])/2.,"lat_0":(coords[2]+coords[3])/2.})
clevs = range(480,600,4)
cmap = plt.cm.nipy_spectral
factor = .1
extend = "both"
level = 50000
layer1 = ContourFilledLayer('Geopotential Height', level, clevs, cmap, factor, extend)
clevs = [480.,552.,600.]
linewidth = 2.
fontsize = 14
fmt = "%d"
inline = 0
labelcolor = 'k'
layer2 = ContourLayer('Geopotential Height', level, clevs, 'k', factor, linewidth, fontsize, fmt, inline, labelcolor)
level = 0
clevs = range(800,1100,5)
factor = .01
linewidth = 1.5
inline = 0
labelcolor = 'k'
layer3 = ContourLayer('Pressure reduced to MSL', level, clevs, 'w', factor, linewidth, fontsize, fmt, inline, labelcolor)
plot1.addLayer(BackgroundLayer('borders', coords))
plot1.addLayer(layer1)
plot1.addLayer(layer2)
plot1.addLayer(layer3)
plot1.plot(data)
2个月后我自己解决了：
如果您的经度范围是 0 到 359.75，Matplotlib 不会填充该区域，因为从 matplotlib 的角度来看，它在那里结束。我通过划分数据然后堆叠来解决它。
selecteddata_all = data.select(name = "Temperature")[0]
selecteddata1, lats1, lons1 = selecteddata_all.data(lat1=20,lat2=60,lon1=335,lon2=360)
selecteddata2, lats2, lons2 = selecteddata_all.data(lat1=20,lat2=60,lon1=0,lon2=30)
lons = np.hstack((lons1,lons2))
lats = np.hstack((lats1,lats2))
selecteddata = np.hstack((selecteddata1,selecteddata2))
0° 左侧不再有白色区域。
如果你想绘制整个半球（0 到 359.75 度），我不知道是否有解决方法。
我自己 运行 也参与过几次，底图模块的 addcyclic 函数实际上工作得很好。 The basemap docs 把语法排好，用得很好。
就代码中的变量而言，您可以在 GribDataLayer class:
中乘以 self.factor 之前或之后添加循环点layerdata, lons = addcyclic(layerdata, lons)
您也可以使用 np.append 并编写您自己的函数来完成同样的任务。它看起来像这样：
layerdata = np.append(layerdata,layerdata[...,0,None],axis=-1)
如果您的输入数据是二维的，那么上面的语法相当于选择第一个经度带中的所有数据（即 layerdata[:,0]）
layerdata = np.append(layerdata,layerdata[:,0,None],axis=-1)
希望对您有所帮助！