Add intelligent overview selection for tif to pmtiles conversion (#601)

* load overviews instead of full detail to avoid memory bomb

* fix overview level indexing & use None when no overview for clarity

* tweak nodata value for RGBA

python/rio-pmtiles/rio_pmtiles/scripts/cli.py

@@ -335,11 +335,12 @@ def pmtiles(
                 warp_options["cutline"] = shapely.wkt.dumps(cutline_rev)
 
         # Resolve the minimum and maximum zoom levels for export.
+        maxzoom_in_file = guess_maxzoom(src.crs, src.bounds, src.width, src.height, tile_size)
         if zoom_levels:
             minzoom, maxzoom = map(int, zoom_levels.split(".."))
         else:
             minzoom = 0
-            maxzoom = guess_maxzoom(src.crs, src.bounds, src.width, src.height, tile_size)
+            maxzoom = maxzoom_in_file
 
         log.debug("Zoom range: %d..%d", minzoom, maxzoom)
 
@@ -358,7 +359,7 @@ def pmtiles(
             {
                 "driver": img_format.upper(),
                 "dtype": "uint8",
-                "nodata": 0,
+                "nodata": 255 if rgba else 0,
                 "height": tile_size,
                 "width": tile_size,
                 "count": count,
@@ -442,6 +443,7 @@ def pmtiles(
                 warp_options,
                 creation_options,
                 exclude_empty_tiles,
+                maxzoom_in_file,
             ),
         ) as executor:
             for tile, contents in executor.map(process_tile, unwrap_tiles(tiles)):

python/rio-pmtiles/rio_pmtiles/worker.py

@@ -25,8 +25,9 @@ def init_worker(
     warp_opts=None,
     creation_opts=None,
     exclude_empties=True,
+    max_zoom=None,
 ):
-    global base_kwds, filename, resampling, open_options, warp_options, creation_options, exclude_empty_tiles
+    global base_kwds, filename, resampling, open_options, warp_options, creation_options, exclude_empty_tiles, max_zoom_level
     resampling = Resampling[resampling_method]
     base_kwds = profile.copy()
     filename = path
@@ -34,6 +35,7 @@ def init_worker(
     warp_options = warp_opts.copy() if warp_opts is not None else {}
     creation_options = creation_opts.copy() if creation_opts is not None else {}
     exclude_empty_tiles = exclude_empties
+    max_zoom_level = max_zoom
 
 
 def process_tile(tile):
@@ -54,7 +56,28 @@ def process_tile(tile):
         Image bytes corresponding to the tile.
 
     """
-    global base_kwds, resampling, filename, open_options, warp_options, creation_options, exclude_empty_tiles
+    global base_kwds, resampling, filename, open_options, warp_options, creation_options, exclude_empty_tiles, max_zoom_level
+
+    # Determine overview level to use
+    temp_src = rasterio.open(filename)
+    overviews = temp_src.overviews(1)
+    temp_src.close()
+    overview_level = None
+    if overviews and tile.z < max_zoom_level:
+        OVERSAMPLING_FACTOR = 4  # oversampling factor to ensure sufficient pixels for resampling operations
+        target_factor = 2 ** (max_zoom_level - tile.z) / OVERSAMPLING_FACTOR
+        best_overview = overview_level
+        best_score = float("inf")
+        for i_overview, factor in enumerate(overviews):
+            if factor <= target_factor:
+                score = abs(factor - target_factor)
+                if score < best_score:
+                    best_score = score
+                    best_overview = i_overview
+        overview_level = best_overview
+
+    if overview_level is not None:
+        open_options["OVERVIEW_LEVEL"] = overview_level
 
     with rasterio.open(filename, **open_options) as src: