Repository navigation
MultiZarrToZarr uses first file's x/y coordinates for all inputs, causing spatial misalignment in concatenated output #584
Description
When using MultiZarrToZarr to combine multiple (time, y, x) NetCDF files along the time dimension, the resulting Zarr store inherits the x and y coordinate arrays exclusively from the first input file. If subsequent input files have different (but overlapping) spatial extents — e.g., different bounding boxes at the same grid resolution — their data seems to be written at the original chunk offsets without remapping to the correct coordinate positions. This causes spatial displacement of data from all files except the first.
All 3D variables in the input NetCDF files have identical (1, 512, 512) chunking along the (time, y, x) dimensions.
Is there a way to merge the data and consolidate x and y extents across such input files?
Minimal reproducible example:
import fsspec
import s3fs
import os
from kerchunk.hdf import SingleHdf5ToZarr
from kerchunk.combine import MultiZarrToZarr
import xarray as xr
s3 = s3fs.S3FileSystem()
bucket = 's3://its-live-data'
granules = [
"velocity_image_pair/landsatOLI/v02/S80W170/LC08_L1GT_020121_20231013_20231102_02_T2_X_LC09_L1GT_020121_20231106_20231106_02_T2_G0120V02_P084.nc",
"velocity_image_pair/landsatOLI/v02/S80W170/LC08_L1GT_020120_20201121_20210315_02_T2_X_LC08_L1GT_020120_20210124_20210305_02_T2_G0120V02_P051.nc",
]
so = dict(
anon=True, default_fill_cache=False, default_cache_type='first'
)
# Read input granules in
singles = []
for each_granule in granules:
s3_path = os.path.join(bucket, each_granule)
print(f'Reading {s3_path}...')
with fsspec.open(s3_path, **so) as inf:
h5chunks = SingleHdf5ToZarr(inf, s3_path, inline_threshold=300)
singles.append(h5chunks.translate())
# Combine granules
mzz = MultiZarrToZarr(
singles,
remote_protocol="s3",
remote_options={'anon': True},
concat_dims=["time"],
identical_dims=['mapping']
)
out = mzz.translate()
# Open merged cube as xarray.Dataset for exploration (plotting)
fs = fsspec.filesystem(
"reference",
fo=out, # or pass mzz_ref dict directly
remote_protocol="s3",
remote_options={
"anon": True,
"asynchronous": False, # must match ReferenceFileSystem's mode
},
skip_instance_cache=True,
)
mapper = fs.get_mapper("")
ds = xr.open_zarr(mapper, consolidated=False, decode_timedelta=False, chunks={})
This is (time, y, x) v data variable from original input data (please note the x extents):
And the same input data as it appears in the virtual cube with x extent of original first input file:
Examining x and y chunking shows reference to the first input file:
for key in out['refs'].keys():
if ( key.startswith('x/') or key.startswith('y/') ) and "/." not in key:
print(f'{key=}: {out['refs'][key]}')
Activity
martindurant commentedon Feb 27, 2026
I'm afraid this is by design - MZZ assumes that your data already conforms to zarr's strictly rectilinear coordinate system.
I believe that virtualizarr, which uses xarray's own coordinate matching system, may be able to do what you are after, but I am not certain (cc @TomNicholas , please link if such a workflow/example already exists).
mliukis commentedon Feb 27, 2026
@martindurant Thank you for your response. That's unfortunate...
I already tried
virtualizarr, filed github issue, and got response fromvirtualizarrthat it's not supported.martindurant commentedon Feb 27, 2026
Sorry ...
I guess reprojecting on-the-fly is too complex. If it's more a case of slice (without reproduction), I can imagine it should not be too complex really - the equivalent of a codec - but there is nowhere to store per-chunk metadata to make it happen.
TomNicholas commentedon Feb 27, 2026
Unless these bounding boxes happen to be perfectly aligned with internal chunking in the netCDF4 files, then I'm afraid what you're trying to do almost certainly can never be done virtually, with Kerchunk or with VirtualiZarr. Your final Zarr store has to have chunks placed in a regular grid. I think what you're trying to do would require chunks being placed at arbitrary locations, not necessarily forming a rectilinear grid.
mliukis commentedon Feb 27, 2026
@TomNicholas All the chunks that exist are grid aligned and the same size, but the resulting layer for each granule will be sparse (only have data chunks for a portion of the whole data cube).
mliukis commentedon Feb 27, 2026
I guess, we could re-adjust all the chunking offsets per our common x, y grid for all layers...
TomNicholas commentedon Feb 27, 2026
In that case it's an advanced use case, but should be possible...
Maybe you could try creating a skeletal zarr store then using the upcoming
regionkwarg (zarr-developers/VirtualiZarr#873) to write in each tile?