revert to merged statistics for MultiBaseTilerFactory (#437)

vincentsarago · web-flow · commit 63822d89f25f · 2022-02-22T13:39:05.000+01:00
diff --git a/CHANGES.md b/CHANGES.md
@@ -3,6 +3,7 @@
 ## 0.5.0 (TBD)
 
 * update rio-tiler/morecantile/rio-cogeo/cogeo-mosaic versions
+* add MultiBaseTilerFactory `/asset_statistics` which will return *per asset* statistics. Returns response in form of `Dict[{asset name}, Dict[{band name}, BandStatistics]]`
 
 **breaking change**
 
@@ -16,6 +17,44 @@ expression = "b1+b2,b2"
 expression = "b1+b2;b2"
 ```
 
+* MultiBaseTilerFactory `/statistics` now returns *merged* statistics in form of `Dict[{asset_band or expression}, BandStatistics]` (instead of `Dict[{asset name}, Dict[{band name}, BandStatistics]]`)
+
+```python
+# before
+response = httpx.get(f"/stac/statistics?url=item.json").json()
+print(response)
+>>> {
+    "asset1": {
+        "1": {
+            "min": ...,
+            "max": ...,
+            ...
+        },
+        "2": {
+            "min": ...,
+            "max": ...,
+            ...
+        }
+    }
+}
+
+# now
+response = httpx.get(f"/stac/statistics?url=item.json").json()
+print(response)
+>>> {
+    "asset1_1": {
+        "min": ...,
+        "max": ...,
+        ...
+    },
+    "asset1_2": {
+        "min": ...,
+        "max": ...,
+        ...
+    },
+}
+```
+
 ## 0.4.3 (2022-02-08)
 
 * add tile `buffer` option to match rio-tiler tile options (https://github.com/developmentseed/titiler/pull/427)
diff --git a/docs/advanced/tiler_factories.md b/docs/advanced/tiler_factories.md
@@ -49,8 +49,9 @@ app.include_router(cog.router, tags=["STAC"])
 | `GET`  | `/assets`                                                       | JSON                                             | return the list of available assets
 | `GET`  | `/info`                                                         | JSON ([Info][multiinfo_model])                   | return assets basic info
 | `GET`  | `/info.geojson`                                                 | GeoJSON ([InfoGeoJSON][multiinfo_geojson_model]) | return assets basic info as a GeoJSON feature
-| `GET`  | `/statistics`                                                   | JSON ([Statistics][multistats_model])            | return assets statistics
-| `POST` | `/statistics`                                                   | GeoJSON ([Statistics][multistats_geojson_model]) | return assets statistics for a GeoJSON
+| `GET`  | `/asset_statistics`                                             | JSON ([Statistics][multistats_model])            | return per asset statistics
+| `GET`  | `/statistics`                                                   | JSON ([Statistics][stats_model])                 | return assets statistics (merged)
+| `POST` | `/statistics`                                                   | GeoJSON ([Statistics][multistats_geojson_model]) | return assets statistics for a GeoJSON (merged)
 | `GET`  | `/tiles/[{TileMatrixSetId}]/{z}/{x}/{y}[@{scale}x][.{format}]`  | image/bin                                        | create a web map tile image from assets
 | `GET`  | `/[{TileMatrixSetId}]/tilejson.json`                            | JSON ([TileJSON][tilejson_model])                | return a Mapbox TileJSON document
 | `GET`  | `/{TileMatrixSetId}/WMTSCapabilities.xml`                       | XML                                              | return OGC WMTS Get Capabilities
diff --git a/docs/endpoints/stac.md b/docs/endpoints/stac.md
@@ -21,6 +21,7 @@ app.include_router(stac.router, prefix="/stac", tags=["SpatioTemporal Asset Cata
 | `GET`  | `/stac/bounds`                                                       | JSON      | return STAC item bounds
 | `GET`  | `/stac/info`                                                         | JSON      | return asset's basic info
 | `GET`  | `/stac/info.geojson`                                                 | GeoJSON   | return asset's basic info as a GeoJSON feature
+| `GET`  | `/stac/asset_statistics`                                             | JSON      | return per asset statistics
 | `GET`  | `/stac/statistics`                                                   | JSON      | return asset's statistics
 | `POST` | `/stac/statistics`                                                   | GeoJSON   | return asset's statistics for a GeoJSON
 | `GET`  | `/stac/tiles/[{TileMatrixSetId}]/{z}/{x}/{y}[@{scale}x][.{format}]`  | image/bin | create a web map tile image from assets
@@ -285,11 +286,36 @@ Example:
 
 ### Statistics
 
+`:endpoint:/stac/asset_statistics - [GET]`
+
+- QueryParams:
+    - **url** (str): STAC Item URL. **Required**
+    - **assets** (array[str]): asset names. Default to all available assets.
+    - **asset_bidx** (array[str]): Per asset band math expression (e.g `Asset1|1;2;3`).
+    - **asset_expression** (array[str]): Per asset band math expression (e.g `Asset1|b1\*b2`).
+    - **max_size** (int): Max image size from which to calculate statistics, default is 1024.
+    - **height** (int): Force image height from which to calculate statistics.
+    - **width** (int): Force image width from which to calculate statistics.
+    - **nodata** (str, int, float): Overwrite internal Nodata value.
+    - **unscale** (bool): Apply dataset internal Scale/Offset.
+    - **resampling** (str): rasterio resampling method. Default is `nearest`.
+    - **categorical** (bool): Return statistics for categorical dataset, default is false.
+    - **c** (array[float]): Pixels values for categories.
+    - **p** (array[int]): Percentile values.
+    - **histogram_bins** (str): Histogram bins.
+    - **histogram_range** (str): Comma (',') delimited Min,Max histogram bounds
+
+Example:
+
+- `https://myendpoint/stac/statistics?url=https://somewhere.com/item.json&assets=B01&categorical=true&c=1&c=2&c=3&p=2&p98`
+
+
 `:endpoint:/stac/statistics - [GET]`
 
 - QueryParams:
     - **url** (str): STAC Item URL. **Required**
     - **assets** (array[str]): asset names. Default to all available assets.
+    - **expression** (str): rio-tiler's math expression with asset names (e.g `Asset1/Asset2`).
     - **asset_bidx** (array[str]): Per asset band math expression (e.g `Asset1|1;2;3`).
     - **asset_expression** (array[str]): Per asset band math expression (e.g `Asset1|b1\*b2`).
     - **max_size** (int): Max image size from which to calculate statistics, default is 1024.
diff --git a/src/titiler/application/titiler/application/routers/stac.py b/src/titiler/application/titiler/application/routers/stac.py
@@ -25,7 +25,7 @@ def stac_demo(request: Request):
             "request": request,
             "tilejson_endpoint": stac.url_for(request, "tilejson"),
             "info_endpoint": stac.url_for(request, "info"),
-            "statistics_endpoint": stac.url_for(request, "statistics"),
+            "statistics_endpoint": stac.url_for(request, "asset_statistics"),
         },
         media_type="text/html",
     )
diff --git a/src/titiler/core/tests/test_factories.py b/src/titiler/core/tests/test_factories.py
@@ -609,7 +609,7 @@ def test_MultiBaseTilerFactory(rio):
     rio.open = mock_rasterio_open
 
     stac = MultiBaseTilerFactory(reader=STACReader)
-    assert len(stac.router.routes) == 26
+    assert len(stac.router.routes) == 27
 
     app = FastAPI()
     app.include_router(stac.router)
@@ -691,7 +691,9 @@ def test_MultiBaseTilerFactory(rio):
     assert meta["count"] == 3
 
     # GET - statistics
-    response = client.get(f"/statistics?url={DATA_DIR}/item.json&assets=B01&assets=B09")
+    response = client.get(
+        f"/asset_statistics?url={DATA_DIR}/item.json&assets=B01&assets=B09"
+    )
     assert response.status_code == 200
     assert response.headers["content-type"] == "application/json"
     resp = response.json()
@@ -714,9 +716,8 @@ def test_MultiBaseTilerFactory(rio):
         "percentile_2",
         "percentile_98",
     }
-
     response = client.get(
-        f"/statistics?url={DATA_DIR}/item.json&assets=B01&assets=B09&asset_bidx=B01|1&asset_bidx=B09|1"
+        f"/asset_statistics?url={DATA_DIR}/item.json&assets=B01&assets=B09&asset_bidx=B01|1&asset_bidx=B09|1"
     )
     assert response.status_code == 200
     assert response.headers["content-type"] == "application/json"
@@ -725,6 +726,40 @@ def test_MultiBaseTilerFactory(rio):
     assert resp["B01"]["1"]
     assert resp["B09"]["1"]
 
+    response = client.get(f"/statistics?url={DATA_DIR}/item.json&assets=B01&assets=B09")
+    assert response.status_code == 200
+    assert response.headers["content-type"] == "application/json"
+    resp = response.json()
+    assert list(resp) == ["B01_1", "B09_1"]
+    assert set(resp["B01_1"].keys()) == {
+        "min",
+        "max",
+        "mean",
+        "count",
+        "sum",
+        "std",
+        "median",
+        "majority",
+        "minority",
+        "unique",
+        "histogram",
+        "valid_percent",
+        "masked_pixels",
+        "valid_pixels",
+        "percentile_2",
+        "percentile_98",
+    }
+
+    response = client.get(
+        f"/statistics?url={DATA_DIR}/item.json&assets=B01&assets=B09&asset_bidx=B01|1&asset_bidx=B09|1"
+    )
+    assert response.status_code == 200
+    assert response.headers["content-type"] == "application/json"
+    resp = response.json()
+    assert len(resp) == 2
+    assert resp["B01_1"]
+    assert resp["B09_1"]
+
     stac_feature = {
         "type": "FeatureCollection",
         "features": [
diff --git a/src/titiler/core/titiler/core/dependencies.py b/src/titiler/core/titiler/core/dependencies.py
@@ -237,6 +237,24 @@ def __post_init__(self):
             }
 
 
+@dataclass
+class AssetsBidxExprParamsOptional(AssetsBidxExprParams):
+    """Assets, Band Indexes and Expression parameters but with no requirement."""
+
+    def __post_init__(self):
+        """Post Init."""
+        if self.asset_indexes:
+            self.asset_indexes: Dict[str, Sequence[int]] = {  # type: ignore
+                idx.split("|")[0]: list(map(int, idx.split("|")[1].split(",")))
+                for idx in self.asset_indexes
+            }
+
+        if self.asset_expression:
+            self.asset_expression: Dict[str, str] = {  # type: ignore
+                idx.split("|")[0]: idx.split("|")[1] for idx in self.asset_expression
+            }
+
+
 @dataclass
 class AssetsBidxParams(AssetsParams):
     """asset and extra."""
diff --git a/src/titiler/core/titiler/core/factory.py b/src/titiler/core/titiler/core/factory.py
@@ -16,6 +16,7 @@
 
 from titiler.core.dependencies import (
     AssetsBidxExprParams,
+    AssetsBidxExprParamsOptional,
     AssetsBidxParams,
     AssetsParams,
     BandsExprParams,
@@ -1055,7 +1056,7 @@ def statistics(self):  # noqa: C901
 
         # GET endpoint
         @self.router.get(
-            "/statistics",
+            "/asset_statistics",
             response_class=JSONResponse,
             response_model=MultiBaseStatistics,
             responses={
@@ -1065,15 +1066,15 @@ def statistics(self):  # noqa: C901
                 }
             },
         )
-        def statistics(
+        def asset_statistics(
             src_path=Depends(self.path_dependency),
             asset_params=Depends(AssetsBidxParams),
             dataset_params=Depends(self.dataset_dependency),
             image_params=Depends(self.img_dependency),
             stats_params=Depends(self.stats_dependency),
             histogram_params=Depends(self.histogram_dependency),
         ):
-            """Create image from a geojson feature."""
+            """Per Asset statistics"""
             with rasterio.Env(**self.gdal_config):
                 with self.reader(src_path) as src_dst:
                     return src_dst.statistics(
@@ -1084,6 +1085,43 @@ def statistics(
                         hist_options={**histogram_params},
                     )
 
+        # MultiBaseReader merged statistics
+        # https://github.com/cogeotiff/rio-tiler/blob/master/rio_tiler/io/base.py#L455-L468
+        # GET endpoint
+        @self.router.get(
+            "/statistics",
+            response_class=JSONResponse,
+            response_model=Statistics,
+            responses={
+                200: {
+                    "content": {"application/json": {}},
+                    "description": "Return dataset's statistics.",
+                }
+            },
+        )
+        def statistics(
+            src_path=Depends(self.path_dependency),
+            layer_params=Depends(AssetsBidxExprParamsOptional),
+            dataset_params=Depends(self.dataset_dependency),
+            image_params=Depends(self.img_dependency),
+            stats_params=Depends(self.stats_dependency),
+            histogram_params=Depends(self.histogram_dependency),
+        ):
+            """Merged assets statistics."""
+            with rasterio.Env(**self.gdal_config):
+                with self.reader(src_path) as src_dst:
+                    # Default to all available assets
+                    if not layer_params.assets and not layer_params.expression:
+                        layer_params.assets = src_dst.assets
+
+                    return src_dst.merged_statistics(
+                        **layer_params,
+                        **image_params,
+                        **dataset_params,
+                        **stats_params,
+                        hist_options={**histogram_params},
+                    )
+
         # POST endpoint
         @self.router.post(
             "/statistics",
@@ -1102,7 +1140,7 @@ def geojson_statistics(
                 ..., description="GeoJSON Feature or FeatureCollection."
             ),
             src_path=Depends(self.path_dependency),
-            asset_params=Depends(AssetsBidxParams),
+            layer_params=Depends(AssetsBidxExprParamsOptional),
             dataset_params=Depends(self.dataset_dependency),
             image_params=Depends(self.img_dependency),
             stats_params=Depends(self.stats_dependency),
@@ -1112,15 +1150,15 @@ def geojson_statistics(
             with rasterio.Env(**self.gdal_config):
                 with self.reader(src_path) as src_dst:
                     # Default to all available assets
-                    if not asset_params.assets:
-                        asset_params.assets = src_dst.assets
+                    if not layer_params.assets and not layer_params.expression:
+                        layer_params.assets = src_dst.assets
 
                     # TODO: stream features for FeatureCollection
                     if isinstance(geojson, FeatureCollection):
                         for feature in geojson:
                             data = src_dst.feature(
                                 feature.dict(exclude_none=True),
-                                **asset_params,
+                                **layer_params,
                                 **image_params,
                                 **dataset_params,
                             )
@@ -1148,7 +1186,7 @@ def geojson_statistics(
                     else:  # simple feature
                         data = src_dst.feature(
                             geojson.dict(exclude_none=True),
-                            **asset_params,
+                            **layer_params,
                             **image_params,
                             **dataset_params,
                         )

Original file line number	Diff line number	Diff line change
`@@ -25,7 +25,7 @@ def stac_demo(request: Request):`
`25`	`25`	`"request": request,`
`26`	`26`	`"tilejson_endpoint": stac.url_for(request, "tilejson"),`
`27`	`27`	`"info_endpoint": stac.url_for(request, "info"),`
`28`		`- "statistics_endpoint": stac.url_for(request, "statistics"),`
	`28`	`+ "statistics_endpoint": stac.url_for(request, "asset_statistics"),`
`29`	`29`	`},`
`30`	`30`	`media_type="text/html",`
`31`	`31`	`)`