Update read raster intensity / fraction data

climada/entity/exposures/base.py

@@ -498,7 +498,7 @@ def from_raster(cls, file_name, band=1, src_crs=None, window=None,
         --------
         Exposures
         """
-        exp = cls()
+
         meta, value = u_coord.read_raster(file_name, [band], src_crs, window,
                                           geometry, dst_crs, transform, width,
                                           height, resampling)
@@ -507,14 +507,11 @@ def from_raster(cls, file_name, band=1, src_crs=None, window=None,
         lry = uly + meta['height'] * yres
         x_grid, y_grid = np.meshgrid(np.arange(ulx + xres / 2, lrx, xres),
                                      np.arange(uly + yres / 2, lry, yres))
+        return cls(
+            {'longitude': x_grid.flatten(), 'latitude': y_grid.flatten(), 'value': value.reshape(-1)}
+            , meta=meta, crs= meta['crs']
+            )
 
-        if exp.crs is None:
-            exp.set_crs()
-        exp.gdf['longitude'] = x_grid.flatten()
-        exp.gdf['latitude'] = y_grid.flatten()
-        exp.gdf['value'] = value.reshape(-1)
-        exp.meta = meta
-        return exp
 
     def plot_scatter(self, mask=None, ignore_zero=False, pop_name=True,
                      buffer=0.0, extend='neither', axis=None, figsize=(9, 13),

climada/hazard/base.py

@@ -340,14 +340,15 @@ def from_raster(cls, files_intensity, files_fraction=None, attrs=None,
         if haz_type is not None:
             hazard_kwargs["haz_type"] = haz_type
 
-        centroids = Centroids.from_raster_file(
+        centroids, meta = Centroids.from_raster_file(
             files_intensity[0], src_crs=src_crs, window=window, geometry=geometry, dst_crs=dst_crs,
-            transform=transform, width=width, height=height, resampling=resampling)
+            transform=transform, width=width, height=height, resampling=resampling, return_meta=True)
+
         if pool:
             chunksize = max(min(len(files_intensity) // pool.ncpus, 1000), 1)
             inten_list = pool.map(
-                centroids.values_from_raster_files,
-                [[f] for f in files_intensity],
+                _values_from_raster_files,
+                [[f] for f in files_intensity], itertools.repeat(meta),
                 itertools.repeat(band), itertools.repeat(src_crs),
                 itertools.repeat(window), itertools.repeat(geometry),
                 itertools.repeat(dst_crs), itertools.repeat(transform),
@@ -356,22 +357,22 @@ def from_raster(cls, files_intensity, files_fraction=None, attrs=None,
             intensity = sparse.vstack(inten_list, format='csr')
             if files_fraction is not None:
                 fract_list = pool.map(
-                    centroids.values_from_raster_files,
-                    [[f] for f in files_fraction],
+                    _values_from_raster_files,
+                    [[f] for f in files_fraction], itertools.repeat(meta),
                     itertools.repeat(band), itertools.repeat(src_crs),
                     itertools.repeat(window), itertools.repeat(geometry),
                     itertools.repeat(dst_crs), itertools.repeat(transform),
                     itertools.repeat(width), itertools.repeat(height),
                     itertools.repeat(resampling), chunksize=chunksize)
                 fraction = sparse.vstack(fract_list, format='csr')
         else:
-            intensity = centroids.values_from_raster_files(
-                files_intensity, band=band, src_crs=src_crs, window=window, geometry=geometry,
+            intensity = _values_from_raster_files(
+                files_intensity, meta=meta, band=band, src_crs=src_crs, window=window, geometry=geometry,
                 dst_crs=dst_crs, transform=transform, width=width, height=height,
                 resampling=resampling)
             if files_fraction is not None:
-                fraction = centroids.values_from_raster_files(
-                    files_fraction, band=band, src_crs=src_crs, window=window, geometry=geometry,
+                fraction = _values_from_raster_files(
+                    files_fraction, meta=meta, band=band, src_crs=src_crs, window=window, geometry=geometry,
                     dst_crs=dst_crs, transform=transform, width=width, height=height,
                     resampling=resampling)
 
@@ -2491,3 +2492,65 @@ def _get_fraction(self, cent_idx=None):
         if cent_idx is None:
             return self.fraction
         return self.fraction[:, cent_idx]
+
+
+#Function to read intensity/fraction values from raster files
+#No a method to allow for parallel computing
+def _values_from_raster_files(file_names, meta, band=None, src_crs=None, window=None,
+                                 geometry=None, dst_crs=None, transform=None, width=None,
+                                 height=None, resampling=Resampling.nearest):
+        """Read raster of bands and set 0 values to the masked ones.
+
+        Each band is an event. Select region using window or geometry. Reproject input by proving
+        dst_crs and/or (transform, width, height).
+
+        Parameters
+        ----------
+        file_names : str
+            path of the file
+        meta : dict
+            description of the centroids raster
+        band : list(int), optional
+            band number to read. Default: [1]
+        src_crs : crs, optional
+            source CRS. Provide it if error without it.
+        window : rasterio.windows.Window, optional
+            window to read
+        geometry : list of shapely.geometry, optional
+            consider pixels only within these shapes
+        dst_crs : crs, optional
+            reproject to given crs
+        transform : rasterio.Affine
+            affine transformation to apply
+        wdith : float
+            number of lons for transform
+        height : float
+            number of lats for transform
+        resampling : rasterio.warp,.Resampling optional
+            resampling function used for reprojection to dst_crs
+
+        Raises
+        ------
+        ValueError
+
+        Returns
+        -------
+        inten : scipy.sparse.csr_matrix
+            Each row is an event.
+        """
+        if band is None:
+            band = [1]
+
+        values = []
+        for file_name in file_names:
+            tmp_meta, data = u_coord.read_raster(
+                file_name, band, src_crs, window, geometry, dst_crs,
+                transform, width, height, resampling)
+            if (tmp_meta['crs'] != meta['crs']
+                    or tmp_meta['transform'] != meta['transform']
+                    or tmp_meta['height'] != meta['height']
+                    or tmp_meta['width'] != meta['width']):
+                raise ValueError('Raster data is inconsistent with contained raster.')
+            values.append(sparse.csr_matrix(data))
+
+        return sparse.vstack(values, format='csr')

climada/hazard/centroids/centr.py

@@ -512,11 +512,10 @@ def get_area_pixel(self, min_resol=1.0e-8):
     I/O methods
     '''
 
-
     @classmethod
     def from_raster_file(cls, file_name, src_crs=None, window=None,
                          geometry=None, dst_crs=None, transform=None, width=None,
-                         height=None, resampling=Resampling.nearest):
+                         height=None, resampling=Resampling.nearest, return_meta=False):
         """Create a new Centroids object from a raster file
 
         Select region using window or geometry. Reproject input by providing
@@ -552,7 +551,9 @@ def from_raster_file(cls, file_name, src_crs=None, window=None,
             file_name, [1], src_crs, window, geometry, dst_crs,
             transform, width, height, resampling)
         lat, lon = _meta_to_lat_lon(meta)
-        return cls(longitude=lon, latitude=lat, crs=dst_crs)
+        if return_meta:
+            return cls(longitude=lon, latitude=lat, crs=meta['crs']), meta
+        return cls(longitude=lon, latitude=lat, crs=meta['crs'])
 
     @classmethod
     def from_meta(cls, meta):
@@ -572,7 +573,6 @@ def from_meta(cls, meta):
         lat, lon = _meta_to_lat_lon(meta)
         return cls(longitude=lon, latitude=lat, crs=crs)
 
-
     @classmethod
     def from_vector_file(cls, file_name, dst_crs=None):
         """Create Centroids object from vector file (any format supported by fiona).