diff --git a/tests/test_collections.py b/tests/test_collections.py index 08615a7..a487ffe 100644 --- a/tests/test_collections.py +++ b/tests/test_collections.py @@ -468,3 +468,43 @@ def test_bbox_collection(rio, app): assert response.headers["content-type"] == "image/tiff; application=geotiff" meta = parse_img(response.content) assert meta["crs"] == "epsg:3857" + + +def test_query_point_collections(app): + """Get values for a Point.""" + response = app.get( + f"/collections/{collection_id}/-85.5,36.1624/values", params={"assets": "cog"} + ) + + assert response.status_code == 200 + resp = response.json() + values = resp["values"] + assert len(values) == 2 + assert values[0][0] == "noaa-emergency-response/20200307aC0853130w361030" + assert values[0][2] == ["cog_b1", "cog_b2", "cog_b3"] + assert values[1][0] == "noaa-emergency-response/20200307aC0853000w361030" + + # with coord-crs + response = app.get( + f"/collections/{collection_id}/-9517816.46282489,4322990.432036275/values", + params={"assets": "cog", "coord_crs": "epsg:3857"}, + ) + assert response.status_code == 200 + resp = response.json() + assert len(resp["values"]) == 2 + + # CollectionId not found + response = app.get( + "/collections/aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa/-85.5,36.1624/values", + params={"assets": "cog"}, + ) + assert response.status_code == 404 + resp = response.json() + assert resp["detail"] == "CollectionId `aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa` not found" + + # at a point with no assets + response = app.get( + f"/collections/{collection_id}/-86.0,-35.0/values", params={"assets": "cog"} + ) + + assert response.status_code == 204 # (no content) diff --git a/tests/test_searches.py b/tests/test_searches.py index bc32cad..aef4bed 100644 --- a/tests/test_searches.py +++ b/tests/test_searches.py @@ -969,3 +969,43 @@ def test_bbox(rio, app, search_no_bbox): assert response.headers["content-type"] == "image/tiff; application=geotiff" meta = parse_img(response.content) assert meta["crs"] == "epsg:3857" + + +def test_query_point_searches(app, search_no_bbox, search_bbox): + """Test getting values for a Point.""" + response = app.get( + f"/searches/{search_no_bbox}/-85.5,36.1624/values", params={"assets": "cog"} + ) + + assert response.status_code == 200 + resp = response.json() + + values = resp["values"] + assert len(values) == 2 + assert values[0][0] == "noaa-emergency-response/20200307aC0853130w361030" + assert values[0][2] == ["cog_b1", "cog_b2", "cog_b3"] + + # with coord-crs + response = app.get( + f"/searches/{search_no_bbox}/-9517816.46282489,4322990.432036275/values", + params={"assets": "cog", "coord_crs": "epsg:3857"}, + ) + assert response.status_code == 200 + resp = response.json() + assert len(resp["values"]) == 2 + + # SearchId not found + response = app.get( + "/searches/aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa/-85.5,36.1624/values", + params={"assets": "cog"}, + ) + assert response.status_code == 404 + resp = response.json() + assert resp["detail"] == "SearchId `aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa` not found" + + # outside of searchid bbox + response = app.get( + f"/searches/{search_bbox}/-86.0,35.0/values", params={"assets": "cog"} + ) + + assert response.status_code == 204 # (no content) diff --git a/titiler/pgstac/factory.py b/titiler/pgstac/factory.py index ffff269..8e41a0b 100644 --- a/titiler/pgstac/factory.py +++ b/titiler/pgstac/factory.py @@ -51,9 +51,10 @@ from titiler.core.models.mapbox import TileJSON from titiler.core.models.responses import MultiBaseStatisticsGeoJSON from titiler.core.resources.enums import ImageType, MediaType, OptionalHeader -from titiler.core.resources.responses import GeoJSONResponse, XMLResponse +from titiler.core.resources.responses import GeoJSONResponse, JSONResponse, XMLResponse from titiler.core.utils import render_image from titiler.mosaic.factory import PixelSelectionParams +from titiler.mosaic.models.responses import Point from titiler.pgstac import model from titiler.pgstac.dependencies import ( BackendParams, @@ -150,6 +151,7 @@ def register_routes(self) -> None: self._tiles_routes() self._tilejson_routes() self._wmts_routes() + self._point_routes() if self.add_part: self._part_routes() @@ -220,7 +222,6 @@ def tile( reader_options={**reader_params}, **backend_params, ) as src_dst: - if MOSAIC_STRICT_ZOOM and ( tile.z < src_dst.minzoom or tile.z > src_dst.maxzoom ): @@ -901,6 +902,53 @@ def feature_image( return Response(content, media_type=media_type, headers=headers) + def _point_routes(self): + """Register point values endpoint.""" + + @self.router.get( + "/{lon},{lat}/values", + response_model=Point, + response_class=JSONResponse, + responses={200: {"description": "Return a value for a point"}}, + ) + def point( + lon: Annotated[float, Path(description="Longitude")], + lat: Annotated[float, Path(description="Latitude")], + search_id=Depends(self.path_dependency), + coord_crs=Depends(CoordCRSParams), + layer_params=Depends(self.layer_dependency), + dataset_params=Depends(self.dataset_dependency), + pgstac_params=Depends(self.pgstac_dependency), + backend_params=Depends(self.backend_dependency), + reader_params=Depends(self.reader_dependency), + env=Depends(self.environment_dependency), + ): + """Get Point value for a Mosaic.""" + threads = int(os.getenv("MOSAIC_CONCURRENCY", MAX_THREADS)) + + with rasterio.Env(**env): + with self.reader( + search_id, + reader_options={**reader_params}, + **backend_params, + ) as src_dst: + values = src_dst.point( + lon, + lat, + coord_crs=coord_crs or WGS84_CRS, + threads=threads, + **layer_params, + **dataset_params, + **pgstac_params, + ) + + return { + "coordinates": [lon, lat], + "values": [ + (src, pts.data.tolist(), pts.band_names) for src, pts in values + ], + } + def add_search_register_route( app: FastAPI, diff --git a/titiler/pgstac/mosaic.py b/titiler/pgstac/mosaic.py index d1112a9..7f5a18c 100644 --- a/titiler/pgstac/mosaic.py +++ b/titiler/pgstac/mosaic.py @@ -1,7 +1,7 @@ """TiTiler.PgSTAC custom Mosaic Backend and Custom STACReader.""" import json -from typing import Any, Dict, List, Optional, Tuple, Type +from typing import Any, Callable, Dict, List, Optional, Sequence, Tuple, Type import attr import rasterio @@ -17,13 +17,13 @@ from psycopg_pool import ConnectionPool from rasterio.crs import CRS from rasterio.warp import transform, transform_bounds, transform_geom -from rio_tiler.constants import WEB_MERCATOR_TMS, WGS84_CRS +from rio_tiler.constants import MAX_THREADS, WEB_MERCATOR_TMS, WGS84_CRS from rio_tiler.errors import InvalidAssetName, PointOutsideBounds from rio_tiler.io import Reader from rio_tiler.io.base import BaseReader, MultiBaseReader -from rio_tiler.models import ImageData +from rio_tiler.models import ImageData, PointData from rio_tiler.mosaic import mosaic_reader -from rio_tiler.tasks import multi_values +from rio_tiler.tasks import create_tasks, filter_tasks from rio_tiler.types import AssetInfo, BBox from titiler.pgstac.settings import CacheSettings, RetrySettings @@ -33,6 +33,30 @@ retry_config = RetrySettings() +def multi_points_pgstac( + asset_list: Sequence[Dict[str, Any]], + reader: Callable[..., PointData], + *args: Any, + threads: int = MAX_THREADS, + allowed_exceptions: Optional[Tuple] = None, + **kwargs: Any, +) -> Dict: + """Merge values returned from tasks. + + Custom version of `rio_tiler.task.multi_values` which + use constructed `item_id` as dict key. + + """ + tasks = create_tasks(reader, asset_list, threads, *args, **kwargs) + + out: Dict[str, Any] = {} + for val, asset in filter_tasks(tasks, allowed_exceptions=allowed_exceptions): + item_id = f"{asset['collection']}/{asset['id']}" + out[item_id] = val + + return out + + @attr.s class CustomSTACReader(MultiBaseReader): """Simplified STAC Reader. @@ -355,16 +379,18 @@ def _reader( item: Dict[str, Any], lon: float, lat: float, - coord_crs=coord_crs, + coord_crs: CRS = coord_crs, **kwargs: Any, - ) -> Dict: + ) -> PointData: with self.reader(item, **self.reader_options) as src_dst: return src_dst.point(lon, lat, coord_crs=coord_crs, **kwargs) if "allowed_exceptions" not in kwargs: kwargs.update({"allowed_exceptions": (PointOutsideBounds,)}) - return list(multi_values(mosaic_assets, _reader, lon, lat, **kwargs).items()) + return list( + multi_points_pgstac(mosaic_assets, _reader, lon, lat, **kwargs).items() + ) def part( self,