如何在 grib 中迭代？

Question

我正在尝试遍历 grib 文件以获取所有数据，但我不知道该怎么做。我是用 Python 做的。

我正在使用此代码

import xarray as xr

ds = xr.open_dataset("myFile.grib2", engine="pynio")

for v in ds:
    print(ds[v])

我得到的结果是这样的一些输出：

<xarray.DataArray 'PERPW_P0_L1_GLL0' (lat_0: 336, lon_0: 720)>
[241920 values with dtype=float32]
Coordinates:
  * lat_0    (lat_0) float32 90.0 89.5 89.0 88.5 ... -76.0 -76.5 -77.0 -77.5
  * lon_0    (lon_0) float32 0.0 0.5 1.0 1.5 2.0 ... 358.0 358.5 359.0 359.5
Attributes:
    center:                                         US National Weather Servi...
    production_status:                              Operational products
    long_name:                                      Primary wave mean period
    units:                                          s
    grid_type:                                      Latitude/longitude
    parameter_discipline_and_category:              Oceanographic products, W...
    parameter_template_discipline_category_number:  [ 0 10  0 11]
    level_type:                                     Ground or water surface
    level:                                          [1.]
    forecast_time:                                  [219]
    forecast_time_units:                            hours
    initial_time:                                   01/18/2021 (00:00)

我正在使用来自 NOAA 的 file。

据我所知，这就像表格的列。例如，如何获取所有 level 值的数据？

Answer 1

xarray 提供了一个网格化数据的接口（如您所见in the documentation）。

我使用包cfgrib 来访问GRIB 数据。它是由 ECMWF (https://github.com/ecmwf/cfgrib) 开发的。

所以，如果你这样做：

import xarray as xr

d = xr.open_dataset('gefs.wave.t00z.c00.global.0p25.f000.grib2', engine = 'cfgrib')
print(d)

你明白了：

<xarray.Dataset>
Dimensions:     (latitude: 721, level: 3, longitude: 1440)
Coordinates:
    time        datetime64[ns] 2020-12-01
    step        timedelta64[ns] 00:00:00
    surface     int32 1
  * latitude    (latitude) float64 90.0 89.75 89.5 89.25 ... -89.5 -89.75 -90.0
  * longitude   (longitude) float64 0.0 0.25 0.5 0.75 ... 359.3 359.5 359.8
    valid_time  datetime64[ns] 2020-12-01
  * level       (level) int32 1 2 3
Data variables:
    ws          (latitude, longitude) float32 nan nan nan nan ... nan nan nan
    wdir        (latitude, longitude) float32 ...
    u           (latitude, longitude) float32 ...
    v           (latitude, longitude) float32 ...
    siconc      (latitude, longitude) float32 ...
    swh         (latitude, longitude) float32 ...
    paramId_0   (latitude, longitude) float32 ...
    mwp         (latitude, longitude) float32 ...
    perpw       (latitude, longitude) float32 ...
    mwd         (latitude, longitude) float32 ...
    dirpw       (latitude, longitude) float32 ...
    shww        (latitude, longitude) float32 ...
    swell       (level, latitude, longitude) float32 ...
    mpww        (latitude, longitude) float32 ...
    swper       (level, latitude, longitude) float32 ...
    wvdir       (latitude, longitude) float32 ...
    swdir       (level, latitude, longitude) float32 ...
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:                 2021-02-23T16:14:32 GRIB to CDM+CF via cfgrib-0....

每个数据变量都有自己的坐标，可以通过xarray访问：

print(d['swh'])
<xarray.DataArray 'swh' (latitude: 721, longitude: 1440)>
array([[nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       ...,
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan]], dtype=float32)
Coordinates:
    time        datetime64[ns] 2020-12-01
    step        timedelta64[ns] 00:00:00
    surface     int32 1
  * latitude    (latitude) float64 90.0 89.75 89.5 89.25 ... -89.5 -89.75 -90.0
  * longitude   (longitude) float64 0.0 0.25 0.5 0.75 ... 359.3 359.5 359.8
    valid_time  datetime64[ns] 2020-12-01
Attributes:
    GRIB_paramId:                             140229
    GRIB_shortName:                           swh
    GRIB_units:                               m
    GRIB_name:                                Significant height of combined ...
    GRIB_cfVarName:                           swh
    GRIB_dataType:                            an
    GRIB_missingValue:                        9999
    GRIB_numberOfPoints:                      1038240
    GRIB_typeOfLevel:                         surface
    GRIB_NV:                                  0
    GRIB_stepUnits:                           1
    GRIB_stepType:                            instant
    GRIB_gridType:                            regular_ll
    GRIB_gridDefinitionDescription:           Latitude/longitude. Also called...
    GRIB_Nx:                                  1440
    GRIB_iDirectionIncrementInDegrees:        0.25
    GRIB_iScansNegatively:                    0
    GRIB_longitudeOfFirstGridPointInDegrees:  0.0
    GRIB_longitudeOfLastGridPointInDegrees:   359.750016
    GRIB_Ny:                                  721
    GRIB_jDirectionIncrementInDegrees:        0.25
    GRIB_jPointsAreConsecutive:               0
    GRIB_jScansPositively:                    0
    GRIB_latitudeOfFirstGridPointInDegrees:   90.0
    GRIB_latitudeOfLastGridPointInDegrees:    -90.0
    long_name:                                Significant height of combined ...
    units:                                    m

变量是DataArray，但您可以使用values 访问底层网格数据，因此：

d['swh'].values
array([[nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       ...,
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan],
       [nan, nan, nan, ..., nan, nan, nan]], dtype=float32)

或使用to_* 函数导出，例如导出到 Pandas DataFrame：

print(d['swh'].to_dataframe().head())
                         time   step  surface valid_time  swh
latitude longitude                                           
90.0     0.00      2020-12-01 0 days        1 2020-12-01  NaN
         0.25      2020-12-01 0 days        1 2020-12-01  NaN
         0.50      2020-12-01 0 days        1 2020-12-01  NaN
         0.75      2020-12-01 0 days        1 2020-12-01  NaN
         1.00      2020-12-01 0 days        1 2020-12-01  NaN

Answer 2

我通常将 pygrib 用于 GRIB 文件，但我认为 xarray 也不错。

G = pygrib.open('file.grb'):
for g in G:
    print(g)

例如获取这样的坐标和值...

latlons = g.latlons()
values = g.values

查看所有其他可用键的数据.. 需要根据数据选择有效时间和预测步距。

print(g.keys())

Answer 3

我不确定您所说的 level 值是什么意思，但对于您链接的文件，所有变量都是随时间变化的二维字段。它们实际上并不像表格的列，而是二维网格的时间序列。