Remapping a dataset

In this tutorial we will remap a dataset from the native ICON grid to HEALPix. The dataset used as an example is from the EERIE project and is available online. It consists of five 2D variables and covers a period of 23 years, with a time step of six hours.

Note

Please note that all of the methods used in this example are generic and can be applied to any combination of source and target grids.

import healpix as hp
import numpy as np
import xarray as xr

import easygems.healpix as egh
import easygems.remap as egr


ds = xr.open_dataset(
    "https://eerie.cloud.dkrz.de/datasets/icon-esm-er.eerie-control-1950.v20240618.atmos.native.2d_6h_inst/kerchunk",
    engine="zarr",
    chunks={},
).squeeze()
ds

<xarray.Dataset> Size: 12TB
Dimensions:             (time: 112496, ncells: 5242880)
Coordinates:
    cell_sea_land_mask  (ncells) int32 21MB dask.array<chunksize=(5242880,), meta=np.ndarray>
    height              float64 8B 10.0
    height_2            float64 8B 2.0
    lat                 (ncells) float64 42MB dask.array<chunksize=(5242880,), meta=np.ndarray>
    lon                 (ncells) float64 42MB dask.array<chunksize=(5242880,), meta=np.ndarray>
  * time                (time) datetime64[ns] 900kB 1991-01-01T05:59:59 ... 2...
Dimensions without coordinates: ncells
Data variables:
    hus2m               (time, ncells) float32 2TB dask.array<chunksize=(1, 5242880), meta=np.ndarray>
    psl                 (time, ncells) float32 2TB dask.array<chunksize=(1, 5242880), meta=np.ndarray>
    ts                  (time, ncells) float32 2TB dask.array<chunksize=(1, 5242880), meta=np.ndarray>
    uas                 (time, ncells) float32 2TB dask.array<chunksize=(1, 5242880), meta=np.ndarray>
    vas                 (time, ncells) float32 2TB dask.array<chunksize=(1, 5242880), meta=np.ndarray>
Attributes: (12/31)
    Conventions:           CF-1.7 CMIP-6.2
    activity_id:           EERIE
    data_specs_version:    01.00.32
    forcing_index:         1
    initialization_index:  1
    license:               EERIE model data produced by MPI-M is licensed und...
    ...                    ...
    parent_activity_id:    EERIE
    sub_experiment_id:     none
    experiment:            coupled control with fixed 1950's forcing (HighRes...
    source:                ICON-ESM-ER (2023): \naerosol: none, prescribed MA...
    institution:           Max Planck Institute for Meteorology, Hamburg 2014...
    sub_experiment:        none

xarray.Dataset

• Dimensions:
  • time: 112496
  • ncells: 5242880
• Coordinates: (6)
  • cell_sea_land_mask
    (ncells)
    int32
    dask.array<chunksize=(5242880,), meta=np.ndarray>

    long_name :

    sea (-2 inner, -1 boundary) land (2 inner, 1 boundary) mask for the cell

    units :

    2,1,-1,-

    grid_type :

    unstructured

    number_of_grid_in_reference :

    1

    dtype :

    int32

    zlib :

    False

    szip :

    False

    zstd :

    False

    bzip2 :

    False

    blosc :

    False

    shuffle :

    False

    complevel :

    0

    fletcher32 :

    False

    contiguous :

    True

    chunksizes :

    None

    source :

    /work/bm1344/DKRZ/kerchunks_batched/erc2002/icon_grid_0033_R02B08_G.nc

    original_shape :

    [5242880]

    Array Chunk Bytes 20.00 MiB 20.00 MiB Shape (5242880,) (5242880,) Dask graph 1 chunks in 2 graph layers Data type int32 numpy.ndarray

  • height
    ()
    float64
    10.0

    axis :

    Z

    long_name :

    height

    positive :

    up

    standard_name :

    height

    units :

    m

    array(10.)

  • height_2
    ()
    float64
    2.0

    axis :

    Z

    long_name :

    height

    positive :

    up

    standard_name :

    height

    units :

    m

    array(2.)

  • lat
    (ncells)
    float64
    dask.array<chunksize=(5242880,), meta=np.ndarray>

    long_name :

    center latitude

    units :

    radian

    standard_name :

    grid_latitude

    bounds :

    clat_vertices

    dtype :

    float64

    zlib :

    False

    szip :

    False

    zstd :

    False

    bzip2 :

    False

    blosc :

    False

    shuffle :

    False

    complevel :

    0

    fletcher32 :

    False

    contiguous :

    True

    chunksizes :

    None

    source :

    /work/bm1344/DKRZ/kerchunks_batched/erc2002/icon_grid_0033_R02B08_G.nc

    original_shape :

    [5242880]

    Array Chunk Bytes 40.00 MiB 40.00 MiB Shape (5242880,) (5242880,) Dask graph 1 chunks in 2 graph layers Data type float64 numpy.ndarray

  • lon
    (ncells)
    float64
    dask.array<chunksize=(5242880,), meta=np.ndarray>

    long_name :

    center longitude

    units :

    radian

    standard_name :

    grid_longitude

    bounds :

    clon_vertices

    dtype :

    float64

    zlib :

    False

    szip :

    False

    zstd :

    False

    bzip2 :

    False

    blosc :

    False

    shuffle :

    False

    complevel :

    0

    fletcher32 :

    False

    contiguous :

    True

    chunksizes :

    None

    source :

    /work/bm1344/DKRZ/kerchunks_batched/erc2002/icon_grid_0033_R02B08_G.nc

    original_shape :

    [5242880]

    Array Chunk Bytes 40.00 MiB 40.00 MiB Shape (5242880,) (5242880,) Dask graph 1 chunks in 2 graph layers Data type float64 numpy.ndarray

  • time
    (time)
    datetime64[ns]
    1991-01-01T05:59:59 ... 2067-12-...

    array(['1991-01-01T05:59:59.000000000', '1991-01-01T11:59:59.000000000',
           '1991-01-01T17:59:59.000000000', ..., '2067-12-31T11:59:59.000000000',
           '2067-12-31T17:59:59.000000000', '2067-12-31T23:59:59.000000000'],
          shape=(112496,), dtype='datetime64[ns]')

• Data variables: (5)
  • hus2m
    (time, ncells)
    float32
    dask.array<chunksize=(1, 5242880), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    specific humidity in 2m

    number_of_grid_in_reference :

    1

    param :

    6.0.0

    standard_name :

    qv2m

    units :

    kg kg-1

    Array Chunk Bytes 2.15 TiB 20.00 MiB Shape (112496, 5242880) (1, 5242880) Dask graph 112496 chunks in 3 graph layers Data type float32 numpy.ndarray

  • psl
    (time, ncells)
    float32
    dask.array<chunksize=(1, 5242880), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    mean sea level pressure

    number_of_grid_in_reference :

    1

    param :

    1.3.0

    standard_name :

    mean sea level pressure

    units :

    Pa

    Array Chunk Bytes 2.15 TiB 20.00 MiB Shape (112496, 5242880) (1, 5242880) Dask graph 112496 chunks in 2 graph layers Data type float32 numpy.ndarray

  • ts
    (time, ncells)
    float32
    dask.array<chunksize=(1, 5242880), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    surface temperature

    number_of_grid_in_reference :

    1

    param :

    0.0.0

    standard_name :

    surface_temperature

    units :

    K

    Array Chunk Bytes 2.15 TiB 20.00 MiB Shape (112496, 5242880) (1, 5242880) Dask graph 112496 chunks in 2 graph layers Data type float32 numpy.ndarray

  • uas
    (time, ncells)
    float32
    dask.array<chunksize=(1, 5242880), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    zonal wind in 10m

    number_of_grid_in_reference :

    1

    param :

    2.2.0

    standard_name :

    uas

    units :

    m s-1

    Array Chunk Bytes 2.15 TiB 20.00 MiB Shape (112496, 5242880) (1, 5242880) Dask graph 112496 chunks in 3 graph layers Data type float32 numpy.ndarray

  • vas
    (time, ncells)
    float32
    dask.array<chunksize=(1, 5242880), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    meridional wind in 10m

    number_of_grid_in_reference :

    1

    param :

    3.2.0

    standard_name :

    vas

    units :

    m s-1

    Array Chunk Bytes 2.15 TiB 20.00 MiB Shape (112496, 5242880) (1, 5242880) Dask graph 112496 chunks in 3 graph layers Data type float32 numpy.ndarray

• Indexes: (1)
  • time
    PandasIndex

    PandasIndex(DatetimeIndex(['1991-01-01 05:59:59', '1991-01-01 11:59:59',
                   '1991-01-01 17:59:59', '1991-01-01 23:59:59',
                   '1991-01-02 05:59:59', '1991-01-02 11:59:59',
                   '1991-01-02 17:59:59', '1991-01-02 23:59:59',
                   '1991-01-03 05:59:59', '1991-01-03 11:59:59',
                   ...
                   '2067-12-29 17:59:59', '2067-12-29 23:59:59',
                   '2067-12-30 05:59:59', '2067-12-30 11:59:59',
                   '2067-12-30 17:59:59', '2067-12-30 23:59:59',
                   '2067-12-31 05:59:59', '2067-12-31 11:59:59',
                   '2067-12-31 17:59:59', '2067-12-31 23:59:59'],
                  dtype='datetime64[ns]', name='time', length=112496, freq=None))

• Attributes: (31)

  Conventions :

  CF-1.7 CMIP-6.2

  activity_id :

  EERIE

  data_specs_version :

  01.00.32

  forcing_index :

  1

  initialization_index :

  1

  license :

  EERIE model data produced by MPI-M is licensed under a Creative Commons Attribution 4.0 International License (https://creativecommons.org/licenses). The data producers and data providers make no warranty, either express or implied, including, but not limited to, warranties of merchantability and fitness for a particular purpose. All liabilities arising from the supply of the information (including any liability arising in negligence) are excluded to the fullest extent permitted by law.

  mip_era :

  EERIE

  project_id :

  EERIE

  physics_index :

  1

  product :

  model-output

  realization_index :

  1

  source_type :

  AOGCM

  variant_label :

  r1i1p1f1

  creation_date :

  2025-07-08T14:09:13

  institution_id :

  MPI-M

  source_id :

  ICON-ESM-ER

  experiment_id :

  eerie-control-1950

  contact :

  juergen.kroeger@mpimet.mpg.de

  references :

  Hohenegger et al., ICON-Sapphire: simulating the components of the Earth system and their interactions at kilometer and subkilometer scales. Geosci. Model Dev., 16, 779–811, 2023, https://doi.org/10.5194/gmd-16-779-2023

  nominal_resolution :

  10 km

  version_id :

  v20240618

  realm :

  atmos

  grid :

  native

  grid_label :

  gn

  parent_experiment_id :

  eerie-spinup-1950

  parent_activity_id :

  EERIE

  sub_experiment_id :

  none

  experiment :

  coupled control with fixed 1950's forcing (HighResMIP equivalent of pre-industrial control)

  source :

  ICON-ESM-ER (2023): aerosol: none, prescribed MACv2-SP atmos: ICON-A (icosahedral/triangles; 10 km; 90 levels; top level 80 km) atmosChem: none land: JSBACH4.20 landIce: none/prescribed ocean: ICON-O (icosahedral/triangles; 5 km; 72 levels; top grid cell 0-2 m) ocnBgchem: none seaIce: unnamed (thermodynamic (Semtner zero-layer) dynamic (Hibler 79) sea ice model)

  institution :

  Max Planck Institute for Meteorology, Hamburg 20146, Germany

  sub_experiment :

  none

The first step is to create a HEALPix grid that is close to the resolution of our source grid. Here we will choose an order of 9 (also known as “zoom”) and “nest” ordering.

order = zoom = 9
nside = hp.order2nside(order)
pixels = np.arange(hp.nside2npix(nside))

hp_lon, hp_lat = hp.pix2ang(nside=nside, ipix=pixels, lonlat=True, nest=True)

Next, we can use our defined source and target grids to compute interpolation weights. The easygems package provides a function to compute these weights using the Delaunay triangulation method.

weights = egr.compute_weights_delaunay(
    points=(ds.lon, ds.lat),
    xi=(hp_lon, hp_lat)
)

Tip

For larger grids, computing the remapping weights can be quite resource-intensive. However, you can save the calculated weights for future use:

weights.to_netcdf("healpix_weights.nc")

These weights can be applied to single fields directly:

egh.healpix_show(egr.apply_weights(ds.ts.isel(time=0), **weights));

/builds/easy/gems/.venv/lib/python3.13/site-packages/cartopy/io/__init__.py:241: DownloadWarning: Downloading: https://naturalearth.s3.amazonaws.com/110m_physical/ne_110m_coastline.zip
  warnings.warn(f'Downloading: {url}', DownloadWarning)

<Figure size 640x480 with 0 Axes>

Alternatively, we can use xarray’s apply_ufunc() function to lift the function onto a full dataset. This requires a coupled of additional information from the user, e.g. the input dimension along which the function should be applied, and the resulting output dimensions name and size.

ds_remap = xr.apply_ufunc(
    egr.apply_weights,
    ds,
    kwargs=weights,
    keep_attrs=True,
    input_core_dims=[["ncells"]],
    output_core_dims=[["cell"]],
    output_dtypes=["f4"],
    vectorize=True,
    dask="parallelized",
    dask_gufunc_kwargs={
        "output_sizes": {"cell": weights.sizes["tgt_idx"]},
    },
)
ds_remap

<xarray.Dataset> Size: 7TB
Dimensions:   (time: 112496, cell: 3145728)
Coordinates:
    height    float64 8B 10.0
    height_2  float64 8B 2.0
  * time      (time) datetime64[ns] 900kB 1991-01-01T05:59:59 ... 2067-12-31T...
Dimensions without coordinates: cell
Data variables:
    hus2m     (time, cell) float32 1TB dask.array<chunksize=(1, 3145728), meta=np.ndarray>
    psl       (time, cell) float32 1TB dask.array<chunksize=(1, 3145728), meta=np.ndarray>
    ts        (time, cell) float32 1TB dask.array<chunksize=(1, 3145728), meta=np.ndarray>
    uas       (time, cell) float32 1TB dask.array<chunksize=(1, 3145728), meta=np.ndarray>
    vas       (time, cell) float32 1TB dask.array<chunksize=(1, 3145728), meta=np.ndarray>
Attributes: (12/31)
    Conventions:           CF-1.7 CMIP-6.2
    activity_id:           EERIE
    data_specs_version:    01.00.32
    forcing_index:         1
    initialization_index:  1
    license:               EERIE model data produced by MPI-M is licensed und...
    ...                    ...
    parent_activity_id:    EERIE
    sub_experiment_id:     none
    experiment:            coupled control with fixed 1950's forcing (HighRes...
    source:                ICON-ESM-ER (2023): \naerosol: none, prescribed MA...
    institution:           Max Planck Institute for Meteorology, Hamburg 2014...
    sub_experiment:        none

xarray.Dataset

• Dimensions:
  • time: 112496
  • cell: 3145728
• Coordinates: (3)
  • height
    ()
    float64
    10.0

    axis :

    Z

    long_name :

    height

    positive :

    up

    standard_name :

    height

    units :

    m

    array(10.)

  • height_2
    ()
    float64
    2.0

    axis :

    Z

    long_name :

    height

    positive :

    up

    standard_name :

    height

    units :

    m

    array(2.)

  • time
    (time)
    datetime64[ns]
    1991-01-01T05:59:59 ... 2067-12-...

    array(['1991-01-01T05:59:59.000000000', '1991-01-01T11:59:59.000000000',
           '1991-01-01T17:59:59.000000000', ..., '2067-12-31T11:59:59.000000000',
           '2067-12-31T17:59:59.000000000', '2067-12-31T23:59:59.000000000'],
          shape=(112496,), dtype='datetime64[ns]')

• Data variables: (5)
  • hus2m
    (time, cell)
    float32
    dask.array<chunksize=(1, 3145728), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    specific humidity in 2m

    number_of_grid_in_reference :

    1

    param :

    6.0.0

    standard_name :

    qv2m

    units :

    kg kg-1

    Array Chunk Bytes 1.29 TiB 12.00 MiB Shape (112496, 3145728) (1, 3145728) Dask graph 112496 chunks in 6 graph layers Data type float32 numpy.ndarray

  • psl
    (time, cell)
    float32
    dask.array<chunksize=(1, 3145728), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    mean sea level pressure

    number_of_grid_in_reference :

    1

    param :

    1.3.0

    standard_name :

    mean sea level pressure

    units :

    Pa

    Array Chunk Bytes 1.29 TiB 12.00 MiB Shape (112496, 3145728) (1, 3145728) Dask graph 112496 chunks in 5 graph layers Data type float32 numpy.ndarray

  • ts
    (time, cell)
    float32
    dask.array<chunksize=(1, 3145728), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    surface temperature

    number_of_grid_in_reference :

    1

    param :

    0.0.0

    standard_name :

    surface_temperature

    units :

    K

    Array Chunk Bytes 1.29 TiB 12.00 MiB Shape (112496, 3145728) (1, 3145728) Dask graph 112496 chunks in 5 graph layers Data type float32 numpy.ndarray

  • uas
    (time, cell)
    float32
    dask.array<chunksize=(1, 3145728), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    zonal wind in 10m

    number_of_grid_in_reference :

    1

    param :

    2.2.0

    standard_name :

    uas

    units :

    m s-1

    Array Chunk Bytes 1.29 TiB 12.00 MiB Shape (112496, 3145728) (1, 3145728) Dask graph 112496 chunks in 6 graph layers Data type float32 numpy.ndarray

  • vas
    (time, cell)
    float32
    dask.array<chunksize=(1, 3145728), meta=np.ndarray>

    CDI_grid_type :

    unstructured

    long_name :

    meridional wind in 10m

    number_of_grid_in_reference :

    1

    param :

    3.2.0

    standard_name :

    vas

    units :

    m s-1

    Array Chunk Bytes 1.29 TiB 12.00 MiB Shape (112496, 3145728) (1, 3145728) Dask graph 112496 chunks in 6 graph layers Data type float32 numpy.ndarray

• Indexes: (1)
  • time
    PandasIndex

    PandasIndex(DatetimeIndex(['1991-01-01 05:59:59', '1991-01-01 11:59:59',
                   '1991-01-01 17:59:59', '1991-01-01 23:59:59',
                   '1991-01-02 05:59:59', '1991-01-02 11:59:59',
                   '1991-01-02 17:59:59', '1991-01-02 23:59:59',
                   '1991-01-03 05:59:59', '1991-01-03 11:59:59',
                   ...
                   '2067-12-29 17:59:59', '2067-12-29 23:59:59',
                   '2067-12-30 05:59:59', '2067-12-30 11:59:59',
                   '2067-12-30 17:59:59', '2067-12-30 23:59:59',
                   '2067-12-31 05:59:59', '2067-12-31 11:59:59',
                   '2067-12-31 17:59:59', '2067-12-31 23:59:59'],
                  dtype='datetime64[ns]', name='time', length=112496, freq=None))

• Attributes: (31)

  Conventions :

  CF-1.7 CMIP-6.2

  activity_id :

  EERIE

  data_specs_version :

  01.00.32

  forcing_index :

  1

  initialization_index :

  1

  license :

  EERIE model data produced by MPI-M is licensed under a Creative Commons Attribution 4.0 International License (https://creativecommons.org/licenses). The data producers and data providers make no warranty, either express or implied, including, but not limited to, warranties of merchantability and fitness for a particular purpose. All liabilities arising from the supply of the information (including any liability arising in negligence) are excluded to the fullest extent permitted by law.

  mip_era :

  EERIE

  project_id :

  EERIE

  physics_index :

  1

  product :

  model-output

  realization_index :

  1

  source_type :

  AOGCM

  variant_label :

  r1i1p1f1

  creation_date :

  2025-07-08T14:09:13

  institution_id :

  MPI-M

  source_id :

  ICON-ESM-ER

  experiment_id :

  eerie-control-1950

  contact :

  juergen.kroeger@mpimet.mpg.de

  references :

  Hohenegger et al., ICON-Sapphire: simulating the components of the Earth system and their interactions at kilometer and subkilometer scales. Geosci. Model Dev., 16, 779–811, 2023, https://doi.org/10.5194/gmd-16-779-2023

  nominal_resolution :

  10 km

  version_id :

  v20240618

  realm :

  atmos

  grid :

  native

  grid_label :

  gn

  parent_experiment_id :

  eerie-spinup-1950

  parent_activity_id :

  EERIE

  sub_experiment_id :

  none

  experiment :

  coupled control with fixed 1950's forcing (HighResMIP equivalent of pre-industrial control)

  source :

  ICON-ESM-ER (2023): aerosol: none, prescribed MACv2-SP atmos: ICON-A (icosahedral/triangles; 10 km; 90 levels; top level 80 km) atmosChem: none land: JSBACH4.20 landIce: none/prescribed ocean: ICON-O (icosahedral/triangles; 5 km; 72 levels; top grid cell 0-2 m) ocnBgchem: none seaIce: unnamed (thermodynamic (Semtner zero-layer) dynamic (Hibler 79) sea ice model)

  institution :

  Max Planck Institute for Meteorology, Hamburg 20146, Germany

  sub_experiment :

  none

The resulting dataset can then be used as usual and the remapping is performed on demand.

egh.healpix_show(ds_remap.ts.sel(time="1991-01-01").mean("time"));

<Figure size 640x480 with 0 Axes>

Attaching a coordinate reference system

It is good practice to store map projection information in the Coordinate Reference Systems (CRS). By making the crs coordinate non-dimensional it will “stick” to the dataset even if individual variables are subselected.

crs = xr.DataArray(
    name="crs",
    attrs={
        "grid_mapping_name": "healpix",
        "healpix_nside": nside,
        "healpix_order": "nest",
    },
)
ds_remap = ds_remap.assign_coords(crs=crs)