matplotlib图上的白色区域,pygrib数据在359.5到360度之间

2024-09-28 20:50:28 发布

您现在位置:Python中文网/ 问答频道 /正文

我尝试用matplotlib绘制gfs天气模型的输出,使用pygrib保存数据,这些数据保存在grib文件中。几乎一切正常,输出如下:

enter image description here

程序似乎没有用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)

Tags: 数据importselfdataplotinitdefinline
2条回答

我自己也运行过几次,basemap模块的addcyclic函数实际上运行得很好。The basemap docs列出语法并使用得很好。在

对于代码中的变量,可以在乘以之前或之后添加循环点自我因素在你的GribDataLayer类中:

layerdata, lons = addcyclic(layerdata, lons)

你也可以使用np.追加编写自己的函数来完成同样的任务。它看起来像这样:

^{pr2}$

如果您的输入数据是二维的,那么上面的语法相当于选择第一个经度带中的所有数据(即layerdata[:,0])

layerdata = np.append(layerdata,layerdata[:,0,None],axis=-1)

希望这有帮助!在

两个月后我自己解决了:

如果经度范围从0到359.75,Matplotlib不会填充该区域,因为从matplotlibs的角度来看,它结束于此。我先把数据分开,然后再叠加,解决了这个问题。在

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度),我不知道是否有修复方法。在

相关问题 更多 >