赞
踩
最近使用的EC和GFS预报数据给的都是grib2格式的,之前用惯nc格式的,python读取grib2数据的时候还走了些弯路,看到很多博客上给的教程其实不能满足我的需求,现在搞明白了分享一下
第一个问题就是我电脑上pygrib安装都折腾了一阵子
我电脑上直接pip是装不上去的
pip3 install pygrib
使用conda的时候最好使用的是conda-forge源的,用这个装成功了
conda install -c conda-forge pygrib
以一个GFS的预报数据为例,我们使用pygrib读取之后,得到的是一个可迭代对象,可以拿来循环
import pygrib
grbs = pygrib.open('gfs.t00z.pgrb2.0p25.f024')
for grb in grbs:
print(grb)
我们来看一下遍历出的grbmessage对象打印出来是什么
1:Pressure reduced to MSL:Pa (instant):regular_ll:meanSea:level 0:fcst time 24 hrs:from 202407210000
2:Cloud mixing ratio:kg kg**-1 (instant):regular_ll:hybrid:level 1:fcst time 24 hrs:from 202407210000
3:Ice water mixing ratio:kg kg**-1 (instant):regular_ll:hybrid:level 1:fcst time 24 hrs:from 202407210000
可以发现,每一个grbMessage对象可以是二维格点数据,同时单位,气压层等维度信息也是有的
以“510:Temperature:K (instant):regular_ll:isobaricInhPa:level 95000 Pa:fcst time 24 hrs:from 202407210000”为例,表示的意思如下
510: 这是记录的索引或序号。它标识这是文件中的第 510 个记录。
Temperature: 这是变量的名称。在这里,它表示“温度”。
K (instant): 这是变量的单位和时间类型。
K
表示开尔文,是温度的单位。instant
表示这是瞬时值。regular_ll: 这是网格类型,regular_ll
表示规则的纬度/经度网格。
isobaricInhPa: 这是变量所属的层。在这里,它表示等压面。
level 95000 Pa: 这是变量的层次信息。这里表示95000帕斯卡 (Pa) 等压面(950 hPa 等压面)。
fcst time 24 hrs: 这是预报时间。这里表示从初始时间起 24 小时后的预报。
from 202407210000: 这是初始时间,表示 2024 年 7 月 21 日 00:00。
我下载的明明有41个高度层的信息,但是这样便利出来每一条数据都只是一个二维变量,我最终想获得的温度数据应该是三维数据,这意味着我需要将遍历的数据中为Temperature的调出来将他们拼成三维数组,这样才获取了一个变量完整的数据
可以看到,遍历出的数据以Temperature开头的就有56个,这是将非等压面的一些奇怪的Temperature也算上的总数
因此我们需要从中获取grbs中真正是气压层的Temperature
import pygrib import numpy as np def select_variable(filepath, variable, typeOfLevel): # 打开GRIB2文件 grbs = pygrib.open(filepath) for grb in grbs: print(grb) temperature_records = grbs.select(name='Temperature', typeOfLevel='isobaricInhPa') # 获取记录数量 num_records = len(temperature_records) # 假设所有记录的维度相同,获取第一个记录的维度 # 这里使用第一个记录来确定数据的维度 first_record = temperature_records[0].values data_shape = first_record.shape # 创建一个三维数组来存储温度数据 select_data = np.zeros((num_records, *data_shape)) # 提取每个记录的数据并存储到三维数组中 for i, record in enumerate(temperature_records): select_data[i, :, :] = record.values return select_data filepath = '../download_res/ncep/20240721/gfs.t00z.pgrb2.0p25.f024' variable = 'Temperature' typeOfLevel='isobaricInhPa' T = select_variable(filepath, variable, typeOfLevel) print(np.shape(T))
这个函数和案例就展示了如何从grib数据中获取需要的变量
这里需要解释的是typeOfLevel=‘isobaricInhPa’,这是因为这个数据中只要是有不同高度层的数据,他的typeOfLevel就是isobaricInhPa,因为除了高度层的Temperature可以看到这种不是描述等压面温度的Temperature,如果不加上typeOfLevel='isobaricInhPa就会把他们也筛选出来
Panoply直接搜索Panoply去官网下载之后直接打开即可,现在需要Java11装好,也很简单,去orical官网下载Java11的非安装版本,之后配置一下环境变量,在cmd输入java -version看到是11就可以了,无脑使用
可以看到,与pygrib不一样,Panoply直接会读出一个叫Temperature isobaric的变量,表示的就是有不同气压层的温度变量,这样可以很直观地看到数据的维度,但需要注意的是panoply和pygrib读出的变量可能长得不一样,需要留意
在Panoply中也可以一目了然总共有什么变量,再回去python中遍历一下找到pygrib读出来的变量名,获取变量
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。