diff --git a/docs/docs/Advanced.md b/docs/docs/Advanced.md index 6c2fd18..b2f64ff 100644 --- a/docs/docs/Advanced.md +++ b/docs/docs/Advanced.md @@ -41,6 +41,20 @@ print(overviews) [2, 4, 8] ``` +### Decimation Base + +As described above, a decimation base of 2 is used by default. However you can provide a custom base, N > 1, with *--decimation-base N*. Optimal overviews are computed assuming a base 2 is used, so using *--decimation-base* also requires that *--overview-level* is provided. Similar to the default example, here are the overviews for base 3: + +```python +overview_level = 3 +decimation_base = 3 +overviews = [decimation_base ** j for j in range(1, overview_level + 1)] +print(overviews) +[3, 9, 27] +``` + +This is primarily useful when working with [custom TileMatrixSets](https://developmentseed.org/morecantile/usage/#define-custom-grid) that also use a non-default decimation base. + ## Band metadata By default rio cogeo DO NOT forward **band** metadata (e.g statistics) to the output dataset. diff --git a/rio_cogeo/cogeo.py b/rio_cogeo/cogeo.py index 3744f9b..9dc12d6 100644 --- a/rio_cogeo/cogeo.py +++ b/rio_cogeo/cogeo.py @@ -98,6 +98,7 @@ def cog_translate( # noqa: C901 colormap: Optional[Dict] = None, additional_cog_metadata: Optional[Dict] = None, use_cog_driver: bool = False, + decimation_base: int = 2, ): """ Create Cloud Optimized Geotiff. @@ -171,6 +172,10 @@ def cog_translate( # noqa: C901 Additional dataset metadata to add to the COG. use_cog_driver: bool, optional (default: False) Use GDAL COG driver if set to True. COG driver is available starting with GDAL 3.1. + decimation_base: int, default: 2 + How overviews are divided at each zoom level (default is 2). Must be greater than 1. + Also requires that `overview_level` is provided for `decimation_base` values greater + than 2. .. deprecated:: 5.1.2 `web_optimized` is deprecated in favor of `tms`. @@ -187,6 +192,15 @@ def cog_translate( # noqa: C901 ) tms = tms or morecantile.tms.get("WebMercatorQuad") + if decimation_base <= 1: + raise ValueError( + "Decimation base must be greater than 1 for building overviews." + ) + elif decimation_base > 2 and overview_level is None: + raise ValueError( + "Decimation base values greater than 2 require that overview_level is defined." + ) + dst_kwargs = dst_kwargs.copy() if isinstance(indexes, int): @@ -343,7 +357,7 @@ def cog_translate( # noqa: C901 if not quiet and overview_level: click.echo("Adding overviews...", err=True) - overviews = [2**j for j in range(1, overview_level + 1)] + overviews = [decimation_base**j for j in range(1, overview_level + 1)] tmp_dst.build_overviews(overviews, ResamplingEnums[overview_resampling]) if not quiet: diff --git a/rio_cogeo/scripts/cli.py b/rio_cogeo/scripts/cli.py index 927bf25..2b45a95 100644 --- a/rio_cogeo/scripts/cli.py +++ b/rio_cogeo/scripts/cli.py @@ -122,6 +122,14 @@ def cogeo(): "selected until the smallest overview is smaller than the value of the " "internal blocksize)", ) +@click.option( + "--decimation-base", + type=int, + help="Decimation base, how to sub-divide a raster for each overview level. " + "Also requires that --overview-level is provided.", + default=2, + show_default=True, +) @click.option( "--overview-resampling", help="Overview creation resampling algorithm.", @@ -230,6 +238,7 @@ def create( add_mask, blocksize, overview_level, + decimation_base, overview_resampling, overview_blocksize, web_optimized, @@ -318,6 +327,7 @@ def create( tms=tilematrixset, use_cog_driver=use_cog_driver, quiet=quiet, + decimation_base=decimation_base, ) diff --git a/tests/test_cli.py b/tests/test_cli.py index 0bdd90d..d7f3ae1 100644 --- a/tests/test_cli.py +++ b/tests/test_cli.py @@ -210,6 +210,28 @@ def test_cogeo_validOvrOption(runner): assert src.overviews(1) == [2, 4] +def test_cogeo_decimation_base_option(runner): + """Should work as expected.""" + with runner.isolated_filesystem(): + result = runner.invoke( + cogeo, + [ + "create", + raster_path_rgb, + "output.tif", + "--decimation-base", + 3, + "--overview-level", + 6, + ], + ) + assert not result.exception + assert result.exit_code == 0 + with rasterio.open("output.tif") as src: + assert len(src.overviews(1)) == 6 + assert src.overviews(1)[0] == 3 + + def test_cogeo_overviewTilesize(monkeypatch, runner): """Should work as expected.""" with runner.isolated_filesystem(): diff --git a/tests/test_cogeo.py b/tests/test_cogeo.py index 2e34e17..d683b17 100644 --- a/tests/test_cogeo.py +++ b/tests/test_cogeo.py @@ -724,3 +724,23 @@ def test_cog_translate_forward_ns_metadata(runner): with rasterio.open("cogeo.tif") as src: assert src.tags(ns="IMD") assert src.tags(ns="RPC") + + +def test_cog_translate_decimation_base(runner): + """Should create proper overviews when using custom decimation base.""" + with runner.isolated_filesystem(): + + base_level_pairs = [(3, 6), (4, 5), (5, 4)] + + for decimation_base, overview_level in base_level_pairs: + cog_translate( + raster_path_rgb, + "cogeo.tif", + cog_profiles.get("deflate"), + decimation_base=decimation_base, + overview_level=overview_level, + quiet=True, + ) + + with rasterio.open("cogeo.tif") as src: + assert src.overviews(1)[0] == decimation_base