added metadata property to GeoArray

io.raster.GeoArray.GeoArray():
- added property metadata: returns GeoDataFrame of all available metadata
- set_gdalDataset_meta(): now also reads metadata for each band from file on disk
- save(): now also writes metadata to file on disk

ptds/io/raster/GeoArray.py

@@ -21,6 +21,7 @@ try:
 except ImportError:
     import gdal
     import gdalnumeric
+from geopandas import GeoDataFrame, GeoSeries
 
 
 from ...geo.coord_calc                  import get_corner_coordinates, calc_FullDataset_corner_positions
@@ -88,6 +89,7 @@ class GeoArray(object):
             self._mask_nodata    = None
             self._footprint_poly = None
             self._gdalDataset_meta_already_set = False
+            self._metadata       = None
 
             if isinstance(self.arg, str):
                 assert ' ' not in self.arg, "The given path contains whitespaces. This is not supported by GDAL."
@@ -105,9 +107,6 @@ class GeoArray(object):
                 self._dtype = self.arr.dtype
                 bands = self.arr.shape[2] if len(self.arr.shape) == 3 else 1
 
-            if self.filePath:
-                self.set_gdalDataset_meta()
-
             if bandnames:
                 assert len(bandnames) == bands, \
                     'Number of given bandnames does not match number of bands in array.'
@@ -123,6 +122,9 @@ class GeoArray(object):
             if projection:
                 self.projection   = projection # use property in order to validate given value
 
+            if self.filePath:
+                self.set_gdalDataset_meta()
+
 
     @property
     def is_inmem(self):
@@ -354,6 +356,42 @@ class GeoArray(object):
             raise ValueError("'footprint_poly' can only be set from a shapely polygon or a WKT string.")
 
 
+    @property
+    def metadata(self):
+        """
+        Returns a GeoDataFrame containing all available metadata (read from file if available).
+        Use 'metadata[band_index].to_dict()' to get a metadata dictionary for a specific band.
+        Use 'metadata.loc[row_name].to_dict()' to get all metadata values of the same key for all bands as dictionary.
+        Use 'metadata.loc[row_name, band_index] = value' to set a new value.
+
+        :return:  geopandas.GeoDataFrame
+        """
+
+        if self._metadata is not None:
+            return self._metadata
+        else:
+            default = GeoDataFrame(columns=range(self.bands))
+            #for bn,idx in self.bandnames.items():
+            #    default.loc['band_index',bn] = idx
+            self._metadata = default
+            if not self.is_inmem:
+                self.set_gdalDataset_meta()
+                return self._metadata
+            else:
+                return self._metadata
+
+
+    @metadata.setter
+    def metadata(self, GDF):
+        assert isinstance(GDF, GeoDataFrame) and len(GDF.columns)==self.bands, \
+            "%s.metadata can only be set with an instance of geopandas.GeoDataFrame of which the column number " \
+            "corresponds to the band number of %s." %(self.__class__.__name__, self.__class__.__name__)
+        self._metadata = GDF
+
+
+    meta = _alias_property('metadata')
+
+
     def __getitem__(self, given):
         # TODO check if array cache contains the needed slice and return data from there
 
@@ -473,6 +511,12 @@ class GeoArray(object):
             if not 'nodata' in self._initParams or self._initParams['nodata'] is None:
                 band           = ds.GetRasterBand(1)
                 self._nodata   = band.GetNoDataValue()   # FIXME this does not support different nodata values within the same file
+
+            # read metadata
+            for b in range(self.bands):
+                band     = ds.GetRasterBand(b+1)
+                self.metadata[b] = GeoSeries(band.GetMetadata())
+
             ds = band = None
 
         self._gdalDataset_meta_already_set = True
@@ -603,10 +647,19 @@ class GeoArray(object):
         if not os.path.isdir(os.path.dirname(out_path)): os.makedirs(os.path.dirname(out_path))
 
         if self.is_inmem:
-            ds     = get_GDAL_ds_inmem(self.arr,self.geotransform, self.projection, self.nodata) # expects rows,columns,bands
             driver = gdal.GetDriverByName(fmt)
             if driver is None:
                 raise Exception("'%s' is not a supported GDAL driver." %fmt)
+
+            ds     = get_GDAL_ds_inmem(self.arr,self.geotransform, self.projection, self.nodata) # expects rows,columns,bands
+
+            # set metadata
+            if not self.metadata.empty:
+                for bn, bidx in self.bandnames.items():
+                    band = ds.GetRasterBand(bidx+1)
+                    band.SetMetadata(self.metadata[bn].to_dict())
+                    band = None
+
             driver.CreateCopy(out_path, ds, options=creationOptions if creationOptions else [])
 
             #out_arr = self.arr if self.ndim == 2 else np.swapaxes(np.swapaxes(self.arr, 0, 2), 1, 2)  # rows, columns, bands => bands, rows, columns