Repository: zarr-developers/VirtualiZarr

Issue:
MetadataError from ValueError: Could not convert object to NumPy datetime
I'm trying to debug @thodson-usgs's example from https://github.com/cubed-dev/cubed/pull/520 (and originally https://github.com/zarr-developers/VirtualiZarr/pull/197).

He is doing a whole serverless reduction of virtual references to multiple files (!!! - relevant to #123), but there seem to be some more basic errors to be fixed first.

Specifically, if I try to use virtualizarr on just one of his files this happens:

```python
import xarray as xr
from virtualizarr import open_virtual_dataset

vds = open_virtual_dataset(
    's3://wrf-se-ak-ar5/ccsm/rcp85/daily/2060/WRFDS_2060-01-01.nc',
    indexes={},
    loadable_variables=['Time'],
    cftime_variables=['Time'],
)
vds
```
```
<xarray.Dataset> Size: 31MB
Dimensions:        (Time: 1, south_north: 250, west_east: 320,
                    interp_levels: 9, soil_layers_stag: 4)
Coordinates:
    interp_levels  (interp_levels) float32 36B ManifestArray<shape=(9,), dtyp...
    Time           (Time) datetime64[ns] 8B 2060-01-01
Dimensions without coordinates: south_north, west_east, soil_layers_stag
Data variables: (12/39)
    SNOWH          (Time, south_north, west_east) float32 320kB ManifestArray...
    ACSNOW         (Time, south_north, west_east) float32 320kB ManifestArray...
    TSK            (Time, south_north, west_east) float32 320kB ManifestArray...
    XLONG          (south_north, west_east) float32 320kB ManifestArray<shape...
    T              (Time, interp_levels, south_north, west_east) float32 3MB ...
    XLAT           (south_north, west_east) float32 320kB ManifestArray<shape...
    ...             ...
    PSFC           (Time, south_north, west_east) float32 320kB ManifestArray...
    ALBEDO         (Time, south_north, west_east) float32 320kB ManifestArray...
    CLDFRA         (Time, interp_levels, south_north, west_east) float32 3MB ...
    SWDNB          (Time, south_north, west_east) float32 320kB ManifestArray...
    PW             (Time, south_north, west_east) float32 320kB ManifestArray...
    SH2O           (Time, soil_layers_stag, south_north, west_east) float32 1MB ManifestArray<shape=(1, 4, 250, 320), dtype=float32, chunks=(1, 4, 250, 32...
Attributes:
    contact:  rtladerjr@alaska.edu
    data:     Downscaled CCSM4
    date:     Mon Oct 21 11:37:23 AKDT 2019
    format:   version 2
    info:     Alaska CASC
```
```python
ds = xr.open_dataset('combined.json', engine="kerchunk")
```
```python
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/meta.py:127, in Metadata2.decode_array_metadata(cls, s)
    126 dimension_separator = meta.get("dimension_separator", None)
--> 127 fill_value = cls.decode_fill_value(meta["fill_value"], dtype, object_codec)
    128 meta = dict(
    129     zarr_format=meta["zarr_format"],
    130     shape=tuple(meta["shape"]),
   (...)
    136     filters=meta["filters"],
    137 )

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/meta.py:260, in Metadata2.decode_fill_value(cls, v, dtype, object_codec)
    259 else:
--> 260     return np.array(v, dtype=dtype)[()]

ValueError: Could not convert object to NumPy datetime

The above exception was the direct cause of the following exception:

MetadataError                             Traceback (most recent call last)
Cell In[8], line 1
----> 1 ds = xr.open_dataset('combined.json', engine="kerchunk")

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/backends/api.py:571, in open_dataset(filename_or_obj, engine, chunks, cache, decode_cf, mask_and_scale, decode_times, decode_timedelta, use_cftime, concat_characters, decode_coords, drop_variables, inline_array, chunked_array_type, from_array_kwargs, backend_kwargs, **kwargs)
    559 decoders = _resolve_decoders_kwargs(
    560     decode_cf,
    561     open_backend_dataset_parameters=backend.open_dataset_parameters,
   (...)
    567     decode_coords=decode_coords,
    568 )
    570 overwrite_encoded_chunks = kwargs.pop("overwrite_encoded_chunks", None)
--> 571 backend_ds = backend.open_dataset(
    572     filename_or_obj,
    573     drop_variables=drop_variables,
    574     **decoders,
    575     **kwargs,
    576 )
    577 ds = _dataset_from_backend_dataset(
    578     backend_ds,
    579     filename_or_obj,
   (...)
    589     **kwargs,
    590 )
    591 return ds

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/kerchunk/xarray_backend.py:12, in KerchunkBackend.open_dataset(self, filename_or_obj, storage_options, open_dataset_options, **kw)
      8 def open_dataset(
      9     self, filename_or_obj, *, storage_options=None, open_dataset_options=None, **kw
     10 ):
     11     open_dataset_options = (open_dataset_options or {}) | kw
---> 12     ref_ds = open_reference_dataset(
     13         filename_or_obj,
     14         storage_options=storage_options,
     15         open_dataset_options=open_dataset_options,
     16     )
     17     return ref_ds

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/kerchunk/xarray_backend.py:46, in open_reference_dataset(filename_or_obj, storage_options, open_dataset_options)
     42     open_dataset_options = {}
     44 m = fsspec.get_mapper("reference://", fo=filename_or_obj, **storage_options)
---> 46 return xr.open_dataset(m, engine="zarr", consolidated=False, **open_dataset_options)

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/backends/api.py:571, in open_dataset(filename_or_obj, engine, chunks, cache, decode_cf, mask_and_scale, decode_times, decode_timedelta, use_cftime, concat_characters, decode_coords, drop_variables, inline_array, chunked_array_type, from_array_kwargs, backend_kwargs, **kwargs)
    559 decoders = _resolve_decoders_kwargs(
    560     decode_cf,
    561     open_backend_dataset_parameters=backend.open_dataset_parameters,
   (...)
    567     decode_coords=decode_coords,
    568 )
    570 overwrite_encoded_chunks = kwargs.pop("overwrite_encoded_chunks", None)
--> 571 backend_ds = backend.open_dataset(
    572     filename_or_obj,
    573     drop_variables=drop_variables,
    574     **decoders,
    575     **kwargs,
    576 )
    577 ds = _dataset_from_backend_dataset(
    578     backend_ds,
    579     filename_or_obj,
   (...)
    589     **kwargs,
    590 )
    591 return ds

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/backends/zarr.py:1182, in ZarrBackendEntrypoint.open_dataset(self, filename_or_obj, mask_and_scale, decode_times, concat_characters, decode_coords, drop_variables, use_cftime, decode_timedelta, group, mode, synchronizer, consolidated, chunk_store, storage_options, stacklevel, zarr_version, store, engine)
   1180 store_entrypoint = StoreBackendEntrypoint()
   1181 with close_on_error(store):
-> 1182     ds = store_entrypoint.open_dataset(
   1183         store,
   1184         mask_and_scale=mask_and_scale,
   1185         decode_times=decode_times,
   1186         concat_characters=concat_characters,
   1187         decode_coords=decode_coords,
   1188         drop_variables=drop_variables,
   1189         use_cftime=use_cftime,
   1190         decode_timedelta=decode_timedelta,
   1191     )
   1192 return ds

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/backends/store.py:43, in StoreBackendEntrypoint.open_dataset(self, filename_or_obj, mask_and_scale, decode_times, concat_characters, decode_coords, drop_variables, use_cftime, decode_timedelta)
     29 def open_dataset(  # type: ignore[override]  # allow LSP violation, not supporting **kwargs
     30     self,
     31     filename_or_obj: str | os.PathLike[Any] | BufferedIOBase | AbstractDataStore,
   (...)
     39     decode_timedelta=None,
     40 ) -> Dataset:
     41     assert isinstance(filename_or_obj, AbstractDataStore)
---> 43     vars, attrs = filename_or_obj.load()
     44     encoding = filename_or_obj.get_encoding()
     46     vars, attrs, coord_names = conventions.decode_cf_variables(
     47         vars,
     48         attrs,
   (...)
     55         decode_timedelta=decode_timedelta,
     56     )

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/backends/common.py:221, in AbstractDataStore.load(self)
    199 def load(self):
    200     """
    201     This loads the variables and attributes simultaneously.
    202     A centralized loading function makes it easier to create
   (...)
    218     are requested, so care should be taken to make sure its fast.
    219     """
    220     variables = FrozenDict(
--> 221         (_decode_variable_name(k), v) for k, v in self.get_variables().items()
    222     )
    223     attributes = FrozenDict(self.get_attrs())
    224     return variables, attributes

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/backends/zarr.py:563, in ZarrStore.get_variables(self)
    562 def get_variables(self):
--> 563     return FrozenDict(
    564         (k, self.open_store_variable(k, v)) for k, v in self.zarr_group.arrays()
    565     )

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/core/utils.py:443, in FrozenDict(*args, **kwargs)
    442 def FrozenDict(*args, **kwargs) -> Frozen:
--> 443     return Frozen(dict(*args, **kwargs))

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/xarray/backends/zarr.py:563, in <genexpr>(.0)
    562 def get_variables(self):
--> 563     return FrozenDict(
    564         (k, self.open_store_variable(k, v)) for k, v in self.zarr_group.arrays()
    565     )

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/hierarchy.py:691, in Group._array_iter(self, keys_only, method, recurse)
    689 if contains_array(self._store, path):
    690     _key = key.rstrip("/")
--> 691     yield _key if keys_only else (_key, self[key])
    692 elif recurse and contains_group(self._store, path):
    693     group = self[key]

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/hierarchy.py:467, in Group.__getitem__(self, item)
    465 path = self._item_path(item)
    466 try:
--> 467     return Array(
    468         self._store,
    469         read_only=self._read_only,
    470         path=path,
    471         chunk_store=self._chunk_store,
    472         synchronizer=self._synchronizer,
    473         cache_attrs=self.attrs.cache,
    474         zarr_version=self._version,
    475         meta_array=self._meta_array,
    476     )
    477 except ArrayNotFoundError:
    478     pass

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/core.py:170, in Array.__init__(self, store, path, read_only, chunk_store, synchronizer, cache_metadata, cache_attrs, partial_decompress, write_empty_chunks, zarr_version, meta_array)
    167     self._metadata_key_suffix = self._hierarchy_metadata["metadata_key_suffix"]
    169 # initialize metadata
--> 170 self._load_metadata()
    172 # initialize attributes
    173 akey = _prefix_to_attrs_key(self._store, self._key_prefix)

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/core.py:193, in Array._load_metadata(self)
    191 """(Re)load metadata from store."""
    192 if self._synchronizer is None:
--> 193     self._load_metadata_nosync()
    194 else:
    195     mkey = _prefix_to_array_key(self._store, self._key_prefix)

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/core.py:207, in Array._load_metadata_nosync(self)
    204     raise ArrayNotFoundError(self._path) from e
    205 else:
    206     # decode and store metadata as instance members
--> 207     meta = self._store._metadata_class.decode_array_metadata(meta_bytes)
    208     self._meta = meta
    209     self._shape = meta["shape"]

File ~/miniconda3/envs/numpy2.0_released/lib/python3.11/site-packages/zarr/meta.py:141, in Metadata2.decode_array_metadata(cls, s)
    139         meta["dimension_separator"] = dimension_separator
    140 except Exception as e:
--> 141     raise MetadataError("error decoding metadata") from e
    142 else:
    143     return meta

MetadataError: error decoding metadata
```

At first I assumed there was something wrong with our handling of the loaded `cftime_variables`, but actually even if I drop the `'Time'` variable I still get exactly the same error:

```python
vds = open_virtual_dataset(
    's3://wrf-se-ak-ar5/ccsm/rcp85/daily/2060/WRFDS_2060-01-01.nc',
    indexes={},
    drop_variables=['Time'],
)
```

I don't know why it's even trying to convert anything to a datetime - none of the other variables have units of time.

What's also weird is that this is raised from within `meta.py:260, in Metadata2.decode_fill_value(cls, v, dtype, object_codec)`, which suggests a problem with the `fill_value`. But I checked and all of the variables in this virtual dataset have a fill_value of either a float or `nan` in their `.encoding`, again nothing about a datetime.



Discussion:
thodson-usgs: @TomAugspurger, referencing back to #201. Pulling you over here, because we might be crossing issues.

I would test this with:
```bash
aws s3api get-object --bucket wrf-se-ak-ar5 --key ccsm/rcp85/daily/2060/WRFDS_2060-01-01.nc WRFDS_2060-01-01.nc
aws s3api get-object --bucket wrf-se-ak-ar5 --key ccsm/rcp85/daily/2060/WRFDS_2060-01-01.nc WRFDS_2060-01-02.nc
```

then run this `example.py`
```python
import xarray as xr
import glob
from virtualizarr import open_virtual_dataset

virtual_datasets = [
    open_virtual_dataset(filepath,
                         indexes={},
                         loadable_variables=['Time'],
                         cftime_variables=['Time'])
    for filepath in glob.glob('*.nc')
]

virtual_ds = xr.combine_nested(virtual_datasets, concat_dim=["Time"], coords="minimal", compat="override")
virtual_ds.virtualize.to_kerchunk('combined.json', format='json')

ds = xr.open_dataset('combined.json', engine="kerchunk")

```

This should fail on `main` with
```
  raise MetadataError("error decoding metadata") from e
zarr.errors.MetadataError: error decoding metadata
```
but passes with this PR.
TomAugspurger: Thanks. I had to modify that slightly to work with my env:

```
aws s3api get-object --region us-west-2 --no-sign-request --bucket wrf-se-ak-ar5 --key ccsm/rcp85/daily/2060/WRFDS_2060-01-01.nc WRFDS_2060-01-01.nc
  aws s3api get-object --region us-west-2 --no-sign-request --bucket wrf-se-ak-ar5 --key ccsm/rcp85/daily/2060/WRFDS_2060-01-01.nc WRFDS_2060-01-02.nc
```

Taking a look.
TomAugspurger: What version of Zarr do you have?
thodson-usgs: > What version of Zarr do you have?

`2.18.2`
TomAugspurger: Thanks, I was able to reproduce it.

I would suggest something like this:

```diff
diff --git a/virtualizarr/zarr.py b/virtualizarr/zarr.py
index e5015b3..3647db5 100644
--- a/virtualizarr/zarr.py
+++ b/virtualizarr/zarr.py
@@ -32,13 +32,14 @@ ZAttrs = NewType(
 )  # just the .zattrs (for one array or for the whole store/group)
 FillValueT = bool | str | float | int | list | None
 
-ZARR_DEFAULT_FILL_VALUE: dict[np.dtype, FillValueT] = {
+ZARR_DEFAULT_FILL_VALUE: dict[str, FillValueT] = {
     # numpy dtypes's hierarchy lets us avoid checking for all the widths
     # https://numpy.org/doc/stable/reference/arrays.scalars.html
-    np.dtype("bool"): False,
-    np.dtype("int"): 0,
-    np.dtype("float"): 0.0,
-    np.dtype("complex"): [0.0, 0.0],
+    np.dtype("bool").kind: False,
+    np.dtype("int").kind: 0,
+    np.dtype("float").kind: 0.0,
+    np.dtype("complex").kind: [0.0, 0.0],
+    np.dtype("datetime64").kind: np.datetime64("NaT").view("i8").item(),
 }
 """
 The value and format of the fill_value depend on the `data_type` of the array.
@@ -67,7 +68,7 @@ class ZArray(BaseModel):
     chunks: tuple[int, ...]
     compressor: dict | None = None
     dtype: np.dtype
-    fill_value: FillValueT = Field(default=0.0, validate_default=True)
+    fill_value: FillValueT = Field(None, validate_default=True)
     filters: list[dict] | None = None
     order: Literal["C", "F"]
     shape: tuple[int, ...]
@@ -90,7 +91,7 @@ class ZArray(BaseModel):
     @model_validator(mode="after")
     def _check_fill_value(self) -> Self:
         if self.fill_value is None:
-            self.fill_value = ZARR_DEFAULT_FILL_VALUE.get(self.dtype, 0.0)
+            self.fill_value = ZARR_DEFAULT_FILL_VALUE.get(self.dtype.kind, 0.0)
         return self
 
     @property
```

That adds a default fill value for datetimes that ends up as `NaT`. Unfortunately, something strange was going on with the lookup of that dtype, so I had to change that mapping to use `.dtype.kind` instead. Dunno if that'll cause other issues.


I think that once virtualizarr is able to use `zarr.Array` internally then things will be smoother.
thodson-usgs: Great suggestion, @TomAugspurger!
I revised the PR, and it passed my example, as well as the test suite. However, I have no sense of who might be affected by this change and whether it's worth the risk.
TomNicholas: Thanks for rooting this out guys! Let's a regression test to be sure the fix is working.

Patch:
diff --git a/docs/releases.rst b/docs/releases.rst
index f472db0..e7f6df2 100644
--- a/docs/releases.rst
+++ b/docs/releases.rst
@@ -9,6 +9,9 @@ v1.0.1 (unreleased)
 New Features
 ~~~~~~~~~~~~
 
+- Load scalar variables by default. (:pull:`205`)
+  By `Gustavo Hidalgo <https://github.com/ghidalgo3>`_.
+
 Breaking changes
 ~~~~~~~~~~~~~~~~
 
@@ -23,6 +26,9 @@ Bug fixes
 
 - Exclude empty chunks during `ChunkDict` construction. (:pull:`198`)
   By `Gustavo Hidalgo <https://github.com/ghidalgo3>`_.
+- Fixed regression in `fill_value` handling for datetime dtypes making virtual
+  Zarr stores unreadable (:pr:`206`)
+  By `Timothy Hodson <https://github.com/thodson-usgs>`_
 
 Documentation
 ~~~~~~~~~~~~~
diff --git a/pyproject.toml b/pyproject.toml
index 9fe0468..4496116 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -45,6 +45,7 @@ test = [
     "fsspec",
     "s3fs",
     "fastparquet",
+    "h5py"
 ]
 
 
diff --git a/virtualizarr/xarray.py b/virtualizarr/xarray.py
index 3a6dd02..528c252 100644
--- a/virtualizarr/xarray.py
+++ b/virtualizarr/xarray.py
@@ -333,13 +333,16 @@ def variable_from_kerchunk_refs(
 
     arr_refs = kerchunk.extract_array_refs(refs, var_name)
     chunk_dict, zarray, zattrs = kerchunk.parse_array_refs(arr_refs)
-
-    manifest = ChunkManifest._from_kerchunk_chunk_dict(chunk_dict)
-
     # we want to remove the _ARRAY_DIMENSIONS from the final variables' .attrs
     dims = zattrs.pop("_ARRAY_DIMENSIONS")
-
-    varr = virtual_array_class(zarray=zarray, chunkmanifest=manifest)
+    if chunk_dict:
+        manifest = ChunkManifest._from_kerchunk_chunk_dict(chunk_dict)
+        varr = virtual_array_class(zarray=zarray, chunkmanifest=manifest)
+    else:
+        # This means we encountered a scalar variable of dimension 0,
+        # very likely that it actually has no numeric value and its only purpose
+        # is to communicate dataset attributes.
+        varr = zarray.fill_value
 
     return xr.Variable(data=varr, dims=dims, attrs=zattrs)
 
diff --git a/virtualizarr/zarr.py b/virtualizarr/zarr.py
index e5015b3..43735b1 100644
--- a/virtualizarr/zarr.py
+++ b/virtualizarr/zarr.py
@@ -32,13 +32,14 @@ ZAttrs = NewType(
 )  # just the .zattrs (for one array or for the whole store/group)
 FillValueT = bool | str | float | int | list | None
 
-ZARR_DEFAULT_FILL_VALUE: dict[np.dtype, FillValueT] = {
+ZARR_DEFAULT_FILL_VALUE: dict[str, FillValueT] = {
     # numpy dtypes's hierarchy lets us avoid checking for all the widths
     # https://numpy.org/doc/stable/reference/arrays.scalars.html
-    np.dtype("bool"): False,
-    np.dtype("int"): 0,
-    np.dtype("float"): 0.0,
-    np.dtype("complex"): [0.0, 0.0],
+    np.dtype("bool").kind: False,
+    np.dtype("int").kind: 0,
+    np.dtype("float").kind: 0.0,
+    np.dtype("complex").kind: [0.0, 0.0],
+    np.dtype("datetime64").kind: 0,
 }
 """
 The value and format of the fill_value depend on the `data_type` of the array.
@@ -67,7 +68,7 @@ class ZArray(BaseModel):
     chunks: tuple[int, ...]
     compressor: dict | None = None
     dtype: np.dtype
-    fill_value: FillValueT = Field(default=0.0, validate_default=True)
+    fill_value: FillValueT = Field(None, validate_default=True)
     filters: list[dict] | None = None
     order: Literal["C", "F"]
     shape: tuple[int, ...]
@@ -90,7 +91,7 @@ class ZArray(BaseModel):
     @model_validator(mode="after")
     def _check_fill_value(self) -> Self:
         if self.fill_value is None:
-            self.fill_value = ZARR_DEFAULT_FILL_VALUE.get(self.dtype, 0.0)
+            self.fill_value = ZARR_DEFAULT_FILL_VALUE.get(self.dtype.kind, 0.0)
         return self
 
     @property


Test patch:
diff --git a/conftest.py b/conftest.py
index 8d5e351..b558abf 100644
--- a/conftest.py
+++ b/conftest.py
@@ -1,3 +1,4 @@
+import h5py
 import pytest
 import xarray as xr
 
@@ -50,3 +51,21 @@ def netcdf4_files(tmpdir):
     ds2.close()
 
     return filepath1, filepath2
+
+
+@pytest.fixture
+def hdf5_empty(tmpdir):
+    filepath = f"{tmpdir}/empty.nc"
+    f = h5py.File(filepath, "w")
+    dataset = f.create_dataset("empty", shape=(), dtype="float32")
+    dataset.attrs["empty"] = "true"
+    return filepath
+
+
+@pytest.fixture
+def hdf5_scalar(tmpdir):
+    filepath = f"{tmpdir}/scalar.nc"
+    f = h5py.File(filepath, "w")
+    dataset = f.create_dataset("scalar", data=0.1, dtype="float32")
+    dataset.attrs["scalar"] = "true"
+    return filepath
diff --git a/virtualizarr/tests/test_integration.py b/virtualizarr/tests/test_integration.py
index 239316a..7210b3f 100644
--- a/virtualizarr/tests/test_integration.py
+++ b/virtualizarr/tests/test_integration.py
@@ -4,6 +4,9 @@ import xarray as xr
 import xarray.testing as xrt
 
 from virtualizarr import open_virtual_dataset
+from virtualizarr.manifests.array import ManifestArray
+from virtualizarr.manifests.manifest import ChunkManifest
+from virtualizarr.zarr import ZArray
 
 
 @pytest.mark.parametrize(
@@ -166,6 +169,50 @@ class TestKerchunkRoundtrip:
         # assert equal to original dataset
         xrt.assert_identical(roundtrip, ds)
 
+    def test_datetime64_dtype_fill_value(self, tmpdir, format):
+        chunks_dict = {
+            "0.0.0": {"path": "foo.nc", "offset": 100, "length": 100},
+        }
+        manifest = ChunkManifest(entries=chunks_dict)
+        chunks = (1, 1, 1)
+        shape = (1, 1, 1)
+        zarray = ZArray(
+            chunks=chunks,
+            compressor={"id": "zlib", "level": 1},
+            dtype=np.dtype("<M8[ns]"),
+            # fill_value=0.0,
+            filters=None,
+            order="C",
+            shape=shape,
+            zarr_format=2,
+        )
+        marr1 = ManifestArray(zarray=zarray, chunkmanifest=manifest)
+        ds = xr.Dataset(
+            {
+                "a": xr.DataArray(
+                    marr1,
+                    attrs={
+                        "_FillValue": np.datetime64("1970-01-01T00:00:00.000000000")
+                    },
+                )
+            }
+        )
+
+        if format == "dict":
+            # write those references to an in-memory kerchunk-formatted references dictionary
+            ds_refs = ds.virtualize.to_kerchunk(format=format)
+
+            # use fsspec to read the dataset from the kerchunk references dict
+            roundtrip = xr.open_dataset(ds_refs, engine="kerchunk")
+        else:
+            # write those references to disk as kerchunk references format
+            ds.virtualize.to_kerchunk(f"{tmpdir}/refs.{format}", format=format)
+
+            # use fsspec to read the dataset from disk via the kerchunk references
+            roundtrip = xr.open_dataset(f"{tmpdir}/refs.{format}", engine="kerchunk")
+
+        assert roundtrip.a.attrs == ds.a.attrs
+
 
 def test_open_scalar_variable(tmpdir):
     # regression test for GH issue #100
diff --git a/virtualizarr/tests/test_xarray.py b/virtualizarr/tests/test_xarray.py
index 31ef832..00140a1 100644
--- a/virtualizarr/tests/test_xarray.py
+++ b/virtualizarr/tests/test_xarray.py
@@ -419,6 +419,16 @@ class TestLoadVirtualDataset:
         }
         mock_read_kerchunk.assert_called_once_with(**args)
 
+    def test_open_dataset_with_empty(self, hdf5_empty, tmpdir):
+        vds = open_virtual_dataset(hdf5_empty)
+        assert vds.empty.dims == ()
+        assert vds.empty.attrs == {"empty": "true"}
+
+    def test_open_dataset_with_scalar(self, hdf5_scalar, tmpdir):
+        vds = open_virtual_dataset(hdf5_scalar)
+        assert vds.scalar.dims == ()
+        assert vds.scalar.attrs == {"scalar": "true"}
+
 
 class TestRenamePaths:
     def test_rename_to_str(self, netcdf4_file):