Merge branch 'develop' into feature/copernicus_forecast

.github/workflows/testing.yml

@@ -25,6 +25,7 @@ jobs:
   build-and-test:
     name: Unit Test Pipeline
     runs-on: ubuntu-latest
+    timeout-minutes: 10
     permissions:
       # For publishing results
       checks: write

Makefile

@@ -16,15 +16,15 @@ help:  ## Use one of the following instructions:
 
 .PHONY : lint
 lint : ## Static code analysis with Pylint
-	pylint -ry climada_petals > pylint.log || true
+	python -m pylint -ry climada_petals > pylint.log || true
 
 .PHONY : unit_test
 unit_test : ## Unit tests execution with coverage and xml reports
-	pytest $(PYTEST_ARGS) --ignore=climada_petals/test climada_petals/
+	python -m pytest $(PYTEST_ARGS) --ignore=climada_petals/test climada_petals/
 
 .PHONY : install_test
 install_test : ## Test installation was successful
-	pytest $(PYTEST_JUNIT_ARGS) --pyargs climada.engine.test.test_cost_benefit \
+	python -m pytest $(PYTEST_JUNIT_ARGS) --pyargs climada.engine.test.test_cost_benefit \
 	climada.engine.test.test_impact
 
 .PHONY : data_test
@@ -37,11 +37,11 @@ notebook_test : ## Test notebooks in doc/tutorial
 
 .PHONY : integ_test
 integ_test : ## Integration tests execution with xml reports
-	pytest $(PYTEST_ARGS) climada_petals/test/
+	python -m pytest $(PYTEST_ARGS) climada_petals/test/
 
 .PHONY : test
 test : ## Unit and integration tests execution with coverage and xml reports
-	pytest $(PYTEST_ARGS) climada_petals/
+	python -m pytest $(PYTEST_ARGS) climada_petals/
 
 .PHONY : ci-clean
 ci-clean :

climada_petals/engine/supplychain.py

@@ -566,12 +566,12 @@ def calc_shock_to_sectors(self,
         self.events_date = impact.date[events_w_imp_bool]
 
         self.secs_exp = pd.DataFrame(
-            0,
+            0.0,
             index=["total_value"],
             columns=self.mriot.Z.columns
         )
         self.secs_imp = pd.DataFrame(
-            0,
+            0.0,
             index=impact.event_id[events_w_imp_bool],
             columns=self.mriot.Z.columns
         )

climada_petals/engine/test/test_supplychain.py

@@ -179,11 +179,16 @@ def dummy_exp_imp():
     " Generate dummy exposure and impacts "
     lat = np.array([1, 3])
     lon = np.array([1.5, 3.5])
-    exp = Exposures(crs=DEF_CRS)
-    exp.gdf['longitude'] = lon
-    exp.gdf['latitude'] = lat
-    exp.gdf['value'] = np.array([150., 80.])
-    exp.gdf["region_id"] = [840, 608] # USA, PHL (ROW)
+
+    exp = Exposures(
+        crs=DEF_CRS,
+        lon = lon,
+        lat = lat,
+        data = dict(
+            value = np.array([150., 80.]),
+            region_id = [840, 608], # USA, PHL (ROW)
+        )
+    )
 
     imp = Impact(
         event_id=np.arange(2) + 10,

climada_petals/entity/exposures/black_marble.py

@@ -61,6 +61,11 @@ class BlackMarble(Exposures):
 
     _metadata = Exposures._metadata + ['nightlight_file']
 
+    def __init__(self, *args, meta=None, nightlight_file=None, **kwargs):
+        super().__init__(*args, meta=meta, **kwargs)
+        meta = meta or {}
+        self.nightlight_file = Exposures._consolidate(meta, 'nightlight_file', nightlight_file)
+
     def set_countries(self, countries, ref_year=2016, res_km=None, from_hr=None,
                       admin_file='admin_0_countries', **kwargs):
         """ Model countries using values at reference year. If GDP or income
@@ -124,13 +129,6 @@ def set_countries(self, countries, ref_year=2016, res_km=None, from_hr=None,
             nightlight_file = Path(fn_nl).name,
         )
 
-        rows, cols, ras_trans = u_coord.pts_to_raster_meta(
-            (self.gdf.longitude.min(), self.gdf.latitude.min(),
-             self.gdf.longitude.max(), self.gdf.latitude.max()),
-            (coord_nl[0, 1], -coord_nl[0, 1])
-        )
-        self.meta = {'width': cols, 'height': rows, 'crs': self.crs, 'transform': ras_trans}
-
     @staticmethod
     def _set_one_country(cntry_info, nightlight, coord_nl, res_fact,
                          res_km, admin1_geom, **kwargs):

climada_petals/entity/exposures/crop_production.py

@@ -148,6 +148,11 @@ class CropProduction(Exposures):
 
     _metadata = Exposures._metadata + ['crop']
 
+    def __init__(self, *args, meta=None, crop=None, **kwargs):
+        super().__init__(*args, meta=meta, **kwargs)
+        meta = meta or {}
+        self.crop = Exposures._consolidate(meta, 'crop', crop)
+
     def set_from_isimip_netcdf(self, *args, **kwargs):
         """This function is deprecated, use LitPop.from_isimip_netcdf instead."""
         LOGGER.warning("The use of LitPop.set_from_isimip_netcdf is deprecated."
@@ -269,13 +274,27 @@ def from_isimip_netcdf(cls, input_dir=None, filename=None, hist_mean=None,
 
         # The latitude and longitude are set; the region_id is determined
         lon, lat = np.meshgrid(data.lon.values, data.lat.values)
-        exp = cls()
-        exp.gdf['latitude'] = lat.flatten()
-        exp.gdf['longitude'] = lon.flatten()
-        exp.gdf['region_id'] = u_coord.get_country_code(exp.gdf.latitude,
-                                                        exp.gdf.longitude,
-                                                        gridded=True)
-
+        latitude = lat.flatten()
+        longitude = lon.flatten()
+        region_id = u_coord.get_country_code(
+            lat=latitude,
+            lon=longitude,
+            gridded=True
+        )
+        exp = cls(
+            lat = latitude,
+            lon = longitude,
+            data = dict(
+                region_id=region_id,
+            ),
+            description = "Crop production exposure from ISIMIP"
+                         f" {CROP_NAME[crop]['print']} {irr}"
+                         f" {yearrange[0]} {yearrange[-1]}",
+            value_unit = 't/y', # input unit, will be reset below if required by user
+            crop = crop,
+            ref_year = yearrange,
+        )
+
         # The indeces of the yearrange to be extracted are determined
         time_idx = (int(yearrange[0] - yearchunk['startyear']),
                     int(yearrange[1] - yearchunk['startyear']))
@@ -310,10 +329,10 @@ def from_isimip_netcdf(cls, input_dir=None, filename=None, hist_mean=None,
                 # if irr=='combined', both 'firr' and 'noirr' are required.
                 raise ValueError(f'Invalid hist_mean provided: {hist_mean}')
             hist_mean_dict = hist_mean
-            lat_mean = exp.gdf.latitude.values
+            lat_mean = exp.latitude
         elif isinstance(hist_mean, np.ndarray) or isinstance(hist_mean, list):
             hist_mean_dict[irr_types[0]] = np.array(hist_mean)
-            lat_mean = exp.gdf.latitude.values
+            lat_mean = exp.latitude
         elif Path(hist_mean).is_dir(): # else if hist_mean is given as path to directory
         # The adequate file from the directory (depending on crop and irrigation) is extracted
         # and the variables hist_mean, lat_mean and lon_mean are set accordingly
@@ -343,12 +362,12 @@ def from_isimip_netcdf(cls, input_dir=None, filename=None, hist_mean=None,
             raise ValueError(f"Invalid hist_mean provided: {hist_mean}")
 
         # The bbox is cut out of the hist_mean data file if needed
-        if len(lat_mean) != len(exp.gdf.latitude.values):
-            idx_mean = np.zeros(len(exp.gdf.latitude.values), dtype=int)
-            for i in range(len(exp.gdf.latitude.values)):
+        if len(lat_mean) != len(exp.latitude):
+            idx_mean = np.zeros(len(exp.latitude), dtype=int)
+            for i in range(len(exp.latitude)):
                 idx_mean[i] = np.where(
-                    (lat_mean == exp.gdf.latitude.values[i])
-                    & (lon_mean == exp.gdf.longitude.values[i])
+                    (lat_mean == exp.latitude[i])
+                    & (lon_mean == exp.longitude[i])
                 )[0][0]
         else:
             idx_mean = np.arange(0, len(lat_mean))
@@ -362,27 +381,6 @@ def from_isimip_netcdf(cls, input_dir=None, filename=None, hist_mean=None,
             value_tmp = np.squeeze(area_crop[irr_val]*hist_mean_dict[irr_val][idx_mean])
             value_tmp = np.nan_to_num(value_tmp) # replace NaN by 0.0
             exp.gdf['value'] += value_tmp
-        exp.description=("Crop production exposure from ISIMIP"
-                        f" {CROP_NAME[crop]['print']} {irr}"
-                        f" {yearrange[0]} {yearrange[-1]}")
-        exp.value_unit = 't/y' # input unit, will be reset below if required by user
-        exp.crop = crop
-        exp.ref_year = yearrange
-        try:
-            rows, cols, ras_trans = u_coord.pts_to_raster_meta(
-                (exp.gdf.longitude.min(), exp.gdf.latitude.min(),
-                 exp.gdf.longitude.max(), exp.gdf.latitude.max()),
-                u_coord.get_resolution(exp.gdf.longitude, exp.gdf.latitude))
-            exp.meta = {
-                'width': cols,
-                'height': rows,
-                'crs': exp.crs,
-                'transform': ras_trans,
-            }
-        except ValueError:
-            LOGGER.warning('Could not write attribute meta, because exposure'
-                           ' has only 1 data point')
-            exp.meta = {}
 
         if 'USD' in unit:
             # set_value_to_usd() is called to compute the exposure in USD/y (country specific)
@@ -471,36 +469,32 @@ def from_area_and_yield_nc4(cls, crop_type, layer_yield, layer_area,
         # The latitude and longitude are set; region_id is determined
         lon, lat = np.meshgrid(data_area.lon.values, data_area.lat.values)
 
-        exp = cls()
         # initiate coordinates and values in GeoDatFrame:
-        exp.gdf['latitude'] = lat.flatten()
-        exp.gdf['longitude'] = lon.flatten()
-        exp.gdf['region_id'] = u_coord.get_country_code(exp.gdf.latitude,
-                                                        exp.gdf.longitude, gridded=True)
+        latitude = lat.flatten()
+        longitude = lon.flatten()
+        region_id = u_coord.get_country_code(
+            lat=latitude,
+            lon=longitude,
+            gridded=True)
+        # calc annual crop production, [t/y] = [ha] * [t/ha/y]:
+        value = np.multiply(data_area.values, data_yield.values).flatten()
+
+        exp = cls(
+            lat=latitude,
+            lon=longitude,
+            data=dict(
+                region_id=region_id,
+                value=value,
+            ),
+            crop=crop_type,
+            description=f"Annual crop production from {var_area} and {var_yield} for"
+                        f" {crop_type} from files {filename_area} and {filename_yield}",
+            value_unit = 't/y'
+        )
+
         exp.gdf[INDICATOR_IMPF + DEF_HAZ_TYPE] = 1
         exp.gdf[INDICATOR_IMPF] = 1
-        # calc annual crop production, [t/y] = [ha] * [t/ha/y]:
-        exp.gdf['value'] = np.multiply(data_area.values, data_yield.values).flatten()
-
-        exp.crop = crop_type
-        exp.description=(f"Annual crop production from {var_area} and {var_yield} for"
-                         f" {exp.crop} from files {filename_area} and {filename_yield}")
-        exp.value_unit = 't/y'
-        try:
-            rows, cols, ras_trans = u_coord.pts_to_raster_meta(
-                (exp.gdf.longitude.min(), exp.gdf.latitude.min(),
-                 exp.gdf.longitude.max(), exp.gdf.latitude.max()),
-                u_coord.get_resolution(exp.gdf.longitude, exp.gdf.latitude))
-            exp.meta = {
-                'width': cols,
-                'height': rows,
-                'crs': exp.crs,
-                'transform': ras_trans,
-            }
-        except ValueError:
-            LOGGER.warning('Could not write attribute meta, because exposure'
-                           ' has only 1 data point')
-            exp.meta = {}
+
         return exp
 
     def set_from_spam_ray_mirca(self, *args, **kwargs):

climada_petals/entity/exposures/gdp_asset.py

@@ -92,8 +92,8 @@ def set_countries(self, countries=[], reg=[], ref_year=2000,
         res = 0.0416666
 
         rows, cols, ras_trans = u_coord.pts_to_raster_meta(
-            (self.gdf.longitude.min(), self.gdf.latitude.min(),
-             self.gdf.longitude.max(), self.gdf.latitude.max()), res)
+            (self.longitude.min(), self.latitude.min(),
+             self.longitude.max(), self.latitude.max()), res)
         self.meta = {'width': cols, 'height': rows, 'crs': self.crs,
                      'transform': ras_trans}
 
@@ -130,13 +130,16 @@ def _set_one_country(countryISO, ref_year, path=None):
         reg_id_info = np.full((len(assets),), reg_id)
         impf_rf_info = np.full((len(assets),), impf_rf)
 
-        exp_gdpasset = GDP2Asset()
-        exp_gdpasset.gdf['value'] = assets
-        exp_gdpasset.gdf['latitude'] = coord[:, 0]
-        exp_gdpasset.gdf['longitude'] = coord[:, 1]
-        exp_gdpasset.gdf[INDICATOR_IMPF + DEF_HAZ_TYPE] = impf_rf_info
-        exp_gdpasset.gdf['region_id'] = reg_id_info
-        return exp_gdpasset
+        return GDP2Asset(
+            ref_year = ref_year,
+            value = assets,
+            lat = coord[:, 0],
+            lon = coord[:, 1],
+            data = {
+                INDICATOR_IMPF + DEF_HAZ_TYPE: impf_rf_info,
+                'region_id': reg_id_info,
+            }
+        )
 
 
 def _read_GDP(shp_exposures, ref_year, path=None):

climada_petals/entity/exposures/spam_agrar.py

@@ -24,6 +24,7 @@
 from pathlib import Path
 import pandas as pd
 import numpy as np
+from shapely.geometry import Point
 
 from climada import CONFIG
 from climada.entity.exposures.base import Exposures, INDICATOR_IMPF
@@ -153,12 +154,11 @@ def init_spam_agrar(self, **parameters):
             i_1 = 7  # get sum over all crops (columns 7 to 48)
             i_2 = 49
         self.gdf['value'] = data.iloc[:, i_1:i_2].sum(axis=1).values
-        self.gdf['latitude'] = lat.values
-        self.gdf['longitude'] = lon.values
+        self.gdf['geometry'] = [Point(x,y) for x,y in zip(lon.values, lat.values)]
         LOGGER.info('Lat. range: {:+.3f} to {:+.3f}.'.format(
-            np.min(self.gdf.latitude), np.max(self.gdf.latitude)))
+            np.min(self.gdf.geometry.y), np.max(self.gdf.geometry.y)))
         LOGGER.info('Lon. range: {:+.3f} to {:+.3f}.'.format(
-            np.min(self.gdf.longitude), np.max(self.gdf.longitude)))
+            np.min(self.gdf.geometry.x), np.max(self.gdf.geometry.x)))
 
         # set region_id (numeric ISO3):
         country_id = data.loc[:, 'iso3']
@@ -189,7 +189,7 @@ def init_spam_agrar(self, **parameters):
             self.value_unit = 'USD'
 
         LOGGER.info('Total {} {} {}: {:.1f} {}.'.format(
-            spam_v, spam_t, region, self.gdf.value.sum(), self.value_unit))
+            spam_v, spam_t, region, self.value.sum(), self.value_unit))
         self.check()
 
     def _set_impf(self, spam_t, haz_type):
@@ -218,7 +218,6 @@ def _set_impf(self, spam_t, haz_type):
             self.description += "\nrainfed portion of crop (= TA - TI)"
         else:
             self.gdf[INDICATOR_IMPF + haz_type] = 1
-        self.set_geometry_points()
 
     def _read_spam_file(self, **parameters):
         """Reads data from SPAM CSV file and cuts out the data for the