Merge pull request #57 from pblottiere/optimize_histogram

Optimize histogram computation
pblottiere · Jul 22, 2024 · cd9913f · cd9913f
2 parents 2f89d6c + 48baa63
commit cd9913f
Show file tree

Hide file tree

Showing 3 changed files with 57 additions and 22 deletions.
diff --git a/docs/src/qsa-api/endpoints.md b/docs/src/qsa-api/endpoints.md
@@ -245,7 +245,7 @@ $ curl "http://localhost:5000/api/projects/my_project/styles" \
 
 | Method  |                      URL                              |         Description                                                  |
 |---------|-------------------------------------------------------|----------------------------------------------------------------------|
-| GET     | `/api/processing/raster/histogram/{project}/{layer}`  | Return an histogram in JSON                                          |
+| POST    | `/api/processing/raster/histogram/{project}/{layer}`  | Return an histogram in JSON                                          |
 | POST    | `/api/processing/raster/calculator/{project}`         | Create a raster based on an `expression` and an `output` filename    |
 
 

diff --git a/qsa-api/qsa_api/api/processing.py b/qsa-api/qsa_api/api/processing.py
@@ -50,17 +50,46 @@ def raster_calculator(project: str):
     return jsonify(rc), 201
 
 
-@processing.get("/raster/histogram/<project>/<layer>")
+@processing.post("/raster/histogram/<project>/<layer>")
 def raster_histogram(project: str, layer: str):
-   proj = QSAProject(project)
-   if proj.exists():
-      layer_infos = proj.layer(layer)
-      if layer_infos:
-         if "type" in layer_infos and layer_infos["type"] != "raster":
-            return {"error": "Histogram is available for raster layer only"}
-         histo = Histogram(proj._qgis_project_uri, layer)
-         return jsonify(histo.process()), 201
-      else:
-         return {"error": "Layer does not exist"}, 415
-   else:
-      return {"error": "Project does not exist"}, 415
+    schema = {
+        "type": "object",
+        "properties": {
+            "min": {"type": "number"},
+            "max": {"type": "number"},
+            "count": {"type": "number"},
+        },
+    }
+
+    data = request.get_json()
+    try:
+        validate(data, schema)
+    except ValidationError as e:
+        return {"error": e.message}, 415
+
+    mini = None
+    if "min" in data:
+        mini = data["min"]
+
+    maxi = None
+    if "max" in data:
+        maxi = data["max"]
+
+    count = 1000
+    if "count" in data:
+        count = data["count"]
+
+    proj = QSAProject(project)
+    if proj.exists():
+        layer_infos = proj.layer(layer)
+        if layer_infos:
+            if "type" in layer_infos and layer_infos["type"] != "raster":
+                return {
+                    "error": "Histogram is available for raster layer only"
+                }
+            histo = Histogram(proj._qgis_project_uri, layer)
+            return jsonify(histo.process(mini, maxi, count)), 201
+        else:
+            return {"error": "Layer does not exist"}, 415
+    else:
+        return {"error": "Project does not exist"}, 415
diff --git a/qsa-api/qsa_api/processing/histogram.py b/qsa-api/qsa_api/processing/histogram.py
@@ -2,15 +2,15 @@
 
 from multiprocessing import Process, Manager
 
-from qgis.core import QgsProject
+from qgis.core import QgsProject, QgsRectangle
 
 
 class Histogram:
     def __init__(self, project_uri: str, layer: str) -> None:
         self.layer = layer
         self.project_uri = project_uri
 
-    def process(self) -> (bool, dict):
+    def process(self, mini, maxi, count) -> (bool, dict):
         # Some kind of cache is bothering us because when a raster layer is
         # added on S3, we cannot open it with GDAL provider later. The
         # QgsApplication needs to be restarted... why???
@@ -19,27 +19,33 @@ def process(self) -> (bool, dict):
 
         p = Process(
             target=Histogram._process,
-            args=(self.project_uri, self.layer, out),
+            args=(self.project_uri, self.layer, mini, maxi, count, out),
         )
         p.start()
         p.join()
 
         if "histo" in out:
-            return out["histo"]
+            return out["histo"].copy()
 
         return {}
 
     @staticmethod
-    def _process(
-        project_uri: str, layer: str, out: dict
-    ) -> None:
+    def _process(project_uri: str, layer: str, mini, maxi, count, out: dict) -> None:
 
         project = QgsProject.instance()
         project.read(project_uri)
         lyr = project.mapLayersByName(layer)[0]
 
         histo = {}
         for band in range(lyr.bandCount()):
-           histo[band+1] = lyr.dataProvider().histogram(band+1).histogramVector
+            h = (
+                lyr.dataProvider()
+                .histogram(band + 1, count, mini, maxi, QgsRectangle(), 250000)
+            )
+
+            histo[band + 1] = {}
+            histo[band + 1]["min"] = h.minimum
+            histo[band + 1]["max"] = h.maximum
+            histo[band + 1]["values"] = h.histogramVector
 
         out["histo"] = histo