Repository navigation
Merging CRSs #4
Description
A useful feature (regardless of the single vs. multi CRS model discussed in #2) would be to merge all the CRSs found in a dataset or dataarray.
Let's take the example below where two datasets differ only by the name of the CRS-indexed coordinate:
>>> ds = xr.Dataset(coords={"lat": [1, 2, 3], "lon": [1, 2, 3]})
>>> ds_wgs84 = ds.proj.assign_crs(crs_wgs84="epsg:4326")
>>> ds_spatial_ref = ds.proj.assign_crs(spatial_ref="epsg:4326")>>> ds_wgs84
<xarray.Dataset> Size: 56B
Dimensions: (lat: 3, lon: 3)
Coordinates:
* lat (lat) int64 24B 1 2 3
* lon (lon) int64 24B 1 2 3
* crs_wgs84 int64 8B 0
Data variables:
*empty*
Indexes:
crs_wgs84 CRSIndex (crs=EPSG:4326)>>> ds_spatial_ref
<xarray.Dataset> Size: 56B
Dimensions: (lat: 3, lon: 3)
Coordinates:
* lat (lat) int64 24B 1 2 3
* lon (lon) int64 24B 1 2 3
* spatial_ref int64 8B 0
Data variables:
*empty*
Indexes:
spatial_ref CRSIndex (crs=EPSG:4326)Merging these two datasets with xarray.merge() will keep both CRS coordinates:
>>> merged = xr.merge([ds_wgs84, ds_spatial_ref])
>>> merged
<xarray.Dataset> Size: 64B
Dimensions: (lat: 3, lon: 3)
Coordinates:
* lat (lat) int64 24B 1 2 3
* lon (lon) int64 24B 1 2 3
* crs_wgs84 int64 8B 0
* spatial_ref int64 8B 0
Data variables:
*empty*
Indexes:
crs_wgs84 CRSIndex (crs=EPSG:4326)
spatial_ref CRSIndex (crs=EPSG:4326)Getting the unique CRS from the merged dataset won't work despite all CRS are the same:
>>> merged.proj.crs
ValueError: found more than one coordinate with a CRSIndex in Dataset or DataArrayFor convenience we could provide a proj.merge_crs() method such that:
>>> merged = merged.proj.merge_crs("spatial_ref")
>>> merged
<xarray.Dataset> Size: 56B
Dimensions: (lat: 3, lon: 3)
Coordinates:
* lat (lat) int64 24B 1 2 3
* lon (lon) int64 24B 1 2 3
* spatial_ref int64 8B 0
Data variables:
*empty*
Indexes:
spatial_ref CRSIndex (crs=EPSG:4326)
>>> merged.proj.crs
<Geographic 2D CRS: EPSG:4326>
...We could first restrict .proj.merge_crs("coord_name") to checking that all CRS indexes found in the dataset are equal and if yes merge all CRS coordinates into one (with the name given as argument). Later on we could expand the functionality with optional arguments (e.g., trigger re-projection).
Activity