Skip to content

xr.concat of virtual datasets silently drops mismatched CF encoding (scale_factor, add_offset, _FillValue) #1004

Description

@TomNicholas

Summary

When two virtual datasets are produced from source files whose data variables share the same dtype/codecs/chunks but differ in their CF decoding attributes (scale_factor, add_offset, _FillValue, missing_value), xr.concat (and the corresponding np.concatenate on ManifestArray) succeeds silently and keeps only the first array's attrs. Once the concatenated dataset is written and later read back with xr.open_zarr(..., decode_cf=True), the surviving attrs are applied to every chunk — including chunks that were packed with a different scale/offset — silently corrupting decoded values.

Reproducer

import h5py, numpy as np, xarray as xr
from virtualizarr import open_virtual_dataset
from virtualizarr.parsers import HDFParser
from obspec_utils.registry import ObjectStoreRegistry
from obstore.store import LocalStore

def write(p, scale, offset, x_start, n=4):
    with h5py.File(p, "w") as f:
        d = f.create_dataset("foo",
            data=np.arange(x_start, x_start+n, dtype="i2").reshape(n,1),
            chunks=(n,1))
        d.attrs["scale_factor"] = np.float64(scale)
        d.attrs["add_offset"]  = np.float64(offset)
        # ... attach dim scales x, y

write("a.nc", scale=0.1,  offset=0.0,   x_start=0)
write("b.nc", scale=0.01, offset=100.0, x_start=4)

registry = ObjectStoreRegistry({"file://": LocalStore()})
parser = HDFParser()
vds1 = open_virtual_dataset(url="file://.../a.nc", parser=parser, registry=registry)
vds2 = open_virtual_dataset(url="file://.../b.nc", parser=parser, registry=registry)

combined = xr.concat([vds1, vds2], dim="x")   # silently succeeds
# combined["foo"].attrs is {'scale_factor': 0.1, 'add_offset': 0.0, ...}
# vds2's scale_factor=0.01 / add_offset=100.0 have been silently dropped.
# Writing this dataset + re-reading with decode_cf=True will mis-decode
# every value originally from b.nc.

Why it happens

ManifestArray.__array_function__ dispatches np.concatenate to virtualizarr.manifests.array_api.concatenate, which calls check_combinable_zarr_arrays (manifests/utils.py). That check validates dtype, codecs, and chunk shape — but not per-array attributes. The result's metadata is then built via copy_and_replace_metadata(first_arr.metadata, new_shape=...), which inherits the first array's attrs.

xarray's combine_attrs setting doesn't save users either: the underlying numpy dispatch happens before xarray's variable-level attr merge, and combine_attrs="override" (the default) also silently picks one.

Affected operations

Any call that funnels through np.concatenate on virtual datasets:

  • xr.concat([vds1, vds2], dim=...)
  • xr.combine_by_coords([...]) / xr.combine_nested([...])
  • virtualizarr.open_virtual_mfdataset(...) style multi-file aggregations

vds.vz.to_icechunk(..., append_dim=...) is not affected — the icechunk-append path calls check_compatible_encodings separately.

Proposed fix

Add an attribute-equality check to check_combinable_zarr_arrays for the four CF decoding keys (scale_factor, add_offset, _FillValue, missing_value), raising ValueError on any mismatch (presence/absence or different values).

This requires the CF attrs to be inspectable at the ManifestArray layer when np.concatenate dispatches. Currently ManifestArray.to_virtual_variable strips all attrs off the array and places them on the wrapping xr.Variable, so the array sees {}. The fix routes only the four CF decoding attrs onto the inner ManifestArray's metadata.attributes (where the check can see them) and leaves arbitrary attrs on xr.Variable.attrs. Writers gain a small helper (extract_cf_encoding_attrs(var)) to keep round-trip behavior identical.

Related

Kerchunk's MultiZarrToZarr.second_pass has the same bug pattern — see combine.py line ~593: "other attributes copied as-is from first occurrence of this array". Worth flagging upstream too.

Metadata

Metadata

Assignees

No one assigned

    Type

    No type

    Fields

    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions