work on supporting writing to geotiff

Author: TonioF

test/core/store/fs/impl/test_geotiff.py

@@ -147,11 +147,11 @@ def test_read_geotiff(self):
 
 class DatasetGeoTiffFsDataAccessorTest(unittest.TestCase):
     """
-    A Test class to test opening a GeoTIFF as multilevel dataset or
+    A Test class to test opening and writing a GeoTIFF as multilevel dataset or
     as normal dataset
     """
 
-    def test_ml_to_dataset(self):
+    def test_open_ml_to_dataset(self):
         fs, cog_path = GeoTIFFMultiLevelDatasetTest.get_params("sample-cog.tif")
         data_accessor = DatasetGeoTiffFsDataAccessor()
         self.assertIsInstance(data_accessor, DatasetGeoTiffFsDataAccessor)
@@ -164,7 +164,7 @@ def test_ml_to_dataset(self):
             data_accessor.get_open_data_params_schema(cog_path), JsonSchema
         )
 
-    def test_dataset(self):
+    def test_open_dataset(self):
         fs, tiff_path = GeoTIFFMultiLevelDatasetTest.get_params("sample-geotiff.tif")
         data_accessor = DatasetGeoTiffFsDataAccessor()
         self.assertIsInstance(data_accessor, DatasetGeoTiffFsDataAccessor)
@@ -178,6 +178,24 @@ def test_dataset(self):
         )
         self.assertIsInstance(dataset, xarray.Dataset)
 
+    def test_write_dataset(self):
+        from xcube.core.new import new_cube
+        fs, _ = GeoTIFFMultiLevelDatasetTest.get_params("sample-geotiff.tif")
+        ds = new_cube(variables=dict(xxx=1, yyy=2))
+        data_accessor = DatasetGeoTiffFsDataAccessor()
+        tiff_path = "geotiff_write_test.tif"
+        data_accessor.write_data(
+            ds, tiff_path, fs=fs
+        )
+        dataset = data_accessor.open_data(
+            data_id=tiff_path,
+            fs=fs,
+            root=None,
+            tile_size=[256, 256],
+            overview_level=None,
+        )
+        self.assertIsNotNone(dataset)
+
 
 class ObjectStorageMultiLevelDatasetTest(S3Test):
     """

xcube/core/store/fs/impl/dataset.py

@@ -13,6 +13,7 @@
 import zarr
 from rasterio.session import AWSSession
 
+from xcube.core.gridmapping import GridMapping
 # Note, we need the following reference to register the
 # xarray property accessor
 # noinspection PyUnresolvedReferences
@@ -413,7 +414,33 @@ def get_write_data_params_schema(self) -> JsonObjectSchema:
     def write_data(
         self, data: xr.Dataset, data_id: str, replace=False, **write_params
     ) -> str:
-        raise NotImplementedError("Writing of GeoTIFF not yet supported")
+        assert_instance(data, xr.Dataset, name="data")
+        assert_instance(data_id, str, name="data_id")
+        fs, root, write_params = self.load_fs(write_params)
+
+        gm = GridMapping.from_dataset(data)
+
+        try:
+            with fs.open(data_id, "wb") as fobj:
+                with rasterio.open(
+                    fobj,
+                    "w",
+                    driver="GTiff",
+                    height=gm.height,
+                    width=gm.height,
+                    count=len(data.data_vars) * len(data.time),
+                    dtype=data[list(data.data_vars.keys())[0]].dtype,
+                    crs=data.rio.crs,
+                    transform=data.rio.transform(),
+                ) as dst:
+                    band_index = 1
+                    for var in data.data_vars:
+                        for time in data.time.values:
+                            dst.write(data[var].sel(time=time).values, band_index)
+                            band_index += 1
+        except ValueError as e:
+            raise DataStoreError(f"Failed to write dataset {data_id!r}: {e}") from e
+        return data_id
 
     @classmethod
     def _sanitize_dataset_attrs(cls, dataset):