I am a little discouraged that we have not managed to align better across projects so far (e.g. having this conversation before the GDAL Zarr CRS convention was implemented)

just seeing this conversation. I wasn't aware of GeoZarr when I implemented the _CRS attribute in the GDAL Zarr driver. There is a practical difficulty in reusing the CF reading & writing code from the GDAL netCDF driver as it is tied to that driver, but I'd guess with sufficient effort it could be made agnostic of the carrier. No opposition from me if someone wants to tackle this. That said, the way CRS is encoding in CF conventions is not super elegant, because you need to create a 0-dim variable, define a set of attributes for each projection method, and have an attribute in the georeferenced variable to point to the variable that holds the CRS attribute. A _CRS attribute on the georeferenced variable is much more simple, and by reusing the WKT or PROJJSON encodings, this save the implementor the business of knowing how to define exactly the CRS (obviously when they have a library that does that job for them)

I am guilty of sidetracking this issue into the politics of CRS encoding. That discussion is important. But in the meantime, @wankoelias's original issue reveals is narrower technical issue with Xarray's Zarr writer: Xarray won't let you serialize a dictionary attribute to zarr, even though zarr has no problem with this. That is a problem we can fix pretty easily.

The _validate_attrs helper function was just borrowed from to_netcdf:

https://github.com/pydata/xarray/blob/586992e8d2998751cb97b1cab4d3caa9dca116e0/xarray/backends/api.py#L133-L135

We could refactor this function to be more flexible to account for zarr's broader range of allowed attribute types (as we have evidently already done for h5netcdf). Or we could just bypass it completely in the to_zarr method. That is the only real decision we need to make here right now.

@wankoelias - you seem to understand the issue pretty well. Would you be game for making a PR? We would be glad to support you along the way.

GDAL doesn't actually need the CRS stored with the _CRS convention as documented in their ZARR driver specification?

GDAL may need to update the ZARR driver. The NetCDF driver does support CF. It would be good for them to be consistent.

so just to make it clear... GDAL doesn't actually need the CRS stored with the _CRS convention as documented in their ZARR driver specification?

Top reasons for using CF over a CRS attribute:

• Supported by geospatial & netCDF software (GDAL/Esri/Unidata)
• Attributes are easily lost ref
• Encourage more usage of standard mechanisms of CRS storage

This document is also a useful reference for storing CRS in xarray: https://corteva.github.io/rioxarray/stable/getting_started/crs_management.html

Instead, it recommends to use CF conventions for encoding CRS. This is more compatible with NetCDF, but won't be parsed correctly by GDAL.

GDAL does support the CF conventions for storing the CRS.

https://gdal.org/drivers/raster/netcdf.html#georeference

"The driver first tries to follow the CF-1 Convention from UNIDATA looking for the Metadata named “grid_mapping”."

One of the main motivations behind the the rioxarray extension is GDAL compatibility. It looks like @snowman2 and @TomAugspurger have discussed saving many geotiffs loaded into xarray as GDAL-compatible Zarr for example https://github.com/corteva/rioxarray/issues/433#issuecomment-967685356.

While it seems that the ultimate solution is agreeing on a format standard, here is another small example using the rioxarray extension where format conversion doesn't currently work as you might expect:

```python

https://github.com/pydata/xarray-data

ds = xr.open_dataset('xarray-data/air_temperature.nc', engine='rasterio')

TooManyDimensions: Only 2D and 3D data arrays supported.

ds.rio.to_raster('test.zarr', driver='ZARR')

Does not error, but output not equivalent to gdal_translate -of ZARR xarray-data/air_temperature.nc gdal_air_temp.zarr

for example, gdalinfo xarray-tutorial-airtemp.zarr gives

Warning 1: Too many samples along the > 2D dimensions of /air.

ds.to_zarr('xarray-tutorial-airtemp.zarr')

```

Some people will prefer Cloud-Optimised GeoTiff, others (geo)Zarr.

By the way, I use rasterio (based on GDAL), xarray, GDAL and all tools without problems with CF conventions.

The impact of GDAL using CRS attribute are only: - when converting a non-netcdf file to Zarr with GDAL, the CRS is encoded in the CRS attribute - when converting a Zarr file to a non-netcdf file which encode CRS in grid-mapping, this is not converted to the CRS attribute.

@christophenoel - I share your perspective. But there is a huge swath of the geospatial world who basically hate NetCDF and avoid it like the plague. These communities prefer to use geotiff and GDAL. We need to reach for interoperability.

Very interesting topic. I assume that GDAL proposed something as Zarr specification has no provision for spatial reference system encoding.

From my point of view, by making it close to NetCDF (which is so widely used in Earth Observation missions) and providing supporting librarires, xArray made the success of Zarr in the EO world.

Indeed, my dream would be GDAL align to the GeoZarr spec which is mostly aligned with NetCDF/xArray.

I think the core problem here is that Zarr itself supports arbitrary json data structures as attributes, but netCDF does not. The Zarr serialization in Xarray is designed to emulate netCDF, but we could make that optional, for example, with a flag to bypass attribute encoding / decoding and just pass the python data directly through to Zarr.

However, my concern would be that netCDF4 C library would not be able to read those files (nczarr). What happens if you try to open up a GDAL-created Zarr with netCDF4?

FWIW, the new GeoZarr Spec by @christophenoel does not use the GDAL convention for CRS. Instead, it recommends to use CF conventions for encoding CRS. This is more compatible with NetCDF, but won't be parsed correctly by GDAL.

I am a little discouraged that we have not managed to align better across projects so far (e.g. having this conversation before the GDAL Zarr CRS convention was implemented). 😞 For example, either of these two GDAL PRs: - https://github.com/OSGeo/gdal/pull/3896 - https://github.com/OSGeo/gdal/pull/4521

However, it is not too late! Let's try to reach for a standard way of encoding CRS in Zarr that can be used across languages and implementations of Zarr.

My own preference would be to try to get GDAL to support the GeoZarr Spec and thus the CF-convention CRS attribute, rather than trying to get Xarray to be able to write the GDAL CRS convention.