Field arithmetic

earthkit-data supports standard arithmetic operators directly on Field and FieldList objects. Operations produce new objects and never modify their inputs.

The table below lists the supported operators.

Operator

Example

Description

+

f + g, f + 2

Addition

-

f - g, f - 2

Subtraction

*

f * g, f * 2

Multiplication

/

f / g, f / 2

True division

//

f // g, f // 2

Floor division

%

f % g, f % 2

Modulo

**

f ** 2

Power

-f

-f

Unary negation

>, <, >=, <=, !=

f > 270

Comparison (returns boolean array field)

[1]:
import numpy as np

import earthkit.data as ekd

Load data

Load a multi-field GRIB file so we have several fields to work with.

[2]:
ds = ekd.from_source("sample", "tuv_pl.grib").to_fieldlist()
len(ds)
[2]:
18
[3]:
# Extract individual fields
t_fields = ds.sel({"parameter.variable": "t"})  # temperature at all levels
u_fields = ds.sel({"parameter.variable": "u"})  # u-component of wind
v_fields = ds.sel({"parameter.variable": "v"})  # v-component of wind

t = t_fields[0]  # single temperature field
u = u_fields[0]  # single u-wind field
v = v_fields[0]  # single v-wind field

print("Shape:", t.geography.shape())
print("Level:", t.vertical.level())
print("Param:", t.parameter.variable())
Shape: (7, 12)
Level: 1000
Param: t

Scalar arithmetic

Arithmetic with a scalar operates element-wise on every grid point. The result is a new Field whose metadata is inherited from the original.

[4]:
# Convert temperature from Kelvin to Celsius
t_celsius = t - 273.15
print(f"T(K)  mean: {t.values.mean():.2f}")
print(f"T(°C) mean: {t_celsius.values.mean():.2f}")
T(K)  mean: 279.71
T(°C) mean: 6.56
[5]:
# Scaling and offset
t_scaled = t * 0.5 + 100
print(f"scaled mean: {t_scaled.values.mean():.2f}")

# Field as a right hand side operand
t_scaled = 0.5 * t + 100
print(f"scaled mean: {t_scaled.values.mean():.2f}")

# Power
t_sq = t**2
print(f"squared mean: {t_sq.values.mean():.2f}")

# Unary negation
t_neg = -t
np.testing.assert_allclose(t_neg.values, -t.values)
print(f"negative mean: {t_neg.values.mean():.2f}")
scaled mean: 239.85
scaled mean: 239.85
squared mean: 78623.10
negative mean: -279.71

Field–field arithmetic

When both operands are Field objects with the same grid, operations are applied point-wise. A common example is computing wind speed from its components.

[6]:
# Wind speed: sqrt(u^2 + v^2)
wind_speed = (u**2 + v**2) ** 0.5
print(f"u mean:          {u.values.mean():.2f} m/s")
print(f"v mean:          {v.values.mean():.2f} m/s")
print(f"wind speed mean: {wind_speed.values.mean():.2f} m/s")
u mean:          -0.38 m/s
v mean:          -0.07 m/s
wind speed mean: 7.45 m/s
[7]:
# Difference between two fields (e.g. temperature anomaly relative to first level)
t_850 = t_fields.sel({"vertical.level": 850})[0]
t_500 = t_fields.sel({"vertical.level": 500})[0]
thickness = t_850 - t_500
print(f"T@850 mean: {t_850.values.mean():.2f} K")
print(f"T@500 mean: {t_500.values.mean():.2f} K")
print(f"diff  mean: {thickness.values.mean():.2f} K")
T@850 mean: 272.73 K
T@500 mean: 252.89 K
diff  mean: 19.84 K

FieldList arithmetic

Arithmetic between two FieldList objects applies the operation element-wise across corresponding fields (i.e. field 0 op field 0, field 1 op field 1, …). The two fieldlists must have the same length.

[8]:
# Convert all temperature fields from K to °C at once
t_celsius_all = t_fields - 273.15
print(f"Original length: {len(t_fields)}, result length: {len(t_celsius_all)}")
for i, (orig, conv) in enumerate(zip(t_fields, t_celsius_all)):
    print(f"  level {orig.vertical.level():4d}:  {orig.values.mean():.2f} K  ->  {conv.values.mean():.2f} °C")
Original length: 6, result length: 6
  level 1000:  279.71 K  ->  6.56 °C
  level  850:  272.73 K  ->  -0.42 °C
  level  700:  265.17 K  ->  -7.98 °C
  level  500:  252.89 K  ->  -20.26 °C
  level  400:  241.53 K  ->  -31.62 °C
  level  300:  227.84 K  ->  -45.31 °C
[9]:
# Element-wise wind speed across all levels
wind_speed_all = (u_fields**2 + v_fields**2) ** 0.5
print("Wind speed at each level (mean m/s):")
for f in wind_speed_all:
    print(f"  level {f.vertical.level():4d}: {f.values.mean():.2f}")
Wind speed at each level (mean m/s):
  level 1000: 7.45
  level  850: 9.14
  level  700: 9.22
  level  500: 12.39
  level  400: 15.83
  level  300: 19.13

Comparison operators

Comparison operators return a new Field containing boolean (0/1) values, which is useful as a mask.

[10]:
# Find grid points above freezing
above_freezing = t > 273.15
fraction = above_freezing.values.mean()
print(f"Fraction of grid points above 0°C: {fraction:.2%}")

# Use as a mask: zero out sub-freezing temperatures
t_warm = t * (t > 273.15)
print(f"Masked field min: {t_warm.values.min():.2f} K  (should be 0 where frozen)")
Fraction of grid points above 0°C: 58.33%
Masked field min: 0.00 K  (should be 0 where frozen)

Using numpy ufuncs

numpy universal functions work transparently on Field objects via earthkit.data.utils.compute.apply_ufunc(). You can also extract values with field.values, apply any numpy operation, and reconstruct a field with field.set(values=...).

[11]:
# Using values directly with numpy, then setting back
log_t = t.set(values=np.log(t.values))
print(f"log(T) mean: {log_t.values.mean():.4f}")

# Verify the metadata is preserved
print(f"param: {log_t.parameter.variable()}, level: {log_t.vertical.level()}")
log(T) mean: 5.6312
param: t, level: 1000

apply_ufunc

apply_ufunc is earthkit-data’s mechanism for applying arbitrary callables — including numpy ufuncs and custom functions — directly to Field and FieldList objects without manually extracting and re-inserting values.

from earthkit.data.utils.compute import apply_ufunc

result = apply_ufunc(func, *args)
  • func receives plain numpy arrays and must return a numpy array.

  • args can be Field, FieldList, scalar, or numpy array objects.

  • earthkit-data iterates over matching fields, calls func on their value arrays, and wraps each result back into a Field, copying metadata from the first FieldList argument.

  • If one FieldList has length 1, it is broadcast across all fields of the other.

[12]:
from earthkit.data.utils.compute import apply_ufunc

Unary functions

Pass a single Field or FieldList to apply any unary callable across every field.

[13]:
# np.sin applied to all temperature fields in one call
t_sin = apply_ufunc(np.sin, t_fields)
print(f"sin(T) fields: {len(t_sin)}")
for orig, result in zip(t_fields, t_sin):
    print(f"  level {orig.vertical.level():4d}: sin(T) mean = {result.values.mean():.6f}")
sin(T) fields: 6
  level 1000: sin(T) mean = 0.398520
  level  850: sin(T) mean = -0.029065
  level  700: sin(T) mean = 0.161524
  level  500: sin(T) mean = -0.067683
  level  400: sin(T) mean = -0.096225
  level  300: sin(T) mean = -0.189231
[14]:
# Custom unary function — clip temperatures to [250, 310] K
def clip_temp(arr):
    return np.clip(arr, 250.0, 310.0)


t_clipped = apply_ufunc(clip_temp, t_fields)
for orig, clipped in zip(t_fields, t_clipped):
    print(
        f"  level {orig.vertical.level():4d}: "
        f"original [{orig.values.min():.1f}, {orig.values.max():.1f}] K  ->  "
        f"clipped  [{clipped.values.min():.1f}, {clipped.values.max():.1f}] K"
    )
  level 1000: original [240.6, 320.6] K  ->  clipped  [250.0, 310.0] K
  level  850: original [232.5, 304.5] K  ->  clipped  [250.0, 304.5] K
  level  700: original [223.3, 287.3] K  ->  clipped  [250.0, 287.3] K
  level  500: original [227.8, 271.8] K  ->  clipped  [250.0, 271.8] K
  level  400: original [216.0, 264.0] K  ->  clipped  [250.0, 264.0] K
  level  300: original [206.7, 250.7] K  ->  clipped  [250.0, 250.7] K

Binary functions

Pass two Field / FieldList arguments for element-wise binary operations. np.hypot(u, v) computes wind speed more precisely than (u**2 + v**2)**0.5 because it avoids intermediate overflow.

[15]:
# Wind speed via np.hypot — avoids overflow for large values
wind_hypot = apply_ufunc(np.hypot, u_fields, v_fields)
wind_manual = (u_fields**2 + v_fields**2) ** 0.5

print("Wind speed comparison (hypot vs manual):")
for h, m in zip(wind_hypot, wind_manual):
    print(f"  level {h.vertical.level():4d}: hypot={h.values.mean():.4f}  manual={m.values.mean():.4f}")
Wind speed comparison (hypot vs manual):
  level 1000: hypot=7.4502  manual=7.4502
  level  850: hypot=9.1436  manual=9.1436
  level  700: hypot=9.2167  manual=9.2167
  level  500: hypot=12.3901  manual=12.3901
  level  400: hypot=15.8330  manual=15.8330
  level  300: hypot=19.1302  manual=19.1302
[16]:
# Custom binary function — wind direction (meteorological convention, degrees from N)
def wind_direction(u_arr, v_arr):
    # atan2 gives mathematical angle; convert to met convention (0=N, clockwise)
    return (270.0 - np.degrees(np.arctan2(v_arr, u_arr))) % 360.0


wdir = apply_ufunc(wind_direction, u_fields, v_fields)
print("Wind direction at each level (mean degrees from N):")
for f in wdir:
    print(f"  level {f.vertical.level():4d}: {f.values.mean():.1f}°")
Wind direction at each level (mean degrees from N):
  level 1000: 162.6°
  level  850: 166.3°
  level  700: 184.6°
  level  500: 241.0°
  level  400: 189.3°
  level  300: 250.3°

Broadcasting a single field across a FieldList

When one argument has length 1 it is broadcast against every field in the other argument. Here we subtract the 1000 hPa temperature from every level to produce a temperature anomaly profile.

[17]:
t_1000 = t_fields.sel({"vertical.level": 1000})  # single-field FieldList (length 1)
print(f"Reference field count: {len(t_1000)}, all-level count: {len(t_fields)}")

# Broadcast: subtract t_1000 from every level
t_anomaly = apply_ufunc(lambda a, b: a - b, t_fields, t_1000)

print("\nTemperature anomaly relative to 1000 hPa (mean K):")
for orig, anom in zip(t_fields, t_anomaly):
    print(f"  level {orig.vertical.level():4d}: {anom.values.mean():+.2f} K")
Reference field count: 1, all-level count: 6

Temperature anomaly relative to 1000 hPa (mean K):
  level 1000: +0.00 K
  level  850: -6.98 K
  level  700: -14.54 K
  level  500: -26.82 K
  level  400: -38.18 K
  level  300: -51.86 K
[ ]: