Writing GRIB data to Zarr

# get input GRIB data
import earthkit.data as ekd
ds = ekd.from_source("sample", "pl.grib")

This data contains 32 fields: several forecasts on pressure levels for 2 parameters. We can check its content with describe().

ds.describe()

		level	date	time	step	paramId	class	stream	type	experimentVersionNumber
shortName	typeOfLevel
r	isobaricInhPa	700,500	20240603,20240604	0,1200	0,6	157	od	oper	fc	0001
t	isobaricInhPa	700,500	20240603,20240604	0,1200	0,6	130	od	oper	fc	0001

Using to_target() on the data object

We use to_target() to write the GRIB fieldlist/field into a zarr store. First, the data is converted to Xarray then xarray.Dataset.to_zarr() is called to generate the zarr store. We need to set the kwargs accordingly.

# with these options each field will be a separate chunk
ds.to_target("zarr",
             earthkit_to_xarray_kwargs={"chunks": {"forecast_reference_time": 1,
                                                   "step": 1,
                                                   "level": 1}},
             xarray_to_zarr_kwargs={"store": "_pl.zarr", "mode": "w"})

/opt/homebrew/Caskroom/miniforge/base/envs/dev/lib/python3.11/site-packages/zarr/api/asynchronous.py:205: UserWarning: Consolidated metadata is currently not part in the Zarr format 3 specification. It may not be supported by other zarr implementations and may change in the future.
  warnings.warn(

import zarr
root = zarr.group("_pl.zarr")
root.tree()

/
├── forecast_reference_time (4,) int64
├── latitude (19,) float64
├── level (2,) int64
├── longitude (36,) float64
├── r (4, 2, 2, 19, 36) float64
├── step (2,) int64
└── t (4, 2, 2, 19, 36) float64

root["t"].info

Type               : Array
Zarr format        : 3
Data type          : DataType.float64
Shape              : (4, 2, 2, 19, 36)
Chunk shape        : (1, 1, 1, 19, 36)
Order              : C
Read-only          : False
Store type         : LocalStore
Filters            : ()
Serializer         : BytesCodec(endian=<Endian.little: 'little'>)
Compressors        : (ZstdCodec(level=0, checksum=False),)
No. bytes          : 87552 (85.5K)

The zarr store can be loaded to Xarray to check its content.

import xarray
xarray.open_dataset("_pl.zarr")

/var/folders/93/w0p869rx17q98wxk83gn9ys40000gn/T/ipykernel_45349/754541422.py:2: FutureWarning: In a future version of xarray decode_timedelta will default to False rather than None. To silence this warning, set decode_timedelta to True, False, or a 'CFTimedeltaCoder' instance.
  xarray.open_dataset("_pl.zarr")

<xarray.Dataset> Size: 176kB
Dimensions:                  (step: 2, longitude: 36,
                              forecast_reference_time: 4, latitude: 19, level: 2)
Coordinates:
  * step                     (step) timedelta64[ns] 16B 00:00:00 06:00:00
  * longitude                (longitude) float64 288B 0.0 10.0 ... 340.0 350.0
  * forecast_reference_time  (forecast_reference_time) datetime64[ns] 32B 202...
  * latitude                 (latitude) float64 152B 90.0 80.0 ... -80.0 -90.0
  * level                    (level) int64 16B 500 700
Data variables:
    r                        (forecast_reference_time, step, level, latitude, longitude) float64 88kB ...
    t                        (forecast_reference_time, step, level, latitude, longitude) float64 88kB ...
Attributes:
    class:        od
    stream:       oper
    levtype:      pl
    type:         fc
    expver:       0001
    date:         20240603
    time:         0
    domain:       g
    number:       0
    Conventions:  CF-1.8
    institution:  ECMWF

xarray.Dataset

• Dimensions:
  • step: 2
  • longitude: 36
  • forecast_reference_time: 4
  • latitude: 19
  • level: 2
• Coordinates: (5)
  • step
    (step)
    timedelta64[ns]
    00:00:00 06:00:00

    array([             0, 21600000000000], dtype='timedelta64[ns]')

  • longitude
    (longitude)
    float64
    0.0 10.0 20.0 ... 330.0 340.0 350.0

    units :

    degrees_east

    standard_name :

    longitude

    long_name :

    longitude

    array([  0.,  10.,  20.,  30.,  40.,  50.,  60.,  70.,  80.,  90., 100., 110.,
           120., 130., 140., 150., 160., 170., 180., 190., 200., 210., 220., 230.,
           240., 250., 260., 270., 280., 290., 300., 310., 320., 330., 340., 350.])

  • forecast_reference_time
    (forecast_reference_time)
    datetime64[ns]
    2024-06-03 ... 2024-06-04T12:00:00

    standard_name :

    forecast_reference_time

    long_name :

    initial time of forecast

    array(['2024-06-03T00:00:00.000000000', '2024-06-03T12:00:00.000000000',
           '2024-06-04T00:00:00.000000000', '2024-06-04T12:00:00.000000000'],
          dtype='datetime64[ns]')

  • latitude
    (latitude)
    float64
    90.0 80.0 70.0 ... -80.0 -90.0

    units :

    degrees_north

    standard_name :

    latitude

    long_name :

    latitude

    array([ 90.,  80.,  70.,  60.,  50.,  40.,  30.,  20.,  10.,   0., -10., -20.,
           -30., -40., -50., -60., -70., -80., -90.])

  • level
    (level)
    int64
    500 700

    units :

    hPa

    positive :

    down

    stored_direction :

    decreasing

    standard_name :

    air_pressure

    long_name :

    pressure

    array([500, 700])

• Data variables: (2)
  • r
    (forecast_reference_time, step, level, latitude, longitude)
    float64
    ...

    param :

    r

    standard_name :

    relative_humidity

    long_name :

    Relative humidity

    paramId :

    157

    units :

    %

    [10944 values with dtype=float64]

  • t
    (forecast_reference_time, step, level, latitude, longitude)
    float64
    ...

    param :

    t

    standard_name :

    air_temperature

    long_name :

    Temperature

    paramId :

    130

    units :

    K

    [10944 values with dtype=float64]

• Indexes: (5)
  • step
    PandasIndex

    PandasIndex(TimedeltaIndex(['0 days 00:00:00', '0 days 06:00:00'], dtype='timedelta64[ns]', name='step', freq=None))

  • longitude
    PandasIndex

    PandasIndex(Index([  0.0,  10.0,  20.0,  30.0,  40.0,  50.0,  60.0,  70.0,  80.0,  90.0,
           100.0, 110.0, 120.0, 130.0, 140.0, 150.0, 160.0, 170.0, 180.0, 190.0,
           200.0, 210.0, 220.0, 230.0, 240.0, 250.0, 260.0, 270.0, 280.0, 290.0,
           300.0, 310.0, 320.0, 330.0, 340.0, 350.0],
          dtype='float64', name='longitude'))

  • forecast_reference_time
    PandasIndex

    PandasIndex(DatetimeIndex(['2024-06-03 00:00:00', '2024-06-03 12:00:00',
                   '2024-06-04 00:00:00', '2024-06-04 12:00:00'],
                  dtype='datetime64[ns]', name='forecast_reference_time', freq=None))

  • latitude
    PandasIndex

    PandasIndex(Index([ 90.0,  80.0,  70.0,  60.0,  50.0,  40.0,  30.0,  20.0,  10.0,   0.0,
           -10.0, -20.0, -30.0, -40.0, -50.0, -60.0, -70.0, -80.0, -90.0],
          dtype='float64', name='latitude'))

  • level
    PandasIndex

    PandasIndex(Index([500, 700], dtype='int64', name='level'))

• Attributes: (11)

  class :

  od

  stream :

  oper

  levtype :

  pl

  type :

  fc

  expver :

  0001

  date :

  20240603

  time :

  0

  domain :

  g

  number :

  0

  Conventions :

  CF-1.8

  institution :

  ECMWF