From 7b01fed939f4600f2612d9bb59ba915430e504e3 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Sun, 5 Feb 2023 20:23:05 +0100
Subject: [PATCH 01/21] add at_reg_id_event

---
 climada/engine/impact.py | 24 ++++++++++++++++++++++--
 1 file changed, 22 insertions(+), 2 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 35241b2bd..57759511f 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -99,6 +99,7 @@ def __init__(self,
                  crs=DEF_CRS,
                  eai_exp=None,
                  at_event=None,
+                 at_reg_id_event=None,
                  tot_value=0,
                  aai_agg=0,
                  unit='',
@@ -129,6 +130,8 @@ def __init__(self,
             expected impact for each exposure within a period of 1/frequency_unit
         at_event : np.array, optional
             impact for each hazard event
+        at_reg_id_event : np.array, optional
+            impact for each hazard event at each exposure region id
         tot_value : float, optional
             total exposure value affected
         aai_agg : float, optional
@@ -150,6 +153,7 @@ def __init__(self,
         self.crs = crs
         self.eai_exp = np.array([], float) if eai_exp is None else eai_exp
         self.at_event = np.array([], float) if at_event is None else at_event
+        self.at_reg_id_event = np.array([], float) if at_reg_id_event is None else at_reg_id_event
         self.frequency = np.array([],float) if frequency is None else frequency
         self.frequency_unit = frequency_unit
         self.tot_value = tot_value
@@ -191,8 +195,6 @@ def __init__(self,
         else:
             self.imp_mat = sparse.csr_matrix(np.empty((0, 0)))
 
-
-
     def calc(self, exposures, impact_funcs, hazard, save_mat=False, assign_centroids=True):
         """This function is deprecated, use ``ImpactCalc.impact`` instead.
         """
@@ -390,6 +392,24 @@ def impact_per_year(self, all_years=True, year_range=None):
             year_set[year] = sum(self.at_event[orig_year == year])
         return year_set
 
+    def impact_per_exp_reg(self, exposures):
+        """Aggregate impact matrix at the regional level, based on the
+        regions specified in the exposure.
+
+        Parameters
+        ----------
+        exposures : climada.entity.Exposures
+            exposure used to compute imp_mat
+        Returns
+        -------
+        np.matrix
+        """
+        self.at_reg_id_event = np.hstack([
+                       self.imp_mat[:, np.where(exposures.gdf.region_id == reg_id)[0]].sum(1)
+                       for reg_id in exposures.gdf.region_id.unique()
+                       ])
+        return self.at_reg_id_event
+
     def calc_impact_year_set(self,all_years=True, year_range=None):
         """This function is deprecated, use Impact.impact_per_year instead."""
         LOGGER.warning("The use of Impact.calc_impact_year_set is deprecated."

From f0096d754bad966b084daaf807f307cf420674b7 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Fri, 17 Feb 2023 10:44:05 +0100
Subject: [PATCH 02/21] remove as attr, add admin_0 option, add unit test

---
 climada/engine/impact.py           | 40 +++++++++++++++++---------
 climada/engine/test/test_impact.py | 45 ++++++++++++++++++++++++++++++
 2 files changed, 72 insertions(+), 13 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 57759511f..cbed6587b 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -99,7 +99,6 @@ def __init__(self,
                  crs=DEF_CRS,
                  eai_exp=None,
                  at_event=None,
-                 at_reg_id_event=None,
                  tot_value=0,
                  aai_agg=0,
                  unit='',
@@ -130,8 +129,6 @@ def __init__(self,
             expected impact for each exposure within a period of 1/frequency_unit
         at_event : np.array, optional
             impact for each hazard event
-        at_reg_id_event : np.array, optional
-            impact for each hazard event at each exposure region id
         tot_value : float, optional
             total exposure value affected
         aai_agg : float, optional
@@ -153,7 +150,6 @@ def __init__(self,
         self.crs = crs
         self.eai_exp = np.array([], float) if eai_exp is None else eai_exp
         self.at_event = np.array([], float) if at_event is None else at_event
-        self.at_reg_id_event = np.array([], float) if at_reg_id_event is None else at_reg_id_event
         self.frequency = np.array([],float) if frequency is None else frequency
         self.frequency_unit = frequency_unit
         self.tot_value = tot_value
@@ -392,9 +388,10 @@ def impact_per_year(self, all_years=True, year_range=None):
             year_set[year] = sum(self.at_event[orig_year == year])
         return year_set
 
-    def impact_per_exp_reg(self, exposures):
-        """Aggregate impact matrix at the regional level, based on the
-        regions specified in the exposure.
+    def impact_at_reg(self, exposures):
+        """Aggregate impact at the region_id level (if these are specified in exposure)
+        or at the admin 0 level (if they are not). This method works only if Impact.imp_mat
+        was stored during the impact calculation.
 
         Parameters
         ----------
@@ -402,13 +399,30 @@ def impact_per_exp_reg(self, exposures):
             exposure used to compute imp_mat
         Returns
         -------
-        np.matrix
+        pd.DataFrame
         """
-        self.at_reg_id_event = np.hstack([
-                       self.imp_mat[:, np.where(exposures.gdf.region_id == reg_id)[0]].sum(1)
-                       for reg_id in exposures.gdf.region_id.unique()
-                       ])
-        return self.at_reg_id_event
+        if self.imp_mat.nnz == 0:
+            LOGGER.warning("No Impact.imp_mat was stored during the impact calculation,"
+                           "thus the impact per region cannot be computed")
+            return None
+
+        if 'region_id' in exposures.gdf.columns:
+            regs = exposures.gdf.region_id
+
+        else:
+            LOGGER.warning("region_id in Exposure not found,"
+                           "impact is calculated at the admin_0 level")
+            regs = pd.Series(u_coord.get_country_code(exposures.gdf.latitude,
+                                                      exposures.gdf.longitude))
+
+        at_reg_event = np.hstack([
+                self.imp_mat[:, np.where(regs == reg)[0]].sum(1)
+                for reg in regs.unique()
+                ])
+
+        at_reg_event = pd.DataFrame(at_reg_event, columns=regs.unique())
+
+        return at_reg_event
 
     def calc_impact_year_set(self,all_years=True, year_range=None):
         """This function is deprecated, use Impact.impact_per_year instead."""
diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index fa8ace1f2..8cf2188d1 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -467,6 +467,50 @@ def test_local_exceedance_imp_pass(self):
         self.assertAlmostEqual(np.max(impact_rp), 2916964966.388219, places=5)
         self.assertAlmostEqual(np.min(impact_rp), 444457580.131494, places=5)
 
+class TestImpactReg(unittest.TestCase):
+    """Test impact aggregation per exposure region id or admin0"""
+    def test_impact_at_reg(self):
+        """Test calc local impacts per region"""
+
+        # Read default hazard file
+        hazard = Hazard.from_hdf5(HAZ_TEST_TC)
+
+        # Read an exposure with no region id
+        ent = Entity.from_excel(ENT_DEMO_TODAY)
+        ent.check()
+
+        # Calculate impact
+        impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=True)
+        at_reg_event = impact.impact_at_reg(ent.exposures)
+
+        self.assertEqual(at_reg_event.sum().sum(), impact.at_event.sum())
+        self.assertEqual(at_reg_event.shape[0], impact.at_event.shape[0])
+
+        # Add dummies region ids
+        region_ids = np.hstack([np.repeat(i, 10) for i in range(5)])
+        ent.exposures.gdf['region_id'] = region_ids
+
+        # Calculate impact
+        impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=True)
+        at_reg_event = impact.impact_at_reg(ent.exposures)
+
+        self.assertAlmostEqual(at_reg_event.sum().sum(), impact.at_event.sum(), places=2)
+        self.assertEqual(at_reg_event.shape[0], impact.at_event.shape[0])
+        self.assertListEqual(at_reg_event.columns.tolist(), ent.exposures.gdf.region_id.unique().tolist())
+
+        self.assertEqual(at_reg_event[0].sum(), 2071193014030.06)
+        self.assertEqual(at_reg_event[1].sum(), 2163005244090.206)
+        self.assertEqual(at_reg_event[2].sum(), 2219969197097.4907)
+        self.assertEqual(at_reg_event[3].sum(), 2203363516026.8047)
+        self.assertEqual(at_reg_event[4].sum(), 1827112892434.1558)
+
+        # Do not save Impact.imp_mat in the impact calculation
+        impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=False)
+        at_reg_event = impact.impact_at_reg(ent.exposures)
+
+        self.assertIsNone(at_reg_event)
+
+
 class TestRiskTrans(unittest.TestCase):
     """Test risk transfer methods"""
     def test_risk_trans_pass(self):
@@ -816,6 +860,7 @@ def test__exp_build_event(self):
     TESTS.addTests(unittest.TestLoader().loadTestsFromTestCase(TestImpactPerYear))
     TESTS.addTests(unittest.TestLoader().loadTestsFromTestCase(TestIO))
     TESTS.addTests(unittest.TestLoader().loadTestsFromTestCase(TestRPmatrix))
+    TESTS.addTests(unittest.TestLoader().loadTestsFromTestCase(TestImpactReg))
     TESTS.addTests(unittest.TestLoader().loadTestsFromTestCase(TestRiskTrans))
     TESTS.addTests(unittest.TestLoader().loadTestsFromTestCase(TestSelect))
     TESTS.addTests(unittest.TestLoader().loadTestsFromTestCase(TestConvertExp))

From a93474fea069be429c2f8fcdbb479f9b1ce05e78 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Fri, 17 Feb 2023 14:10:06 +0100
Subject: [PATCH 03/21] work with agg regs and not exposure

---
 climada/engine/impact.py           | 36 ++++++++++++++++--------------
 climada/engine/test/test_impact.py | 20 ++++++++---------
 2 files changed, 28 insertions(+), 28 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index cbed6587b..70288d091 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -388,15 +388,18 @@ def impact_per_year(self, all_years=True, year_range=None):
             year_set[year] = sum(self.at_event[orig_year == year])
         return year_set
 
-    def impact_at_reg(self, exposures):
-        """Aggregate impact at the region_id level (if these are specified in exposure)
-        or at the admin 0 level (if they are not). This method works only if Impact.imp_mat
-        was stored during the impact calculation.
+    def impact_at_reg(self, agg_regions=None):
+        """Aggregate impact on given aggregation regions. This method works 
+        only if Impact.imp_mat was stored during the impact calculation.
 
         Parameters
         ----------
-        exposures : climada.entity.Exposures
-            exposure used to compute imp_mat
+        agg_regions : np.array, list, pd.Series
+            It has the same lenght as the number of coordinates points in exposure 
+            and reports what macro-regions these points belong to. If no aggregation
+            regions are passed, the method aggregate impact at the admin_0 level.
+            Default is None.
+
         Returns
         -------
         pd.DataFrame
@@ -406,21 +409,20 @@ def impact_at_reg(self, exposures):
                            "thus the impact per region cannot be computed")
             return None
 
-        if 'region_id' in exposures.gdf.columns:
-            regs = exposures.gdf.region_id
-
-        else:
-            LOGGER.warning("region_id in Exposure not found,"
-                           "impact is calculated at the admin_0 level")
-            regs = pd.Series(u_coord.get_country_code(exposures.gdf.latitude,
-                                                      exposures.gdf.longitude))
+        if agg_regions is None:
+            LOGGER.warning("Aggregation regions were not specified,"
+                           "impact is aggregated at the admin_0 level")
+            agg_regions = pd.Series(u_coord.get_country_code(self.coord_exp[:,0],
+                                                             self.coord_exp[:,1]))
+        elif not isinstance(agg_regions, pd.Series):
+            agg_regions = pd.Series(agg_regions)
 
         at_reg_event = np.hstack([
-                self.imp_mat[:, np.where(regs == reg)[0]].sum(1)
-                for reg in regs.unique()
+                self.imp_mat[:, np.where(agg_regions == reg)[0]].sum(1)
+                for reg in agg_regions.unique()
                 ])
 
-        at_reg_event = pd.DataFrame(at_reg_event, columns=regs.unique())
+        at_reg_event = pd.DataFrame(at_reg_event, columns=agg_regions.unique())
 
         return at_reg_event
 
diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index 8cf2188d1..f1387faff 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -475,28 +475,26 @@ def test_impact_at_reg(self):
         # Read default hazard file
         hazard = Hazard.from_hdf5(HAZ_TEST_TC)
 
-        # Read an exposure with no region id
+        # Read an exposure
         ent = Entity.from_excel(ENT_DEMO_TODAY)
         ent.check()
 
         # Calculate impact
         impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=True)
-        at_reg_event = impact.impact_at_reg(ent.exposures)
+
+        # Aggregate impact at the admin 0 level
+        at_reg_event = impact.impact_at_reg()
 
         self.assertEqual(at_reg_event.sum().sum(), impact.at_event.sum())
         self.assertEqual(at_reg_event.shape[0], impact.at_event.shape[0])
 
-        # Add dummies region ids
+        # Aggregate impact at user-defined aggregation regions
         region_ids = np.hstack([np.repeat(i, 10) for i in range(5)])
-        ent.exposures.gdf['region_id'] = region_ids
-
-        # Calculate impact
-        impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=True)
-        at_reg_event = impact.impact_at_reg(ent.exposures)
+        at_reg_event = impact.impact_at_reg(region_ids)
 
         self.assertAlmostEqual(at_reg_event.sum().sum(), impact.at_event.sum(), places=2)
         self.assertEqual(at_reg_event.shape[0], impact.at_event.shape[0])
-        self.assertListEqual(at_reg_event.columns.tolist(), ent.exposures.gdf.region_id.unique().tolist())
+        self.assertListEqual(at_reg_event.columns.tolist(), np.unique(region_ids).tolist())
 
         self.assertEqual(at_reg_event[0].sum(), 2071193014030.06)
         self.assertEqual(at_reg_event[1].sum(), 2163005244090.206)
@@ -504,9 +502,9 @@ def test_impact_at_reg(self):
         self.assertEqual(at_reg_event[3].sum(), 2203363516026.8047)
         self.assertEqual(at_reg_event[4].sum(), 1827112892434.1558)
 
-        # Do not save Impact.imp_mat in the impact calculation
+        # Do not save Impact.imp_mat in the impact calculation and does do not aggregate
         impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=False)
-        at_reg_event = impact.impact_at_reg(ent.exposures)
+        at_reg_event = impact.impact_at_reg()
 
         self.assertIsNone(at_reg_event)
 

From 09254d12e056ab9a284cb37175b93bd54f5d4b02 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Fri, 17 Feb 2023 14:13:01 +0100
Subject: [PATCH 04/21] grammar

---
 climada/engine/impact.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 70288d091..14ec1a26a 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -397,7 +397,7 @@ def impact_at_reg(self, agg_regions=None):
         agg_regions : np.array, list, pd.Series
             It has the same lenght as the number of coordinates points in exposure 
             and reports what macro-regions these points belong to. If no aggregation
-            regions are passed, the method aggregate impact at the admin_0 level.
+            regions are passed, the method aggregates impact at the admin_0 level.
             Default is None.
 
         Returns

From 42360255668d1d2138f891cc64ba09477a800436 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Mon, 27 Feb 2023 16:21:39 +0100
Subject: [PATCH 05/21] Update climada/engine/impact.py

Co-authored-by: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
---
 climada/engine/impact.py | 2 --
 1 file changed, 2 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 14ec1a26a..b8738e6e7 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -410,8 +410,6 @@ def impact_at_reg(self, agg_regions=None):
             return None
 
         if agg_regions is None:
-            LOGGER.warning("Aggregation regions were not specified,"
-                           "impact is aggregated at the admin_0 level")
             agg_regions = pd.Series(u_coord.get_country_code(self.coord_exp[:,0],
                                                              self.coord_exp[:,1]))
         elif not isinstance(agg_regions, pd.Series):

From af28fefd35b9c515d870608de82cbc1f1c8c9d25 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Mon, 27 Feb 2023 16:21:51 +0100
Subject: [PATCH 06/21] Update climada/engine/impact.py

Co-authored-by: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
---
 climada/engine/impact.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index b8738e6e7..d10b94c8d 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -397,7 +397,7 @@ def impact_at_reg(self, agg_regions=None):
         agg_regions : np.array, list, pd.Series
             It has the same lenght as the number of coordinates points in exposure 
             and reports what macro-regions these points belong to. If no aggregation
-            regions are passed, the method aggregates impact at the admin_0 level.
+            regions are passed, the method aggregates impact at the country (admin_0) level.
             Default is None.
 
         Returns

From cff3318932ac1e2170c801e57811b15d6d8ab0d9 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Mon, 27 Feb 2023 16:22:01 +0100
Subject: [PATCH 07/21] Update climada/engine/impact.py

Co-authored-by: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
---
 climada/engine/impact.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index d10b94c8d..f8ff78200 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -394,7 +394,7 @@ def impact_at_reg(self, agg_regions=None):
 
         Parameters
         ----------
-        agg_regions : np.array, list, pd.Series
+        agg_regions : np.array, list, pd.Series (optional)
             It has the same lenght as the number of coordinates points in exposure 
             and reports what macro-regions these points belong to. If no aggregation
             regions are passed, the method aggregates impact at the country (admin_0) level.

From f96ae47addda493b2bed269e6c5f56052beb4914 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Mon, 27 Feb 2023 16:43:05 +0100
Subject: [PATCH 08/21] add haz indexes as df indexes

---
 climada/engine/impact.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index f8ff78200..5586747e7 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -420,7 +420,7 @@ def impact_at_reg(self, agg_regions=None):
                 for reg in agg_regions.unique()
                 ])
 
-        at_reg_event = pd.DataFrame(at_reg_event, columns=agg_regions.unique())
+        at_reg_event = pd.DataFrame(at_reg_event, columns=agg_regions.unique(), index=self.event_id)
 
         return at_reg_event
 

From bdcf0d40c06499a4361256bac0c2299d255f49c2 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Mon, 27 Feb 2023 16:47:43 +0100
Subject: [PATCH 09/21] raise error when no imp_mat

---
 climada/engine/impact.py | 5 ++---
 1 file changed, 2 insertions(+), 3 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 5586747e7..a066fa955 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -405,9 +405,8 @@ def impact_at_reg(self, agg_regions=None):
         pd.DataFrame
         """
         if self.imp_mat.nnz == 0:
-            LOGGER.warning("No Impact.imp_mat was stored during the impact calculation,"
-                           "thus the impact per region cannot be computed")
-            return None
+            raise ValueError("The aggregated impact cannot be computed as no"
+                             "Impact.imp_mat was stored during the impact calculation")
 
         if agg_regions is None:
             agg_regions = pd.Series(u_coord.get_country_code(self.coord_exp[:,0],

From c99ce83a8623ce84efe3ba7db2d8e69bfff91767 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Mon, 27 Feb 2023 16:48:27 +0100
Subject: [PATCH 10/21] Update climada/engine/impact.py

Co-authored-by: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
---
 climada/engine/impact.py | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index a066fa955..a3b7d6621 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -403,6 +403,8 @@ def impact_at_reg(self, agg_regions=None):
         Returns
         -------
         pd.DataFrame
+            Contains the aggregated data per event.
+            Rows: Events. Columns: Unique aggregate regions.
         """
         if self.imp_mat.nnz == 0:
             raise ValueError("The aggregated impact cannot be computed as no"

From 8cc47694cf44b0b78ece5734b7a7d0b03abace02 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Mon, 27 Feb 2023 17:02:43 +0100
Subject: [PATCH 11/21] update docstring

---
 climada/engine/impact.py | 14 +++++++++-----
 1 file changed, 9 insertions(+), 5 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index a3b7d6621..77ab7aa35 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -389,22 +389,26 @@ def impact_per_year(self, all_years=True, year_range=None):
         return year_set
 
     def impact_at_reg(self, agg_regions=None):
-        """Aggregate impact on given aggregation regions. This method works 
+        """Aggregate impact on given aggregation regions. This method works
         only if Impact.imp_mat was stored during the impact calculation.
 
         Parameters
         ----------
         agg_regions : np.array, list, pd.Series (optional)
-            It has the same lenght as the number of coordinates points in exposure 
-            and reports what macro-regions these points belong to. If no aggregation
-            regions are passed, the method aggregates impact at the country (admin_0) level.
+            It has a lenght equal to the number of coordinates points in exposure
+            and it reports what macro-regions these points belong to. For example,
+            asuming there are three centroids and agg_regions = ['A', 'A', 'B']
+            then impact of the first and second centroids will be assigned to
+            region A, whereas impact from the second centroid will be assigned
+            to area B. If no aggregation regions are passed, the method aggregates
+            impact at the country (admin_0) level.
             Default is None.
 
         Returns
         -------
         pd.DataFrame
             Contains the aggregated data per event.
-            Rows: Events. Columns: Unique aggregate regions.
+            Rows: Hzard events. Columns: Aggregation regions.
         """
         if self.imp_mat.nnz == 0:
             raise ValueError("The aggregated impact cannot be computed as no"

From 17d6d5eeb60a1cdaee66a17945a69eb60b336465 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Tue, 28 Feb 2023 11:44:31 +0100
Subject: [PATCH 12/21] update test

---
 climada/engine/test/test_impact.py | 90 +++++++++++++++++++++---------
 1 file changed, 64 insertions(+), 26 deletions(-)

diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index f1387faff..f36e18483 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -472,42 +472,80 @@ class TestImpactReg(unittest.TestCase):
     def test_impact_at_reg(self):
         """Test calc local impacts per region"""
 
-        # Read default hazard file
-        hazard = Hazard.from_hdf5(HAZ_TEST_TC)
+        # # Read default hazard file
+        # hazard = Hazard.from_hdf5(HAZ_TEST_TC)
 
-        # Read an exposure
-        ent = Entity.from_excel(ENT_DEMO_TODAY)
-        ent.check()
+        # # Read an exposure
+        # ent = Entity.from_excel(ENT_DEMO_TODAY)
+        # ent.check()
 
-        # Calculate impact
-        impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=True)
+        # # Calculate impact
+        # impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=True)
+
+        # # Aggregate impact at the admin 0 level
+        # at_reg_event = impact.impact_at_reg()
+
+        imp = dummy_impact()
+
+        # Aggregate over a single region
+        region_ids = ['A', 'A']
+        at_reg_event = imp.impact_at_reg(region_ids)
 
-        # Aggregate impact at the admin 0 level
-        at_reg_event = impact.impact_at_reg()
+        self.assertEqual(at_reg_event.sum().sum(), imp.at_event.sum())
+        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
+        self.assertEqual(at_reg_event.shape[1], np.unique(region_ids).shape[0])
 
-        self.assertEqual(at_reg_event.sum().sum(), impact.at_event.sum())
-        self.assertEqual(at_reg_event.shape[0], impact.at_event.shape[0])
+        # Aggregate over two different regions
+        region_ids = ['A', 'B']
 
-        # Aggregate impact at user-defined aggregation regions
-        region_ids = np.hstack([np.repeat(i, 10) for i in range(5)])
-        at_reg_event = impact.impact_at_reg(region_ids)
+        at_reg_event = imp.impact_at_reg(region_ids)
 
-        self.assertAlmostEqual(at_reg_event.sum().sum(), impact.at_event.sum(), places=2)
-        self.assertEqual(at_reg_event.shape[0], impact.at_event.shape[0])
-        self.assertListEqual(at_reg_event.columns.tolist(), np.unique(region_ids).tolist())
+        self.assertEqual(at_reg_event['A'].sum(), imp.imp_mat[:,0].sum())
+        self.assertEqual(at_reg_event['B'].sum(), imp.imp_mat[:,1].sum())
 
-        self.assertEqual(at_reg_event[0].sum(), 2071193014030.06)
-        self.assertEqual(at_reg_event[1].sum(), 2163005244090.206)
-        self.assertEqual(at_reg_event[2].sum(), 2219969197097.4907)
-        self.assertEqual(at_reg_event[3].sum(), 2203363516026.8047)
-        self.assertEqual(at_reg_event[4].sum(), 1827112892434.1558)
+        self.assertEqual(at_reg_event.sum().sum(), imp.at_event.sum())
+        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
+        self.assertEqual(at_reg_event.shape[1], np.unique(region_ids).shape[0])
 
-        # Do not save Impact.imp_mat in the impact calculation and does do not aggregate
-        impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=False)
-        at_reg_event = impact.impact_at_reg()
+        # Let's specify sample cities' coords and countries' code
+        CHE_num_code = 756
+        ITA_num_code = 380
+        zurich_lat, zurich_lon = 47.37, 8.55
+        bern_lat, bern_lon = 46.94, 7.44
+        rome_lat, rome_lon = 41.89, 12.51
 
-        self.assertIsNone(at_reg_event)
+        # Test admin 0 with one country
+        imp.coord_exp = np.array([[zurich_lat, zurich_lon], [bern_lat, bern_lon]])
 
+        at_reg_event = imp.impact_at_reg()
+
+        self.assertEqual(len(at_reg_event.columns), 1)
+        self.assertEqual(at_reg_event.columns[0], CHE_num_code)
+
+        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
+        self.assertEqual(at_reg_event[CHE_num_code].sum(),
+                         at_reg_event.sum().sum(),
+                         imp.at_event.sum())
+
+        # Test admin 0 with two countries
+        imp.coord_exp = np.array([[rome_lat, rome_lon], [bern_lat, bern_lon]])
+
+        at_reg_event = imp.impact_at_reg()
+
+        self.assertEqual(len(at_reg_event.columns), 2)
+        self.assertEqual(at_reg_event.columns[0], ITA_num_code)
+        self.assertEqual(at_reg_event.columns[1], CHE_num_code)
+
+        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
+        self.assertEqual(at_reg_event[ITA_num_code].sum(), imp.imp_mat[:,0].sum())
+        self.assertEqual(at_reg_event[CHE_num_code].sum(), imp.imp_mat[:,1].sum())
+        self.assertEqual(at_reg_event.sum().sum(), imp.at_event.sum())
+
+        # Test error when no imp_mat is stored
+        imp.imp_mat = sparse.csr_matrix((0, 0))
+
+        with self.assertRaises(ValueError):
+                imp.impact_at_reg()
 
 class TestRiskTrans(unittest.TestCase):
     """Test risk transfer methods"""

From f4ea128651442f86ccd44d168bd0c586cc8324c1 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Tue, 28 Feb 2023 11:53:27 +0100
Subject: [PATCH 13/21] Update climada/engine/impact.py

Co-authored-by: Chahan M. Kropf <chahan.kropf@usys.ethz.ch>
---
 climada/engine/impact.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 77ab7aa35..35fb76b41 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -408,7 +408,7 @@ def impact_at_reg(self, agg_regions=None):
         -------
         pd.DataFrame
             Contains the aggregated data per event.
-            Rows: Hzard events. Columns: Aggregation regions.
+            Rows: Hazard events. Columns: Aggregation regions.
         """
         if self.imp_mat.nnz == 0:
             raise ValueError("The aggregated impact cannot be computed as no"

From afb5f3b0b685eec1cbcd7aafd90c3566c46aa1fa Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Tue, 28 Feb 2023 11:54:16 +0100
Subject: [PATCH 14/21] delete commented code

---
 climada/engine/test/test_impact.py | 15 +--------------
 1 file changed, 1 insertion(+), 14 deletions(-)

diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index f36e18483..a3c4c40c8 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -468,23 +468,10 @@ def test_local_exceedance_imp_pass(self):
         self.assertAlmostEqual(np.min(impact_rp), 444457580.131494, places=5)
 
 class TestImpactReg(unittest.TestCase):
-    """Test impact aggregation per exposure region id or admin0"""
+    """Test impact aggregation per aggregation region or admin 0"""
     def test_impact_at_reg(self):
         """Test calc local impacts per region"""
 
-        # # Read default hazard file
-        # hazard = Hazard.from_hdf5(HAZ_TEST_TC)
-
-        # # Read an exposure
-        # ent = Entity.from_excel(ENT_DEMO_TODAY)
-        # ent.check()
-
-        # # Calculate impact
-        # impact = ImpactCalc(ent.exposures, ent.impact_funcs, hazard).impact(save_mat=True)
-
-        # # Aggregate impact at the admin 0 level
-        # at_reg_event = impact.impact_at_reg()
-
         imp = dummy_impact()
 
         # Aggregate over a single region

From f262bb154cc07b4a71b0989c384025302ffbe96d Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Tue, 28 Feb 2023 13:43:41 +0100
Subject: [PATCH 15/21] use np.array and cntries iso_codes

---
 climada/engine/impact.py           | 16 +++++++++-------
 climada/engine/test/test_impact.py | 17 +++++++++--------
 2 files changed, 18 insertions(+), 15 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 35fb76b41..3d7d56bd3 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -394,7 +394,7 @@ def impact_at_reg(self, agg_regions=None):
 
         Parameters
         ----------
-        agg_regions : np.array, list, pd.Series (optional)
+        agg_regions : np.array, list (optional)
             It has a lenght equal to the number of coordinates points in exposure
             and it reports what macro-regions these points belong to. For example,
             asuming there are three centroids and agg_regions = ['A', 'A', 'B']
@@ -415,17 +415,19 @@ def impact_at_reg(self, agg_regions=None):
                              "Impact.imp_mat was stored during the impact calculation")
 
         if agg_regions is None:
-            agg_regions = pd.Series(u_coord.get_country_code(self.coord_exp[:,0],
-                                                             self.coord_exp[:,1]))
-        elif not isinstance(agg_regions, pd.Series):
-            agg_regions = pd.Series(agg_regions)
+            agg_regions = u_coord.country_to_iso(
+                u_coord.get_country_code(self.coord_exp[:,0], self.coord_exp[:,1])
+            )
+
+        if isinstance(agg_regions, list):
+            agg_regions = np.array(agg_regions)
 
         at_reg_event = np.hstack([
                 self.imp_mat[:, np.where(agg_regions == reg)[0]].sum(1)
-                for reg in agg_regions.unique()
+                for reg in np.unique(agg_regions)
                 ])
 
-        at_reg_event = pd.DataFrame(at_reg_event, columns=agg_regions.unique(), index=self.event_id)
+        at_reg_event = pd.DataFrame(at_reg_event, columns=np.unique(agg_regions), index=self.event_id)
 
         return at_reg_event
 
diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index a3c4c40c8..2c5a84c88 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -495,8 +495,9 @@ def test_impact_at_reg(self):
         self.assertEqual(at_reg_event.shape[1], np.unique(region_ids).shape[0])
 
         # Let's specify sample cities' coords and countries' code
-        CHE_num_code = 756
-        ITA_num_code = 380
+        # CHE_num_code = 756
+        # ITA_num_code = 380
+
         zurich_lat, zurich_lon = 47.37, 8.55
         bern_lat, bern_lon = 46.94, 7.44
         rome_lat, rome_lon = 41.89, 12.51
@@ -507,10 +508,10 @@ def test_impact_at_reg(self):
         at_reg_event = imp.impact_at_reg()
 
         self.assertEqual(len(at_reg_event.columns), 1)
-        self.assertEqual(at_reg_event.columns[0], CHE_num_code)
+        self.assertEqual(at_reg_event.columns[0], 'CHE')
 
         self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
-        self.assertEqual(at_reg_event[CHE_num_code].sum(),
+        self.assertEqual(at_reg_event['CHE'].sum(),
                          at_reg_event.sum().sum(),
                          imp.at_event.sum())
 
@@ -520,12 +521,12 @@ def test_impact_at_reg(self):
         at_reg_event = imp.impact_at_reg()
 
         self.assertEqual(len(at_reg_event.columns), 2)
-        self.assertEqual(at_reg_event.columns[0], ITA_num_code)
-        self.assertEqual(at_reg_event.columns[1], CHE_num_code)
+        self.assertEqual(at_reg_event.columns[0], 'CHE')
+        self.assertEqual(at_reg_event.columns[1], 'ITA')
 
         self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
-        self.assertEqual(at_reg_event[ITA_num_code].sum(), imp.imp_mat[:,0].sum())
-        self.assertEqual(at_reg_event[CHE_num_code].sum(), imp.imp_mat[:,1].sum())
+        self.assertEqual(at_reg_event['CHE'].sum(), imp.imp_mat[:,0].sum())
+        self.assertEqual(at_reg_event['ITA'].sum(), imp.imp_mat[:,1].sum())
         self.assertEqual(at_reg_event.sum().sum(), imp.at_event.sum())
 
         # Test error when no imp_mat is stored

From 5fcf72b484265df0fd882fa07c81c4687e166fa4 Mon Sep 17 00:00:00 2001
From: Alessio Ciullo <29018864+aleeciu@users.noreply.github.com>
Date: Tue, 28 Feb 2023 13:46:06 +0100
Subject: [PATCH 16/21] remove commented code

---
 climada/engine/test/test_impact.py | 5 +----
 1 file changed, 1 insertion(+), 4 deletions(-)

diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index 2c5a84c88..e49d7f642 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -494,10 +494,7 @@ def test_impact_at_reg(self):
         self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
         self.assertEqual(at_reg_event.shape[1], np.unique(region_ids).shape[0])
 
-        # Let's specify sample cities' coords and countries' code
-        # CHE_num_code = 756
-        # ITA_num_code = 380
-
+        # Let's specify sample cities' coords
         zurich_lat, zurich_lon = 47.37, 8.55
         bern_lat, bern_lon = 46.94, 7.44
         rome_lat, rome_lon = 41.89, 12.51

From 6cce7ce835dcb4503dd3c2dbc987c6661bdfcfdc Mon Sep 17 00:00:00 2001
From: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
Date: Tue, 28 Feb 2023 17:06:52 +0100
Subject: [PATCH 17/21] Only compute unique aggregation regions once

---
 climada/engine/impact.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 3d7d56bd3..445a3f699 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -419,15 +419,15 @@ def impact_at_reg(self, agg_regions=None):
                 u_coord.get_country_code(self.coord_exp[:,0], self.coord_exp[:,1])
             )
 
-        if isinstance(agg_regions, list):
-            agg_regions = np.array(agg_regions)
+        agg_regions = np.asanyarray(agg_regions)
+        agg_reg_unique = np.unique(agg_regions)
 
         at_reg_event = np.hstack([
                 self.imp_mat[:, np.where(agg_regions == reg)[0]].sum(1)
-                for reg in np.unique(agg_regions)
+                for reg in np.unique(agg_reg_unique)
                 ])
 
-        at_reg_event = pd.DataFrame(at_reg_event, columns=np.unique(agg_regions), index=self.event_id)
+        at_reg_event = pd.DataFrame(at_reg_event, columns=np.unique(agg_reg_unique), index=self.event_id)
 
         return at_reg_event
 

From cf9c79cf215d1314d6540d9021bd11f046911244 Mon Sep 17 00:00:00 2001
From: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
Date: Tue, 28 Feb 2023 17:12:24 +0100
Subject: [PATCH 18/21] Split TestImpactReg into three test cases

---
 climada/engine/test/test_impact.py | 58 ++++++++++++++++--------------
 1 file changed, 32 insertions(+), 26 deletions(-)

diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index e49d7f642..7558bb9be 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -469,68 +469,74 @@ def test_local_exceedance_imp_pass(self):
 
 class TestImpactReg(unittest.TestCase):
     """Test impact aggregation per aggregation region or admin 0"""
-    def test_impact_at_reg(self):
-        """Test calc local impacts per region"""
 
-        imp = dummy_impact()
+    def setUp(self):
+        """Build the impact object for testing"""
+        self.imp = dummy_impact()
 
+    def test_agg_regions(self):
+        """Test calc local impacts per region"""
         # Aggregate over a single region
         region_ids = ['A', 'A']
-        at_reg_event = imp.impact_at_reg(region_ids)
+        at_reg_event = self.imp.impact_at_reg(region_ids)
 
-        self.assertEqual(at_reg_event.sum().sum(), imp.at_event.sum())
-        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
+        self.assertEqual(at_reg_event.sum().sum(), self.imp.at_event.sum())
+        self.assertEqual(at_reg_event.shape[0], self.imp.at_event.shape[0])
         self.assertEqual(at_reg_event.shape[1], np.unique(region_ids).shape[0])
 
         # Aggregate over two different regions
         region_ids = ['A', 'B']
+        at_reg_event = self.imp.impact_at_reg(region_ids)
 
-        at_reg_event = imp.impact_at_reg(region_ids)
-
-        self.assertEqual(at_reg_event['A'].sum(), imp.imp_mat[:,0].sum())
-        self.assertEqual(at_reg_event['B'].sum(), imp.imp_mat[:,1].sum())
+        self.assertEqual(at_reg_event['A'].sum(), self.imp.imp_mat[:,0].sum())
+        self.assertEqual(at_reg_event['B'].sum(), self.imp.imp_mat[:,1].sum())
 
-        self.assertEqual(at_reg_event.sum().sum(), imp.at_event.sum())
-        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
+        self.assertEqual(at_reg_event.sum().sum(), self.imp.at_event.sum())
+        self.assertEqual(at_reg_event.shape[0], self.imp.at_event.shape[0])
         self.assertEqual(at_reg_event.shape[1], np.unique(region_ids).shape[0])
 
+    def test_admin0(self):
+        """Test with aggregation to countries"""
         # Let's specify sample cities' coords
         zurich_lat, zurich_lon = 47.37, 8.55
         bern_lat, bern_lon = 46.94, 7.44
         rome_lat, rome_lon = 41.89, 12.51
 
         # Test admin 0 with one country
-        imp.coord_exp = np.array([[zurich_lat, zurich_lon], [bern_lat, bern_lon]])
+        self.imp.coord_exp = np.array([[zurich_lat, zurich_lon], [bern_lat, bern_lon]])
 
-        at_reg_event = imp.impact_at_reg()
+        at_reg_event = self.imp.impact_at_reg()
 
         self.assertEqual(len(at_reg_event.columns), 1)
         self.assertEqual(at_reg_event.columns[0], 'CHE')
 
-        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
+        self.assertEqual(at_reg_event.shape[0], self.imp.at_event.shape[0])
         self.assertEqual(at_reg_event['CHE'].sum(),
                          at_reg_event.sum().sum(),
-                         imp.at_event.sum())
+                         self.imp.at_event.sum())
 
         # Test admin 0 with two countries
-        imp.coord_exp = np.array([[rome_lat, rome_lon], [bern_lat, bern_lon]])
-
-        at_reg_event = imp.impact_at_reg()
+        self.imp.coord_exp = np.array([[rome_lat, rome_lon], [bern_lat, bern_lon]])
+        at_reg_event = self.imp.impact_at_reg()
 
         self.assertEqual(len(at_reg_event.columns), 2)
         self.assertEqual(at_reg_event.columns[0], 'CHE')
         self.assertEqual(at_reg_event.columns[1], 'ITA')
 
-        self.assertEqual(at_reg_event.shape[0], imp.at_event.shape[0])
-        self.assertEqual(at_reg_event['CHE'].sum(), imp.imp_mat[:,0].sum())
-        self.assertEqual(at_reg_event['ITA'].sum(), imp.imp_mat[:,1].sum())
-        self.assertEqual(at_reg_event.sum().sum(), imp.at_event.sum())
+        self.assertEqual(at_reg_event.shape[0], self.imp.at_event.shape[0])
+        self.assertEqual(at_reg_event['CHE'].sum(), self.imp.imp_mat[:,0].sum())
+        self.assertEqual(at_reg_event['ITA'].sum(), self.imp.imp_mat[:,1].sum())
+        self.assertEqual(at_reg_event.sum().sum(), self.imp.at_event.sum())
+
 
+    def test_no_imp_mat(self):
+        """Check error if no impact matrix is stored"""
         # Test error when no imp_mat is stored
-        imp.imp_mat = sparse.csr_matrix((0, 0))
+        self.imp.imp_mat = sparse.csr_matrix((0, 0))
 
-        with self.assertRaises(ValueError):
-                imp.impact_at_reg()
+        with self.assertRaises(ValueError) as cm:
+            self.imp.impact_at_reg()
+        self.assertIn("no Impact.imp_mat was stored", str(cm.exception))
 
 class TestRiskTrans(unittest.TestCase):
     """Test risk transfer methods"""

From c84deb0a3867d85225ecfa0b8e623b06c58bedd5 Mon Sep 17 00:00:00 2001
From: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
Date: Tue, 28 Feb 2023 17:14:00 +0100
Subject: [PATCH 19/21] Format new code with black

---
 climada/engine/impact.py           | 18 ++++++++++++------
 climada/engine/test/test_impact.py | 27 ++++++++++++++-------------
 2 files changed, 26 insertions(+), 19 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 445a3f699..93ba3c55b 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -411,23 +411,29 @@ def impact_at_reg(self, agg_regions=None):
             Rows: Hazard events. Columns: Aggregation regions.
         """
         if self.imp_mat.nnz == 0:
-            raise ValueError("The aggregated impact cannot be computed as no"
-                             "Impact.imp_mat was stored during the impact calculation")
+            raise ValueError(
+                "The aggregated impact cannot be computed as no Impact.imp_mat was "
+                "stored during the impact calculation"
+            )
 
         if agg_regions is None:
             agg_regions = u_coord.country_to_iso(
-                u_coord.get_country_code(self.coord_exp[:,0], self.coord_exp[:,1])
+                u_coord.get_country_code(self.coord_exp[:, 0], self.coord_exp[:, 1])
             )
 
         agg_regions = np.asanyarray(agg_regions)
         agg_reg_unique = np.unique(agg_regions)
 
-        at_reg_event = np.hstack([
+        at_reg_event = np.hstack(
+            [
                 self.imp_mat[:, np.where(agg_regions == reg)[0]].sum(1)
                 for reg in np.unique(agg_reg_unique)
-                ])
+            ]
+        )
 
-        at_reg_event = pd.DataFrame(at_reg_event, columns=np.unique(agg_reg_unique), index=self.event_id)
+        at_reg_event = pd.DataFrame(
+            at_reg_event, columns=np.unique(agg_reg_unique), index=self.event_id
+        )
 
         return at_reg_event
 
diff --git a/climada/engine/test/test_impact.py b/climada/engine/test/test_impact.py
index 7558bb9be..4d05feedd 100644
--- a/climada/engine/test/test_impact.py
+++ b/climada/engine/test/test_impact.py
@@ -467,6 +467,7 @@ def test_local_exceedance_imp_pass(self):
         self.assertAlmostEqual(np.max(impact_rp), 2916964966.388219, places=5)
         self.assertAlmostEqual(np.min(impact_rp), 444457580.131494, places=5)
 
+
 class TestImpactReg(unittest.TestCase):
     """Test impact aggregation per aggregation region or admin 0"""
 
@@ -477,7 +478,7 @@ def setUp(self):
     def test_agg_regions(self):
         """Test calc local impacts per region"""
         # Aggregate over a single region
-        region_ids = ['A', 'A']
+        region_ids = ["A", "A"]
         at_reg_event = self.imp.impact_at_reg(region_ids)
 
         self.assertEqual(at_reg_event.sum().sum(), self.imp.at_event.sum())
@@ -485,11 +486,11 @@ def test_agg_regions(self):
         self.assertEqual(at_reg_event.shape[1], np.unique(region_ids).shape[0])
 
         # Aggregate over two different regions
-        region_ids = ['A', 'B']
+        region_ids = ["A", "B"]
         at_reg_event = self.imp.impact_at_reg(region_ids)
 
-        self.assertEqual(at_reg_event['A'].sum(), self.imp.imp_mat[:,0].sum())
-        self.assertEqual(at_reg_event['B'].sum(), self.imp.imp_mat[:,1].sum())
+        self.assertEqual(at_reg_event["A"].sum(), self.imp.imp_mat[:, 0].sum())
+        self.assertEqual(at_reg_event["B"].sum(), self.imp.imp_mat[:, 1].sum())
 
         self.assertEqual(at_reg_event.sum().sum(), self.imp.at_event.sum())
         self.assertEqual(at_reg_event.shape[0], self.imp.at_event.shape[0])
@@ -508,27 +509,26 @@ def test_admin0(self):
         at_reg_event = self.imp.impact_at_reg()
 
         self.assertEqual(len(at_reg_event.columns), 1)
-        self.assertEqual(at_reg_event.columns[0], 'CHE')
+        self.assertEqual(at_reg_event.columns[0], "CHE")
 
         self.assertEqual(at_reg_event.shape[0], self.imp.at_event.shape[0])
-        self.assertEqual(at_reg_event['CHE'].sum(),
-                         at_reg_event.sum().sum(),
-                         self.imp.at_event.sum())
+        self.assertEqual(
+            at_reg_event["CHE"].sum(), at_reg_event.sum().sum(), self.imp.at_event.sum()
+        )
 
         # Test admin 0 with two countries
         self.imp.coord_exp = np.array([[rome_lat, rome_lon], [bern_lat, bern_lon]])
         at_reg_event = self.imp.impact_at_reg()
 
         self.assertEqual(len(at_reg_event.columns), 2)
-        self.assertEqual(at_reg_event.columns[0], 'CHE')
-        self.assertEqual(at_reg_event.columns[1], 'ITA')
+        self.assertEqual(at_reg_event.columns[0], "CHE")
+        self.assertEqual(at_reg_event.columns[1], "ITA")
 
         self.assertEqual(at_reg_event.shape[0], self.imp.at_event.shape[0])
-        self.assertEqual(at_reg_event['CHE'].sum(), self.imp.imp_mat[:,0].sum())
-        self.assertEqual(at_reg_event['ITA'].sum(), self.imp.imp_mat[:,1].sum())
+        self.assertEqual(at_reg_event["CHE"].sum(), self.imp.imp_mat[:, 0].sum())
+        self.assertEqual(at_reg_event["ITA"].sum(), self.imp.imp_mat[:, 1].sum())
         self.assertEqual(at_reg_event.sum().sum(), self.imp.at_event.sum())
 
-
     def test_no_imp_mat(self):
         """Check error if no impact matrix is stored"""
         # Test error when no imp_mat is stored
@@ -538,6 +538,7 @@ def test_no_imp_mat(self):
             self.imp.impact_at_reg()
         self.assertIn("no Impact.imp_mat was stored", str(cm.exception))
 
+
 class TestRiskTrans(unittest.TestCase):
     """Test risk transfer methods"""
     def test_risk_trans_pass(self):

From 1188c4501b90fd6072223924663428e9e261c6a8 Mon Sep 17 00:00:00 2001
From: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
Date: Wed, 1 Mar 2023 09:53:31 +0100
Subject: [PATCH 20/21] Apply suggestions for docstring

Co-authored-by: Emanuel Schmid <51439563+emanuel-schmid@users.noreply.github.com>
---
 climada/engine/impact.py | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/climada/engine/impact.py b/climada/engine/impact.py
index 93ba3c55b..5aeec134d 100755
--- a/climada/engine/impact.py
+++ b/climada/engine/impact.py
@@ -395,11 +395,11 @@ def impact_at_reg(self, agg_regions=None):
         Parameters
         ----------
         agg_regions : np.array, list (optional)
-            It has a lenght equal to the number of coordinates points in exposure
-            and it reports what macro-regions these points belong to. For example,
+            The length of the array must equal the number of centroids in exposures.
+            It reports what macro-regions these centroids belong to. For example,
             asuming there are three centroids and agg_regions = ['A', 'A', 'B']
             then impact of the first and second centroids will be assigned to
-            region A, whereas impact from the second centroid will be assigned
+            region A, whereas impact from the third centroid will be assigned
             to area B. If no aggregation regions are passed, the method aggregates
             impact at the country (admin_0) level.
             Default is None.

From 57053be705ca4b3e11d089c577d9664aa6f00009 Mon Sep 17 00:00:00 2001
From: Lukas Riedel <34276446+peanutfun@users.noreply.github.com>
Date: Wed, 1 Mar 2023 09:58:27 +0100
Subject: [PATCH 21/21] Update CHANGELOG.md

---
 CHANGELOG.md | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index 95082dcfc..1def61605 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -18,6 +18,8 @@ Removed:
 
 ### Added
 
+- `Impact.impact_at_reg` method for aggregating impacts per country or custom region [#642](https://github.com/CLIMADA-project/climada_python/pull/642)
+
 ### Changed
 
 ### Fixed