Automated Static Map Generation from GeoJSON

Producing publication-ready cartography from GeoJSON in headless environments is harder than it looks: basemaps misalign when CRS is left unmanaged, legends shift when attribute ranges vary between runs, and export resolution clips silently under memory pressure. This guide builds a deterministic pipeline — from raw GeoJSON ingestion through high-resolution raster export — that runs reliably in scheduled CI/CD jobs, Docker containers, and PDF report generation workflows without manual GIS intervention.

Pipeline Architecture

The diagram below shows the end-to-end data flow from a GeoJSON source file to a document-ready raster asset.

Prerequisites

Install the full stack in an isolated virtual environment:

pip install geopandas matplotlib contextily pyproj shapely

On Debian/Ubuntu, install system-level geospatial libraries before pip: sudo apt-get install libgdal-dev proj-bin. On macOS with Homebrew: brew install gdal proj. Always pin versions in requirements.txt or pyproject.toml — silent breaking changes in PROJ string parsing and tile provider APIs are common across minor releases.

This pipeline fits within the broader Dynamic Map & Data Embedding Workflows context, where spatial data must be transformed, styled, and embedded without manual intervention at any stage.

Step-by-step Implementation

Step 1 — Ingest and Validate GeoJSON

Load the source file with geopandas.read_file(). Mixed-geometry FeatureCollection inputs (Points and Polygons in the same file) will fail downstream rendering unless filtered first.

import geopandas as gpd
from pathlib import Path

def load_and_validate(path: str) -> gpd.GeoDataFrame:
    gdf = gpd.read_file(Path(path))
    if gdf.empty:
        raise ValueError(f"No features found in {path}")
    # Drop null and invalid geometries
    gdf = gdf[gdf.geometry.notna() & gdf.geometry.is_valid]
    # Enforce a single geometry type (keep Polygon/MultiPolygon for area maps)
    dominant = gdf.geometry.geom_type.value_counts().idxmax()
    gdf = gdf[gdf.geometry.geom_type.isin([dominant, f"Multi{dominant}"])]
    return gdf

Automated pipelines fail silently on non-compliant JSON structures. Validate coordinate ordering against the RFC 7946 specification ([longitude, latitude]) and confirm type and geometry fields are present before any downstream operation. Use gdf.is_valid to catch self-intersecting rings that pass JSON parsing but break rasterization.

Step 2 — Normalize Coordinate Reference Systems

Raw GeoJSON always arrives in unprojected geographic coordinates (EPSG:4326). Web basemap tiles use Web Mercator (EPSG:3857). Mixing them produces the most common spatial offset artifact in automated pipelines. The diagram below illustrates why the reprojection step must happen before any basemap fetch or extent calculation.

def normalize_crs(
    gdf: gpd.GeoDataFrame,
    target_crs: str = "EPSG:3857"
) -> gpd.GeoDataFrame:
    if gdf.crs is None:
        gdf = gdf.set_crs("EPSG:4326")
    return gdf.to_crs(target_crs)

If your reporting template requires metric scaling for legends, scale bars, or buffer zones, chain a secondary transformation to a local projected CRS after basemap alignment — for example, a UTM zone or a national grid system. Always verify CRS alignment between vector layers and any raster basemap before compositing; mismatched transforms do not raise exceptions, they silently shift features.

Step 3 — Apply Programmatic Styling and Basemap Alignment

Initialize a Matplotlib figure with explicit dimensions matching your target output resolution. Use contextily.add_basemap() to fetch and align web tiles to your axis extent. The zoom="auto" parameter lets contextily calculate an appropriate zoom level from the bounding box, but set it explicitly for batch jobs where bounding box size varies significantly between runs to prevent over- or under-sampling tiles.

import matplotlib.pyplot as plt
import contextily as ctx

def render_map(
    gdf: gpd.GeoDataFrame,
    figsize: tuple[float, float] = (10.0, 8.0),
    dpi: int = 300
) -> tuple[plt.Figure, plt.Axes]:
    fig, ax = plt.subplots(figsize=figsize, dpi=dpi)
    gdf.plot(
        ax=ax,
        color="#2a9d8f",
        edgecolor="#264653",
        linewidth=0.8,
        alpha=0.85
    )
    ax.set_axis_off()
    try:
        ctx.add_basemap(
            ax,
            source=ctx.providers.OpenStreetMap.Mapnik,
            crs=gdf.crs.to_string(),
            zoom="auto"
        )
    except Exception as exc:
        # Log and continue — basemap is supplementary, not structural
        import logging
        logging.warning("Basemap unavailable (offline mode): %s", exc)
    return fig, ax

For thematic styling on datasets with variable attribute ranges — common when a report covers multiple administrative regions — pair this step with Dynamic Legend Injection for Variable Datasets to automate classification breaks, color ramp normalization, and legend positioning without manual figure tweaking between runs.

Step 4 — Configure Layout and Export

Static map generation demands precise control over margins, DPI, and background transparency. Remove axis spines, ticks, and gridlines with ax.axis('off') to achieve a clean cartographic aesthetic. Use bbox_inches='tight' to prevent clipping and pad_inches=0.1 to retain a minimal visual margin.

def export_map(
    fig: plt.Figure,
    output_path: str,
    dpi: int = 300,
    transparent: bool = True
) -> None:
    out = Path(output_path)
    out.parent.mkdir(parents=True, exist_ok=True)
    fig.savefig(
        out,
        dpi=dpi,
        bbox_inches="tight",
        pad_inches=0.1,
        transparent=transparent
    )
    plt.close(fig)

High-resolution exports frequently exhaust memory in constrained CI environments. Render at 150 DPI for draft validation, then scale to 300 DPI for the production artifact. If your pipeline must embed these exports into complex document generators, see Embedding Interactive Mapbox Exports into WeasyPrint PDFs for raster-to-vector fallback strategies and print-optimized DPI handling.

Step 5 — Integrate with Legends, Tables, and Multi-page Reports

Static maps rarely exist in isolation. Production reporting pipelines combine cartographic outputs with attribute summaries, statistical charts, and metadata blocks. Use Matplotlib’s GridSpec or constrained_layout to allocate dedicated subplots for legends, north arrows, and scale bars.

import matplotlib.gridspec as gridspec

def build_composite_figure(
    gdf: gpd.GeoDataFrame,
    title: str,
    figsize: tuple[float, float] = (14.0, 10.0),
    dpi: int = 300
) -> plt.Figure:
    fig = plt.figure(figsize=figsize, dpi=dpi)
    gs = gridspec.GridSpec(2, 2, figure=fig, height_ratios=[4, 1], hspace=0.05)
    ax_map = fig.add_subplot(gs[0, :])
    ax_legend = fig.add_subplot(gs[1, 0])
    ax_meta = fig.add_subplot(gs[1, 1])

    gdf.plot(ax=ax_map, color="#2a9d8f", edgecolor="#264653", linewidth=0.8)
    ax_map.set_title(title, fontsize=12, pad=8)
    ax_map.axis("off")

    for spine in ax_legend.spines.values():
        spine.set_visible(False)
    ax_legend.axis("off")
    ax_meta.axis("off")

    return fig

When attribute tables exceed page boundaries, implement programmatic row splitting and header repetition. For detailed pagination patterns, see Table Pagination Strategies for Large Attribute Tables. Always cache rendered assets during iterative development to avoid redundant tile downloads.

Production-Ready Script

The following script combines all five steps into a single, headless-ready function with structured logging and graceful error handling.

import geopandas as gpd
import matplotlib.pyplot as plt
import contextily as ctx
import logging
from pathlib import Path

logging.basicConfig(level=logging.INFO, format="%(levelname)s: %(message)s")
logger = logging.getLogger(__name__)


def generate_static_map(
    geojson_path: str,
    output_path: str,
    crs_target: str = "EPSG:3857",
    dpi: int = 300,
    figsize: tuple[float, float] = (10.0, 8.0),
    draft_mode: bool = False
) -> Path:
    """
    Generate a publication-ready static map from a GeoJSON file.

    Args:
        geojson_path: Path to source GeoJSON.
        output_path:  Destination path for the raster export (PNG or PDF).
        crs_target:   Target projected CRS for basemap alignment.
        dpi:          Output resolution; use 150 for drafts, 300 for production.
        figsize:      Matplotlib figure dimensions in inches (width, height).
        draft_mode:   When True, overrides dpi to 150 and skips basemap fetching.

    Returns:
        Path to the written output file.
    """
    input_file = Path(geojson_path)
    if not input_file.exists():
        raise FileNotFoundError(f"GeoJSON not found: {input_file}")

    # 1. Ingest & validate
    logger.info("Loading GeoJSON from %s", input_file)
    gdf = gpd.read_file(input_file)
    if gdf.empty:
        raise ValueError("GeoJSON contains no features.")
    initial_count = len(gdf)
    gdf = gdf[gdf.geometry.notna() & gdf.geometry.is_valid]
    dropped = initial_count - len(gdf)
    if dropped:
        logger.warning("Dropped %d invalid/null geometries.", dropped)
    logger.info("Validated %d features.", len(gdf))

    # 2. Normalize CRS
    if gdf.crs is None:
        logger.warning("No CRS detected — assuming EPSG:4326.")
        gdf = gdf.set_crs("EPSG:4326")
    gdf = gdf.to_crs(crs_target)
    logger.info("Reprojected to %s.", crs_target)

    # 3. Render
    render_dpi = 150 if draft_mode else dpi
    fig, ax = plt.subplots(figsize=figsize, dpi=render_dpi)
    gdf.plot(ax=ax, color="#2a9d8f", edgecolor="#264653", linewidth=0.8, alpha=0.85)
    ax.set_axis_off()

    if not draft_mode:
        try:
            ctx.add_basemap(
                ax,
                source=ctx.providers.OpenStreetMap.Mapnik,
                crs=gdf.crs.to_string(),
                zoom="auto"
            )
            logger.info("Basemap aligned successfully.")
        except Exception as exc:
            logger.warning("Basemap fetch failed (continuing without): %s", exc)

    # 4. Export
    out = Path(output_path)
    out.parent.mkdir(parents=True, exist_ok=True)
    fig.savefig(out, dpi=render_dpi, bbox_inches="tight", pad_inches=0.1, transparent=True)
    plt.close(fig)
    logger.info("Map exported to %s (DPI=%d).", out, render_dpi)
    return out


if __name__ == "__main__":
    generate_static_map(
        geojson_path="data/parcels.geojson",
        output_path="output/parcels_map.png",
        draft_mode=False
    )

Edge Cases & Advanced Configuration

Null and Sparse Data

Empty GeoJSON — files with a valid FeatureCollection wrapper but zero features — pass gpd.read_file() without error and silently produce blank exports. Always check gdf.empty immediately after loading and raise before entering the render stage. Sparse data (valid feature count below a reporting threshold) should trigger a warning logged to your artifact system so downstream consumers know the map is statistically incomplete.

Scale and Performance Tuning

For large datasets (more than 500,000 polygons), apply spatial indexing via gdf.sindex before any bbox-filter operations, and consider pre-aggregating features at coarser administrative levels before plotting. Alternatively, use rasterio to rasterize vector data directly to a pixel grid before handing off to matplotlib — this decouples geometry complexity from render time.

Memory consumption scales with geometry complexity and DPI. Profile peak RSS during render using resource.getrusage(resource.RUSAGE_SELF) and set CI memory limits accordingly. A 10-by-8-inch figure at 300 DPI allocates roughly 3000×2400 pixels per channel — approximately 86 MB for an RGBA output before compression.

Multi-Format Outputs

The same pipeline function handles PNG, JPEG, and PDF outputs by matching the file extension passed to output_path. PDF output from matplotlib.savefig() produces a vector-embedded raster at the specified DPI — suitable for single-page deliverables but not for multi-page report generation. For multi-page spatial reports, hand raster exports to Chart-to-PDF Sync with Matplotlib or a WeasyPrint document assembly step.

Headless Environments

In Docker and headless CI runners, matplotlib requires a non-interactive backend. Set it before importing pyplot:

import matplotlib
matplotlib.use("Agg")  # Must be set before importing pyplot
import matplotlib.pyplot as plt

Include this in your container’s entrypoint script or pass it via the MPLBACKEND=Agg environment variable to avoid import-order bugs across modules. Web tile providers enforce rate limits in concurrent batch jobs; mount a persistent cache directory and configure contextily’s cache_path parameter, or pre-download tiles for the expected bounding box and serve them via a local tile source.

Troubleshooting

Detailed Guides in This Section

• Embedding Interactive Mapbox Exports into WeasyPrint PDFs — strategies for raster-to-vector fallback, print-optimized DPI, and bleed handling when embedding map images into WeasyPrint page layouts.

Related

• Dynamic Legend Injection for Variable Datasets — automate classification breaks and color ramp normalization when attribute ranges vary across report runs.
• Table Pagination Strategies for Large Attribute Tables — programmatic row splitting and header repetition for attribute tables that accompany map exports.
• Chart-to-PDF Sync with Matplotlib — synchronize chart rendering with map exports for multi-panel spatial reports.
• Document Architecture & Layout Rules for Spatial Reports — margin, bleed, and grid constraints that govern how map assets are placed in print-ready documents.
• Dynamic Map & Data Embedding Workflows ↑ parent section

Automated static map generation from GeoJSON resolves the core reproducibility gap in spatial reporting: maps that look correct on one machine but shift, clip, or go blank on another. By enforcing CRS normalization, designing for headless execution, and exporting at locked DPI, this pipeline slots directly into the broader Dynamic Map & Data Embedding Workflows stack — delivering consistent, document-ready cartography across every automated publishing run.