This guide is part of Raster Mosaicking & Tiling for Pipeline Outputs, within the broader Automated Vector & Raster Cleaning Workflows reference.

Writing Cloud-Optimized GeoTIFFs with rio-cogeo

A cloud-optimized GeoTIFF is an ordinary GeoTIFF with two properties: it is internally tiled, and its metadata and overviews are arranged so that a reader can find what it needs in the first range request. Both are cheap to produce and expensive to retrofit.

Why the Layout Decides the Read Cost

• Striped files force whole-row reads. A 200-pixel window in a striped file pulls every row it touches, across the full width of the raster.
• External overviews are invisible to a remote reader. A .ovr sidecar is a second object most clients never request, so zoomed-out reads fall back to full resolution.
• Header at the end means two round trips minimum. A reader must locate the directory before it can request anything useful.
• Wrong block size wastes bandwidth in both directions. Blocks that are too large pull unnecessary pixels; too small and the file’s block index becomes large enough to matter.

Version and Environment Compatibility

pip install "rio-cogeo>=5.3" "rasterio>=1.3" "numpy>=1.26"

The Internal Layout, and Why It Matters

The write_cog Recipe

from __future__ import annotations

import logging
from pathlib import Path

import rasterio
from rio_cogeo.cogeo import cog_translate, cog_validate
from rio_cogeo.profiles import cog_profiles

logger = logging.getLogger(__name__)

CATEGORICAL_DTYPES = {"uint8", "int8", "uint16", "int16"}


def write_cog(
    source: Path,
    destination: Path,
    categorical: bool = False,
    blocksize: int = 512,
    overview_levels: int = 5,
    compression: str = "deflate",
) -> Path:
    """Write a validated COG, choosing predictor and overview resampling by data kind.

    `categorical` is explicit rather than inferred: an integer dtype alone does
    not distinguish a class code from an elevation in centimetres.
    """
    with rasterio.open(source) as src:
        dtype = src.dtypes[0]
        nodata = src.nodata

    profile = cog_profiles.get(compression)
    profile.update(
        blockxsize=blocksize,
        blockysize=blocksize,
        # Horizontal differencing helps continuous data and hurts categorical data.
        predictor=1 if categorical else (3 if dtype.startswith("float") else 2),
        BIGTIFF="IF_SAFER",
    )

    config = {
        "GDAL_NUM_THREADS": "ALL_CPUS",
        "GDAL_TIFF_INTERNAL_MASK": True,
        "GDAL_TIFF_OVR_BLOCKSIZE": str(blocksize),
    }

    temporary = destination.with_suffix(".part.tif")
    cog_translate(
        source,
        temporary,
        profile,
        config=config,
        nodata=nodata,
        overview_level=overview_levels,
        overview_resampling="nearest" if categorical else "average",
        web_optimized=False,
        quiet=True,
    )

    valid, errors, warnings = cog_validate(temporary)
    if not valid:
        temporary.unlink(missing_ok=True)
        raise ValueError(f"{destination.name} failed COG validation: {errors}")
    if warnings:
        logger.warning("%s: %s", destination.name, warnings)

    with rasterio.open(temporary) as check:
        if not check.overviews(1):
            temporary.unlink(missing_ok=True)
            raise ValueError(f"{destination.name}: no overviews were written")

    temporary.replace(destination)   # atomic within a filesystem
    logger.info("wrote %s (%s, %s, %d overview levels)",
                destination.name, dtype, compression, overview_levels)
    return destination

Key Implementation Notes

• categorical is a parameter, not an inference. An int16 band can be a land-cover code or an elevation in centimetres, and the two want opposite predictor and resampling settings — the distinction argued in choosing resampling methods for categorical and continuous rasters.
• Predictor 2 for integers, 3 for floats, 1 for categorical. Horizontal differencing compresses smoothly varying data well and actively harms class codes, where neighbouring values have no numeric relationship.
• Overview resampling follows the same rule. Averaging class codes produces values that exist in no legend, and the damage is invisible until someone zooms out.
• The write goes to a .part object and is renamed after validation. A half-written COG that a consumer can open is worse than no file, and the rename is the same atomic-publish pattern used across this section.
• BIGTIFF="IF_SAFER" avoids the four-gigabyte ceiling without unconditionally paying the BigTIFF header cost on small tiles.
• Validation is not optional. cog_translate can produce a file that opens perfectly and is not cloud-optimized at all — for example when the source’s mask forces a layout change — and only the validator notices.

Choosing the Overview Depth

Each overview level halves both dimensions, so five levels reduce a 10 000-pixel scene to about 312 pixels — small enough that a client showing the whole tile reads one block. Stopping short of that leaves zoomed-out clients pulling full-resolution data; going much further adds levels that no viewport ever requests.

The storage cost is bounded and small: the full overview pyramid adds roughly a third to the file size, because each level is a quarter of the one below. That is almost always worth paying, and the exception — a product read only at full resolution by machines, never displayed — is rarer than it sounds, because someone always ends up looking at it.

Build overviews as part of the COG creation rather than afterwards. An overview added later either lands in an external sidecar or forces a full rewrite, and the rewrite is exactly the operation the tiled product was designed to avoid.

Troubleshooting COG Output

Integration Note

write_cog is the publication step for each tile produced in raster mosaicking and tiling for pipeline outputs, and it belongs inside the per-tile task so that a failed validation fails only that tile. Upload the validated object rather than writing directly to object storage: cog_translate performs many small writes, and doing those over a network filesystem is both slow and a common source of corrupt output.

Masks, Alpha and Nodata

Three mechanisms can mark a pixel as absent, and mixing them produces files that behave differently in different readers.

A nodata value is a number declared in the metadata that readers treat as absent. It is the simplest and the most widely honoured, and its weakness is that the value must not occur legitimately — a nodata of zero on an elevation raster silently erases sea level.

An internal mask band stores validity per pixel independently of the data values, which removes the collision problem entirely. GDAL writes it inside the COG when GDAL_TIFF_INTERNAL_MASK is set, and every modern reader honours it. It costs roughly one extra bit per pixel after compression.

An alpha band is a full extra band with per-pixel opacity, meaningful for visual products and wasteful for analytical ones, since it doubles the band count for RGB imagery to carry information a mask expresses more cheaply.

Pick one per product and state it in the metadata. The failure to avoid is a file that declares a nodata value and carries a mask that disagrees with it, because which one wins then depends on the reader rather than on the data.

Related

• Raster Mosaicking & Tiling for Pipeline Outputs
• Building Virtual Raster Mosaics with GDAL VRT
• Choosing Resampling Methods for Categorical and Continuous Rasters
• Generating XYZ Tiles from a Raster