A DataArray within a Dataset has different CRS

Hello!

I'm accessing the Dataset MODIS/061/MCD64A1 (MODIS Burned Area Monthly Global 500m) and am noticing the Dataset and DataArray report different projections. From what I understand, the Dataset appears to be EPSG:4326, but the DataArray contained within it is reported as SR-ORG:6974. This seems strange as the default appears to be EPSG:4326.

The output seems to be getting closer to usable after the commit Add 'm' as another synonym for 'meter'

Is this intentional behaviour or am I misunderstanding or doing something silly here?

Here's a minimum example:

import xarray as xr
import rioxarray
import ee

ee.Initialize(project="FILL ME")

ic = ee.ImageCollection('MODIS/061/MCD64A1') \
    .filterDate('2019-09-01', '2020-03-01') \
    .select('BurnDate')

ds = xr.open_dataset(ic, engine='ee', scale=0.25)

print(f"Dataset:\n{ds}")

print(f"Dataset CRS: {ds.crs}")
print(f"DataArray CRS: {ds['BurnDate'].crs}")

Standard out:

Dataset:
<xarray.Dataset>
Dimensions:   (time: 6, lon: 1440, lat: 720)
Coordinates:
  * time      (time) datetime64[ns] 2019-09-01 2019-10-01 ... 2020-02-01
  * lon       (lon) float32 -179.9 -179.6 -179.4 -179.1 ... 179.4 179.6 179.9
  * lat       (lat) float32 89.88 89.62 89.38 89.12 ... -89.38 -89.62 -89.88
Data variables:
    BurnDate  (time, lon, lat) int32 ...
Attributes:
    crs:      EPSG:4326
Dataset CRS: EPSG:4326
DataArray CRS: SR-ORG:6974

When written using ds.to_netcdf(..), and examining it with gdalinfo:

Driver: netCDF/Network Common Data Format
Files: via_xarray.nc
Size is 720, 1440
Metadata:
  BurnDate#bounds={-180,-90,180,90}
  BurnDate#crs=SR-ORG:6974
  BurnDate#crs_transform={463.3127165279165,0,-20015109.354,0,-463.3127165279167,7783653.637664001}
  BurnDate#data_type={'type': 'PixelType', 'precision': 'int', 'min': -32768, 'max': 32767}
  BurnDate#dimensions={86400,31200}
  BurnDate#id=BurnDate
  BurnDate#scale_factor=0.25
  lat#_FillValue=nan
  lon#_FillValue=nan
  NC_GLOBAL#crs=EPSG:4326
  NETCDF_DIM_EXTRA={time}
  NETCDF_DIM_time_DEF={6,10}
  NETCDF_DIM_time_VALUES={0,30,61,91,122,153}
  time#calendar=proleptic_gregorian
  time#units=days since 2019-09-01 00:00:00
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 1440.0)
Upper Right (  720.0,    0.0)
Lower Right (  720.0, 1440.0)
Center      (  360.0,  720.0)
Band 1 Block=720x1 Type=Int32, ColorInterp=Undefined
  NoData Value=-2147483647
  Offset: 0,   Scale:0.25
  Metadata:
    bounds={-180,-90,180,90}
    crs=SR-ORG:6974
    crs_transform={463.3127165279165,0,-20015109.354,0,-463.3127165279167,7783653.637664001}
    data_type={'type': 'PixelType', 'precision': 'int', 'min': -32768, 'max': 32767}
    dimensions={86400,31200}
    id=BurnDate
    NETCDF_DIM_time=0
    NETCDF_VARNAME=BurnDate
    scale_factor=0.25

** TRUNCATED, REMAINING BANDS ARE SIMILAR**

Displaying it in QGIS, which reports "Layer has no coordinate reference set", so forcing EPSG:4326:

Versions:

xee @ git+https://github.com/google/Xee.git@49d82f8e80ffa3db09a204d7dd67544e5846f271
xarray==2023.1.0

Any hints for me?

google / Xee

A DataArray within a Dataset has different CRS #108