Remove implicit use of CRS in code base

climada/util/coordinates.py

@@ -1174,13 +1174,15 @@ def match_coordinates(
 
                 check_if_geo_coords(coords[:, 0], coords[:, 1])
                 check_if_geo_coords(coords_to_assign[:, 0], coords_to_assign[:, 1])
-            elif unit == "m":
-                coords *= 1e-3
-                coords_to_assign *= 1e-3
-            elif unit == "km":
-                pass
-            else:
-                raise ValueError("Unit must be one of 'degree', 'm' or 'km'.")
+            elif unit != "km":
+                LOGGER.warning(
+                    "You are using coordinates systems defined in %s. "
+                    "The following distance threshold will be used for coordinates "
+                    "matching: %i %s. Please adapt the distance threshold if needed. ",
+                    unit,
+                    threshold,
+                    unit,
+                )
 
             not_assigned_idx_mask = assigned_idx == -1
             assigned_idx[not_assigned_idx_mask] = nearest_neighbor_funcs[distance](
@@ -1245,21 +1247,17 @@ def match_centroids(
 
     try:
         if not equal_crs(coord_gdf.crs, centroids.crs):
-            raise ValueError("Set hazard and GeoDataFrame to same CRS first!")
-        if coord_gdf.crs is None:
             raise ValueError(
-                "Please provide a coord_gdf (coordinate GeoDataFrame)"
-                "object with a valid crs attribute."
+                "Please provide coord_gdf and centroids defined in a common CRS."
             )
-        if centroids.crs is None :
+        if coord_gdf.crs is None or centroids.crs is None:
             raise ValueError(
-                "Please provide centroids "
-                "object with a valid crs attribute."
+                "Please provide coord_gdf and centroids with valid crs attributes."
             )
     except AttributeError as exc:
         # a crs attribute is needed for unit inference
         raise ValueError(
-            "Please provide coordinate GeoDataFrame and Hazard object with a valid crs attribute."
+            "Please provide coord_gdf and centroids with valid crs attributes."
         ) from exc
 
     # get unit of coordinate systems from axis of crs
@@ -1288,8 +1286,8 @@ def _nearest_neighbor_approx(
     centroids, coordinates, unit, threshold, check_antimeridian=True
 ):
     """Compute the nearest centroid for each coordinate using the
-    euclidean distance d = ((dlon)cos(lat))^2+(dlat)^2. For distant points
-    (e.g. more than 100km apart) use the haversine distance.
+    squared equirectangular approximation distance d = ((dlon)cos(lat))^2+(dlat)^2.
+    For distant points (e.g. more than 100km apart) use the haversine distance.
 
     Parameters
     ----------
@@ -1318,8 +1316,10 @@ def _nearest_neighbor_approx(
     # first check that unit is in degree
     if unit != "degree":
         raise ValueError(
-            "Only degree unit is supported for nearest neighbor approximation."
-            "Please use euclidean distance for non-degree units."
+            "Only degree unit is supported for nearest neighbor matching using"
+            "'approx' distance. Please use euclidean distance by setting "
+            "'distance=euclidean' in the match_coordinates or match_centroids "
+            "function when using units that are not degrees."
         )
     # Compute only for the unique coordinates. Copy the results for the
     # not unique coordinates
@@ -1386,8 +1386,10 @@ def _nearest_neighbor_haversine(centroids, coordinates, unit, threshold):
     # first check that unit is in degree
     if unit != "degree":
         raise ValueError(
-            "Only degree unit is supported for nearest neighbor approximation."
-            "Please use euclidean distance for non-degree units."
+            "Only degree unit is supported for nearest neighbor matching using"
+            "'haversine' distance. Please use euclidean distance by setting "
+            "'distance=euclidean' in the match_coordinates or match_centroids "
+            "function when using units that are not degrees."
         )
     # Construct tree from centroids
     tree = BallTree(np.radians(centroids), metric="haversine")
@@ -1471,17 +1473,28 @@ def _nearest_neighbor_euclidean(
         centroids = np.rad2deg(centroids)
         coordinates = np.rad2deg(coordinates)
         dist = dist * EARTH_RADIUS_KM
+        threshold_unit = "km"
     else:
-        # if unit is not in degree, check_antimeridian is forced to False
-        check_antimeridian = False
+        threshold_unit = (
+            unit  # the unit of the threshold is considered to be in input unit
+        )
+        if check_antimeridian:
+            # if unit is not in degree, check_antimeridian is forced to False
+            check_antimeridian = False
+            LOGGER.warning(
+                "Handling of antimeridian crossing is not implemented for non-degree"
+                " coordinates ('check_antimeridian' has been set to False). Please use"
+                " degree-based coordinates system if you want to enable antimeridian crossing."
+            )
 
     # Raise a warning if the minimum distance is greater than the
     # threshold and set an unvalid index -1
     num_warn = np.sum(dist > threshold)
     if num_warn:
         LOGGER.warning(
-            "Distance to closest centroid is greater than %i km for %i coordinates.",
+            "Distance to closest centroid is greater than %i %s for %i coordinates.",
             threshold,
+            threshold_unit,
             num_warn,
         )
         assigned[dist > threshold] = -1

climada/util/test/test_coordinates.py

@@ -954,7 +954,8 @@ def test_match_centroids(self):
         with self.assertRaises(ValueError) as cm:
             u_coord.match_centroids(gdf, centroids)
         self.assertIn(
-            "Set hazard and GeoDataFrame to same CRS first!", str(cm.exception)
+            "Please provide coord_gdf and centroids defined in a common CRS.",
+            str(cm.exception),
         )
 
         # Test 4: no crs
@@ -971,8 +972,7 @@ def test_match_centroids(self):
         with self.assertRaises(ValueError) as cm:
             u_coord.match_centroids(gdf, centroids)
         self.assertIn(
-            "Please provide coordinate GeoDataFrame and "
-            "Hazard object with a valid crs attribute.",
+            "Please provide coord_gdf and centroids with valid crs attributes.",
             str(cm.exception),
         )
 
@@ -1256,7 +1256,7 @@ def repeat_coord_pass(self, dist):
         # Check copied coordinates have same neighbors
         self.assertEqual(neighbors[2], neighbors[0])
 
-    def antimeridian_warning(self, dist):
+    def distance_threshold_warning_antimeridian(self, dist):
         """Check that a warning is raised when minimum distance greater than threshold"""
         # Load input
         centroids, exposures = self.data_antimeridian_values()
@@ -1285,7 +1285,7 @@ def test_approx_repeat_coord_pass(self):
 
     def test_approx_antimeridian_warning(self):
         """Call normal_warning test for approximate distance"""
-        self.antimeridian_warning("approx")
+        self.distance_threshold_warning_antimeridian("approx")
 
     def test_haver_normal_pass(self):
         """Call normal_pass test for haversine distance"""
@@ -1301,7 +1301,7 @@ def test_haver_repeat_coord_pass(self):
 
     def test_haver_antimeridian_warning(self):
         """Call normal_warning test for haversine distance"""
-        self.antimeridian_warning("haversine")
+        self.distance_threshold_warning_antimeridian("haversine")
 
     def test_euc_normal_pass(self):
         """Call normal_pass test for euclidean distance"""
@@ -1317,7 +1317,7 @@ def test_euc_repeat_coord_pass(self):
 
     def test_euc_antimeridian_warning(self):
         """Call normal_warning test for euclidean distance"""
-        self.antimeridian_warning("euclidean")
+        self.distance_threshold_warning_antimeridian("euclidean")
 
     def test_diff_outcomes(self):
         """Different NN interpolation outcomes"""
@@ -1360,15 +1360,15 @@ def setUp_match_coordinates(self):
         # test with different unit
         self.expected_results = {
             "degree": [
-                # test with different thresholds (in degree)
+                # test with different thresholds (in degree (converted to km))
                 (100, [2, 1, 2, 0, 3, -1, 4]),
                 (20, [-1, 1, 2, 0, 3, -1, -1]),
                 (0, [-1, 1, 2, 0, -1, -1, -1]),
             ],
             "m": [
                 # test with different thresholds (in m)
-                (100, [2, 1, 2, 0, 3, -1, -1]),
-                (20, [-1, 1, 2, 0, 3, -1, -1]),
+                (100e3, [2, 1, 2, 0, 3, -1, -1]),
+                (20e3, [-1, 1, 2, 0, 3, -1, -1]),
                 (0, [-1, 1, 2, 0, -1, -1, -1]),
             ],
             "km": [
@@ -1502,6 +1502,27 @@ def test_nearest_neighbor_haversine_invalid_unit(self):
                 u_coord.NEAREST_NEIGHBOR_THRESHOLD,
             )
 
+    def antimeridian_warning_invalid_unit(self):
+        """Check that a warning is raised when coords are
+        non-degree and check_antimeridian is True"""
+
+        self.setUp_match_coordinates()
+        coords_to_assign = np.deg2rad(self.coords_to_assign) * EARTH_RADIUS_KM  # to km
+        coords = np.deg2rad(self.coords) * EARTH_RADIUS_KM
+
+        with self.assertLogs("climada.util.coordinates", level="WARNING") as cm:
+            neighbors = u_coord.match_coordinates(
+                coords_to_assign, coords, check_antimeridian=True
+            )
+        self.assertIn(
+            "Handling of antimeridian crossing is not implemented for non-degree"
+            " coordinates ('check_antimeridian' has been set to False). Please use"
+            " degree-based coordinates system if you want to enable antimeridian crossing.",
+            cm.output[1],
+        )
+
+        np.testing.assert_array_equal(neighbors, self.data_ref_antimeridian())
+
 
 class TestGetGeodata(unittest.TestCase):
     def test_nat_earth_resolution_pass(self):