Rasterizing Labels into Segmentation Masks

To burn polygon labels into a mask that aligns exactly with an image chip, pass the chip’s own transform and shape to rasterio.features.rasterize:

from rasterio.features import rasterize

mask = rasterize(
    shapes=[(geom, cls) for geom, cls in zip(gdf.geometry, gdf.class_id)],
    out_shape=(256, 256),
    transform=chip_transform,    # the CHIP transform, not the scene transform
    fill=0,                      # 0 = unlabelled
    dtype="uint8",
)

Everything that makes this hard is in the detail: which class wins where polygons overlap, what happens to pixels the boundary crosses, and how to tell the model “I do not know what this pixel is” rather than lying to it. This page sits under building training datasets from satellite imagery in Raster Machine Learning & Model Inference.

Why Burn Order Is a Design Decision

rasterize paints geometries in the order it receives them, so a later geometry overwrites an earlier one wherever they overlap. That is not a quirk to work around — it is the mechanism you use to express priority.

Last geometry burned wins A large orchard polygon overlaps a small pond polygon. Burning the orchard last covers the pond entirely and the water class disappears from the mask. Burning by descending area, so the small pond is painted after the large orchard, preserves both. The rule generalises: the more specific class must be burned later. file order — orchard last sorted by descending area orchard — class 2 the pond is gone; 340 water pixels lost orchard — class 2 pond — 4 both classes survive; the pond is burned second Sorting by descending area is the safe default; an explicit priority column is better when you have one.

Sorting by descending area works because containment is the common case: ponds inside fields, buildings inside plots, clearings inside forest. When priority is not a function of size — a road crossing several land cover polygons should always win — add an explicit priority column and sort on that instead.

Environment & Setup

Package Version pin Used for
rasterio >=1.3.0 features.rasterize and the transform handling
geopandas >=0.14 Reading, reprojecting and sorting the label layer
shapely >=2.0 Boundary buffers for the ignore region
numpy >=1.23 Mask arithmetic and validity checks 
pip install "rasterio>=1.3.0" "geopandas>=0.14" "shapely>=2.0" "numpy>=1.23"

Complete Working Example

This function produces a two-layer target: the class mask, and a boolean array marking pixels the loss should ignore.

import geopandas as gpd
import numpy as np
from rasterio.features import rasterize

IGNORE = 255


def burn_mask(labels: gpd.GeoDataFrame, transform, shape: tuple[int, int],
              *, class_field: str = "class_id",
              priority_field: str | None = None,
              boundary_buffer_m: float = 10.0,
              exhaustive: bool = False) -> tuple[np.ndarray, np.ndarray]:
    """Burn a multi-class mask plus an ignore mask for uncertain pixels."""
    if labels.empty:
        empty = np.full(shape, 0 if exhaustive else IGNORE, dtype="uint8")
        return empty, (empty == IGNORE)

    # Priority: explicit column if given, else largest first so small wins
    if priority_field is not None:
        ordered = labels.sort_values(priority_field)
    else:
        ordered = labels.assign(_a=labels.geometry.area).sort_values("_a", ascending=False)

    mask = rasterize(
        shapes=[(g, int(v)) for g, v in zip(ordered.geometry, ordered[class_field])],
        out_shape=shape,
        transform=transform,
        fill=0,
        all_touched=False,
        dtype="uint8",
    )

    # Boundaries are uncertain: burn a buffered ring as the ignore value
    rings = ordered.geometry.boundary.buffer(boundary_buffer_m)
    ignore = rasterize(
        shapes=[(g, 1) for g in rings if not g.is_empty],
        out_shape=shape, transform=transform, fill=0, dtype="uint8",
    ).astype(bool)

    if not exhaustive:
        # Unlabelled means unknown, not background
        ignore |= (mask == 0)

    return mask, ignore


if __name__ == "__main__":
    import rasterio

    with rasterio.open("s2_stack_36NYF.tif") as src:
        win = rasterio.windows.Window(2048, 4096, 256, 256)
        chip = src.read(window=win)
        tr = src.window_transform(win)
        gdf = gpd.read_file("field_labels.gpkg").to_crs(src.crs)

    m, ig = burn_mask(gdf.cx[tr.c:tr.c + 2560, tr.f - 2560:tr.f], tr, (256, 256))
    print("classes present:", np.unique(m), "ignored:", f"{ig.mean():.1%}")
Unlabelled is not the same as background A chip contains four fields but only two were digitised. Treating the rest as background labels two real fields as negative examples, and the model learns to suppress exactly the pattern it should detect. Marking unlabelled pixels as ignore removes them from the loss instead, so the two digitised fields still teach and the undigitised ones teach nothing. exhaustive = True (wrong here) exhaustive = False two real fields taught as negatives undigitised pixels excluded from the loss labelled field background ignored

The exhaustive flag encodes the single most consequential assumption in the whole file. If the label layer covers everything in the chip, an unlabelled pixel is background and the model should learn it. If only fields of interest were digitised, an unlabelled pixel might be anything, and training the model to call it background actively teaches it to miss real targets.

Variant Patterns

1. A boundary-aware three-class target

Some architectures learn better from an explicit boundary class than from an ignore region: interior, boundary, background. The boundary class gives the network something to predict rather than something to skip.

Interior, boundary, background The same polygon is rendered three ways. As a plain binary mask the model must place the edge exactly. With an ignore ring the edge pixels are excluded from the loss. With an explicit boundary class the model is asked to predict the edge itself, which is what separates touching fields at inference time. binary mask with ignore ring boundary as a class touching fields merge at inference edge errors not penalised predicted edges separate instances interior boundary ignored 
import numpy as np
from rasterio.features import rasterize


def interior_boundary_mask(labels, transform, shape, *, edge_m: float = 15.0):
    """0 = background, 1 = interior, 2 = boundary."""
    interior = rasterize([(g, 1) for g in labels.geometry],
                         out_shape=shape, transform=transform, fill=0, dtype="uint8")
    edges = rasterize([(g.boundary.buffer(edge_m), 1) for g in labels.geometry],
                      out_shape=shape, transform=transform, fill=0, dtype="uint8")
    return np.where(edges == 1, 2, interior).astype("uint8")

2. Burning lines and points

Linear features fall between pixel centres and vanish under the default rule. Buffer them to a realistic width instead of reaching for all_touched, which produces a one-pixel line regardless of whether the road is 3 m or 30 m wide.

import geopandas as gpd
from rasterio.features import rasterize

roads = gpd.read_file("roads.gpkg").to_crs(chip_crs)
roads["geometry"] = roads.buffer(roads["width_m"] / 2)   # real width, not one pixel
road_mask = rasterize([(g, 1) for g in roads.geometry],
                      out_shape=(256, 256), transform=chip_transform,
                      fill=0, dtype="uint8")

3. Masks for a whole scene rather than a chip

The same call scales to a full scene by passing the scene transform and shape. This is how a reference raster for accuracy assessment is made, and the output feeds directly into the comparison in computing confusion matrices from raster pairs.

import rasterio
from rasterio.features import rasterize

with rasterio.open("s2_stack_36NYF.tif") as src:
    profile = src.profile | {"count": 1, "dtype": "uint8", "nodata": 0,
                             "compress": "deflate", "tiled": True}
    scene_mask = rasterize(
        [(g, int(v)) for g, v in zip(gdf.geometry, gdf.class_id)],
        out_shape=(src.height, src.width), transform=src.transform,
        fill=0, dtype="uint8",
    )

with rasterio.open("reference.tif", "w", **profile) as dst:
    dst.write(scene_mask, 1)

Common Errors

The mask is entirely zero

Either no geometry intersects the chip, or the labels are in a different CRS from the transform. rasterize returns the fill value silently in both cases. Reproject first, and assert that a known-overlapping chip produces a non-empty mask.

ValueError: Invalid value for dtype

A class code exceeds the dtype’s range — commonly a class id of 300 in a uint8 mask, or a float class column. Cast the values to int explicitly and keep codes inside 1–254 so 0 and 255 stay available for background and ignore.

Thin features disappear entirely

The default pixel-centre rule skips any polygon that does not cover a pixel centre. Buffer linear features to their real width, or set all_touched=True if a one-pixel-wide representation is genuinely what you want.

Frequently Asked Questions

Q: What happens when two label polygons overlap? The last geometry burned wins, so the iteration order is the priority order. Sort the pairs deliberately — usually smallest area or most specific class last — rather than relying on whatever order the file happened to store.

Q: How do I stop the model being punished for ambiguous boundaries? Burn a thin buffer around every polygon boundary as a dedicated ignore value, then exclude that value from the loss. Field boundaries are rarely digitised to the pixel, and forcing the model to match them exactly wastes capacity on noise.

Q: Should background be a class or an ignore value? It depends on whether unlabelled really means background. If the layer is exhaustive, background is a real class and should be learned. If only some features were digitised, unlabelled pixels are unknown and must be ignored, otherwise the model is trained to call real targets background.