GRIB: modifying values

In this notebook we will show how to do computations with GRIB data.

First we read some GRIB data containing pressure level fields.

[1]:
import earthkit.data as ekd

fl = ekd.from_source("sample", "tuv_pl.grib").to_fieldlist()

Using fieldlist arithmetic

Select temperature fields.

[2]:
t = fl.sel({"parameter.variable": "t"})

We can use maths operators directly on the fieldlist to modify the values.

[3]:
# add 1 to all the field values in t
t1 = t + 1

The result is a new fieldlist with modified field data values but the same metadata. All the values are stored in memory.

We compare the first field from t and t1 to verify the results.

[4]:
t[0].ls()
[4]:
parameter.variable time.valid_datetime time.base_datetime time.step vertical.level vertical.level_type ensemble.member geography.grid_type
0 t 2018-08-01 12:00:00 2018-08-01 12:00:00 0 days 1000 pressure 0 regular_ll
[5]:
t1[0].ls()
[5]:
parameter.variable time.valid_datetime time.base_datetime time.step vertical.level vertical.level_type ensemble.member geography.grid_type
0 t 2018-08-01 12:00:00 2018-08-01 12:00:00 0 days 1000 pressure 0 regular_ll
[6]:
print(t[0].values.max(), t1[0].values.max())
320.5641784667969 321.5641784667969

The modified fieldlist can be saved into a GRIB file.

[7]:
t1.to_target("file", "_res_1.grib")

# read back the data and compare the values in the first field
t1_w = ekd.from_source("file", "_res_1.grib").to_fieldlist()
print(t[0].values.max(), t1_w[0].values.max())
320.5641784667969 321.5641784667969

Looping over fields

The same computation can be carried out by iterating over the fields and modifying them one by one using set(). Then we can create a new fieldlist from the resulting fields with create_fieldlist(). Just like in the previous example all the resulting fields are stored in memory.

[8]:
res = []
for f in t:
    v = f.values
    v1 = v + 1
    f1 = f.set(values=v1)
    res.append(f1)

t1 = ekd.create_fieldlist(res)

The modified fieldlist can be saved into a GRIB file.

[9]:
t1.to_target("file", "_res_2.grib")

# read back the data and compare the values in the first field
t1_w = ekd.from_source("file", "_res_2.grib").to_fieldlist()
print(t[0].values.max(), t1_w[0].values.max())
320.5641784667969 321.5641784667969

Writing modified fields to disk immediately

This example demonstrates how to perform the computation field by field and write the resulting field immediately to disk. This is the most efficient solution in terms of memory usage.

[10]:
with ekd.create_target("file", "_res_3.grib") as target:
    for f in t:
        v = f.values
        v1 = v + 1
        f1 = f.set(values=v1)
        target.write(f1)
[11]:
# read back the data and compare the values in the first field
t1_w = ekd.from_source("file", "_res_3.grib").to_fieldlist()
print(t[0].values.max(), t1_w[0].values.max())
320.5641784667969 321.5641784667969

Modified fields and the associated GRIB message

When a field is created from GRIB data the associated GRIB message can be accessed via the field using message(). Having modified the field values this GRIB message is still available and the data that message() returns will contain the updated values.

The following example demonstrates this by creating new fields from the messages in the original and modified fields and comparing their values.

[12]:
from earthkit.data.field.grib.create import create_grib_field_from_message

f = fl[0]
f1 = f.set(values=f.values + 1)

f_m = create_grib_field_from_message(f.message())
f1_m = create_grib_field_from_message(f1.message())

f_m.values[0], f1_m.values[1]
[12]:
(np.float64(272.5641784667969), np.float64(273.5641784667969))

If the metadata is also modified the associated GRIB message is not updated and we cannot access it any longer in the new field. The same is true for any raw GRIB metadata.

[13]:
f = fl[0]
f1 = f.set({"values": f.values + 1, "vertical.level": 850})
[14]:
f1.message()
[15]:
f1.get("metadata.shortName")