How to Serve Raster Tiles from a COG with TiTiler

On this page

Problem statement

A GeoTIFF needs to appear on a web map. The traditional route is to cut a tile pyramid — hours of processing, gigabytes of small files, and a rerun every time the raster changes.

A Cloud Optimised GeoTIFF plus a dynamic tile server removes both. The COG is internally tiled and carries overviews, so a tile server can read exactly the bytes one tile needs over a range request and render it on demand. TiTiler is that server: a FastAPI application built on rasterio and rio-tiler that turns any COG URL into an XYZ tile endpoint.

The parts that go wrong are predictable: the source is not actually a valid COG, the rescaling is wrong so every tile is black or white, nothing is cached, and the server reads from slow storage on every request.

Quick answer

pip install "titiler.application"
uvicorn titiler.application.main:app --port 8000

http://localhost:8000/cog/tiles/WebMercatorQuad/{z}/{x}/{y}.png
    ?url=https://example.org/scene.tif
    &rescale=0,3000
    &colormap_name=viridis

and in MapLibre:

map.addSource("dem", {
  type: "raster",
  tiles: ["http://localhost:8000/cog/tiles/WebMercatorQuad/{z}/{x}/{y}.png" +
          "?url=" + encodeURIComponent(COG_URL) + "&rescale=0,3000"],
  tileSize: 256,
});

The two parameters that decide whether anything appears are url and rescale. Everything else is styling.

Validate the source before blaming the server.

A real COG, inspected, rescaled and coloured

1. Make sure the source is a real COG

"Cloud optimised" is a specific internal layout: the image is tiled rather than striped, overviews are present, and the metadata sits at the front of the file. A plain GeoTIFF works and is slow, because reading one tile means reading whole strips.

from rio_cogeo.cogeo import cog_validate

valid, errors, warnings = cog_validate("scene.tif")
print("valid COG:", valid)
for message in errors + warnings:
    print("  ", message)

valid COG: False
   The file is greater than 512xH or 512xW, but is not tiled
   The file is greater than 512xH or 512xW, it is recommended to include internal overviews

Converting is one command:

from rio_cogeo.cogeo import cog_translate
from rio_cogeo.profiles import cog_profiles

cog_translate("scene.tif", "scene_cog.tif", cog_profiles.get("deflate"),
              overview_level=6, web_optimized=True, quiet=True)
print("valid COG:", cog_validate("scene_cog.tif")[0])

valid COG: True

web_optimized=True aligns the internal tiling with the Web Mercator grid, so a map tile maps onto whole internal blocks rather than straddling them.

2. Inspect the raster before styling it

Half the "all my tiles are black" reports are a rescale problem, and the fix starts with knowing the data's actual range:

GET /cog/info?url=…          bands, dtype, nodata, overviews, bounds
GET /cog/statistics?url=…    min, max, percentiles per band

Use the 2nd and 98th percentiles rather than the true minimum and maximum: a single hot pixel otherwise compresses the whole visible range.

3. Set `rescale` explicitly

TiTiler maps the source values onto 0–255 for display. For anything that is not already 8-bit — a DEM in metres, Sentinel-2 reflectance in scaled integers, a float index — the default is almost certainly wrong.

&rescale=0,3000            one range for all bands
&rescale=0,3000&rescale=0,3000&rescale=0,3000    per band, in order

Everything black usually means the range is far too wide; everything white means it is far too narrow. On the Snowdonia DEM (Example 4), rescale=0,100000 left 86.8% of a zoom-10 tile pure black and the rest nearly so, a mean grey level of 0.1 out of 255; rescale=0,1 turned every land pixel, 35.4% of the tile, pure white.

4. Choose bands and a colormap

&bidx=4&bidx=3&bidx=2                     false colour from a multiband scene
&expression=(b8-b4)/(b8%2Bb4)             NDVI, computed per tile
&colormap_name=viridis&rescale=-1,1

expression is evaluated per tile on the fly, which is what makes a dynamic tiler genuinely different from a pyramid: an index can be recomputed with a different formula without regenerating anything. Note that + must be URL-encoded as %2B.

5. Cache, because tile generation is not free

Every tile is a range read plus a decode plus a resample plus an encode. It is fast — 21.4 to 23.9 ms per tile from the local COG in Example 4 — and it is not free, and identical tiles are requested repeatedly. A map view is a dozen or more tiles, each one rendered again for every visitor unless something stores it. Put a CDN or a caching proxy in front, and set Cache-Control on the responses, as Example 1's middleware does.

6. Mind where the COG actually lives

A tile served from a COG on the same machine is a local read. The same COG on object storage in another region is a range request per tile, with that latency on every uncached tile.

Both work. The second needs a cache much more urgently, and benefits from GDAL's block cache being sized sensibly:

export GDAL_CACHEMAX=512                 # MB
export VSI_CACHE=TRUE
export VSI_CACHE_SIZE=536870912
export GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR

That last one prevents GDAL listing the whole bucket prefix every time it opens a file, which on object storage halved the time of a cold open in Example 4.

Everything black means the range is far too wide; everything white, too narrow.

Serving Snowdonia's elevation with TiTiler

The data used on this page

The Copernicus DEM GLO-30 tile N53 W005 — Snowdonia, Anglesey and the Llŷn peninsula at 1 arc-second, published on AWS as a COG — downloaded as dem.tif. The setup also writes a plain, striped GeoTIFF copy without overviews, scene.tif, as the kind of file step 1 is about.

import os
import urllib.request

import rasterio

URL = ("https://copernicus-dem-30m.s3.amazonaws.com/Copernicus_DSM_COG_10_N53_00_W005_00_DEM/"
       "Copernicus_DSM_COG_10_N53_00_W005_00_DEM.tif")
if not os.path.exists("dem.tif"):
    urllib.request.urlretrieve(URL, "dem.tif")
with rasterio.open("dem.tif") as src:
    data, profile = src.read(), src.profile
    print(f"dem.tif: {src.width} x {src.height} {src.dtypes[0]}, blocks {src.block_shapes[0]}, "
          f"overviews {src.overviews(1)}, nodata {src.nodata}")
    print(f"  exactly 0 m (sea, not flagged): {(data == 0).mean():.1%}, highest {data.max():,.1f} m")
profile.update(tiled=False, compress="deflate")
for key in ("blockxsize", "blockysize"):
    profile.pop(key, None)
with rasterio.open("scene.tif", "w", **profile) as dst:
    dst.write(data)
with rasterio.open("scene.tif") as src:
    print(f"scene.tif: blocks {src.block_shapes[0]}, overviews {src.overviews(1)}")

dem.tif: 2400 x 3600 float32, blocks (1024, 1024), overviews [2, 4, 8], nodata None
  exactly 0 m (sea, not flagged): 83.4%, highest 1,074.5 m
scene.tif: blocks (1, 2400), overviews []

Example 1 — a TiTiler application with your own defaults and caching

from fastapi import FastAPI, Request, Response
from titiler.core.factory import TilerFactory
from titiler.core.errors import DEFAULT_STATUS_CODES, add_exception_handlers

app = FastAPI(title="Raster tiles")
cog = TilerFactory()
app.include_router(cog.router, prefix="/cog", tags=["Cloud Optimised GeoTIFF"])
add_exception_handlers(app, DEFAULT_STATUS_CODES)


@app.middleware("http")
async def cache_tiles(request: Request, call_next):
    response = await call_next(request)
    if "/tiles/" in request.url.path and response.status_code == 200:
        response.headers["Cache-Control"] = "public, max-age=86400"
    return response


@app.get("/health")
def health():
    return {"status": "ok"}

To try it, run the app in a background thread of the same Python process — in production it is uvicorn module:app:

import threading
import time

import httpx
import uvicorn

server = uvicorn.Server(uvicorn.Config(app, host="127.0.0.1", port=8765, log_level="warning"))
threading.Thread(target=server.run, daemon=True).start()
while not server.started:
    time.sleep(0.05)
BASE = "http://127.0.0.1:8765"
r = httpx.get(f"{BASE}/cog/tiles/WebMercatorQuad/10/500/333.png",
              params={"url": "dem.tif", "rescale": "0,1000", "colormap_name": "terrain"})
print(httpx.get(f"{BASE}/health").json(), r.status_code, r.headers["content-type"],
      f"{len(r.content):,} bytes", "Cache-Control:", r.headers.get("cache-control"))

{'status': 'ok'} 200 image/png 39,620 bytes Cache-Control: public, max-age=86400

Wrapping the factory rather than running titiler.application directly is what lets you add auth, restrict which URLs may be tiled, and set caching policy — all of which a public deployment needs.

Example 2 — restricting which rasters may be served

from fastapi import HTTPException
from urllib.parse import urlparse

ALLOWED_HOSTS = {"data.example.org", "s3.eu-west-2.amazonaws.com"}
ALLOWED_PREFIXES = ("s3://our-bucket/", "/data/cogs/")


def validate_url(url: str) -> str:
    """A tile server that accepts any URL is an open proxy."""
    if url.startswith(ALLOWED_PREFIXES):
        return url
    parsed = urlparse(url)
    if parsed.scheme in ("http", "https") and parsed.hostname in ALLOWED_HOSTS:
        return url
    raise HTTPException(403, f"rasters from {parsed.hostname or url!r} are not served here")


for url in ("https://data.example.org/dem.tif", "/data/cogs/dem.tif", "http://169.254.169.254/latest/meta-data/"):
    try:
        print("served:", validate_url(url))
    except HTTPException as e:
        print("refused:", e.status_code, e.detail)

served: https://data.example.org/dem.tif
served: /data/cogs/dem.tif
refused: 403 rasters from '169.254.169.254' are not served here

This matters more than it sounds. A ?url= parameter that accepts anything lets a stranger use your server to fetch arbitrary URLs, including internal ones.

Example 3 — picking the rescale range from the data

import httpx


def suggest_rescale(base_url, cog_url, percentiles=(2, 98), **params):
    """Ask the server for statistics and produce the rescale parameter."""
    stats = httpx.get(f"{base_url}/cog/statistics",
                      params={"url": cog_url, **params}, timeout=60).json()

    ranges = []
    for band, values in stats.items():
        low = values.get(f"percentile_{percentiles[0]}", values["min"])
        high = values.get(f"percentile_{percentiles[1]}", values["max"])
        ranges.append((low, high))
        print(f"{band}: min {values['min']:.1f}  "
              f"p{percentiles[0]} {low:.1f}  p{percentiles[1]} {high:.1f}  "
              f"max {values['max']:.1f}")

    return "&".join(f"rescale={low:.0f},{high:.0f}" for low, high in ranges)


print("suggested:", suggest_rescale(BASE, "dem.tif"))
print("suggested:", suggest_rescale(BASE, "dem.tif", nodata=0))

b1: min -0.6  p2 0.0  p98 335.3  max 1056.9
suggested: rescale=0,335
b1: min -0.6  p2 1.1  p98 641.4  max 1056.9
suggested: rescale=1,641

The first suggestion shows what percentiles cannot fix. 83.4% of this DEM is exactly 0 m — the sea, stored as a value rather than flagged as nodata — so the 2nd percentile is sea level and the 98th is 335.3 m, and everything above it, most of Snowdonia's high ground, would render as flat white. Passing nodata=0 gives land-only percentiles and rescale=1,641. The statistics come from a downsampled read, which is why the maximum is 1,056.9 m rather than the file's 1,074.5 m; that is fine for choosing a stretch.

Example 4 — what the page's other claims measure

import statistics
import subprocess
import sys

import numpy as np
from PIL import Image
import io


def tile_ms(path, z, x, y, repeats=5, **params):
    times = []
    for _ in range(repeats):
        start = time.perf_counter()
        r = httpx.get(f"{BASE}/cog/tiles/WebMercatorQuad/{z}/{x}/{y}.png",
                      params={"url": path, "rescale": "0,1000", **params}, timeout=120)
        times.append((time.perf_counter() - start) * 1000)
    return statistics.median(times), r


# 1. COG against the striped copy, at low and high zoom (tiles over Snowdon)
for label, (z, x, y) in {"z7, whole tile": (7, 62, 41), "z12, Snowdon": (12, 2001, 1332)}.items():
    cog, _ = tile_ms("dem.tif", z, x, y)
    plain, _ = tile_ms("scene.tif", z, x, y)
    print(f"{label:15} COG {cog:7.1f} ms   striped, no overviews {plain:7.1f} ms")

# 2. rescale: far too wide, sensible, far too narrow (mean grey level of a z10 tile, 0-255)
for rescale in ("0,100000", "0,1000", "0,1"):
    _, r = tile_ms("dem.tif", 10, 500, 333, repeats=1, rescale=rescale)
    grey = np.asarray(Image.open(io.BytesIO(r.content)).convert("L"))
    print(f"rescale={rescale:9} mean grey {grey.mean():5.1f}, black {(grey == 0).mean():6.1%}, "
          f"white {(grey == 255).mean():6.1%}")

# 3. An unencoded + in an expression
raw = f"{BASE}/cog/tiles/WebMercatorQuad/10/500/333.png?url=dem.tif&rescale=0,1000&expression="
for label, expr in (("b1+100", "b1+100"), ("b1%2B100", "b1%2B100")):
    print(f"expression={label:9} -> HTTP {httpx.get(raw + expr).status_code}")

# 4. GDAL_DISABLE_READDIR_ON_OPEN, one open per fresh process (GDAL caches listings)
probe = """import sys, time, rasterio
env = {"GDAL_DISABLE_READDIR_ON_OPEN": "EMPTY_DIR"} if sys.argv[1] == "empty" else {}
href = "https://sentinel-cogs.s3.us-west-2.amazonaws.com/sentinel-s2-l2a-cogs/30/U/VD/2025/7/S2C_30UVD_20250712_0_L2A/B0%s.tif"
with rasterio.Env(**env):
    t = time.perf_counter()
    with rasterio.open(href % sys.argv[2]) as ds:
        ds.read(1, window=((0, 256), (0, 256)))
    print(f"{time.perf_counter() - t:.2f}")"""
for mode in ("default", "empty"):
    seconds = [float(subprocess.run([sys.executable, "-c", probe, mode, band], capture_output=True, text=True).stdout)
               for band in "234"]
    print(f"open + one block from S3, readdir {mode:7}: {', '.join(f'{t:.2f}' for t in seconds)} s")

z7, whole tile  COG    21.4 ms   striped, no overviews    63.6 ms
z12, Snowdon    COG    23.9 ms   striped, no overviews    15.4 ms
rescale=0,100000  mean grey   0.1, black  86.8%, white   0.0%
rescale=0,1000    mean grey  32.1, black  64.6%, white   0.0%
rescale=0,1       mean grey  90.3, black  64.6%, white  35.4%
expression=b1+100    -> HTTP 400
expression=b1%2B100  -> HTTP 200
open + one block from S3, readdir default: 3.98, 4.19, 4.15 s
open + one block from S3, readdir empty  : 2.17, 2.14, 2.10 s

Why a COG makes dynamic tiling possible

Why a COG makes dynamic raster tiling possible

Measured in Example 4. At high zoom both read a small window; at low zoom only the COG has a smaller copy to read.

An ordinary GeoTIFF is stored in strips, so reading a small square means reading every strip it touches — potentially the full width of the image. Internal tiling stores the file in blocks, so a reader can fetch one block.

Overviews add the second half: a low-zoom tile is rendered from a downsampled overview, not from the full-resolution data. Without overviews, a low-zoom tile reads and downsamples the entire image. On this 2,400 × 3,600 DEM the striped copy took 63.6 ms for the zoom-7 tile that covers all of it, against 21.4 ms for the COG, and was no slower at zoom 12 (15.4 against 23.9 ms). The low-zoom penalty grows with the image: a full Sentinel-2 band is 10,980 pixels square, 14 times as many pixels.

Why `rescale` is the parameter that decides whether anything appears

Display is 8-bit; source data usually is not. A DEM in metres can run from about −430 to 8,849; Sentinel-2 reflectance is scaled integers to 10,000; an NDVI is −1 to 1.

TiTiler cannot know which part of that range is interesting, so it needs telling. Set it from the 2nd and 98th percentiles of the actual data, not from the theoretical range, or one extreme pixel compresses everything else into a few grey levels.

Why `GDAL_DISABLE_READDIR_ON_OPEN` matters so much on object storage

By default, when GDAL opens a file it lists the containing directory to find sidecar files. On a local disk that is instant. On object storage with a prefix containing thousands of objects, it is a listing request per open — and a tile server opens the file constantly.

Setting it to EMPTY_DIR tells GDAL not to look. Opening a Sentinel-2 COG on S3 and reading one block took 3.98–4.19 s from here with the default and 2.10–2.17 s with EMPTY_DIR (Example 4): the listing was half the cost of every cold open.

Why dynamic raster tiles and pre-rendered ones are not the same product

A pyramid is fixed: one styling, one band combination, one stretch, decided when it was cut. A dynamic tiler takes the styling as query parameters, so the same COG serves true colour, false colour, an NDVI and three different stretches without any regeneration.

That flexibility is the reason to accept the per-request cost. If only one fixed rendering is ever needed, cutting a pyramid once is cheaper to operate — and the decision should be made on that basis rather than on which is more modern.

GDAL lists the containing directory on open unless told not to.

Validation, percentiles and URL allowlists

Validate the COG. A plain GeoTIFF works and is slow, especially at low zoom.
web_optimized=True when converting aligns internal blocks with the tile grid.
Percentiles, not min/max, for the rescale range.
+ in an expression must be %2B.
Restrict the url parameter, or the service is an open proxy.
Set GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR for object storage.
GDAL_CACHEMAX is per process, so it multiplies by the worker count.
Nodata must be declared in the file, or it renders as data.

Internal links

Cloud optimised GeoTIFF explained — why the format matters here
How to write a cloud optimised GeoTIFF — converting the source
How to serve raster tiles from a COG — the tiling pipeline
Dynamic tile server or static tiles: what to serve — the underlying decision
Fixing a COG read that is slow — the storage-side problems
How to set cache headers on a spatial API and tile service — the caching layer
Fixing a tile endpoint that returns 404 or blank tiles — when nothing appears
Fixing rasterio output that is black or empty — the same rescale problem offline

FAQ

What is TiTiler?

A FastAPI application built on rasterio and rio-tiler that serves XYZ tiles rendered on demand from a Cloud Optimised GeoTIFF, with the band selection, stretch and colormap as query parameters.

Why are all my tiles black?

Almost always the rescale range. Ask the server for /cog/statistics and set it from the 2nd and 98th percentiles rather than the full data range.

Does the source have to be a COG?

No, but it should be. A plain GeoTIFF has no internal tiling and no overviews, so a low-zoom tile reads and downsamples the whole image.

Can I compute an index on the fly?

Yes — expression=(b8-b4)/(b8%2Bb4) computes NDVI per tile. That flexibility is the main reason to serve dynamically rather than cutting a pyramid.

Is it fast enough for a public map?

With a cache in front, yes. Each tile took 21–24 ms to render from a local COG here, paid again on every request that nothing caches; a CDN answers repeats without rendering at all.

What should I set for object storage?

GDAL_DISABLE_READDIR_ON_OPEN=EMPTY_DIR first — it stops GDAL listing the bucket prefix on every open — then GDAL_CACHEMAX and the VSI cache settings.