除夕快乐 | xarray2024.11.0读取GRIB数据进阶

用户11172986

发布于 2025-02-05 14:08:06

11200

代码可运行

文章被收录于专栏：气python风雨气python风雨

运行总次数：0

代码可运行

黑暗中的我们都没有说话 | xarray2024.11.0读取GRIB数据进阶

前言

上期内容有读者朋友留言如下

grib、bufr这种表驱动格式企图“规范化”所有要素，但是气象海洋等各种学科各种要素茫茫多，数据发布机构茫茫多，搞出各种不公开的local表，事实上成了各是各的规范。而且eccodes解析包本身设计得贼难用。我辈不要再给这俩格式提供热度，让它们像ncl慢慢死掉才是人道的做法。

诚然如此格式形成了各式各样的数据壁垒，但是欧洲中心等庞然大物仍然以grib格式等提供产品

且还需用grib格式进行数值模拟输入，我辈还是学一学

项目概述

GRIB作为气象领域的"老大难"格式，其复杂的层次结构与多样的变量组合常让分析者望而生畏。本文将带你学习几个数据处理小技巧

三维数据读取
xarray数据写入GRIB格式
快速打印变量
三维变量可视化

安装依赖

# 推荐使用conda进行依赖管理  
conda install -c conda-forge eccodes cfgrib xarray=2024.11.0

技巧一：三维数据读取测试

同时读取多变量多高度层数据

import xarray as xr
# 一次性读取所有等压面数据
ds = xr.open_dataset('/home/mw/input/GFS1824/gfs_4_20230902_0000_021.grb2', engine='cfgrib',
                    filter_by_keys={'typeOfLevel': 'isobaricInhPa'})

ds

读取单层数据

# 读取单个层次测试
ds_single = xr.open_dataset(
    '/home/mw/input/GFS1824/gfs_4_20230902_0000_021.grb2',
    engine='cfgrib',
    backend_kwargs={
        'filter_by_keys': {'typeOfLevel': 'isobaricInhPa', 'level': 500}

    }
)

ds_single

多层单变量读取

ds1 = xr.open_dataset('/home/mw/input/GFS1824/gfs_4_20230902_0000_021.grb2', engine='cfgrib',
                     filter_by_keys={'typeOfLevel': 'isobaricInhPa', 'shortName': 't'})
ds1

你可以试试写一个单变量单层读取当然直接使用xarray的sel对多层变量进行裁剪更简单

技巧二：xarray 写入 grib

help(xarray_to_grib)

Help on module cfgrib.xarray_to_grib in cfgrib:

NAME
    cfgrib.xarray_to_grib

DESCRIPTION
    # Copyright 2017-2021 European Centre for Medium-Range Weather Forecasts (ECMWF).
    #
    # Licensed under the Apache License, Version 2.0 (the "License");
    # you may not use this file except in compliance with the License.
    # You may obtain a copy of the License at
    #
    #   http://www.apache.org/licenses/LICENSE-2.0
    #
    # Unless required by applicable law or agreed to in writing, software
    # distributed under the License is distributed on an "AS IS" BASIS,
    # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    # See the License for the specific language governing permissions and
    # limitations under the License.
    #
    # Authors:
    #   Alessandro Amici - B-Open - https://bopen.eu
    #   Aureliana Barghini - B-Open - https://bopen.eu
    #   Leonardo Barcaroli - B-Open - https://bopen.eu
    #

FUNCTIONS
    canonical_dataarray_to_grib(data_var, file, grib_keys={}, default_grib_keys={'centre': 255, 'typeOfLevel': 'surface'}, **kwargs)
        Write a ``xr.DataArray`` in *canonical* form to a GRIB file.
    
    canonical_dataset_to_grib(dataset, path, mode='wb', no_warn=False, grib_keys={}, **kwargs)
        Write a ``xr.Dataset`` in *canonical* form to a GRIB file.
    
    detect_grib_keys(data_var, default_grib_keys, grib_keys={})
    
    detect_regular_ll_grib_keys(lon, lat)
    
    detect_sample_name(grib_keys, sample_name_template='{geography}_{vertical}_grib{edition}')
    
    expand_dims(data_var: xarray.core.dataarray.DataArray) -> Tuple[List[str], xarray.core.dataarray.DataArray]
    
    make_template_message(merged_grib_keys, template_path=None, sample_name=None)
    
    merge_grib_keys(grib_keys, detected_grib_keys, default_grib_keys)
    
    regular_ll_params(values, min_value=-180.0, max_value=360.0)
    
    to_grib = canonical_dataset_to_grib(dataset, path, mode='wb', no_warn=False, grib_keys={}, **kwargs)
        Write a ``xr.Dataset`` in *canonical* form to a GRIB file.

DATA
    ALL_TYPE_OF_LEVELS = ['surface', 'meanSea', 'cloudBase', 'cloudTop', '...
    DEFAULT_GRIB_KEYS = {'centre': 255, 'typeOfLevel': 'surface'}
    GRID_TYPES = ['polar_stereographic', 'reduced_gg', 'reduced_ll', 'regu...
    LOGGER = <Logger cfgrib.xarray_to_grib (WARNING)>
    MESSAGE_DEFINITION_KEYS = ['productDefinitionTemplateNumber', 'units']
    TYPE_OF_LEVELS_ML = ['hybrid']
    TYPE_OF_LEVELS_PL = ['isobaricInhPa', 'isobaricInPa']
    TYPE_OF_LEVELS_SFC = ['surface', 'meanSea', 'cloudBase', 'cloudTop']

FILE
    /opt/conda/lib/python3.9/site-packages/cfgrib/xarray_to_grib.py

将 DataArray 写入 GRIB 文件

import xarray as xr
from cfgrib.xarray_to_grib import canonical_dataarray_to_grib

# 创建 DataArray（例如地表温度）
temperature = xr.DataArray(
    data=[[25.5, 26.0], [24.8, 25.3]],
    dims=["latitude", "longitude"],
    coords={
        "latitude": [30.0, 30.5],
        "longitude": [120.0, 120.5],
    },
    name="temperature",
    attrs={
        "units": "degC",
        "long_name": "Surface Temperature",
        "GRIB_shortName": "t",
    },
)

# 定义 GRIB 编码参数
grib_keys = {
    "centre": 98,       # 编码中心代码（例如中国气象局）
    "typeOfLevel": "surface",
    "level": 0,
}

# 打开文件句柄并以二进制写入模式写入
with open("output_single.grib", "wb") as f:
    canonical_dataarray_to_grib(
        data_var=temperature,
        file=f,  # 传递文件对象，而非路径字符串
        grib_keys=grib_keys,
        default_grib_keys={"centre": 98, "typeOfLevel": "surface"}
    )

技巧三：快速打印变量

场景：未知文件快速解析他山之石：pygrib meteva 变量过多，展示前20个意思一下

pygrib

import pygrib
def quick_scan(grib_path, max_show=20):
    with pygrib.open(grib_path) as grbs:
        for i, grb in enumerate(grbs):
            if i >= max_show:
                break
            print(f"{grb.messagenumber}: {grb.name} @ {grb.level} {grb.units}")
            
quick_scan("/home/mw/input/GFS1824/gfs_4_20230902_0000_021.grb2")

1: Pressure reduced to MSL @ 0 Pa
2: Cloud mixing ratio @ 1 kg kg**-1
3: Ice water mixing ratio @ 1 kg kg**-1
4: Rain mixing ratio @ 1 kg kg**-1
5: Snow mixing ratio @ 1 kg kg**-1
6: Graupel (snow pellets) @ 1 kg kg**-1
7: Derived radar reflectivity @ 1 dB
8: Derived radar reflectivity @ 2 dB
9: Maximum/Composite radar reflectivity @ 0 dB
10: Visibility @ 0 m
11: U component of wind @ 0 m s**-1
12: V component of wind @ 0 m s**-1
13: Ventilation Rate @ 0 m**2 s**-1
14: Wind speed (gust) @ 0 m s**-1
15: Geopotential height @ 1 gpm
16: Temperature @ 1 K
17: Relative humidity @ 1 %
18: Specific humidity @ 1 kg kg**-1
19: Vertical velocity @ 1 Pa s**-1
20: Geometric vertical velocity @ 1 m s**-1

meteva

#看看数据结构
import meteva.base as meb
meb.print_grib_file_info('/home/mw/input/GFS1824/gfs_4_20230902_0000_021.grb2',filter_by_keys = {'typeOfLevel': 'isobaricInhPa'})

<xarray.Dataset> Size: 343MB
Dimensions:        (isobaricInhPa: 33, latitude: 361, longitude: 720)
Coordinates:
    time           datetime64[ns] 8B ...
    step           timedelta64[ns] 8B ...
  * isobaricInhPa  (isobaricInhPa) float64 264B 1e+03 975.0 950.0 ... 2.0 1.0
  * latitude       (latitude) float64 3kB 90.0 89.5 89.0 ... -89.0 -89.5 -90.0
  * longitude      (longitude) float64 6kB 0.0 0.5 1.0 1.5 ... 358.5 359.0 359.5
    valid_time     datetime64[ns] 8B ...
Data variables:
    gh             (isobaricInhPa, latitude, longitude) float32 34MB ...
    t              (isobaricInhPa, latitude, longitude) float32 34MB ...
    r              (isobaricInhPa, latitude, longitude) float32 34MB ...
    q              (isobaricInhPa, latitude, longitude) float32 34MB ...
    w              (isobaricInhPa, latitude, longitude) float32 34MB ...
    wz             (isobaricInhPa, latitude, longitude) float32 34MB ...
    u              (isobaricInhPa, latitude, longitude) float32 34MB ...
    v              (isobaricInhPa, latitude, longitude) float32 34MB ...
    absv           (isobaricInhPa, latitude, longitude) float32 34MB ...
    o3mr           (isobaricInhPa, latitude, longitude) float32 34MB ...
Attributes:
    GRIB_edition:            2
    GRIB_centre:             kwbc
    GRIB_centreDescription:  US National Weather Service - NCEP
    GRIB_subCentre:          0
    Conventions:             CF-1.7
    institution:             US National Weather Service - NCEP
    history:                 2025-01-27T03:34 GRIB to CDM+CF via cfgrib-0.9.1...

技巧四 : 三维变量

总结

参考以下

https://blog.perillaroc.wang/post/2020/03/2020-03-16-grib-notebook-read-grib-using-cfgrib/

https://www.showdoc.com.cn/meteva/3975601833484385

当然还有更多特性，例如同时打开多个grib文件和使用dask加速处理等源自xarray的优势功能

这些留给有需要的同学自行探索吧

另外shortname 可以到以下网址查询

https://codes.ecmwf.int/grib/param-db/

本文参与腾讯云自媒体同步曝光计划，分享自微信公众号。

原始发表：2025-01-27，如有侵权请联系 cloudcommunity@tencent.com 删除

变量

本文分享自气python风雨微信公众号，前往查看

如有侵权，请联系 cloudcommunity@tencent.com 删除。

本文参与腾讯云自媒体同步曝光计划，欢迎热爱写作的你一起参与！

登录后参与评论

0 条评论

热度