GRIB: getting latitudes, longitudes and values (regular LL grid)

In this example we will work with GRIB data on a regular latitude-longitude grid.

import earthkit.data as ekd

We will only use the temperature fields, so we read the file and extract them with sel().

ds_in = ekd.from_source("sample", "tuv_pl.grib").to_fieldlist()
ds = ds_in.sel({"parameter.variable": "t"})
ds.ls()

	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
1	t	2018-08-01 12:00:00	2018-08-01 12:00:00	0 days	850	pressure	0	regular_ll
2	t	2018-08-01 12:00:00	2018-08-01 12:00:00	0 days	700	pressure	0	regular_ll
3	t	2018-08-01 12:00:00	2018-08-01 12:00:00	0 days	500	pressure	0	regular_ll
4	t	2018-08-01 12:00:00	2018-08-01 12:00:00	0 days	400	pressure	0	regular_ll
5	t	2018-08-01 12:00:00	2018-08-01 12:00:00	0 days	300	pressure	0	regular_ll

In the data each field is defined on a regular latitude-longitude grid. These grids can be represented as a 2D array so the field shape is always 2 dimensional:

ds[0].shape

(7, 12)

Using the data() method

The simplest way to access the latitudes, longitudes and values is using the data() methods:

• earthkit.data.core.field.Field.data()

• earthkit.data.core.fieldlist.FieldList.data()

Fields

data() returns the latitude, longitude and values arrays from a field as an ndarray. By default the field’s shape is kept.

llv = ds[0].data()
type(llv), llv.shape

(numpy.ndarray, (3, 7, 12))

Here llv[0] contains the latitudes, llv[1] the longitudes, while llv[2] the values. Each of these arrays has the same shape that of the field (7, 12)).

print(llv[0].shape)
# first latitude
llv[0, 0, 0]

(7, 12)

np.float64(90.0)

print(llv[1].shape)
# first longitude
llv[1, 0, 0]

(7, 12)

np.float64(0.0)

print(llv[2].shape)
# first value
llv[2, 0, 0]

(7, 12)

np.float64(272.5641784667969)

When using the flatten keyword 1D arrays are returned:

llv = ds[0].data(flatten=True)
llv.shape

(3, 84)

FieldLists

data() only works on a fieldlist if all the fields have the same grid. The first two elements of the resulting ndarray are the latitude and longitude arrays (shared between fields), while the rest of the elements are the value arrays per field. Since we have 6 fields in our data the size of the first axis of the resulting ndarray is 2+6=8. Each of these arrays has the same shape that of the field (7, 12)).

llv = ds.data()
type(llv), llv.shape

(numpy.ndarray, (8, 7, 12))

Latitudes can be accessed as llv[0].

print(llv[0].shape)
# first latitude
llv[0, 0, 0]

(7, 12)

np.float64(90.0)

Longitudes can be accessed as llv[1].

print(llv[1].shape)
# first longitude
llv[1, 0, 0]

(7, 12)

np.float64(0.0)

The field values can be accessed as llv[2:].

print(llv[2:].shape)

# first value in first field
print(llv[2, 0, 0])

# first value from all the fields
print(llv[2:, 0, 0])

(6, 7, 12)
272.5641784667969
[272.56417847 272.53916931 271.26531982 255.84306335 244.00323486
 226.65315247]

When using the flatten keyword 1D arrays are returned for latitude and longitude, and the values array will be flattened per field:

llv = ds.data(flatten=True)
llv.shape

(8, 84)

Using the latlons() and to_numpy() methods

We can also get the latitudes, longitudes and values by using the latlons() and to_numpy() methods.

Fields

latlons() returns a dict:

lat, lon = ds[0].geography.latlons()
lat.shape, lon.shape

((7, 12), (7, 12))

The field values can be accessed with to_numpy().

v = ds[0].to_numpy()
v.shape

(7, 12)

By default both methods keep the field’s shape, but we can use the flatten keyword to get 1D arrays:

lat, lon = ds[0].geography.latlons(flatten=True)
v = ds[0].to_numpy(flatten=True)
lat.shape, lon.shape, v.shape

((84,), (84,), (84,))

FieldLists

latlons() only works on a fieldlist if all the fields have the same grid. In this case it returns the same dict that we would get for any of the fields:

lat, lon = ds.geography.latlons()
lat.shape, lon.shape

((7, 12), (7, 12))

to_numpy() only works on a fieldlist if all the fields has the same grid. It returns the array of the field value arrays.

v = ds.to_numpy()
v.shape

(6, 7, 12)

E.g. the first value from each field can be extracted as:

v[:, 0, 0]

array([272.56417847, 272.53916931, 271.26531982, 255.84306335,
       244.00323486, 226.65315247])

Specifying the array type

For all the methods above we can set the array type with the dtype keyword both for fields and fieldlists:

import numpy as np

v = ds.to_numpy(dtype=np.float32)
v[:, 0, 0]

array([272.56418, 272.53918, 271.26532, 255.84306, 244.00323, 226.65315],
      dtype=float32)

llv = ds.data(dtype=np.float32)

llv[:, 0, 0]

array([ 90.     ,   0.     , 272.56418, 272.53918, 271.26532, 255.84306,
       244.00323, 226.65315], dtype=float32)