diff --git a/CHANGES.md b/CHANGES.md index b27c12e2..868db7cd 100644 --- 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 index 5ac6db53..50e4ad46 100644 --- 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 index 9f996d1c..ece1186a 100644 --- 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