Chunked Arrays and Zarr Explained

Problem statement

A satellite time series is naturally a four-dimensional array: time, band, y, x. Nothing in the raster world stores that well. A GeoTIFF is two dimensions plus bands; a directory of GeoTIFFs makes time an afterthought and every read a file open.

Zarr stores the array as it is — n-dimensional, chunked into blocks, with each chunk a separate object and the metadata in a small JSON file.

The consequence that matters: the chunk shape decides which queries are fast. A cube chunked for maps is slow for time series and vice versa, and the difference can be two orders of magnitude in the number of objects read.

Quick answer

import xarray as xr

cube = xr.open_zarr("s3://bucket/cube.zarr")
print(cube)
print(cube["reflectance"].encoding["chunks"])

<xarray.Dataset>
Dimensions:  (time: 120, band: 4, y: 509, x: 543)
Chunks:      (1, 4, 509, 543)

That chunking — one time step per chunk, the whole spatial extent — makes "give me one date's map" a single chunk read and "give me one pixel's history" 120 chunk reads, each of which discards all but one value.

The chunk shape is the whole performance model. Choose it from the access pattern.

Step-by-step solution

1. Understand what a chunk is

A chunk is the unit of storage and the unit of transfer. Each is compressed independently and stored as one object with a name encoding its position — 0.3.1.2 for the chunk at index (0, 3, 1, 2).

Reading any element requires fetching, decompressing and discarding the rest of its chunk. That is the same principle as a COG tile, extended to n dimensions.

2. Choose the chunk shape from the query, not the data

For a cube of shape (time=120, band=4, y=509, x=543):

access pattern	chunk shape	reads per query
one date, whole scene	`(1, 4, 509, 543)`	1
one pixel, all dates	`(1, 4, 509, 543)`	120
one pixel, all dates	`(120, 4, 64, 64)`	1
one date, whole scene	`(120, 4, 64, 64)`	72

Neither is wrong; they answer different questions. If both matter, either store two chunkings, or compromise — (30, 4, 128, 128) costs 4 reads for a map and 4 for a pixel history rather than 1 and 120.

3. Size the chunks sensibly

Aim for a few megabytes uncompressed. Too small and the metadata and per-object overhead dominate — a million tiny objects is slow on any store and expensive on some. Too large and every read pulls data you discard.

chunk_bytes = np.prod(chunk_shape) * dtype_size

(1, 4, 509, 543) of float32 is 4.4 MB — a reasonable single chunk. (120, 4, 64, 64) is 7.9 MB, also fine.

4. Understand that the metadata is tiny and separate

Zarr stores array shape, dtype, chunk shape, compression and attributes in small JSON documents. A client reads those first — a few kilobytes — and then knows exactly which chunk objects to request.

That is why opening a Zarr store is fast regardless of its size, and why xr.open_zarr returns instantly on a terabyte: it has read the metadata and nothing else.

5. Combine with dask for larger-than-memory work

cube = xr.open_zarr(path, chunks={"time": 30})
median = cube["reflectance"].median("time").compute()

chunks= makes the arrays dask arrays, so operations build a graph and execute in chunk-sized pieces. Aligning the dask chunks with the storage chunks avoids reading each stored chunk more than once.

A Zarr store as a small JSON metadata document plus many independently compressed chunk objects named by index. — Metadata first, then only the chunks the query touches. No server, no index scan.

Code examples

Example 1 — writing a cube with a deliberate chunking

import numpy as np
import xarray as xr


def write_cube(arrays, times, bands, transform, crs, path,
               time_chunk=1, spatial_chunk=None):
    """Write a time-band-y-x cube with an explicit chunk shape."""
    stack = np.stack(arrays)                      # (time, band, y, x)
    n_time, n_band, height, width = stack.shape

    a, b, c, d, e, f = transform[:6]
    xs = c + (np.arange(width) + 0.5) * a
    ys = f + (np.arange(height) + 0.5) * e

    cube = xr.Dataset(
        {"reflectance": (("time", "band", "y", "x"), stack.astype("float32"))},
        coords={"time": np.array(times, dtype="datetime64[s]"),
                "band": list(bands), "y": ys, "x": xs},
        attrs={"crs": str(crs)},
    )

    chunks = (time_chunk, n_band,
              spatial_chunk or height, spatial_chunk or width)
    bytes_per_chunk = np.prod(chunks) * 4
    print(f"  cube {stack.shape}, chunks {chunks} "
          f"({bytes_per_chunk / 1e6:.1f} MB each, "
          f"{np.prod([int(np.ceil(s / c)) for s, c in zip(stack.shape, chunks)])} chunks)")
    if bytes_per_chunk < 1e6:
        print("  ! chunks under 1 MB — per-object overhead will dominate")
    if bytes_per_chunk > 100e6:
        print("  ! chunks over 100 MB — every read pulls a lot you discard")

    cube["reflectance"].encoding["chunks"] = chunks
    cube.to_zarr(path, mode="w")
    return cube

Printing the chunk size and count at write time is the cheapest possible guard. A store with a million 4 kB chunks is unusable and looks fine until someone tries to read it.

Example 2 — measuring which access pattern your chunking suits

import numpy as np


def chunk_read_cost(shape, chunks, queries):
    """How many chunks does each query touch?"""
    print(f"  array {shape}, chunks {chunks}")
    for name, selection in queries.items():
        touched = 1
        for dim, (size, chunk) in enumerate(zip(shape, chunks)):
            span = selection.get(dim, size)
            touched *= int(np.ceil(span / chunk)) if span > 1 else 1
        useful = np.prod([selection.get(d, s) for d, s in enumerate(shape)])
        fetched = touched * np.prod(chunks)
        print(f"    {name:28} {touched:6,} chunks, "
              f"{fetched / max(useful, 1):8.1f}x amplification")


chunk_read_cost(
    (120, 4, 509, 543), (1, 4, 509, 543),
    {"one date, whole scene": {0: 1},
     "one pixel, all dates": {1: 1, 2: 1, 3: 1},
     "small window, all dates": {1: 1, 2: 64, 3: 64}},
)

  array (120, 4, 509, 543), chunks (1, 4, 509, 543)
    one date, whole scene             1 chunks,      1.0x amplification
    one pixel, all dates            120 chunks, 276,381.0x amplification
    small window, all dates         120 chunks,     67.5x amplification

The amplification factor — bytes fetched divided by bytes used — is the number to design against. Above about 10 the chunking is wrong for that query.

Example 3 — rechunking an existing store

import xarray as xr


def rechunk(src_path, dst_path, chunks, variable=None):
    """Rewrite a store with a different chunk shape."""
    source = xr.open_zarr(src_path)
    names = [variable] if variable else list(source.data_vars)

    for name in names:
        old = source[name].encoding.get("chunks")
        source[name] = source[name].chunk(chunks)
        source[name].encoding = {}
        source[name].encoding["chunks"] = tuple(
            chunks.get(d, source[name].sizes[d]) for d in source[name].dims)
        print(f"  {name}: {old} -> {source[name].encoding['chunks']}")

    source.to_zarr(dst_path, mode="w")
    print(f"  written to {dst_path}")
    return dst_path

Rechunking reads everything and writes everything, so it is expensive — but it is a one-off, and for a store queried in a pattern its chunking does not suit, it pays back immediately.

For very large stores, a dedicated rechunking tool that stages through an intermediate is more memory-efficient than a straight read-write.

Explanation

Why chunking is the whole design

In a chunked format there is no index scan, no query planner and no server. The client computes which chunks a query touches from the array's shape and chunk shape, and requests exactly those objects.

That makes performance completely predictable — and completely determined by the chunk shape. There is no runtime optimisation to compensate for a bad choice.

The amplification factor above makes the stakes concrete: 276,381× for a pixel time series against a time-chunked cube. The client fetches 4.4 MB to return four bytes, 120 times over.

Why Zarr and not NetCDF or HDF5

All three store chunked n-dimensional arrays with metadata. The difference is the storage model.

NetCDF and HDF5 put everything in one file with an internal index, which suits a filesystem and needs a library that can seek within it. Reading them from object storage means either downloading the file or a complicated protocol over range requests.

Zarr puts each chunk in a separate object and the metadata in separate small documents. On object storage that is exactly right: every chunk is a key, and fetching one is a plain GET.

The trade is a very large number of objects, which some filesystems handle badly. Zarr v3's sharding addresses that by grouping chunks into larger objects with an internal index.

Why chunk size has a sweet spot

Too small: each chunk carries per-object overhead — an HTTP request, a metadata lookup, a compression header. At 4 kB chunks, the overhead exceeds the data.

Too large: every read fetches the whole chunk regardless of how little you wanted. At 500 MB chunks, reading a single pixel costs 500 MB.

A few megabytes balances the two: large enough that per-object overhead is negligible, small enough that a targeted read is targeted. That is also roughly the size at which HTTP transfer becomes bandwidth-bound rather than latency-bound.

Why to align dask chunks with storage chunks

Dask chunks describe the computation; Zarr chunks describe the storage. When they differ, each dask chunk reads parts of several stored chunks, and each stored chunk is read by several dask tasks.

Aligning them — or making the dask chunks whole multiples of the storage chunks — means each stored chunk is fetched exactly once. Misaligned by one, and a chunk can be fetched twice for no benefit.

Amplification is a property of the layout and the query. It is the number to design against.

Edge cases or notes

The chunk shape decides the performance, and there is no runtime compensation.
Aim for a few megabytes per chunk uncompressed.
A million tiny chunks is slow everywhere and expensive on object storage.
Align dask chunks with storage chunks, or read each one twice.
Rechunking rewrites everything. Choose well the first time.
Zarr v3 sharding groups chunks into larger objects with an internal index.
Metadata is separate and tiny, which is why opening a huge store is instant.
Consolidated metadata (zarr.consolidate_metadata) turns many small metadata reads into one.

Internal links

Cloud-native geospatial explained — the pattern Zarr belongs to
Lazy loading explained: why xarray reads nothing until you ask — how a cube stays cheap to open
How to turn a STAC search into an xarray data cube — building one from imagery
How to choose chunk and tile sizes that actually help — sizing chunks in practice
Cloud-optimised GeoTIFF explained — the two-dimensional equivalent
How to scale a GeoPandas job with dask-geopandas — the same chunking idea for tables
How to process a very large GeoPackage in chunks with Python — chunking without a chunked format
How to build a cloud-free composite from many satellite scenes — a natural cube operation

FAQ

What is Zarr?

A format for chunked, compressed, n-dimensional arrays where each chunk is a separate object and the metadata is a small JSON document — designed for object storage.

How do I choose a chunk shape?

From the query you will run most. Time-chunked cubes make maps cheap and pixel histories expensive; spatially chunked cubes do the reverse.

How big should a chunk be?

A few megabytes uncompressed. Much smaller and per-object overhead dominates; much larger and every read pulls data you discard.

Why is my pixel time series so slow?

Because the cube is chunked by time. Reading one pixel across 120 dates touches 120 chunks, each of which is a whole scene.

Should I use Zarr or NetCDF?

Zarr for object storage, because each chunk is an independently fetchable object. NetCDF and HDF5 for local filesystems, where their single-file model is an advantage.

Do I need dask to use Zarr?

No — xarray reads Zarr without it. Use dask when the computation is larger than memory, and align the dask chunks with the storage chunks.

Can I change the chunking later?

Yes, by rewriting the store. It reads and writes everything, so it is expensive but one-off.

Chunked Arrays and Zarr Explained #

Problem statement #

Quick answer #

Step-by-step solution #

1. Understand what a chunk is #

2. Choose the chunk shape from the query, not the data #

3. Size the chunks sensibly #

4. Understand that the metadata is tiny and separate #

5. Combine with dask for larger-than-memory work #

Code examples #

Example 1 — writing a cube with a deliberate chunking #

Example 2 — measuring which access pattern your chunking suits #

Example 3 — rechunking an existing store #

Explanation #

Why chunking is the whole design #

Why Zarr and not NetCDF or HDF5 #

Why chunk size has a sweet spot #

Why to align dask chunks with storage chunks #

Edge cases or notes #

Internal links #

FAQ #

What is Zarr? #

How do I choose a chunk shape? #

How big should a chunk be? #

Why is my pixel time series so slow? #

Should I use Zarr or NetCDF? #

Do I need dask to use Zarr? #

Can I change the chunking later? #