Consider array band count when writing Image

[MWieland]

Currently write_to_file() takes width and height from the Image array but its band number (count) from the dataset.

In case one initializes an Image() from an existing dataset and an array it may happen that the wrong band number is used internally to write to file. Take for example the case where someone computes a simple NVDI from a multi-band Image() and wants to write the result (a single-band image) back to file.

img = Image(path=multiband_file.tif)
ndvi_arr = (img.arr[3] - img.arr[2]) / (img.arr[3] + img.arr[2])
ndvi = Image(dataset=img.dataset, array=ndvi_arr)
ndvi.write_to_file(ndvi_file.tif)

This would result in an error upon writing the array to file, because img band number != ndvi band number.

The following could be a possible solution to account for this:

    def write_to_file(self, path_to_file, dtype, driver="GTiff", compress=None):
        """
        Write a dataset to file.
        :param path_to_file: str, path to new file
        :param dtype: datatype, like np.uint16, 'float32' or 'min' to use the minimum type to represent values
        :param driver: str, optional (default: 'GTiff')
        :param compress: compression, e.g. 'lzw' (default: None)
        """
        if dtype == 'min':
            dtype = get_minimum_dtype(self.__arr)

        profile = self.dataset.meta
        profile.update(
            {
                "driver": driver,
                "height": self.__arr.shape[-2],
                "width": self.__arr.shape[-1],
                "count": self.__arr.shape[0],
                "dtype": dtype,
                "transform": self.transform,
                "crs": self.crs,
            }
        )

        if compress:
            profile.update({"compress": compress})

        with rasterio.open(path_to_file, "w", **profile) as dst:
            dst.write(self.__arr.astype(dtype), indexes=range(1, self.__arr.shape[0] + 1)

[fwfichtner]

OT: Please be aware that you are quoting an old version of write_to_file().

I also don't quite understand what you're doing, img is not a dataset. This method is also not meant to write a new array to file, at least that's not how I intended it, for such purposes I would suggest using rasterio directly with something like this:

        with rasterio.open(
            destination,
            "w",
            driver="GTiff",
            width=ndvi_arr.shape[-1],
            height=ndvi_arr.shape[-2],
            count=1,
            dtype=np.float32,
            transform=self.transform,
        ) as dst:
            dst.write(ndvi_arr, indexes=1)

[MWieland]

Just edited the above comment. Should be img.dataset.

[fwfichtner]

The dataset is not the dataset you're looking for and this all becomes very confusing then. But I agree that it is very inviting to do something like that.
As an alternative we could consider something else:
Instead of initializing Image with a dataset (itself confusing, GDAL datasets will not work I guess) and array we could initalize it with crs, affine_transform and array. Then we could use an adaption of __update_dataset() to create the correct dataset internally. This would also be nice to go from GDAL to ukis-pysat in memory because you just have to convert the transform.

[MWieland]

I totally agree with that. So from a user's perspective we would call something like this:

img = Image(array=myarray, crs=mycrs, transform=mytransform, dimorder="first")

This would be more meaningful, less confusing and more flexible (I like the link to easily go from GDAL to ukis-pysat).