ICON2VTK - ICON netCDF to VTK converter

ICON2VTK is a lightweight Python command-line tool that converts cell-based ICON model output from netCDF into visualization-ready files for direct use in ParaView, with legacy VTK and XDMF/HDF5 field export support.

The script supports a broader workflow for exploratory visualization:

• export one or more 2-D ICON cell fields by selecting time and level indices
• subset the domain by a bounding box or circular region
• filter triangles by exported field value ranges
• coarsen the mesh by one or more ICON refinement levels using parent-child metadata
• apply configurable radius offsets to fields and overlays
• write legacy VTK in ASCII or binary format
• write XDMF with a companion HDF5 file for large HPC-style datasets
• add coastline, river, country-boundary, province-boundary, and longitude-latitude graticule overlays
• print field statistics and optionally write them to CSV

The main script is:

• icon2vtk.py

What the script does

ICON commonly stores atmospheric fields on an unstructured triangular grid. In the bundled sample files in this repository, the mesh geometry is stored separately from the actual model variables:

• the data file contains the variable values, for example ts, pr, or ta
• the grid file contains the triangular connectivity and the vertex coordinates

The script combines both:

1. It reads the ICON mesh from the grid file.
2. It reads one selected variable from the netCDF data file.
3. It resolves that variable to one or more 2-D slices by selecting time and level indices when needed.
4. It writes one or more field files that ParaView can display as unstructured surfaces, either as legacy VTK or as XDMF with an HDF5 sidecar, on the sphere or in plate-carree projection.

Optionally, it can also write additional VTK polyline files containing:

• coastlines
• rivers
• country boundaries
• province boundaries
• graticule lines

These can be loaded into the same ParaView scene as overlays.

Input assumptions

The current implementation is designed for cell-based ICON fields, meaning the exported variable must have an ncells dimension.

Typical supported shapes are:

• (time, ncells) for 2-D time-dependent fields
• (time, height, ncells) for 3-D fields on vertical levels
• (time, singleton_dim, ncells) for variables like 2 m temperature or 10 m wind

The script supports at most one non-singleton non-cell dimension besides time. In practice this means one horizontal field per timestep, optionally with one additional vertical or level-like index.

The script uses the ICON grid file to obtain:

• vertex coordinates
• triangle connectivity
• cell center coordinates for region selection

Main files in this directory

• icon2vtk.py: the converter
• data/icon_grid_*.nc: ICON horizontal grid description
• data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc: example 2-D ICON output
• data/aes_amip_atm_3d_qp_ml_19790101T000000Z.nc: example 3-D ICON output

Requirements

Required Python packages:

• numpy
• netCDF4

Optional packages used only for Natural Earth overlays:

• cartopy
• shapely

If you do not request Natural Earth overlays, the script does not need Cartopy or Shapely.

Installation

For a simple local setup, you can use the bundled setup.sh helper script to create a virtual environment and install the dependencies from requirements.txt:

bash setup.sh

The script creates .venv, activates it, upgrades pip, and installs the required Python packages.

If you prefer to set up the environment manually, run:

python3 -m venv .venv
source .venv/bin/activate
pip install -r requirements.txt

After that, you can run the converter with:

python3 icon2vtk.py --help

If you do not need Natural Earth overlay support, the required packages are only:

• numpy
• netCDF4

If cartopy is difficult to install on your system, you can still use the script without coastline, river, country-boundary, or province-boundary export.

Quick start with example data

The easiest way to try the converter is to run the bundled example script from the example/ directory:

cd example
bash run.sh

This writes a small set of VTK files back into example/:

• one sphere example using a 3-D field, one time index, one vertical level, coarsening, and overlays
• one plate-carree example using a 2-D field, clipped seam handling, and overlays

This is a good first check that:

• your Python environment can read the ICON netCDF files from data/
• the optional overlay dependencies are installed if you requested coastlines
• ParaView can open the generated field and overlay files together

The example directory also contains saved ParaView state files and the corresponding generated files for both scenes:

• example/sphere.pvsm
• example/sphere_ta_t1_l45.vtk
• example/sphere_coastlines.vtk
• example/sphere_graticule.vtk
• example/plate_carree.pvsm
• example/plate_carree_ts_bbox.vtk
• example/plate_carree_coastlines.vtk
• example/plate_carree_graticule.vtk

Example screenshots are included as a quick visual reference:

Sphere example

Plate-carree example

Basic command structure

The minimal command is:

python3 icon2vtk.py DATA.nc GRID.nc VARIABLE

For example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts

This reads the variable ts from the ICON netCDF data file, uses the ICON grid from the grid file, and writes ts.vtk in binary VTK format by default when you run the command from the repository root.

For larger datasets, you can switch the field export to XDMF with a companion HDF5 file:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --field-format xdmf \
  -o example/ts.xdmf

This writes:

• example/ts.xdmf
• example/ts.h5

The XDMF file contains the mesh and array metadata, while the HDF5 sidecar holds the heavy numeric arrays. This is usually the better choice when the legacy VTK output becomes too large for efficient storage or transfer.

If you want multiple XDMF field exports to reuse one mesh file, add --xdmf-shared-grid and point it to a common HDF5 path:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  --field-format xdmf \
  --xdmf-shared-grid example/shared_grid.h5 \
  -o example/ts.xdmf

Then export another variable with the same processed mesh settings:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  pr \
  --time-index 1 \
  --field-format xdmf \
  --xdmf-shared-grid example/shared_grid.h5 \
  -o example/pr.xdmf

In that setup:

• example/shared_grid.h5 stores the shared points and cells once
• example/ts.h5 stores only the ts values
• example/pr.h5 stores only the pr values

This saves disk space when several exported variables share exactly the same processed mesh. The tool validates that the shared grid file matches the requested mesh and raises an error if it does not. That means all geometry-affecting options must stay the same across those exports, including projection, precision, subsetting, coarsening, seam handling, and field radius offset.

--time-index and --level-index accept either a single index or a comma-separated list of indices. When you pass multiple indices, the script exports one field file per selected slice combination and appends suffixes such as _t0, _t1, _l45, or _t1_l45 to the output filename.

For example, to export two timesteps of a 2-D field:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 0,1 \
  -o example/ts_batch.vtk

This writes:

• example/ts_batch_t0.vtk
• example/ts_batch_t1.vtk

If you want the per-slice field statistics as a CSV file, add --stats-output:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 0,1 \
  -o example/ts_batch.vtk \
  --stats-output example/ts_stats.csv

This writes one CSV row per exported field slice with the columns:

• output_file
• variable
• time_index
• level_index
• min
• max
• mean
• count
• nan_count

If a selected dimension is not actually used by the exported variable, the corresponding time_index or level_index value is written as nan in the CSV output.

If the grid file provides the ICON variable parent_cell_index, you can coarsen the exported field by one or more refinement levels:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --coarsen-level 1 \
  -o example/ts_coarse.vtk

This groups complete four-child sibling families using parent_cell_index, reconstructs the parent triangle, and writes the average of the sibling values onto that coarser cell. Higher values such as --coarsen-level 2 or 3 request repeated collapse of the same 4:1 refinement pattern.

Coarsening is applied on the full global ICON mesh first. If you also request --bbox or --circle, that regional filtering is applied afterward on the already coarsened mesh using the coarse cell centers.

Conceptually, each coarsening step replaces every complete four-child family by its parent face, carries the mean of the four child values onto that parent, and then repeats the same process on the next coarser mesh level.

Important details:

• --coarsen-level requires parent_cell_index in the grid file; otherwise the script exits with an error
• the requested level is an upper bound, not a guarantee; the script stops early if no further complete sibling families can be collapsed and reports requested=... applied=... after the export
• if you also subset the domain, incomplete families near the subset boundary are kept at their current resolution instead of being forced to coarsen
• higher coarsening levels should currently be treated as approximate for visualization; the cell count may collapse as expected, but exact parent-face connectivity is not guaranteed at the deepest levels

You can also remove uninteresting triangles directly during export by filtering on the selected field values:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  --field-min 270 \
  --field-max 280 \
  -o example/ts_270_280.vtk

This keeps only cells whose scalar value is finite and lies inside the closed interval [270, 280]. The filter is applied after optional coarsening and spatial subsetting. --field-min and --field-max affect only the exported field triangles; overlays are not filtered by field value and still only follow --bbox and --circle. This lets you combine value-filtered fields with the same contextual overlays to build layered ParaView scenes from different variables or value bands.

Overlay-only mode

You can also generate overlay VTK files without exporting any ICON field.

This is useful when you only want to recreate the overlays for an existing ParaView scene.

For example, to generate only a graticule in plate-carree projection:

python3 icon2vtk.py \
  --projection plate-carree \
  --graticule-output example/graticule_only.vtk \
  --graticule-spacing 60 30

Or to generate only country boundaries on the sphere:

python3 icon2vtk.py \
  --country-output example/countries_only.vtk \
  --country-resolution 10m

In overlay-only mode:

• no DATA.nc GRID.nc VARIABLE positional arguments are required
• --bbox and --circle still work
• --projection, seam handling, and radius offsets still apply
• the radius defaults to 6371229 m unless you pass --radius

The same fixed default radius is also used for field export. If you need a different scale, pass --radius explicitly.

Listing variables in a netCDF file

Before exporting, it is often useful to inspect which variables are available in a netCDF file without using ncdump.

You can do that with:

python3 icon2vtk.py data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc --list-variables

or, for the 3-D example file:

python3 icon2vtk.py data/aes_amip_atm_3d_qp_ml_19790101T000000Z.nc --list-variables

The output lists each variable with its metadata and dimensions, followed by a coordinate summary that maps available time and level indices to their values.

Example:

- ta: long_name="Temperature"; standard_name="air_temperature"; units=[K]
    dims=time, height, ncells; shape=(2, 90, 20480); grid=unstructured; dtype=float32

Coordinate values:

    time[0] = 1979-01-01T00:00:00
    time[1] = 1979-01-02T00:00:00

    height[0] = 1
    ...

This makes it easier to decide:

• which variable to export
• whether the variable is a 2-D field directly or needs time/level selection
• whether it lives on the ICON cell grid
• whether you will need --time-index and --level-index

A first useful example

This is a good starting point for testing with a 2-D field:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  -o example/ts_t1.vtk

This tells the script:

• use the 2-D daily output file
• use the ICON grid file
• export the variable ts
• select the second time record via --time-index 1
• write the output to example/ts_t1.vtk

After the export, the script prints a short summary including:

• number of cells written
• number of finite values
• minimum value
• maximum value
• mean value

This is useful because some fields can be nearly zero, and the summary immediately tells you whether the data range is meaningful for visualization.

Exporting 3-D variables

For variables with one additional non-singleton dimension besides ncells, the script needs:

• a time index
• a level index for that extra dimension

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_3d_qp_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ta \
  --time-index 1 \
  --level-index 45 \
  -o example/ta_t1_l45.vtk

Here:

• ta is air temperature
• --time-index 1 selects the second available timestep
• --level-index 45 selects one vertical model level

Both options also accept comma-separated lists such as --time-index 0,1 or --level-index 10,20,50. In that case, the script exports one field file per selected slice combination.

The output is still a surface mesh file, not a full 3-D volume. In other words, the script writes one horizontal slice over the ICON sphere for the chosen level.

If a variable has more than one non-singleton non-cell dimension besides time, the script rejects it instead of guessing how --level-index should be applied.

ASCII versus binary legacy VTK

By default, the script writes binary legacy VTK files:

--vtk-format binary

You usually do not need to specify this explicitly.

Binary output is recommended for most use cases:

• smaller file size
• faster read/write performance
• better suited for large ICON grids and batch processing
• when ParaView loading speed matters

ASCII output is mainly useful for debugging or inspection:

• human-readable text format
• easier to inspect or diff
• convenient for small test cases

However, ASCII files are significantly larger and slower to read and write.

Floating-point VTK output precision is controlled separately:

--vtk-precision float32

This setting affects VTK point coordinates and field scalar values. The default is float32, which usually reduces file size and write time. Use --vtk-precision float64 if you prefer double-precision geometry and scalar output.

To use ASCII instead of binary:

--vtk-format ascii

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  --vtk-format ascii \
  -o example/ts_ascii.vtk

Natural Earth line overlays

The script can generate separate VTK polyline files from several Natural Earth line datasets via Cartopy:

• coastlines
• rivers
• country boundaries
• province boundaries

This is useful because the ICON grid itself does not carry general-purpose linework for these map references in a clean overlay-ready form.

Coastlines

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  -o example/ts_with_coast.vtk \
  --coastline-output example/coastlines_10m.vtk \
  --coastline-resolution 10m \
  --coastline-radius-offset 10000

This command produces two files:

• example/ts_with_coast.vtk: the field on the ICON sphere
• example/coastlines_10m.vtk: the coastline overlay

Coastline resolution

You can select one of three Natural Earth resolutions:

• 110m: coarse, fast, small files
• 50m: medium detail
• 10m: highest detail, much larger files

Typical choice:

• use 110m for quick checks
• use 10m for more detailed regional figures

Coastline radius offset

Coastlines can be written at a slightly larger radius than the field:

--coastline-radius-offset 10000

In practice, offsets on the order of 10 to 20 km tend to work better than 1 to 2 km for avoiding visual overlap or z-fighting when multiple datasets lie at nearly the same radius.

Rivers

Rivers are exported from Natural Earth rivers_lake_centerlines as polyline overlays.

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  -o example/ts_with_rivers.vtk \
  --river-output example/rivers_10m.vtk \
  --river-resolution 10m \
  --river-radius-offset 13000

Country boundaries

Country boundaries are exported from Natural Earth admin_0_boundary_lines_land.

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  -o example/ts_with_countries.vtk \
  --country-output example/country_boundaries_10m.vtk \
  --country-resolution 10m \
  --country-radius-offset 16000

Province boundaries

Province or state boundaries are exported from Natural Earth admin_1_states_provinces_lines.

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  -o example/ts_with_provinces.vtk \
  --province-output example/province_boundaries_10m.vtk \
  --province-resolution 10m \
  --province-radius-offset 16000

All Natural Earth line overlays support the same 110m, 50m, and 10m resolution choices and the same region filtering, projection, seam handling, and radius-offset behavior.

Graticule overlays

The script can also create a geographic graticule, meaning a set of longitude and latitude lines.

This is often useful when showing regional subsets or when you want geographic orientation in a global view.

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  -o example/ts_with_graticule.vtk \
  --graticule-output example/graticule_30x15.vtk \
  --graticule-spacing 30 15 \
  --graticule-radius-offset 19000

This writes:

• the field file
• a separate polyline VTK file containing longitude lines every 30 degrees and latitude lines every 15 degrees

Why the name “graticule”?

Graticule is the standard cartographic term for the network of meridians and parallels, that is, longitude and latitude grid lines. It is the precise technical term, although “lat/lon grid” would also be understandable.

Graticule spacing

The option:

--graticule-spacing DLON DLAT

controls the spacing in degrees.

Examples:

• --graticule-spacing 30 15
• --graticule-spacing 60 30
• --graticule-spacing 10 10

Smaller spacing gives more lines and a denser visual grid.

Graticule radius offset

Like coastlines, the graticule can be lifted above the sphere:

--graticule-radius-offset 19000

This is often useful when field surfaces, coastlines, and grid lines are shown together.

Field radius offset

The script can also lift the field surface itself:

--field-radius-offset 5000

This is especially useful when you want to display a regional subset above another surface without adding a separate Transform step in ParaView.

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  --circle 10 50 1500 \
  --field-radius-offset 5000 \
  -o example/ts_circle_offset.vtk \
  --coastline-output example/coastlines_circle_offset.vtk \
  --coastline-resolution 10m \
  --coastline-radius-offset 10000

In that example:

• the field is lifted by 5000 m
• the coastlines are lifted by 10000 m

So the coastlines still remain slightly above the field.

Regional subsetting by bounding box

You can restrict the export to a longitude-latitude box:

--bbox lon_min lat_min lon_max lat_max

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  --bbox -15 30 35 72 \
  -o example/ts_europe.vtk \
  --coastline-output example/coastlines_europe.vtk \
  --coastline-resolution 10m \
  --coastline-radius-offset 10000

How bbox selection works

For the ICON field itself:

• each triangular cell has a center longitude and latitude
• a cell is included if its center lies inside the requested box

For overlays:

• the overlay linework is clipped or filtered to the same region

Dateline crossing

If lon_min > lon_max, the script interprets the box as crossing the dateline.

That makes boxes such as:

170 -20 -170 20

possible.

Regional subsetting by circle

You can also select a circular region defined by:

• center longitude
• center latitude
• radius in kilometers

Syntax:

--circle lon_center lat_center radius_km

Example:

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  --circle 10 50 1500 \
  -o example/ts_circle.vtk \
  --coastline-output example/coastlines_circle.vtk \
  --coastline-resolution 10m \
  --coastline-radius-offset 10000

How circle selection works

For ICON cells:

• the script computes the great-circle distance from each cell center to the requested center point
• cells inside the requested radius are kept

For overlays:

• the lines are filtered to the same circular region

This is useful for compact regional views centered on a city, country, or basin.

Mutual exclusivity

You can use either:

• --bbox
• or --circle

but not both in the same command.

Full example with field and multiple overlays

python3 icon2vtk.py \
  data/aes_amip_atm_2d_P1D_ml_19790101T000000Z.nc \
  data/icon_grid_0049_R02B04_G.nc \
  ts \
  --time-index 1 \
  --circle 10 50 1500 \
  --field-radius-offset 5000 \
  -o example/ts_combo.vtk \
  --coastline-output example/coastlines_combo.vtk \
  --coastline-resolution 10m \
  --coastline-radius-offset 10000 \
  --river-output example/rivers_combo.vtk \
  --river-resolution 10m \
  --river-radius-offset 13000 \
  --country-output example/countries_combo.vtk \
  --country-resolution 10m \
  --country-radius-offset 16000 \
  --graticule-output example/graticule_combo.vtk \
  --graticule-spacing 30 15 \
  --graticule-radius-offset 19000

This produces:

• a field surface
• a coastline overlay
• a river overlay
• a country-boundary overlay
• a longitude-latitude graticule

all as separate files that can be opened together in ParaView.

ParaView workflow

A typical workflow in ParaView is:

1. Open the main field file.
2. Open any overlay files such as coastlines, rivers, boundaries, and graticule.
3. Click Apply.
4. Color the field by the exported scalar.
5. Set overlays to solid contrasting colors.
6. Increase Line Width for the line datasets.

Useful settings for overlays:

• Representation: Surface or Wireframe
• Coloring: Solid Color
• Line Width: 2 to 4

If the line datasets lie exactly on the field surface, they can visually interfere with it. That is the main reason the script supports separate radius offsets.

Notes on the output format

The script writes legacy VTK field and overlay files by default, and also supports XDMF field export with a companion HDF5 sidecar.

That means:

• field surfaces are written as UNSTRUCTURED_GRID
• overlays are written as legacy VTK POLYDATA
• XDMF field exports reference triangle topology and cell-centered data stored in HDF5

Legacy VTK remains the default because it is simple and easy to generate directly from Python without additional VTK dependencies in the writer path, while XDMF/HDF5 is available as the more scalable option for large field datasets.

Current limitations

At the moment, the script is focused on one specific use case:

• cell-based ICON variables on the horizontal grid

It does not currently handle:

• edge-based variables
• vertex-based variables
• vector field export as VTK vectors
• XML VTK output formats such as .vtu or .vtp
• XDMF export for overlays

Getting help

To see the full command-line help:

python3 icon2vtk.py --help

Summary of the most important options

• --time-index: select one or more time records
• --level-index: select one or more vertical levels for slice export
• --stats-output: write per-slice field statistics to a CSV file
• --coarsen-level N: coarsen by up to N ICON refinement levels when parent_cell_index is available
• --vtk-format ascii|binary: choose legacy VTK encoding
• --vtk-precision float32|float64: choose floating-point precision for VTK points and scalar arrays
• --bbox ...: select a rectangular region
• --circle ...: select a circular region
• --field-radius-offset: lift the field surface
• --coastline-output: write coastline overlay
• --coastline-resolution: choose Natural Earth detail
• --coastline-radius-offset: lift coastline overlay
• --river-output: write river overlay
• --river-resolution: choose Natural Earth river detail
• --river-radius-offset: lift river overlay
• --country-output: write country-boundary overlay
• --country-resolution: choose Natural Earth country-boundary detail
• --country-radius-offset: lift country-boundary overlay
• --province-output: write province-boundary overlay
• --province-resolution: choose Natural Earth province-boundary detail
• --province-radius-offset: lift province-boundary overlay
• --graticule-output: write longitude-latitude grid overlay
• --graticule-spacing: choose graticule spacing
• --graticule-radius-offset: lift graticule overlay

Contact

If you have questions or need support, please feel free to get in touch.

Dr. Lars Hoffmann

Jülich Supercomputing Centre, Forschungszentrum Jülich

e-mail: l.hoffmann@fz-juelich.de