Orthomosaic Generation and Tiling for Archaeological Site Maps #

An orthomosaic is the flattened, map-accurate photograph that most heritage projects actually publish, print, and dig against — the deliverable a site director points at during excavation and the base layer a planning authority expects in a submission pack. Getting there from a raw photogrammetric block is not a single button-press: it demands orthorectification against an elevation surface, radiometric reconciliation across flights shot in different light, deliberate seamline routing, and a storage format that a web client can stream without downloading gigabytes. This section sits inside the broader Photogrammetry & 3D Site Mapping Pipelines workflow and takes the georeferenced products from automated drone image processing workflows through to tiled, archive-ready rasters. The two guides below — generating Cloud-Optimized GeoTIFF orthomosaics and tiling orthomosaics for web map delivery — carry the practical detail; this overview frames the decisions that bind them together.

Orthorectification, DSM-based ortho, and true-ortho #

Orthorectification removes the perspective and relief displacement of a raw aerial photo so that every pixel sits at its correct planimetric position. The transform needs an elevation model: for most open archaeological landscapes a Digital Surface Model (DSM) derived from the dense cloud is sufficient, and a standard DSM-based orthomosaic drapes imagery over that surface. The distinction that matters on structured sites is true-ortho. A conventional DSM ortho corrects for terrain but not for the lean of tall vertical features — the parapet of a standing wall, a spoil-heap edge, or a scaffold tower still topples outward and occludes the ground behind it. A true-ortho pipeline uses per-pixel visibility from the camera network to select occlusion-free source pixels, so a Roman fort rampart reads as a crisp footprint rather than a smeared skirt. True-ortho costs more processing and denser imagery, so reserve it for sites where standing structures matter and accept DSM ortho for open trench plans.

Radiometric and colour balancing across flights #

Multi-flight surveys almost never share exposure. A morning grid over a barrow cemetery and an afternoon grid over the same field will differ in sun angle, white balance, and haze, and the join between them shows as a visible tonal step. Radiometric balancing normalises these blocks before mosaicking — histogram matching against a reference flight, or a gain-and-offset model fitted in overlap zones — so that the published sheet reads as one continuous surface. This is not cosmetic: uneven radiometry defeats any later index-based analysis (crop-mark enhancement, vegetation masking) and misleads a specialist reading soil colour off the map. Balance before you cut seamlines, never after.

Seamlines, resolution, and storage #

A seamline is the path along which two overlapping source images are cut and joined. Left to a naive nearest-centre rule, seamlines run straight through buildings and produce ghosted doorways. Route them along low-contrast ground — field edges, paths, the base of a wall rather than its top — and the mosaic hides its joins. Resolution is the other governing trade-off. Ground sample distance (GSD) sets how much real ground each pixel covers; a 2 cm GSD survey of a 4-hectare site is visually superb and enormous on disk. The uncompressed footprint of a single-band raster is approximately

$\text{bytes} \approx \frac{A}{g^2} \times b \times d$

where $A$ is the mapped area in square metres, $g$ the GSD in metres, $b$ the band count, and $d$ the bytes per sample. Halving GSD quadruples storage, so choose a delivery resolution matched to the question the map answers, and keep the full-resolution master separately from the tiled derivative.

GeoTIFF, Cloud-Optimized GeoTIFF, and overviews #

The archival container is a GeoTIFF: a TIFF carrying embedded georeferencing (an affine geotransform plus a CRS, here EPSG:27700 — British National Grid; # substitute your site's EPSG for other regions). A plain GeoTIFF, however, forces a client to read the whole file to display any part of it. A Cloud-Optimized GeoTIFF (COG) reorganises the same pixels into internal tiles plus pre-computed reduced-resolution overviews, and lays them out so an HTTP range request can fetch only the bytes for the current viewport and zoom. That single change lets a heritage archive host a multi-gigabyte orthomosaic on object storage and serve it to QGIS or a browser without a tile server. The Cloud-Optimized GeoTIFF guide walks the exact gdal_translate and validation steps.

Tiling: gdal2tiles, XYZ, and WMTS #

Where a COG serves the raster directly, a tile pyramid pre-renders it into millions of small PNG or JPEG images addressed by zoom/column/row. This is what a Leaflet or MapLibre front end consumes as an XYZ layer, and what an OGC WMTS service advertises to desktop GIS. gdal2tiles.py builds either scheme, and the pyramid can be packed into a single MBTiles SQLite container for portable offline delivery to the field. The web map delivery guide covers zoom-range selection and serving.

Prerequisites #

Assumptions: a georeferenced dense cloud and DSM exist for each flight, all blocks share a target CRS, and you have write access to the archive object store. GDAL binaries must be built with GeoTIFF, PNG, and SQLite (MBTiles) drivers — confirm with gdalinfo --formats | grep -E "COG|MBTiles".

Implementation #

The field-to-archive order is: orthorectify each block, balance radiometry, mosaic with routed seamlines, convert to COG with overviews, then tile.

1. Orthorectify each flight against the DSM. gdalwarp resamples a block into the target grid using the elevation surface implicitly baked into the georeferenced ortho export, or explicitly via a DSM-driven pipeline in the photogrammetry package. At the GDAL stage you normalise CRS and pixel grid:

# GDAL 3.8.4
# Reproject and snap each block to a common 0.02 m grid in British National Grid
# substitute your site's EPSG
gdalwarp -t_srs EPSG:27700 -tr 0.02 0.02 -tap \
  -r bilinear -overwrite \
  flight_a_ortho.tif flight_a_bng.tif

2. Balance radiometry across blocks. Fit a per-band gain/offset against a reference flight in the overlap region with rasterio:

# requirements.txt
# rasterio==1.3.9
# numpy==1.26.4

# rasterio==1.3.9, numpy==1.26.4
import numpy as np
import rasterio

def gain_offset(reference, target):
    """Least-squares gain/offset per band from co-located valid pixels."""
    with rasterio.open(reference) as ref, rasterio.open(target) as tgt:
        r = ref.read(masked=True)
        t = tgt.read(masked=True)
    coeffs = []
    for band in range(r.shape[0]):
        m = ~(r[band].mask | t[band].mask)
        x = t[band].data[m].astype("f4")
        y = r[band].data[m].astype("f4")
        g, o = np.polyfit(x, y, 1)  # y = g*x + o
        coeffs.append((float(g), float(o)))
    return coeffs

3. Mosaic with routed seamlines. Build a virtual mosaic and cut lines that avoid structures, then materialise the sheet:

# GDAL 3.8.4
gdalbuildvrt -resolution highest site_mosaic.vrt flight_*_bng.tif
gdal_translate -of GTiff -co TILED=YES site_mosaic.vrt site_ortho.tif

4. Convert to COG with overviews — detailed in the COG guide. In short, gdaladdo builds the pyramid and gdal_translate -of COG writes the tiled layout.

5. Tile for web delivery — gdal2tiles.py against the COG, detailed in the tiling guide.

CRS & Geometry Validation #

Confirm the written CRS and geotransform with gdalinfo site_ortho.tif: the Coordinate System is: block must report your target EPSG and the Pixel Size must match the intended GSD. A common failure is an axis-order surprise — some tools emit EPSG:4326 with latitude first, which shifts a mosaic into the sea. Validate the corner coordinates fall inside the site’s known bounding box before promotion. For balanced blocks, check that overlap-zone residuals are near zero; a persistent tonal step means the gain/offset fit used too few valid pixels.

Automated QA/QC #

Gate promotion on four machine checks: (1) CRS equals the mandated EPSG; (2) no unexpected nodata swathes inside the survey boundary (gdalinfo -stats reporting an all-nodata band signals a failed warp); (3) overview levels present via gdalinfo showing Overviews:; and (4) COG validity from rio cogeo validate. Route failures to a quarantine directory with the failing check name, and record source flights, GSD, CRS, seamline strategy, and validation results to an audit table so any published sheet is reproducible.

Integration with adjacent workflows #

The tiled and COG outputs feed two sibling workflows directly. Georeferenced meshes and the same site datum flow into exporting 3D models to spatial databases, where the orthomosaic becomes the draped base for 3D archive records. Upstream, the alignment and control established in automated drone image processing workflows determines the absolute accuracy this section can only preserve, never improve.

Troubleshooting #

• ERROR 1: Too many points (X out of Y) failed to transform from gdalwarp — the source and target CRS pairing lacks a datum transform or the input georeferencing is corrupt. Check gdalsrsinfo on the source and install the relevant PROJ grid.
• ERROR 6: SetColorTable() only supported for Byte or UInt16 bands — you are mosaicking paletted and RGB tiles together; expand with gdal_translate -expand rgb first.
• Visible seams after balancing — the gain/offset fit ran on too few overlap pixels; widen the overlap mask or fall back to full histogram matching.
• ERROR 1: TIFFReadDirectory:Sanity check on directory count failed — a truncated or partially written GeoTIFF; re-run the translate step and verify disk space.
• Blurry mosaic at native zoom — a lossy JPEG-compressed intermediate was re-warped; keep the master in DEFLATE or LZW and only apply JPEG at the tile stage.
• Mosaic offset by tens of metres — axis-order or a local-grid export; reproject explicitly and verify corners against the site bounding box.

Compliance notes #

Record ISO 19115 lineage for every published sheet: source flights and dates, the CRS as an EPSG code, GSD as spatial resolution, the orthorectification method (DSM vs true-ortho), and the radiometric processing applied. MIDAS Heritage expects an event/date and a monument reference; attach both to the raster metadata or an accompanying sidecar. To meet FAIR reuse, publish the COG on a stable, range-request-capable URL with a licence statement, so the interoperable raster is genuinely retrievable rather than merely archived.

Related #

• Automated Drone Image Processing Workflows — the alignment and georeferencing that feeds orthomosaic generation.
• Exporting 3D Models to Spatial Databases — where orthomosaics become draped base layers for 3D archive records.
• Generating Cloud-Optimized GeoTIFF Orthomosaics — COG creation, compression, and validation.
• Tiling Orthomosaics for Web Map Delivery — XYZ/TMS pyramids, MBTiles, and serving.
• Mesh Generation & Optimization for Ruins — the mesh products alongside the raster deliverable.

Part of Photogrammetry & 3D Site Mapping Pipelines.