issue_comments: 407404131

Here is an example of how it would look on a dataset: <xarray.Dataset> Dimensions: (x: 65, y: 31) Coordinates: * x (x) float64 ... * y (y) float64 ... time datetime64[ns] ... crs int64 ... Data variables: ndvi (y, x) float64 ... Attributes: Here is how the crs or spatial_ref coodinate variable would look: <xarray.DataArray 'crs' ()> array(0) Coordinates: time datetime64[ns] ... crs int64 0 Attributes: spatial_ref: PROJCS["UTM Zone 15, Northern Hemisphere",GEOGCS["WGS 84",D... And here is how it would look on the variables: <xarray.DataArray 'ndvi' (y: 31, x: 65)> array([[ ...]]) Coordinates: * x (x) float64 ... * y (y) float64 ... time datetime64[ns] ... crs int64 0 Attributes: grid_mapping: crs

@djhoese Whether or not we use the CF convention is not what I am concerned about. What I think would benefit the most people is with the file format to be able to do to_netcdf() with the file and be able to have it read in with standard GIS tools such as rasterio, GDAL, and QGIS. With this schema, this is possible.

Another benefit is that it keeps the crs or spatial_ref with it when you do your operations as it is a coordinate of the variable.

Also, as a side note if you use the center pixel coordinates, then GDAL, rasterio, and QGIS are able to read in the file and determine it's affine/transform without a problem.

For the new library, if you have a crs method attached to it, it isn't too difficult to convert it to whatever format you need it to be in. With the rasterio.crs.CRS.from_string you can "Make a CRS from an EPSG, PROJ, or WKT string" and you can also get back the EPSG, PROJ, or WKT string with a simple method call.

For example, using the recommended method to extend xarray, you could add a crs property:

python from rasterio.crs import CRS ........ @property def crs(self): """:obj:`rasterio.crs.CRS`: Projection from `xarray.DataArray` """ if self._crs is not None: return self._crs try: # look in grid_mapping self._crs = CRS.from_string(self._obj.coords[self._obj.grid_mapping].spatial_ref) except AttributeError: raise ValueError("Spatial reference not found.") return self._crs

And if you call your extension geo, all you would need to get the CRS would be: python ds.geo.crs To get proj.4 string: python ds.geo.crs.to_string() To get WKT string: python ds.geo.crs.wkt To get EPSG code: python ds.geo.crs.to_epsg()