[Backport 2.9] Get profile and height from service

doc/integrator/raster.rst

@@ -21,8 +21,12 @@ To configure these web services, you need to set the ``raster`` variable in the
             file: /var/sig/altimetrie/mnt.vrt
             type: gdal
             round: 1
+        DTM2:
+            url: 'https://api3.geo.admin.ch/rest/services'
+            type: external_url
+            elevation_name: 'height'
 
-``raster`` is a list of "DEM layers". There are only two entries in this example, but there could be more.
+``raster`` is a list of "DEM layers".
 
 ``file`` provides the path to the shape index that references the raster files.
 The raster files should be in the Binary Terrain (BT/VTP .bt 1.3) format.
@@ -38,6 +42,13 @@ By default the nodata value is directly read from the source.
 ``round`` specifies how the result values should be rounded.
 For instance '1': round to the unit, '0.01': round to the hundredth, etc.
 
+``url`` specifies the URL of a service to use to get the elevation or profile. The service must accept a request with the this format for profiles
+``{base_url}/profile.json?geom={geom}&nbPoints={nb_points}`` and
+``{base_url}/height?easting={lon}&northing={lat}`` for points, where
+``{geom}`` is the geometry of the line, ``{nb_points}`` is the number of points to get, ``{lon}`` is the longitude of the point and ``{lat}`` is the latitude of the point.
+
+``elevation_name`` specifies the name of the elevation field in the response.
+
 .. note:: gdalbuildvrt usage example:
 
     Set the environment variables (example for Exoscale):

geoportal/c2cgeoportal_geoportal/views/profile.py

@@ -1,4 +1,4 @@
-# Copyright (c) 2012-2023, Camptocamp SA
+# Copyright (c) 2012-2024, Camptocamp SA
 # All rights reserved.
 
 # Redistribution and use in source and binary forms, with or without
@@ -25,20 +25,26 @@
 # of the authors and should not be interpreted as representing official policies,
 # either expressed or implied, of the FreeBSD Project.
 
-
+import json
+import logging
 import math
+import urllib.parse
 from decimal import Decimal
+from json.decoder import JSONDecodeError
 from typing import Any
 
 import geojson
 import pyramid.request
-from pyramid.httpexceptions import HTTPNotFound
+import requests
+from pyramid.httpexceptions import HTTPBadRequest, HTTPInternalServerError, HTTPNotFound
 from pyramid.i18n import TranslationStringFactory
 from pyramid.view import view_config
 
 from c2cgeoportal_geoportal.lib.common_headers import Cache, set_common_headers
 from c2cgeoportal_geoportal.views.raster import Raster
 
+_LOG = logging.getLogger(__name__)
+
 _ = TranslationStringFactory("c2cgeoportal")
 
 
@@ -48,6 +54,41 @@ class Profile(Raster):
     def __init__(self, request: pyramid.request.Request):
         Raster.__init__(self, request)
 
+    @staticmethod
+    def _to_filtered(points: list[dict[str, Any]], layers: list[str]) -> list[dict[str, Any]]:
+        profile = []
+        for point in points:
+            filtered_alts = {key: value for key, value in point["alts"].items() if key in layers}
+            profile.append(
+                {
+                    "dist": point["dist"],
+                    "values": filtered_alts,
+                    "x": point["easting"],
+                    "y": point["northing"],
+                }
+            )
+        return profile
+
+    def _get_profile_service_data(
+        self, layers: list[str], geom: dict[str, Any], rasters: dict[str, Any], nb_points: int
+    ) -> list[dict[str, Any]]:
+        request = f"{rasters[layers[0]]['url']}/profile.json?{urllib.parse.urlencode({'geom': geom, 'nbPoints': nb_points, 'distinct_points': 'true'})}"
+        response = requests.get(request, timeout=10)
+        if not response.ok:
+            _LOG.error("profile request %s failed with status code %s", request, response.status_code)
+            raise HTTPInternalServerError(
+                f"Failed to fetch profile data from internal request: \
+                {response.status_code} {response.reason}"
+            )
+
+        try:
+            points = json.loads(response.content)
+        except (TypeError, JSONDecodeError) as exc:
+            _LOG.exception("profile request %s failed", request)
+            raise HTTPInternalServerError("Failed to decode JSON response from internal request") from exc
+
+        return self._to_filtered(points, layers)
+
     @view_config(route_name="profile.json", renderer="fast_json")  # type: ignore
     def json(self) -> dict[str, Any]:
         """Answer to /profile.json."""
@@ -58,6 +99,9 @@ def json(self) -> dict[str, Any]:
     def _compute_points(self) -> tuple[list[str], list[dict[str, Any]]]:
         """Compute the alt=fct(dist) array."""
         geom = geojson.loads(self.request.params["geom"], object_hook=geojson.GeoJSON.to_instance)
+        nb_points = int(self.request.params["nbPoints"])
+        coords = []
+        service_results: list[dict[str, Any]] = []
 
         layers: list[str]
         if "layers" in self.request.params:
@@ -73,26 +117,47 @@ def _compute_points(self) -> tuple[list[str], list[dict[str, Any]]]:
             layers = list(rasters.keys())
             layers.sort()
 
-        points: list[dict[str, Any]] = []
+        service_layers = [layer for layer in layers if rasters[layer].get("type") == "external_url"]
 
-        dist = 0
-        prev_coord = None
-        coords = self._create_points(geom.coordinates, int(self.request.params["nbPoints"]))
-        for coord in coords:
-            if prev_coord is not None:
-                dist += self._dist(prev_coord, coord)
-
-            values = {}
-            for ref in list(rasters.keys()):
-                value = self._get_raster_value(self.rasters[ref], ref, coord[0], coord[1])
-                values[ref] = value
-
-            # 10cm accuracy is enough for distances
-            rounded_dist = Decimal(str(dist)).quantize(Decimal("0.1"))
-            points.append({"dist": rounded_dist, "values": values, "x": coord[0], "y": coord[1]})
-            prev_coord = coord
-
-        return layers, points
+        if len(service_layers) > 0:
+            urls = [rasters[layer]["url"] for layer in service_layers]
+            if len(set(urls)) != 1:
+                raise HTTPBadRequest("All service layers must have the same URL.")
+            service_results = self._get_profile_service_data(service_layers, geom, rasters, nb_points)
+            if len(service_layers) < len(layers):
+                coords = [(point["x"], point["y"]) for point in service_results]
+            else:
+                return layers, service_results
+
+        if len(service_results) == 0:
+            points: list[dict[str, Any]] = []
+
+            dist = 0
+            prev_coord = None
+            coords = self._create_points(geom.coordinates, nb_points)
+            for coord in coords:
+                if prev_coord is not None:
+                    dist += self._dist(prev_coord, coord)
+                _LOG.info("new dist %s", dist)
+
+                values = {}
+                for ref in list(rasters.keys()):
+                    value = self._get_raster_value(self.rasters[ref], ref, coord[0], coord[1])
+                    values[ref] = value
+                _LOG.info("values %s", values)
+
+                # 10cm accuracy is enough for distances
+                rounded_dist = Decimal(str(dist)).quantize(Decimal("0.1"))
+                points.append({"dist": rounded_dist, "values": values, "x": coord[0], "y": coord[1]})
+                prev_coord = coord
+            return layers, points
+        else:
+            additional_layers = [layer for layer in layers if layer not in service_layers]
+            for point in service_results:
+                for ref in additional_layers:
+                    value = self._get_raster_value(self.rasters[ref], ref, point["x"], point["y"])
+                    point["values"][ref] = value
+            return layers, service_results
 
     @staticmethod
     def _dist(coord1: tuple[float, float], coord2: tuple[float, float]) -> float:

geoportal/c2cgeoportal_geoportal/views/raster.py

@@ -27,16 +27,20 @@
 
 
 import decimal
+import json
 import logging
 import math
 import os
 import traceback
+import urllib.parse
+from json.decoder import JSONDecodeError
 from typing import TYPE_CHECKING, Any
 
 import numpy
 import pyramid.request
+import requests
 import zope.event.classhandler
-from pyramid.httpexceptions import HTTPBadRequest, HTTPNotFound
+from pyramid.httpexceptions import HTTPBadRequest, HTTPInternalServerError, HTTPNotFound
 from pyramid.view import view_config
 from rasterio.io import DatasetReader
 
@@ -95,10 +99,21 @@ def raster(self) -> dict[str, Any]:
                     raise HTTPNotFound(f"Layer {layer} not found")
         else:
             rasters = self.rasters
+            layers = list(rasters.keys())
+            layers.sort()
+
+        result: dict[str, Any] = {}
+
+        service_layers = [layer for layer in layers if rasters[layer].get("type") == "external_url"]
+
+        if len(service_layers) > 0:
+            for layer in service_layers:
+                service_result: dict[str, Any] = self._get_service_data(layer, lat, lon, rasters)
+                result.update(service_result)
 
-        result = {}
         for ref in list(rasters.keys()):
-            result[ref] = self._get_raster_value(rasters[ref], ref, lon, lat)
+            if ref not in service_layers:
+                result[ref] = self._get_raster_value(rasters[ref], ref, lon, lat)
 
         set_common_headers(self.request, "raster", Cache.PUBLIC_NO)
         return result
@@ -181,6 +196,31 @@ def get_index(index_: int) -> tuple[int, int]:
 
         return result
 
+    def _get_service_data(
+        self, layer: str, lat: float, lon: float, rasters: dict[str, Any]
+    ) -> dict[str, Any]:
+        request = (
+            f"{rasters[layer]['url']}/height?{urllib.parse.urlencode({'easting': lon, 'northing': lat})}"
+        )
+        _LOG.info("Doing height request to %s", request)
+        response = requests.get(request, timeout=10)
+        if not response.ok:
+            _LOG.error("Elevation request %s failed with status code %s", request, response.status_code)
+            raise HTTPInternalServerError(
+                f"Failed to fetch elevation data from the internal request: \
+                {response.status_code} {response.reason}"
+            )
+
+        try:
+            result = json.loads(response.content).get(rasters[layer]["elevation_name"])
+        except (TypeError, JSONDecodeError) as exc:
+            _LOG.exception("Height request to %s failed", request)
+            raise HTTPInternalServerError("Failed to decode JSON response from the internal request") from exc
+
+        set_common_headers(self.request, "raster", Cache.PUBLIC_NO)
+
+        return {layer: result}
+
     @staticmethod
     def _round(value: numpy.float32, round_to: float) -> decimal.Decimal | None:
         if value is not None: