[715] Enable handling of inverted latitude in reader.point (#716)

* Update point to handle inverted latitude * add formatting * add warning and test * fix test * update changelog --------- Co-authored-by: vincentsarago <[email protected]>
cogeotiff · Jul 31, 2024 · 8573dbe · 8573dbe
1 parent fbdaff8
commit 8573dbe
Show file tree

Hide file tree

Showing 3 changed files with 28 additions and 3 deletions.
diff --git a/CHANGES.md b/CHANGES.md
@@ -1,6 +1,7 @@
 # unreleased
 
 * better error message for `TileOutsideBounds` errors (author @abarciauskas-bgse, https://github.com/cogeotiff/rio-tiler/pull/712)
+* handle of inverted latitude in `reader.point` (author @georgespill, https://github.com/cogeotiff/rio-tiler/pull/716)
 
 # 6.6.1 (2024-05-17)
 

diff --git a/rio_tiler/reader.py b/rio_tiler/reader.py
@@ -199,7 +199,9 @@ def read(
                 values = dataset.read(
                     indexes=indexes,
                     window=window,
-                    out_shape=(len(indexes), height, width) if height and width else None,
+                    out_shape=(
+                        (len(indexes), height, width) if height and width else None
+                    ),
                     resampling=io_resampling,
                     boundless=boundless,
                 )
@@ -559,9 +561,23 @@ def point(
             xs, ys = transform_coords(coord_crs, dataset.crs, [lon], [lat])
             lon, lat = xs[0], ys[0]
 
+        dataset_min_lon, dataset_min_lat, dataset_max_lon, dataset_max_lat = (
+            dataset.bounds
+        )
+        # check if latitude is inverted
+        if dataset_min_lat > dataset_max_lat:
+            warnings.warn(
+                "BoundingBox of the dataset is inverted (minLat > maxLat).",
+                UserWarning,
+            )
+
+        dataset_min_lat, dataset_max_lat = (
+            min(dataset_min_lat, dataset_max_lat),
+            max(dataset_min_lat, dataset_max_lat),
+        )
         if not (
-            (dataset.bounds[0] < lon < dataset.bounds[2])
-            and (dataset.bounds[1] < lat < dataset.bounds[3])
+            (dataset_min_lon < lon < dataset_max_lon)
+            and (dataset_min_lat < lat < dataset_max_lat)
         ):
             raise PointOutsideBounds("Point is outside dataset bounds")
 

diff --git a/tests/test_reader.py b/tests/test_reader.py
@@ -35,6 +35,7 @@
 COG_NODATA_FLOAT_NAN = os.path.join(
     os.path.dirname(__file__), "fixtures", "cog_nodata_float_nan.tif"
 )
+COG_INVERTED = os.path.join(os.path.dirname(__file__), "fixtures", "inverted_lat.tif")
 
 
 @pytest.fixture(autouse=True)
@@ -850,3 +851,10 @@ def test_tile_read_nodata_float():
         prev = reader.read(src_dst, max_size=100)
         assert prev.mask[0, 0] == 0
         assert not numpy.all(prev.mask)
+
+
+def test_inverted_latitude_point():
+    """Make sure we can read a point from a file with inverted latitude."""
+    with rasterio.open(COG_INVERTED) as src_dst:
+        pt = reader.point(src_dst, [-104.77519499, 38.95367054])
+        assert pt.data[0] == -9999.0