Extruded buildings float above or sink into the ground

Problem statement

The extrusion worked, the heights are right, and the buildings hover ten metres above the terrain โ€” or the ground floor of every one of them is buried. On a slope, half the block is in the air and the other half is underground.

There are four distinct causes and they need different fixes: the base height was never set, the base came from the wrong vertical datum, the base is a single value for a sloping footprint, or the model and the terrain are in different height systems entirely. All four produce a picture that is obviously wrong and gives no hint which.

Quick answer

Diagnose by the size and sign of the offset:

import numpy as np, rasterio
from rasterio.mask import mask

with rasterio.open("dtm.tif") as src:
    terrain = []
    for geom in b.geometry:
        arr, _ = mask(src, [geom], crop=True, filled=False)
        v = arr.compressed()
        terrain.append(np.median(v) if v.size else np.nan)

b["terrain_z"] = terrain
b["offset_m"] = b["base_z"] - b["terrain_z"]
print(b["offset_m"].describe().round(2))
offset cause
equal to the base_z, terrain near 0 base never set โ€” buildings at z = 0
a constant 30โ€“60 m vertical datum mismatch (geoid separation)
varies with slope, mean near 0 one base per footprint on sloping ground
a constant few metres different DTM or a different reference epoch
Triage of four causes of floating or sunken extruded buildings and the fix for each.
The size of the offset identifies the cause before you look at any code.

Step-by-step solution

1. Check whether a base was set at all

trimesh.creation.extrude_polygon produces a prism from z = 0 to z = height. If nothing translated it, every building sits at the vertical origin โ€” sea level for a projected CRS with orthometric heights, or the ellipsoid for a geographic one.

print([round(m.bounds[0][2], 2) for m in meshes[:10]])   # all 0.0 means no base was applied

2. Check the vertical datum

A constant offset of tens of metres is a geoid separation. The Netherlands is about 43 m; Britain about 45โ€“50 m; parts of the Mediterranean and the Baltic differ by similar amounts. A model whose z values are heights above a national datum, drawn against a terrain in ellipsoidal height, sits that far underground.

from pyproj import Transformer
t = Transformer.from_crs("EPSG:7415", "EPSG:4979", always_xy=True)   # NAP โ†’ ellipsoidal
print(t.transform(78_642, 457_940, 0.0)[2])

3. Sample the terrain under each footprint, not at its centroid

A centroid can fall in a courtyard, on the far side of a river, or outside the raster. Masking by the polygon and taking a statistic of the cells inside is both more robust and more honest.

4. Choose the statistic deliberately

  • Minimum โ€” the building never floats; the uphill side is buried.
  • Median โ€” split the difference; both sides are out by half the slope.
  • Per-vertex โ€” a sloping base, which needs a mesh rather than a prism.

For most work the median is right. For rendering, the minimum looks better because a floating building is more obviously wrong than a buried one.

5. Handle the sloping case properly if it matters

On steep ground, a flat base over a 20 m footprint at 5% is a metre out at each end. The honest fix is to extrude from a sloping base โ€” build the prism with per-vertex base heights sampled from the DTM.

6. Check the terrain and the model share a CRS

Including the vertical part. EPSG:7415 is not EPSG:28992: the first carries NAP heights, the second is horizontal only, and reprojecting between them silently drops or keeps z depending on the tool.

7. In a web map, remember the base is relative to the terrain

With a terrain source enabled, fill-extrusion-base is measured from the terrain surface, so 0 is correct. Copying a base from the flat case buries every building.

Two scenes contrasting a flat base from the minimum terrain value with a base from the median on a sloping footprint.
On a slope there is no single right base; there is only a documented choice.

Code examples

Example 1 โ€” robust terrain sampling with a report

import numpy as np, rasterio
from rasterio.mask import mask

def terrain_under(buildings, dtm_path, stat="median", min_cells=4):
    stats, problems = [], []
    with rasterio.open(dtm_path) as src:
        if buildings.crs != src.crs:
            raise ValueError(f"CRS mismatch: buildings {buildings.crs}, dtm {src.crs}")
        for i, geom in enumerate(buildings.geometry):
            try:
                arr, _ = mask(src, [geom], crop=True, filled=False)
            except ValueError:
                problems.append((i, "outside the raster"))
                stats.append(np.nan)
                continue
            v = arr.compressed()
            if v.size < min_cells:
                # fall back to a small buffer so slivers still get a value
                arr, _ = mask(src, [geom.buffer(src.res[0])], crop=True, filled=False)
                v = arr.compressed()
            if v.size == 0:
                problems.append((i, "all nodata"))
                stats.append(np.nan)
            else:
                stats.append(float(np.median(v) if stat == "median" else v.min()))
    return np.array(stats), problems

base, problems = terrain_under(b, "dtm.tif")
print(f"{np.isnan(base).sum()} buildings without a terrain value; "
      f"{len(problems)} reported problems")

Example 2 โ€” detect a datum offset rather than guessing

import numpy as np

diff = b["model_base_z"] - b["terrain_z"]
print(f"offset: median {np.nanmedian(diff):.2f} m, "
      f"std {np.nanstd(diff):.2f} m, range {np.nanmin(diff):.1f} to {np.nanmax(diff):.1f}")

if abs(np.nanmedian(diff)) > 20 and np.nanstd(diff) < 2:
    print("constant large offset โ€” this is a vertical datum mismatch, not a sampling problem")
elif np.nanstd(diff) > 1 and abs(np.nanmedian(diff)) < 1:
    print("offset varies and averages to zero โ€” this is terrain slope under a flat base")

A small standard deviation with a large median is the signature of a datum problem; the reverse is a slope problem. Printing both saves an hour of reading code.

Example 3 โ€” a sloping base

import numpy as np, trimesh, rasterio
from shapely.geometry import Polygon

def extrude_on_terrain(poly: Polygon, height, dtm_src):
    """Base follows the terrain under each footprint vertex."""
    coords = np.array(poly.exterior.coords[:-1])
    base_z = np.array([v[0] for v in dtm_src.sample(coords)], dtype=float)

    n = len(coords)
    verts = np.vstack([np.c_[coords, base_z],
                       np.c_[coords, base_z + height]])
    faces = []
    for i in range(n):
        j = (i + 1) % n
        faces += [[i, j, j + n], [i, j + n, i + n]]
    cap = trimesh.creation.triangulate_polygon(poly)[1]
    faces += [[a, b, c] for a, b, c in cap]                      # bottom
    faces += [[a + n, c + n, b + n] for a, b, c in cap]          # top, reversed
    return trimesh.Trimesh(vertices=verts, faces=np.array(faces), process=True)

The reversed winding on the top cap is what keeps the solid watertight; without it the normals disagree and the volume is wrong.

Explanation

Why the geoid separation is the commonest large offset

National height datums are defined relative to mean sea level, which follows the geoid. GNSS and most global terrain products use heights above the WGS84 ellipsoid. The two differ by the geoid separation โ€” roughly 43 m in the Netherlands, 45โ€“55 m over much of western Europe, and negative in parts of the Indian Ocean. Nothing in a file forces the two to be labelled, so the mismatch is invisible until buildings sink.

Why the centroid is a bad sampling point

A building's centroid can be outside its own footprint when the shape is concave, which for a U-shaped block means sampling the courtyard. It can also be a nodata cell in a DTM that masks buildings out. Masking by the polygon and taking a statistic over the cells inside is barely more code and cannot do either.

Why a flat base is a decision and not a bug

Buildings have level floors; the ground does not. A real building on a slope has a plinth, a stepped foundation or a sunken side, and none of that is in a footprint layer. A flat base is therefore a reasonable model โ€” but which flat base matters, and the minimum, median and mean answer different questions.

Why web maps invert the rule

Renderers with terrain support place extrusions on the terrain surface, so the base is an offset from the ground rather than an absolute height. That is the opposite convention from a desktop 3D scene, and code moved between the two without a change buries or floats everything.

Stack showing ellipsoidal height, the geoid separation of about 43 metres, orthometric height and an unlabelled terrain surface.
A constant tens-of-metres offset is almost always this and nothing else.

Edge cases or notes

  • Basements need a negative base. base = terrain โˆ’ depth.
  • DTM versus DSM. Sampling a DSM under a building gives the roof.
  • Some DTMs mask buildings as nodata. Buffer outwards to find ground.
  • Bridges and buildings over water have no terrain beneath them.
  • Check the CRS on both layers, including the vertical component.
  • Terrain exaggeration in viewers is display only.
  • Report the buildings you could not place. NaN is better than zero.
  • Record the statistic used. Minimum and median give different volumes.

FAQ

Why are my extruded buildings floating above the terrain?

Most often no base height was applied, so the prisms start at z = 0, or the model and the terrain use different vertical datums.

What is a 40 m constant offset?

A geoid separation. National height datums and ellipsoidal heights differ by tens of metres โ€” about 43 m in the Netherlands.

Should I use the minimum or the median terrain height?

The median splits the error on a slope; the minimum guarantees the building does not float and buries the uphill side. Pick one and record it.

Why does sampling the centroid fail?

A concave footprint's centroid can lie outside it, and some DTMs mask buildings as nodata. Mask by the polygon and take a statistic of the cells inside.

How do I handle a steep site?

Extrude from a sloping base, sampling the DTM at each footprint vertex, and reverse the winding on the top cap so the solid stays watertight.

Why do my web map buildings sink into a hill?

With a terrain source enabled, fill-extrusion-base is relative to the terrain. A base copied from the flat case is subtracted twice.