Tiles Module — Web-tile Basemaps

The cleopatra.basemap.tiles module adds an optional, pure-Python web-tile basemap helper: add_tiles fetches XYZ map tiles covering an axes' current extent, stitches them with Pillow, and renders the composite underneath your data. No GDAL is required.

For a service published as OGC WMS or WMTS rather than an XYZ template, pass one of the provider objects from cleopatra.basemap.ogc as the source — the fetch, stitch and compose path below is the same one they use.

It is gated behind the cleopatra[tiles] optional extra (pillow, pyproj, xyzservices):

pipconda

pip install "cleopatra[tiles]"

conda install -c conda-forge cleopatra-tiles

If the extra is not installed, the functions raise a clear ImportError with the install hint.

Usage

import matplotlib
matplotlib.use("Agg")  # any backend; Agg shown for headless rendering
import matplotlib.pyplot as plt

from cleopatra.basemap.tiles import add_tiles

fig, ax = plt.subplots()
# plot something in Web Mercator (EPSG:3857) coordinates ...
ax.plot([1_000_000.0, 1_200_000.0], [6_000_000.0, 6_200_000.0])

# ... and drop an OpenStreetMap basemap underneath it
add_tiles(ax, crs=3857)

# a different provider, a fixed zoom, a custom User-Agent (recommended in production):
add_tiles(
    ax,
    source="CartoDB.Positron",
    crs="EPSG:4326",
    zoom=8,
    user_agent="my-app/1.0 (+https://example.org)",
)

fig.savefig("map.png")

Note

add_tiles reads the axes' current xlim/ylim, so plot your data first. When the data CRS is Web Mercator the tiles are placed in-place; for any other crs= the mosaic's own Web-Mercator coverage is reprojected into the target CRS and used as the image extent (the axis limits stay at the data bounds), so the basemap aligns with the data even though the fetched tiles cover a tile-snapped area larger than it. A residual Mercator-vs-linear nonlinearity remains for very large extents — for pixel-accurate results reproject the source data to EPSG:3857 before plotting. If a coarse mosaic overflows a limited-domain target CRS (e.g. a whole-world mosaic into a UTM zone), the reprojection is skipped with a warning and the basemap falls back to the data bounds (slightly misaligned); use a higher zoom or reproject the data to EPSG:3857 to avoid it. The automatic zoom uses a min_tiles_across floor (default 2) so a mid-range extent is not rendered from one or two coarse tiles; the number of tiles is capped by max_tiles (default MAX_TILES = 256), and the zoom is stepped down if a level would need more.

Module Documentation

cleopatra.basemap.tiles

Web-tile basemap helper for matplotlib axes.

Provides add_tiles -- a single entry point that fetches XYZ web tiles for the current axes extent, stitches them into a composite image with Pillow, and renders the image underneath the existing data layer.

The implementation is a pure-Python port of the pyramids.basemap module (basemap.py + tiles.py). It supports any XYZ provider listed in xyzservices. CRS handling is done with pyproj -- there is no GDAL dependency, so the module is safe to use in environments that only have matplotlib + numpy installed.

Notes

For data in CRSes other than Web Mercator (EPSG:3857) the stitched tile image is placed at the mosaic's own coverage: its Web-Mercator bounds are reprojected (with edge densification) into the target CRS and used as the imshow extent, while the axis limits stay at the data bounds. This aligns the basemap with the data even when the fetched tiles cover a tile-snapped area larger than the data. A residual Mercator-vs-linear nonlinearity remains for very large extents (the Mercator pixels are placed on a linear axis); if pixel-accurate warping is required, reproject the source data to Web Mercator (EPSG:3857) before plotting.

Examples:

Add a default OpenStreetMap basemap to an axes that already has data plotted in Web Mercator coordinates:

>>> import matplotlib
>>> matplotlib.use("Agg")
>>> import matplotlib.pyplot as plt
>>> fig, ax = plt.subplots()
>>> _ = ax.plot([1000000.0, 1200000.0], [6000000.0, 6200000.0])
>>> _ = add_tiles(ax, source=None, crs=3857)

Tile

Bases: NamedTuple

An XYZ web-map tile: column x, row y, at zoom level z.

The standard "slippy map" tile-coordinate triple used by every XYZ tile provider (OpenStreetMap, CartoDB, Esri, ...): at zoom z the world is divided into a 2**z by 2**z grid, x counted west to east and y north to south. Hashable and immutable, so it doubles as a dict key (see fetch_tiles's tile -> PNG bytes mapping).

Attributes:

Name	Type	Description
x	int	Column index, 0 (west edge) to 2**z - 1 (east edge).
y	int	Row index, 0 (north edge) to 2**z - 1 (south edge).
z	int	Zoom level; the grid is 2**z tiles wide and tall.

Examples:

• The single tile covering the whole world at zoom 0:

  >>> from cleopatra.basemap.tiles import Tile
  >>> Tile(0, 0, 0)
  Tile(x=0, y=0, z=0)
  

• Two tiles compare equal by value, so one can look the other up in a {tile: data} mapping (as fetch_tiles's return value does):

  >>> from cleopatra.basemap.tiles import Tile
  >>> tile_data = {Tile(548, 335, 10): b"...png bytes..."}
  >>> tile_data[Tile(x=548, y=335, z=10)]
  b'...png bytes...'
  

• Fields are accessible by name or by position:

  >>> from cleopatra.basemap.tiles import Tile
  >>> tile = Tile(548, 335, 10)
  >>> tile.z
  10
  >>> tile[:2]
  (548, 335)
  

Source code in src/cleopatra/basemap/tiles.py

class Tile(NamedTuple):
    """An XYZ web-map tile: column `x`, row `y`, at zoom level `z`.

    The standard "slippy map" tile-coordinate triple used by every XYZ tile
    provider (OpenStreetMap, CartoDB, Esri, ...): at zoom `z` the world is
    divided into a `2**z` by `2**z` grid, `x` counted west to east and `y`
    north to south. Hashable and immutable, so it doubles as a dict key
    (see `fetch_tiles`'s `tile -> PNG bytes` mapping).

    Attributes:
        x: Column index, `0` (west edge) to `2**z - 1` (east edge).
        y: Row index, `0` (north edge) to `2**z - 1` (south edge).
        z: Zoom level; the grid is `2**z` tiles wide and tall.

    Examples:
        - The single tile covering the whole world at zoom 0:
            ```python
            >>> from cleopatra.basemap.tiles import Tile
            >>> Tile(0, 0, 0)
            Tile(x=0, y=0, z=0)

            ```
        - Two tiles compare equal by value, so one can look the other up
            in a `{tile: data}` mapping (as `fetch_tiles`'s return value does):
            ```python
            >>> from cleopatra.basemap.tiles import Tile
            >>> tile_data = {Tile(548, 335, 10): b"...png bytes..."}
            >>> tile_data[Tile(x=548, y=335, z=10)]
            b'...png bytes...'

            ```
        - Fields are accessible by name or by position:
            ```python
            >>> from cleopatra.basemap.tiles import Tile
            >>> tile = Tile(548, 335, 10)
            >>> tile.z
            10
            >>> tile[:2]
            (548, 335)

            ```
    """

    x: int
    y: int
    z: int

add_tiles(ax, source=None, *, crs=None, zoom='auto', alpha=1.0, attribution=True, zorder=-1, interpolation='bilinear', timeout=10, retries=2, user_agent=None, max_tiles=MAX_TILES, min_tiles_across=2)

Overlay a web-tile basemap on a matplotlib axes.

Fetches XYZ web tiles that cover the axes' current extent, stitches them into a single composite image, and renders the image below the existing data layer. When the data is already in Web Mercator (EPSG:3857) the tiles are placed in-place; for any other CRS the mosaic's own Web-Mercator coverage is reprojected into the target CRS and used as the image extent (the axis limits stay at the data bounds), so the basemap aligns with the data even though the fetched tiles cover a tile-snapped area larger than it.

Parameters:

Name	Type	Description	Default
ax	Any	Matplotlib matplotlib.axes.Axes to add the basemap to. Data must already be plotted so the axis limits define the geographic extent.	required
source	Any | None	Tile provider. None defaults to OpenStreetMap.Mapnik. A dot-separated string such as "CartoDB.Positron" is resolved via get_provider. An xyzservices.TileProvider is used directly, and so is any object exposing build_url(x=, y=, z=) -- which is how cleopatra.basemap.ogc.WMSProvider and WMTSProvider bring OGC services through this same path.	None
crs	int | str | None	CRS of the data on ax. An integer is interpreted as an EPSG code; a string is passed through ("EPSG:XXXX" or WKT). None is treated as EPSG:3857.	None
zoom	int | str	Tile zoom level. "auto" derives a level from the axes extent. Integers must be in [0, 19].	'auto'
alpha	float	Opacity of the basemap (0.0--1.0).	1.0
attribution	str | bool	True adds the provider's attribution text, a string overrides it, False skips it.	True
zorder	int	Matplotlib zorder for the basemap (-1 puts it behind all data).	-1
interpolation	str	Interpolation method passed to ax.imshow.	'bilinear'
timeout	int	Per-tile HTTP timeout in seconds.	10
retries	int	Per-tile retry count.	2
user_agent	str | None	User-Agent header to send on tile requests. None (default) uses USER_AGENT (cleopatra/<version> (+repo-url)). Pass your own string when embedding cleopatra in an application so the traffic is attributed to that app (recommended for production use, and required by some providers' usage policies).	None
max_tiles	int	Cap on how many tiles to fetch. If the chosen zoom would need more than this, the zoom is stepped down until the count fits (or reaches 0). Defaults to MAX_TILES (256). Must be a positive int.	MAX_TILES
min_tiles_across	int	Floor for the automatic zoom, forwarded to auto_zoom when zoom="auto" (ignored for an explicit zoom=). Higher values give a sharper basemap at the cost of more tiles. Defaults to 2. See auto_zoom.	2

Returns:

Type	Description
Any	matplotlib.axes.Axes: The same axes, for chaining.

Raises:

Type	Description
ImportError	If the [tiles] extra is not installed.
TypeError	If ax is not a matplotlib Axes.
ValueError	If the axes have no data extent or zoom is invalid.
ConnectionError	If tiles cannot be fetched from the provider.

Examples:

Add a default OpenStreetMap basemap to an existing plot:

>>> import matplotlib
>>> matplotlib.use("Agg")
>>> import matplotlib.pyplot as plt
>>> fig, ax = plt.subplots()
>>> _ = ax.plot([1000000.0, 1200000.0], [6000000.0, 6200000.0])
>>> _ = add_tiles(ax, crs=3857)

Source code in src/cleopatra/basemap/tiles.py

def add_tiles(
    ax: Any,
    source: Any | None = None,
    *,
    crs: int | str | None = None,
    zoom: int | str = "auto",
    alpha: float = 1.0,
    attribution: str | bool = True,
    zorder: int = -1,
    interpolation: str = "bilinear",
    timeout: int = 10,
    retries: int = 2,
    user_agent: str | None = None,
    max_tiles: int = MAX_TILES,
    min_tiles_across: int = 2,
) -> Any:
    """Overlay a web-tile basemap on a matplotlib axes.

    Fetches XYZ web tiles that cover the axes' current extent, stitches
    them into a single composite image, and renders the image below the
    existing data layer. When the data is already in Web Mercator
    (EPSG:3857) the tiles are placed in-place; for any other CRS the
    mosaic's own Web-Mercator coverage is reprojected into the target CRS
    and used as the image extent (the axis limits stay at the data
    bounds), so the basemap aligns with the data even though the fetched
    tiles cover a tile-snapped area larger than it.

    Args:
        ax: Matplotlib `matplotlib.axes.Axes` to add the
            basemap to. Data must already be plotted so the axis
            limits define the geographic extent.
        source: Tile provider. `None` defaults to
            `OpenStreetMap.Mapnik`. A dot-separated string such as
            `"CartoDB.Positron"` is resolved via
            `get_provider`. An `xyzservices.TileProvider` is
            used directly, and so is any object exposing
            `build_url(x=, y=, z=)` -- which is how
            `cleopatra.basemap.ogc.WMSProvider` and `WMTSProvider` bring OGC
            services through this same path.
        crs: CRS of the data on `ax`. An integer is interpreted as an
            EPSG code; a string is passed through (`"EPSG:XXXX"` or
            WKT). `None` is treated as EPSG:3857.
        zoom: Tile zoom level. `"auto"` derives a level from the
            axes extent. Integers must be in `[0, 19]`.
        alpha: Opacity of the basemap (`0.0`--`1.0`).
        attribution: `True` adds the provider's attribution text, a
            string overrides it, `False` skips it.
        zorder: Matplotlib zorder for the basemap (`-1` puts it
            behind all data).
        interpolation: Interpolation method passed to `ax.imshow`.
        timeout: Per-tile HTTP timeout in seconds.
        retries: Per-tile retry count.
        user_agent: `User-Agent` header to send on tile requests.
            `None` (default) uses `USER_AGENT`
            (`cleopatra/<version> (+repo-url)`). Pass your own string
            when embedding cleopatra in an application so the traffic is
            attributed to *that* app (recommended for production use, and
            required by some providers' usage policies).
        max_tiles: Cap on how many tiles to fetch. If the chosen `zoom`
            would need more than this, the zoom is stepped down until the
            count fits (or reaches 0). Defaults to `MAX_TILES`
            (`256`). Must be a positive int.
        min_tiles_across: Floor for the automatic zoom, forwarded to
            `auto_zoom` when `zoom="auto"` (ignored for an explicit `zoom=`).
            Higher values give a sharper basemap at the cost of more tiles.
            Defaults to 2. See `auto_zoom`.

    Returns:
        matplotlib.axes.Axes: The same axes, for chaining.

    Raises:
        ImportError: If the `[tiles]` extra is not installed.
        TypeError: If `ax` is not a matplotlib Axes.
        ValueError: If the axes have no data extent or `zoom` is
            invalid.
        ConnectionError: If tiles cannot be fetched from the provider.

    Examples:
        Add a default OpenStreetMap basemap to an existing plot:

        >>> import matplotlib
        >>> matplotlib.use("Agg")
        >>> import matplotlib.pyplot as plt
        >>> fig, ax = plt.subplots()
        >>> _ = ax.plot([1000000.0, 1200000.0], [6000000.0, 6200000.0])
        >>> _ = add_tiles(ax, crs=3857)  # doctest: +SKIP
    """
    _require_tiles_extra()
    from pyproj import Transformer

    if not hasattr(ax, "get_xlim") or not hasattr(ax, "get_ylim"):
        raise TypeError(
            f"ax must be a matplotlib.axes.Axes instance, got {type(ax).__name__}"
        )

    if not isinstance(max_tiles, int) or isinstance(max_tiles, bool) or max_tiles < 1:
        raise ValueError(f"max_tiles must be a positive int, got {max_tiles!r}.")

    if (
        not isinstance(min_tiles_across, int)
        or isinstance(min_tiles_across, bool)
        or min_tiles_across < 1
    ):
        raise ValueError(
            f"min_tiles_across must be a positive int, got {min_tiles_across!r}."
        )

    x0, x1 = ax.get_xlim()
    y0, y1 = ax.get_ylim()
    west, east = min(x0, x1), max(x0, x1)
    south, north = min(y0, y1), max(y0, y1)

    if (west, east) == (0.0, 1.0) and (south, north) == (0.0, 1.0):
        raise ValueError("Axes have no data extent. Plot data before adding a basemap.")

    if west == east or south == north:
        raise ValueError(
            f"Axes have zero-area extent (west={west}, east={east}, "
            f"south={south}, north={north}). A basemap requires a "
            f"non-degenerate geographic extent."
        )

    if isinstance(source, str) or source is None:
        provider = get_provider(source)
    else:
        provider = source

    crs_value: int | str = 3857 if crs is None else crs
    crs_str = f"EPSG:{crs_value}" if isinstance(crs_value, int) else str(crs_value)
    is_3857 = crs_str == "EPSG:3857"

    if is_3857:
        w3857, s3857, e3857, n3857 = west, south, east, north
    else:
        try:
            w3857, s3857, e3857, n3857 = _densify_and_reproject_bounds(
                west, south, east, north, crs_str, "EPSG:3857"
            )
        except Exception as e:
            if "CRS" in type(e).__name__ or "Invalid" in str(e):
                raise ValueError(
                    f"Invalid CRS: {crs_str!r}. Provide a valid "
                    f"EPSG code or WKT string."
                ) from e
            raise

    transformer_to_4326 = Transformer.from_crs("EPSG:3857", "EPSG:4326", always_xy=True)
    w4326, s4326 = transformer_to_4326.transform(w3857, s3857)
    e4326, n4326 = transformer_to_4326.transform(e3857, n3857)

    if not all(np.isfinite(v) for v in (w4326, s4326, e4326, n4326)):
        raise ValueError(
            "Reprojection to EPSG:4326 produced infinite or NaN "
            "coordinates. The data extent may be outside the valid "
            "Web Mercator domain."
        )

    s4326 = max(s4326, -85.06)
    n4326 = min(n4326, 85.06)

    bounds_4326 = (w4326, s4326, e4326, n4326)

    if zoom == "auto":
        tile_zoom = auto_zoom(bounds_4326, min_tiles_across=min_tiles_across)
    else:
        try:
            tile_zoom = int(zoom)
        except (ValueError, TypeError) as e:
            raise ValueError(f"zoom must be 'auto' or int 0-19, got {zoom!r}") from e
        if not 0 <= tile_zoom <= 19:
            raise ValueError(f"zoom must be 0-19, got {tile_zoom}")

    original_zoom = tile_zoom
    tiles = _tiles_for_bbox(w4326, s4326, e4326, n4326, zoom=tile_zoom)

    while len(tiles) > max_tiles and tile_zoom > 0:
        tile_zoom -= 1
        tiles = _tiles_for_bbox(w4326, s4326, e4326, n4326, zoom=tile_zoom)

    if tile_zoom != original_zoom:
        logger.warning(
            "Zoom reduced from %d to %d (extent requires > %d tiles).",
            original_zoom,
            tile_zoom,
            max_tiles,
        )

    if not tiles:
        raise ValueError(
            f"No tiles found for bounds {bounds_4326} at zoom "
            f"{tile_zoom}. The extent may be outside valid tile "
            f"coverage."
        )

    tile_data = fetch_tiles(
        tiles,
        provider,
        timeout=timeout,
        retries=retries,
        user_agent=user_agent if user_agent is not None else USER_AGENT,
    )
    image, extent_3857 = stitch_tiles(tile_data, tiles, tile_zoom)

    if is_3857:
        extent = extent_3857
    else:
        try:
            extent = _densify_and_reproject_bounds(
                extent_3857[0],
                extent_3857[1],
                extent_3857[2],
                extent_3857[3],
                "EPSG:3857",
                crs_str,
            )
        except ValueError:
            logger.warning(
                "Tile mosaic bounds could not be reprojected from EPSG:3857 "
                "to %s (they overflow the target CRS domain); falling back to "
                "the data bounds. Use a higher zoom or reproject the data to "
                "EPSG:3857 for accurate tile placement.",
                crs_str,
            )
            extent = (west, south, east, north)

    xlim = ax.get_xlim()
    ylim = ax.get_ylim()

    ax.imshow(
        image,
        extent=[extent[0], extent[2], extent[1], extent[3]],
        interpolation=interpolation,
        alpha=alpha,
        zorder=zorder,
        aspect=ax.get_aspect(),
    )

    ax.set_xlim(xlim)
    ax.set_ylim(ylim)

    if attribution is True:
        raw = getattr(provider, "attribution", None) or getattr(
            provider, "html_attribution", ""
        )
        if raw:
            attr_text = html.unescape(re.sub(r"<[^>]+>", "", raw)).strip() or None
        else:
            attr_text = None
    elif isinstance(attribution, str):
        attr_text = attribution
    else:
        attr_text = None

    if attr_text:
        ax.text(
            0.99,
            0.01,
            attr_text,
            transform=ax.transAxes,
            fontsize=6,
            ha="right",
            va="bottom",
            alpha=0.7,
            bbox=dict(boxstyle="round,pad=0.2", fc="white", alpha=0.5),
        )

    return ax

auto_zoom(bounds_4326, min_tiles_across=2)

Compute a default zoom level for the given bounds in EPSG:4326.

Picks the smallest zoom at which the larger of the two extents spans at least min_tiles_across tiles, i.e. zoom = ceil(log2(min_tiles_across * 360 / max(lon_extent, lat_extent))), clamped to 0--19. The min_tiles_across floor (default 2) stops a mid-range regional extent from collapsing onto a single coarse tile stretched over the whole area (a 6--11 degree window would otherwise fetch just 2 tiles); min_tiles_across=1 reproduces the older one-tile-across heuristic.

This is a coarse heuristic that treats degrees of longitude and latitude as interchangeable; it does not account for Web Mercator's latitude distortion, so the result tends to be conservative (under-zoomed) for extents far from the equator. For high-latitude data, pass an explicit zoom= to add_tiles rather than relying on the auto value. The MAX_TILES cap in add_tiles will still step the zoom back down if the chosen level would require too many tiles.

Parameters:

Name	Type	Description	Default
bounds_4326	tuple[float, float, float, float]	(west, south, east, north) in EPSG:4326 degrees.	required
min_tiles_across	int	Minimum number of tiles the larger extent should span; higher values pick a sharper (higher) zoom. Values below 1 are clamped to 1 (the older one-tile-across heuristic). Defaults to 2.	2

Returns:

Name	Type	Description
int	int	Zoom level between 0 and 19.

Examples:

• Worldwide extent maps to zoom 1 (two tiles across the globe):

  >>> from cleopatra.basemap.tiles import auto_zoom
  >>> auto_zoom((-180, -85, 180, 85))
  1
  

• A 0.6 by 0.2 degree window over Berlin yields zoom 11:

  >>> from cleopatra.basemap.tiles import auto_zoom
  >>> auto_zoom((13.0, 52.4, 13.6, 52.6))
  11
  

• min_tiles_across=1 restores the older, coarser one-tile heuristic (worldwide -> zoom 0):

  >>> from cleopatra.basemap.tiles import auto_zoom
  >>> auto_zoom((-180, -85, 180, 85), min_tiles_across=1)
  0
  

• Tiny extents are clamped to the maximum zoom (19):

  >>> from cleopatra.basemap.tiles import auto_zoom
  >>> auto_zoom((0.0, 0.0, 1e-9, 1e-9))
  19
  

Source code in src/cleopatra/basemap/tiles.py

def auto_zoom(
    bounds_4326: tuple[float, float, float, float],
    min_tiles_across: int = 2,
) -> int:
    """Compute a default zoom level for the given bounds in EPSG:4326.

    Picks the smallest zoom at which the larger of the two extents spans at
    least `min_tiles_across` tiles, i.e.
    `zoom = ceil(log2(min_tiles_across * 360 / max(lon_extent, lat_extent)))`,
    clamped to 0--19. The `min_tiles_across` floor (default 2) stops a
    mid-range regional extent from collapsing onto a single coarse tile
    stretched over the whole area (a 6--11 degree window would otherwise
    fetch just 2 tiles); `min_tiles_across=1` reproduces the older
    one-tile-across heuristic.

    This is a coarse heuristic that treats degrees of longitude and
    latitude as interchangeable; it does **not** account for Web
    Mercator's latitude distortion, so the result tends to be
    conservative (under-zoomed) for extents far from the equator. For
    high-latitude data, pass an explicit `zoom=` to `add_tiles`
    rather than relying on the auto value. The `MAX_TILES` cap in
    `add_tiles` will still step the zoom back down if the chosen
    level would require too many tiles.

    Args:
        bounds_4326: `(west, south, east, north)` in EPSG:4326 degrees.
        min_tiles_across: Minimum number of tiles the larger extent should
            span; higher values pick a sharper (higher) zoom. Values below 1
            are clamped to 1 (the older one-tile-across heuristic). Defaults
            to 2.

    Returns:
        int: Zoom level between 0 and 19.

    Examples:
        - Worldwide extent maps to zoom 1 (two tiles across the globe):
            ```python
            >>> from cleopatra.basemap.tiles import auto_zoom
            >>> auto_zoom((-180, -85, 180, 85))
            1

            ```
        - A 0.6 by 0.2 degree window over Berlin yields zoom 11:
            ```python
            >>> from cleopatra.basemap.tiles import auto_zoom
            >>> auto_zoom((13.0, 52.4, 13.6, 52.6))
            11

            ```
        - `min_tiles_across=1` restores the older, coarser one-tile
            heuristic (worldwide -> zoom 0):
            ```python
            >>> from cleopatra.basemap.tiles import auto_zoom
            >>> auto_zoom((-180, -85, 180, 85), min_tiles_across=1)
            0

            ```
        - Tiny extents are clamped to the maximum zoom (19):
            ```python
            >>> from cleopatra.basemap.tiles import auto_zoom
            >>> auto_zoom((0.0, 0.0, 1e-9, 1e-9))
            19

            ```
    """
    west, south, east, north = bounds_4326
    lon_extent = abs(east - west)
    lat_extent = abs(north - south)
    max_extent = max(lon_extent, lat_extent, 1e-10)
    across = max(1, min_tiles_across)
    zoom = math.ceil(math.log2(across * 360.0 / max_extent))
    result = max(0, min(zoom, 19))
    return result

fetch_single_tile(tile, provider, timeout, retries, user_agent=USER_AGENT)

Fetch a single tile, retrying on transient failures.

Every failed attempt is logged at debug level with a redacted URL. A tile URL is not always safe to write to a log -- an XYZ template can embed an API key, and cleopatra.basemap.ogc documents extra_params as the place to put a token -- and a debug log outlives the session. Every parameter name survives, and so do the values of the OGC parameters that say which tile was asked for, so an OGC line still identifies the failing tile; every other value is replaced with ... (see _redact_url). An XYZ URL carries its tile in the path, which is not masked, so its line stays distinguishable too. The request that goes on the wire is the unredacted URL; only the log line is masked.

Parameters:

Name	Type	Description	Default
tile	Any	Tile to fetch (has x, y, z attributes).	required
provider	Any	xyzservices.TileProvider with a URL template.	required
timeout	int	HTTP request timeout in seconds.	required
retries	int	Number of retry attempts on failure.	required
user_agent	str	User-Agent header sent on every request. Defaults to USER_AGENT (cleopatra/<version> (+repo-url)).	USER_AGENT

Returns:

Type	Description
tuple[Any, bytes]	tuple[Any, bytes]: The original tile and its PNG/JPEG bytes.

Raises:

Type	Description
ValueError	If the provider's URL template is not an http(s) URL.
ConnectionError	If the tile cannot be fetched after all retries are exhausted.

Examples:

• Fetch a single OpenStreetMap tile (network-dependent, hence skipped under doctest):

  >>> from cleopatra.basemap.tiles import Tile, fetch_single_tile, get_provider
  >>> tile = Tile(0, 0, 0)
  >>> provider = get_provider("OpenStreetMap.Mapnik")
  >>> tile_obj, data = fetch_single_tile(  # doctest: +SKIP
  ...     tile, provider, timeout=10, retries=2
  ... )
  >>> from cleopatra.basemap.tiles import _looks_like_image
  >>> _looks_like_image(data)  # doctest: +SKIP
  True
  

• Tile failures raise ConnectionError after retries are exhausted:

  >>> from cleopatra.basemap.tiles import Tile, fetch_single_tile
  >>> from xyzservices import TileProvider
  >>> bad = TileProvider(
  ...     name="bad",
  ...     url="http://127.0.0.1:1/{z}/{x}/{y}.png",
  ...     attribution="",
  ... )
  >>> fetch_single_tile(  # doctest: +SKIP
  ...     Tile(0, 0, 0), bad, timeout=1, retries=0
  ... )
  Traceback (most recent call last):
      ...
  ConnectionError: Failed to fetch tile z=0/x=0/y=0 ...
  

Source code in src/cleopatra/basemap/tiles.py

def fetch_single_tile(
    tile: Any,
    provider: Any,
    timeout: int,
    retries: int,
    user_agent: str = USER_AGENT,
) -> tuple[Any, bytes]:
    """Fetch a single tile, retrying on transient failures.

    Every failed attempt is logged at debug level with a *redacted* URL. A tile
    URL is not always safe to write to a log -- an XYZ template can embed an API
    key, and `cleopatra.basemap.ogc` documents `extra_params` as the place to
    put a token -- and a debug log outlives the session. Every parameter name
    survives, and so do the values of the OGC parameters that say which tile was
    asked for, so an OGC line still identifies the failing tile; every other
    value is replaced with `...` (see `_redact_url`). An XYZ URL carries its
    tile in the path, which is not masked, so its line stays distinguishable
    too. The request that goes on the wire is the unredacted URL; only the log
    line is masked.

    Args:
        tile: Tile to fetch (has `x`, `y`, `z` attributes).
        provider: `xyzservices.TileProvider` with a URL template.
        timeout: HTTP request timeout in seconds.
        retries: Number of retry attempts on failure.
        user_agent: `User-Agent` header sent on every request. Defaults
            to `USER_AGENT` (`cleopatra/<version> (+repo-url)`).

    Returns:
        tuple[Any, bytes]: The original tile and its PNG/JPEG bytes.

    Raises:
        ValueError: If the provider's URL template is not an `http(s)` URL.
        ConnectionError: If the tile cannot be fetched after all
            retries are exhausted.

    Examples:
        - Fetch a single OpenStreetMap tile (network-dependent, hence
            skipped under doctest):
            ```python
            >>> from cleopatra.basemap.tiles import Tile, fetch_single_tile, get_provider
            >>> tile = Tile(0, 0, 0)
            >>> provider = get_provider("OpenStreetMap.Mapnik")
            >>> tile_obj, data = fetch_single_tile(  # doctest: +SKIP
            ...     tile, provider, timeout=10, retries=2
            ... )
            >>> from cleopatra.basemap.tiles import _looks_like_image
            >>> _looks_like_image(data)  # doctest: +SKIP
            True

            ```
        - Tile failures raise `ConnectionError` after retries
            are exhausted:
            ```python
            >>> from cleopatra.basemap.tiles import Tile, fetch_single_tile
            >>> from xyzservices import TileProvider
            >>> bad = TileProvider(
            ...     name="bad",
            ...     url="http://127.0.0.1:1/{z}/{x}/{y}.png",
            ...     attribution="",
            ... )
            >>> fetch_single_tile(  # doctest: +SKIP
            ...     Tile(0, 0, 0), bad, timeout=1, retries=0
            ... )
            Traceback (most recent call last):
                ...
            ConnectionError: Failed to fetch tile z=0/x=0/y=0 ...

            ```
    """
    url = provider.build_url(x=tile.x, y=tile.y, z=tile.z)
    if urllib.parse.urlsplit(url).scheme.lower() not in ("http", "https"):
        raise ValueError(f"Refusing to fetch non-http(s) tile URL: {url!r}")
    last_error: Exception | None = None
    result_bytes: bytes | None = None
    for attempt in range(retries + 1):
        try:
            request = urllib.request.Request(
                url,
                headers={"User-Agent": user_agent},
            )
            response = urlopen_http(request, timeout=timeout)
            data = response.read()
            if not _looks_like_image(data):
                raise OSError(
                    f"Tile response is not a valid image "
                    f"({len(data)} bytes, starts with "
                    f"{data[:8]!r})"
                )
            result_bytes = data
            break
        except (OSError, urllib.error.URLError, ConnectionError) as e:
            last_error = e
            # Redacted: an XYZ template can embed an API key, and an OGC
            # provider's `extra_params` is documented as the place for a token.
            logger.debug(
                "Tile fetch attempt %d/%d failed for %s: %s",
                attempt + 1,
                retries + 1,
                _redact_url(url),
                e,
            )
    if result_bytes is None:
        raise ConnectionError(
            f"Failed to fetch tile z={tile.z}/x={tile.x}/y={tile.y} "
            f"after {retries + 1} attempts: {last_error}"
        )
    return tile, result_bytes

fetch_tiles(tiles, provider, max_workers=8, timeout=10, retries=2, user_agent=USER_AGENT)

Fetch tile images in parallel over HTTP.

Uses concurrent.futures.ThreadPoolExecutor for parallel downloads. Each tile URL is constructed via the provider's build_url(). A User-Agent header (cleopatra/<version> (+repo-url) by default) is sent on every request so tile providers can attribute the traffic — OpenStreetMap's usage policy requires an identifiable agent.

Parameters:

Name	Type	Description	Default
tiles	list	Tiles to fetch (each has x, y, z attributes).	required
provider	Any	xyzservices.TileProvider with a URL template.	required
max_workers	int	Maximum concurrent HTTP connections.	8
timeout	int	Per-tile HTTP request timeout in seconds.	10
retries	int	Per-tile retry count on failure.	2
user_agent	str	User-Agent header sent on every request. Defaults to USER_AGENT.	USER_AGENT

Returns:

Name	Type	Description
dict	dict	Mapping of Tile to PNG/JPEG bytes.

Raises:

Type	Description
ConnectionError	If any tile cannot be fetched after all retries.

Examples:

• Fetch a small tile grid in parallel (network-dependent, hence skipped under doctest):

  >>> from cleopatra.basemap.tiles import _tiles_for_bbox, fetch_tiles, get_provider
  >>> tiles = _tiles_for_bbox(13.0, 52.4, 13.6, 52.6, 10)
  >>> provider = get_provider("OpenStreetMap.Mapnik")
  >>> data = fetch_tiles(tiles, provider, max_workers=4)  # doctest: +SKIP
  >>> len(data) == len(tiles)  # doctest: +SKIP
  True
  

• Pass an empty list to short-circuit and get an empty dict:

  >>> from cleopatra.basemap.tiles import fetch_tiles, get_provider
  >>> provider = get_provider("OpenStreetMap.Mapnik")
  >>> fetch_tiles([], provider)
  {}
  

Source code in src/cleopatra/basemap/tiles.py

def fetch_tiles(
    tiles: list,
    provider: Any,
    max_workers: int = 8,
    timeout: int = 10,
    retries: int = 2,
    user_agent: str = USER_AGENT,
) -> dict:
    """Fetch tile images in parallel over HTTP.

    Uses `concurrent.futures.ThreadPoolExecutor` for parallel
    downloads. Each tile URL is constructed via the provider's
    `build_url()`. A `User-Agent` header (`cleopatra/<version> (+repo-url)`
    by default) is sent on every request so tile providers can attribute
    the traffic — OpenStreetMap's usage policy requires an identifiable
    agent.

    Args:
        tiles: Tiles to fetch (each has `x`, `y`, `z` attributes).
        provider: `xyzservices.TileProvider` with a URL template.
        max_workers: Maximum concurrent HTTP connections.
        timeout: Per-tile HTTP request timeout in seconds.
        retries: Per-tile retry count on failure.
        user_agent: `User-Agent` header sent on every request. Defaults
            to `USER_AGENT`.

    Returns:
        dict: Mapping of Tile to PNG/JPEG bytes.

    Raises:
        ConnectionError: If any tile cannot be fetched after all
            retries.

    Examples:
        - Fetch a small tile grid in parallel (network-dependent, hence
            skipped under doctest):
            ```python
            >>> from cleopatra.basemap.tiles import _tiles_for_bbox, fetch_tiles, get_provider
            >>> tiles = _tiles_for_bbox(13.0, 52.4, 13.6, 52.6, 10)
            >>> provider = get_provider("OpenStreetMap.Mapnik")
            >>> data = fetch_tiles(tiles, provider, max_workers=4)  # doctest: +SKIP
            >>> len(data) == len(tiles)  # doctest: +SKIP
            True

            ```
        - Pass an empty list to short-circuit and get an empty dict:
            ```python
            >>> from cleopatra.basemap.tiles import fetch_tiles, get_provider
            >>> provider = get_provider("OpenStreetMap.Mapnik")
            >>> fetch_tiles([], provider)
            {}

            ```
    """
    tile_data: dict[Any, bytes] = {}
    with ThreadPoolExecutor(max_workers=max_workers) as executor:
        futures = {
            executor.submit(
                fetch_single_tile, tile, provider, timeout, retries, user_agent
            ): tile
            for tile in tiles
        }
        try:
            for future in as_completed(futures):
                tile_obj, png_bytes = future.result()
                tile_data[tile_obj] = png_bytes
        except ConnectionError:
            for f in futures:
                f.cancel()
            raise
        except Exception as e:
            for f in futures:
                f.cancel()
            logger.error("Unexpected error during tile fetching: %s", e)
            raise
    return tile_data

get_provider(name=None)

Resolve an XYZ tile provider by name.

Parameters:

Name	Type	Description	Default
name	str | None	Dot-separated provider name (e.g. "OpenStreetMap.Mapnik", "CartoDB.Positron", "Esri.WorldImagery"). None returns the default (OpenStreetMap.Mapnik).	None

Returns:

Type	Description
Any	xyzservices.TileProvider: The resolved tile provider with
Any	url template and attribution metadata.

Raises:

Type	Description
ImportError	If the [tiles] extra is not installed.
ValueError	If the provider name cannot be resolved.

Examples:

• Resolve the default OpenStreetMap provider and inspect its URL template:

  >>> from cleopatra.basemap.tiles import get_provider
  >>> provider = get_provider()
  >>> provider.name
  'OpenStreetMap.Mapnik'
  >>> "{z}" in provider.url and "{x}" in provider.url
  True
  

• Resolve a named provider via dot-path syntax:

  >>> from cleopatra.basemap.tiles import get_provider
  >>> provider = get_provider("CartoDB.Positron")
  >>> provider.name
  'CartoDB.Positron'
  

• Invalid provider names raise ValueError:

  >>> from cleopatra.basemap.tiles import get_provider
  >>> get_provider("NonExistent.Provider")
  Traceback (most recent call last):
      ...
  ValueError: Unknown tile provider: 'NonExistent.Provider'. Failed at 'NonExistent'. Use xyzservices.providers to list available providers.
  

Source code in src/cleopatra/basemap/tiles.py

def get_provider(name: str | None = None) -> Any:
    """Resolve an XYZ tile provider by name.

    Args:
        name: Dot-separated provider name (e.g. `"OpenStreetMap.Mapnik"`,
            `"CartoDB.Positron"`, `"Esri.WorldImagery"`). `None`
            returns the default (`OpenStreetMap.Mapnik`).

    Returns:
        xyzservices.TileProvider: The resolved tile provider with
        `url` template and `attribution` metadata.

    Raises:
        ImportError: If the `[tiles]` extra is not installed.
        ValueError: If the provider name cannot be resolved.

    Examples:
        - Resolve the default OpenStreetMap provider and inspect its
            URL template:
            ```python
            >>> from cleopatra.basemap.tiles import get_provider
            >>> provider = get_provider()
            >>> provider.name
            'OpenStreetMap.Mapnik'
            >>> "{z}" in provider.url and "{x}" in provider.url
            True

            ```
        - Resolve a named provider via dot-path syntax:
            ```python
            >>> from cleopatra.basemap.tiles import get_provider
            >>> provider = get_provider("CartoDB.Positron")
            >>> provider.name
            'CartoDB.Positron'

            ```
        - Invalid provider names raise `ValueError`:
            ```python
            >>> from cleopatra.basemap.tiles import get_provider
            >>> get_provider("NonExistent.Provider")
            Traceback (most recent call last):
                ...
            ValueError: Unknown tile provider: 'NonExistent.Provider'. Failed at 'NonExistent'. Use xyzservices.providers to list available providers.

            ```
    """
    _require_tiles_extra()
    import xyzservices.providers as xyz

    provider: Any
    if name is None:
        provider = xyz.OpenStreetMap.Mapnik
    else:
        parts = name.split(".")
        provider = xyz
        for part in parts:
            try:
                provider = provider[part]
            except (KeyError, TypeError) as e:
                raise ValueError(
                    f"Unknown tile provider: '{name}'. "
                    f"Failed at '{part}'. Use "
                    f"xyzservices.providers to list "
                    f"available providers."
                ) from e
    return provider

mercator_to_equirectangular(mosaic, bounds, n_lon=2880, n_lat=1440)

Resample a Web Mercator image onto an equirectangular lon/lat grid.

XYZ tiles -- and the mosaic stitch_tiles builds from them -- are in Web Mercator (EPSG:3857), where a degree of latitude occupies more pixels near the poles. A globe, or a plain lon/lat (EPSG:4326) axis, is textured in equirectangular longitude and latitude, so the mosaic must be resampled from its Mercator rows onto evenly spaced latitude rows.

Each output cell area-averages the whole contiguous block of source pixels that falls inside it, found from the cell's edges -- not a fixed number of point samples, which over-reads the compressed equatorial rows and skips source rows near the poles, exactly where Mercator's stretch makes aliasing worst. Latitudes beyond the Mercator limit (+/-85.051 deg) clamp onto the edge rows; a cell thinner than one source pixel keeps that single row.

Parameters:

Name	Type	Description	Default
mosaic	ndarray	The Web Mercator image, (H, W) or (H, W, C), north-up (row 0 = north), as returned by stitch_tiles.	required
bounds	tuple[float, float, float, float]	(west, south, east, north) of mosaic in EPSG:3857 metres, as returned by stitch_tiles. The output always spans the whole globe (-180..180, 90..-90); latitudes and longitudes outside the mosaic's coverage clamp onto its edge rows/columns rather than being cropped, so a partial mosaic is placed within a full frame.	required
n_lon	int	Width of the output grid, spanning -180..180 degrees.	2880
n_lat	int	Height of the output grid, spanning 90..-90 degrees.	1440

Returns:

Type	Description
ndarray	numpy.ndarray: The (n_lat, n_lon) or (n_lat, n_lon, C)
ndarray	equirectangular resample as float32, in the same value scale as
ndarray	mosaic (a uint8 0..255 mosaic yields 0..255 floats -- averaging is
ndarray	inherently fractional; world_texture normalises to 0..1).

Raises:

Type	Description
ValueError	If mosaic is not 2-D or 3-D, n_lon/n_lat are not positive, or bounds is degenerate (east <= west or north <= south).

Examples:

• A west-blue / east-red Mercator mosaic keeps its longitudinal order after resampling (Web Mercator x is linear in longitude):

  >>> import numpy as np
  >>> from cleopatra.basemap.tiles import mercator_to_equirectangular
  >>> m = np.zeros((4, 8, 3), dtype="uint8")
  >>> m[:, :4] = (0, 0, 255)
  >>> m[:, 4:] = (255, 0, 0)
  >>> world = (-20037508.34, -20037508.34, 20037508.34, 20037508.34)
  >>> tex = mercator_to_equirectangular(m, world, n_lon=8, n_lat=4)
  >>> tex.shape
  (4, 8, 3)
  >>> bool(tex[0, 0, 2] > tex[0, 0, 0])  # west stays blue
  True
  >>> bool(tex[0, -1, 0] > tex[0, -1, 2])  # east stays red
  True
  

Source code in src/cleopatra/basemap/tiles.py

def mercator_to_equirectangular(
    mosaic: np.ndarray,
    bounds: tuple[float, float, float, float],
    n_lon: int = 2880,
    n_lat: int = 1440,
) -> np.ndarray:
    """Resample a Web Mercator image onto an equirectangular lon/lat grid.

    XYZ tiles -- and the mosaic `stitch_tiles` builds from them -- are in Web
    Mercator (EPSG:3857), where a degree of latitude occupies more pixels near
    the poles. A globe, or a plain lon/lat (EPSG:4326) axis, is textured in
    equirectangular longitude and latitude, so the mosaic must be resampled
    from its Mercator rows onto evenly spaced latitude rows.

    Each output cell area-averages the whole contiguous block of source pixels
    that falls inside it, found from the cell's **edges** -- not a fixed number
    of point samples, which over-reads the compressed equatorial rows and skips
    source rows near the poles, exactly where Mercator's stretch makes aliasing
    worst. Latitudes beyond the Mercator limit (`+/-85.051` deg) clamp onto the
    edge rows; a cell thinner than one source pixel keeps that single row.

    Args:
        mosaic: The Web Mercator image, `(H, W)` or `(H, W, C)`, north-up
            (row 0 = north), as returned by `stitch_tiles`.
        bounds: `(west, south, east, north)` of `mosaic` in EPSG:3857 metres,
            as returned by `stitch_tiles`. The output always spans the whole
            globe (-180..180, 90..-90); latitudes and longitudes outside the
            mosaic's coverage clamp onto its edge rows/columns rather than
            being cropped, so a partial mosaic is placed within a full frame.
        n_lon: Width of the output grid, spanning -180..180 degrees.
        n_lat: Height of the output grid, spanning 90..-90 degrees.

    Returns:
        numpy.ndarray: The `(n_lat, n_lon)` or `(n_lat, n_lon, C)`
        equirectangular resample as `float32`, in the same value scale as
        `mosaic` (a `uint8` 0..255 mosaic yields 0..255 floats -- averaging is
        inherently fractional; `world_texture` normalises to 0..1).

    Raises:
        ValueError: If `mosaic` is not 2-D or 3-D, `n_lon`/`n_lat` are not
            positive, or `bounds` is degenerate (`east <= west` or
            `north <= south`).

    Examples:
        - A west-blue / east-red Mercator mosaic keeps its longitudinal order
            after resampling (Web Mercator x is linear in longitude):
            ```python
            >>> import numpy as np
            >>> from cleopatra.basemap.tiles import mercator_to_equirectangular
            >>> m = np.zeros((4, 8, 3), dtype="uint8")
            >>> m[:, :4] = (0, 0, 255)
            >>> m[:, 4:] = (255, 0, 0)
            >>> world = (-20037508.34, -20037508.34, 20037508.34, 20037508.34)
            >>> tex = mercator_to_equirectangular(m, world, n_lon=8, n_lat=4)
            >>> tex.shape
            (4, 8, 3)
            >>> bool(tex[0, 0, 2] > tex[0, 0, 0])  # west stays blue
            True
            >>> bool(tex[0, -1, 0] > tex[0, -1, 2])  # east stays red
            True

            ```
    """
    mosaic = np.asarray(mosaic)
    if mosaic.ndim not in (2, 3):
        raise ValueError(f"mosaic must be a 2-D or 3-D array, got {mosaic.ndim}-D")
    if n_lon < 1 or n_lat < 1:
        raise ValueError(
            f"n_lon and n_lat must be positive, got n_lon={n_lon}, n_lat={n_lat}"
        )
    west, south, east, north = (float(v) for v in bounds)
    if not (east > west and north > south):
        raise ValueError(
            "bounds must be (west, south, east, north) with east > west and "
            f"north > south, got {bounds!r}"
        )
    h, w = mosaic.shape[:2]

    lat_edges = np.clip(
        np.linspace(90.0, -90.0, n_lat + 1), -_MAX_LATITUDE, _MAX_LATITUDE
    )
    y_edges = np.log(np.tan(np.pi / 4 + np.deg2rad(lat_edges) / 2)) * _MERC_MAX / np.pi
    row_edges = np.clip(
        np.round((north - y_edges) / (north - south) * h).astype(int), 0, h
    )
    x_edges = np.linspace(-180.0, 180.0, n_lon + 1) / 180.0 * _MERC_MAX
    col_edges = np.clip(
        np.round((x_edges - west) / (east - west) * w).astype(int), 0, w
    )

    # reduceat sums each block band[start[i]:start[i+1]] (the last runs to the
    # array end). The block END is the next *capped* start, not the raw edge, so
    # size the mean's divisor from the actual block widths: otherwise the south
    # clamp (row_edges == h collapses start to h-1) sums one row fewer than it
    # divides by, darkening the row near -85 deg S. max(.,1) keeps a collapsed
    # polar row (next start == start) on its own.
    row_start = np.minimum(row_edges[:-1], h - 1)
    col_start = np.minimum(col_edges[:-1], w - 1)
    row_next = np.concatenate([row_start[1:], [h]])
    col_next = np.concatenate([col_start[1:], [w]])
    row_w = np.maximum(row_next - row_start, 1)[:, None].astype("float32")
    col_w = np.maximum(col_next - col_start, 1)[None, :].astype("float32")

    flat = mosaic.reshape(h, w, -1)  # (H, W, C), C == 1 for a 2-D mosaic
    bands = []
    for band_index in range(flat.shape[2]):
        # One band at a time: the float32 copy of a large mosaic is heavy, and
        # converting every band at once multiplies the transient peak.
        band = flat[..., band_index].astype("float32")
        band = np.add.reduceat(band, row_start, axis=0) / row_w
        bands.append(np.add.reduceat(band, col_start, axis=1) / col_w)
    grid = np.stack(bands, axis=-1)
    return grid[..., 0] if mosaic.ndim == 2 else grid

stitch_tiles(tile_data, tiles, zoom)

Stitch tile images into a single RGBA array.

Arranges tiles in a grid based on their x, y positions. The tile size is read from the first decoded image (typically 256 or 512 px). Computes the geographic extent of the stitched image in EPSG:3857 using _tile_xy_bounds on the corner tiles.

Parameters:

Name	Type	Description	Default
tile_data	dict	Mapping of Tile to PNG bytes (from fetch_tiles).	required
tiles	list	All tiles in the grid, defining grid dimensions.	required
zoom	int	Zoom level of the tiles.	required

Returns:

Type	Description
ndarray	tuple[numpy.ndarray, tuple[float, float, float, float]]: The
tuple[float, float, float, float]	stitched RGBA image with shape (H, W, 4) and dtype
tuple[ndarray, tuple[float, float, float, float]]	uint8, plus (west, south, east, north) in EPSG:3857
tuple[ndarray, tuple[float, float, float, float]]	meters.

Raises:

Type	Description
ValueError	If any tile bytes cannot be decoded as an image.

Examples:

• Stitch a single synthetic tile into a 256x256 RGBA image:

  >>> import io
  >>> from PIL import Image
  >>> from cleopatra.basemap.tiles import Tile, stitch_tiles
  >>> buf = io.BytesIO()
  >>> Image.new("RGBA", (256, 256), (255, 0, 0, 255)).save(buf, "PNG")
  >>> tile = Tile(0, 0, 0)
  >>> image, extent = stitch_tiles({tile: buf.getvalue()}, [tile], 0)
  >>> image.shape
  (256, 256, 4)
  >>> image.dtype.name
  'uint8'
  

• The returned EPSG:3857 extent comes from _tile_xy_bounds on the corner tiles:

  >>> import io
  >>> from PIL import Image
  >>> from cleopatra.basemap.tiles import Tile, stitch_tiles
  >>> buf = io.BytesIO()
  >>> Image.new("RGBA", (256, 256), (0, 255, 0, 255)).save(buf, "PNG")
  >>> tile = Tile(0, 0, 0)
  >>> _, (w, s, e, n) = stitch_tiles({tile: buf.getvalue()}, [tile], 0)
  >>> w < e and s < n
  True
  

• Invalid tile bytes raise ValueError:

  >>> from cleopatra.basemap.tiles import Tile, stitch_tiles
  >>> tile = Tile(0, 0, 0)
  >>> try:
  ...     stitch_tiles({tile: b"not an image"}, [tile], 0)
  ... except ValueError as exc:
  ...     str(exc).startswith("Failed to decode tile image:")
  True
  

Source code in src/cleopatra/basemap/tiles.py

def stitch_tiles(
    tile_data: dict,
    tiles: list,
    zoom: int,
) -> tuple[np.ndarray, tuple[float, float, float, float]]:
    """Stitch tile images into a single RGBA array.

    Arranges tiles in a grid based on their `x`, `y` positions. The
    tile size is read from the first decoded image (typically 256 or
    512 px). Computes the geographic extent of the stitched image in
    EPSG:3857 using `_tile_xy_bounds` on the corner tiles.

    Args:
        tile_data: Mapping of Tile to PNG bytes (from
            `fetch_tiles`).
        tiles: All tiles in the grid, defining grid dimensions.
        zoom: Zoom level of the tiles.

    Returns:
        tuple[numpy.ndarray, tuple[float, float, float, float]]: The
        stitched RGBA image with shape `(H, W, 4)` and dtype
        `uint8`, plus `(west, south, east, north)` in EPSG:3857
        meters.

    Raises:
        ValueError: If any tile bytes cannot be decoded as an image.

    Examples:
        - Stitch a single synthetic tile into a 256x256 RGBA image:
            ```python
            >>> import io
            >>> from PIL import Image
            >>> from cleopatra.basemap.tiles import Tile, stitch_tiles
            >>> buf = io.BytesIO()
            >>> Image.new("RGBA", (256, 256), (255, 0, 0, 255)).save(buf, "PNG")
            >>> tile = Tile(0, 0, 0)
            >>> image, extent = stitch_tiles({tile: buf.getvalue()}, [tile], 0)
            >>> image.shape
            (256, 256, 4)
            >>> image.dtype.name
            'uint8'

            ```
        - The returned EPSG:3857 extent comes from
            `_tile_xy_bounds` on the corner tiles:
            ```python
            >>> import io
            >>> from PIL import Image
            >>> from cleopatra.basemap.tiles import Tile, stitch_tiles
            >>> buf = io.BytesIO()
            >>> Image.new("RGBA", (256, 256), (0, 255, 0, 255)).save(buf, "PNG")
            >>> tile = Tile(0, 0, 0)
            >>> _, (w, s, e, n) = stitch_tiles({tile: buf.getvalue()}, [tile], 0)
            >>> w < e and s < n
            True

            ```
        - Invalid tile bytes raise `ValueError`:
            ```python
            >>> from cleopatra.basemap.tiles import Tile, stitch_tiles
            >>> tile = Tile(0, 0, 0)
            >>> try:
            ...     stitch_tiles({tile: b"not an image"}, [tile], 0)
            ... except ValueError as exc:
            ...     str(exc).startswith("Failed to decode tile image:")
            True

            ```
    """
    _require_tiles_extra()
    from PIL import Image

    try:
        first_img = Image.open(io.BytesIO(next(iter(tile_data.values()))))
    except Exception as e:
        raise ValueError(
            f"Failed to decode tile image: {e}. The tile server "
            f"may have returned invalid data."
        ) from e
    tile_size = first_img.width

    x_indices = sorted({t.x for t in tiles})
    y_indices = sorted({t.y for t in tiles})
    width = len(x_indices) * tile_size
    height = len(y_indices) * tile_size

    merged = Image.new("RGBA", (width, height))
    for tile, png_bytes in tile_data.items():
        try:
            img = Image.open(io.BytesIO(png_bytes)).convert("RGBA")
        except Exception as e:
            raise ValueError(
                f"Failed to decode tile z={tile.z}/x={tile.x}/y={tile.y}: {e}"
            ) from e
        x_offset = (tile.x - x_indices[0]) * tile_size
        y_offset = (tile.y - y_indices[0]) * tile_size
        merged.paste(img, (x_offset, y_offset))

    image = np.array(merged)

    tl_left, _tl_bottom, _tl_right, tl_top = _tile_xy_bounds(
        Tile(x_indices[0], y_indices[0], zoom)
    )
    _br_left, br_bottom, br_right, _br_top = _tile_xy_bounds(
        Tile(x_indices[-1], y_indices[-1], zoom)
    )
    extent_3857 = (tl_left, br_bottom, br_right, tl_top)

    return image, extent_3857

world_texture(provider=None, *, zoom=5, n_lon=2880, n_lat=1440, cache=True, max_workers=8, timeout=10, retries=2, user_agent=USER_AGENT)

Fetch a whole-world XYZ tile basemap as an equirectangular texture.

The XYZ-provider analogue of cleopatra.basemap.reference.relief: it fetches every tile of the zoom-level world grid, stitches them into a Web Mercator mosaic, and resamples that onto an equirectangular lon/lat grid (mercator_to_equirectangular), so the result can texture a globe or back an EPSG:4326 axis. Unlike relief (a fixed re-hosted asset), this pulls live tiles from any xyzservices provider.

The grid is 2**zoom tiles per side (4**zoom total -- 1024 at the default zoom 5, capped at zoom 6 / 4096 tiles), so the texture is cached on disk and fetched once per (provider, zoom, n_lon, n_lat).

Note

A whole-world fetch pulls thousands of tiles. The default provider, OpenStreetMap.Mapnik, prohibits bulk downloading in its usage policy; pass a bulk-permitting imagery provider (e.g. "Esri.WorldImagery") for whole-world textures.

Parameters:

Name	Type	Description	Default
provider	Any	An xyzservices provider name (e.g. "Esri.WorldImagery") or a resolved xyzservices.TileProvider. Unlike add_tiles this needs a Mapping-like provider, because the disk cache is keyed on provider.get("name", ...) -- so the OGC providers in cleopatra.basemap.ogc raise AttributeError here and are supported only through add_tiles; None uses the default (OpenStreetMap.Mapnik). A name is resolved by get_provider.	None
zoom	int	Tile zoom level (0..6); the world grid is 2**zoom tiles per side.	5
n_lon	int	Width of the returned texture, spanning -180..180 degrees.	2880
n_lat	int	Height of the returned texture, spanning 90..-90 degrees.	1440
cache	bool	When True (default), read/write the texture under Config.get_cache_dir()/"world-textures" so the tiles are fetched once. False always refetches.	True
max_workers	int	Maximum concurrent tile HTTP connections.	8
timeout	int	Per-tile HTTP request timeout in seconds.	10
retries	int	Per-tile retry count on failure.	2
user_agent	str	User-Agent header sent on every tile request.	USER_AGENT

Returns:

Type	Description
ndarray	numpy.ndarray: The (n_lat, n_lon, 3) texture as float32 in [0, 1]
ndarray	(RGB; the tiles' alpha is dropped). This differs from
ndarray	reference.relief, which returns uint8 0..255.

Raises:

Type	Description
ImportError	If the [tiles] extra (xyzservices, Pillow, pyproj) is not installed.
ValueError	If zoom is outside 0..6, or n_lon/n_lat are not positive.
ConnectionError	If any tile cannot be fetched after all retries.

Examples:

• Fetch a coarse world texture (network-dependent, hence skipped under doctest):

  >>> from cleopatra.basemap.tiles import world_texture
  >>> tex = world_texture("Esri.WorldImagery", zoom=2)  # doctest: +SKIP
  >>> tex.shape, tex.dtype.name  # doctest: +SKIP
  ((1440, 2880, 3), 'float32')
  

Source code in src/cleopatra/basemap/tiles.py

def world_texture(
    provider: Any = None,
    *,
    zoom: int = 5,
    n_lon: int = 2880,
    n_lat: int = 1440,
    cache: bool = True,
    max_workers: int = 8,
    timeout: int = 10,
    retries: int = 2,
    user_agent: str = USER_AGENT,
) -> np.ndarray:
    """Fetch a whole-world XYZ tile basemap as an equirectangular texture.

    The XYZ-provider analogue of `cleopatra.basemap.reference.relief`: it
    fetches every tile of the `zoom`-level world grid, stitches them into a Web
    Mercator mosaic, and resamples that onto an equirectangular lon/lat grid
    (`mercator_to_equirectangular`), so the result can texture a globe or back
    an EPSG:4326 axis. Unlike `relief` (a fixed re-hosted asset), this pulls
    live tiles from any `xyzservices` provider.

    The grid is `2**zoom` tiles per side (`4**zoom` total -- 1024 at the default
    zoom 5, capped at zoom 6 / 4096 tiles), so the texture is cached on disk and
    fetched once per `(provider, zoom, n_lon, n_lat)`.

    Note:
        A whole-world fetch pulls thousands of tiles. The default provider,
        `OpenStreetMap.Mapnik`, prohibits bulk downloading in its usage policy;
        pass a bulk-permitting imagery provider (e.g. `"Esri.WorldImagery"`) for
        whole-world textures.

    Args:
        provider: An `xyzservices` provider name (e.g. `"Esri.WorldImagery"`) or
            a resolved `xyzservices.TileProvider`. Unlike `add_tiles` this
            needs a Mapping-like provider, because the disk cache is keyed on
            `provider.get("name", ...)` -- so the OGC providers in
            `cleopatra.basemap.ogc` raise `AttributeError` here and are
            supported only through `add_tiles`;
            `None` uses the default (`OpenStreetMap.Mapnik`). A name is resolved
            by `get_provider`.
        zoom: Tile zoom level (0..6); the world grid is `2**zoom` tiles per side.
        n_lon: Width of the returned texture, spanning -180..180 degrees.
        n_lat: Height of the returned texture, spanning 90..-90 degrees.
        cache: When `True` (default), read/write the texture under
            `Config.get_cache_dir()/"world-textures"` so the tiles are fetched
            once. `False` always refetches.
        max_workers: Maximum concurrent tile HTTP connections.
        timeout: Per-tile HTTP request timeout in seconds.
        retries: Per-tile retry count on failure.
        user_agent: `User-Agent` header sent on every tile request.

    Returns:
        numpy.ndarray: The `(n_lat, n_lon, 3)` texture as `float32` in `[0, 1]`
        (RGB; the tiles' alpha is dropped). This differs from
        `reference.relief`, which returns `uint8` 0..255.

    Raises:
        ImportError: If the `[tiles]` extra (xyzservices, Pillow, pyproj) is
            not installed.
        ValueError: If `zoom` is outside `0..6`, or `n_lon`/`n_lat` are not
            positive.
        ConnectionError: If any tile cannot be fetched after all retries.

    Examples:
        - Fetch a coarse world texture (network-dependent, hence skipped under
            doctest):
            ```python
            >>> from cleopatra.basemap.tiles import world_texture
            >>> tex = world_texture("Esri.WorldImagery", zoom=2)  # doctest: +SKIP
            >>> tex.shape, tex.dtype.name  # doctest: +SKIP
            ((1440, 2880, 3), 'float32')

            ```
    """
    _require_tiles_extra()
    if not 0 <= zoom <= _MAX_WORLD_ZOOM:
        raise ValueError(
            f"zoom must be between 0 and {_MAX_WORLD_ZOOM} (a world grid is "
            f"4**zoom tiles; a higher zoom would exhaust memory), got {zoom}"
        )
    if n_lon < 1 or n_lat < 1:
        raise ValueError(
            f"n_lon and n_lat must be positive, got n_lon={n_lon}, n_lat={n_lat}"
        )
    # Accept a name string or a resolved TileProvider, matching add_tiles.
    if isinstance(provider, str) or provider is None:
        resolved = get_provider(provider)
    else:
        resolved = provider

    cache_path = None
    if cache:
        name = re.sub(r"[^0-9A-Za-z]+", "_", str(resolved.get("name", "default")))
        digest = hashlib.blake2s(
            repr((name, zoom, n_lon, n_lat)).encode(), digest_size=8
        ).hexdigest()
        cache_dir = Config.get_cache_dir() / "world-textures"
        cache_dir.mkdir(parents=True, exist_ok=True)
        cache_path = cache_dir / f"{name}_z{zoom}_{digest}.npy"
        if cache_path.exists():
            try:
                return np.load(cache_path)
            except (OSError, ValueError, EOFError):
                # A truncated / foreign / pre-atomic .npy: drop it and rebuild
                # (as reference.relief does with a poisoned cached asset).
                cache_path.unlink(missing_ok=True)

    side = 2**zoom
    tiles = [Tile(x, y, zoom) for x in range(side) for y in range(side)]
    mosaic, bounds = stitch_tiles(
        fetch_tiles(
            tiles,
            resolved,
            max_workers=max_workers,
            timeout=timeout,
            retries=retries,
            user_agent=user_agent,
        ),
        tiles,
        zoom,
    )
    rgb = np.asarray(mosaic)[..., :3]
    texture = mercator_to_equirectangular(rgb, bounds, n_lon=n_lon, n_lat=n_lat)
    texture = np.clip(texture / 255.0, 0.0, 1.0).astype("float32")

    if cache_path is not None:
        # Write to a unique temp file, then rename: np.save straight onto the
        # final path leaves a truncated array under the name the next run trusts
        # if it is interrupted, and a fixed .tmp name lets concurrent same-key
        # writes clobber each other. Remove the temp if the write fails.
        fd, tmp_name = tempfile.mkstemp(dir=cache_path.parent, suffix=".npy.tmp")
        try:
            with os.fdopen(fd, "wb") as fh:
                np.save(fh, texture)
            os.replace(tmp_name, cache_path)
        except BaseException:
            with contextlib.suppress(FileNotFoundError):
                os.remove(tmp_name)
            raise
    return texture