Quantized Mesh Protocol

quantized-mesh:// (lib/quantized-mesh-protocol.ts) consumes Cesium terrain — the quantized-mesh-1.0 format — as an ordinary MapLibre raster-dem source, so 3D terrain and every visualization mode treat it like any other DEM.

It is the most involved of the protocols here, because unlike LERC or float32dem the source is not a grid at all.

The direction the meshing libraries do not go

Worth stating first, because the obvious tools are the wrong ones.

pymartini and pydelatin — and their JavaScript ports @mapbox/martini and delatin — are terrain mesh generation tools. They take a raster heightmap and produce a TIN. That is what a titiler-style service uses to serve quantized mesh from a DEM, and it is why titiler's own algorithm list (hillshade, slope, contours, terrarium, terrainrgb, …) contains no quantized-mesh entry: the meshing lives in a separate extension, and it points the other way.

Consuming a TIN needs the inverse — mesh → raster — which is plain triangle rasterisation. No package ships it because it is a loop, and that loop is rasterise() in this file.

@loaders.gl/terrain's QuantizedMeshLoader does the decoding, which is the part that genuinely benefits from a library.

Two URL shapes

quantized-mesh://ion/<assetId>/{z}/{x}/{y}     a Cesium ion asset
quantized-mesh://<host+path>/{z}/{x}/{y}       any other service

{z}/{x}/{y} are MapLibre's own Web Mercator coordinates. Quantized mesh is tiled geographically, so the mapping happens inside the protocol rather than in the template — see below.

An ion asset needs two calls: api.cesium.com/v1/assets/<id>/endpoint returns a tile base URL plus a short-lived access token. That pair is cached per asset and re-resolved once on a 401, so an expired session token recovers without a reload.

The tile mapping

Quantized mesh uses the geographic tiling scheme: EPSG:4326, TMS row order (row 0 at the south), and level L holds 2^(L+1) × 2^L tiles over the whole −180…180 / −90…90 world. Two consequences:

• Level. A Mercator tile at zoom z spans 360/2^z degrees of longitude; a geographic tile at level L spans 360/2^(L+1). They match when L = z − 1.
• Count. A geographic tile is square in degrees, while a Mercator tile is shorter in latitude than in longitude. So one output tile usually falls inside a single source tile, and sometimes straddles two. Both are fetched and rasterised into the same output buffer.

Rasterising

The trick that keeps this cheap is QuantizedMeshLoader's bounds option. The mesh's vertices are normalised within their tile; passing the tile's own geographic rectangle maps them straight to (lon, lat, height):

QuantizedMeshLoader.parseSync(buf, {
  "quantized-mesh": { bounds: [b.west, b.south, b.east, b.north], skirtHeight: null },
})

The Mercator projection then costs one operation per vertex rather than a resampling pass over every pixel:

px[v] = ((lon - out.west) / lonSpan) * TILE_SIZE
py[v] = ((out.northY - lat2merc(lat)) / ySpan) * TILE_SIZE

Each triangle is then scan-converted over its pixel bounding box with barycentric weights, interpolating height. Forward rasterisation is deliberate: the alternative — asking, per output pixel, which triangle contains it — is a point-in-TIN query over thousands of triangles, and would also have to do the projection per pixel instead of per vertex.

The scan-conversion, step by step

The loop is short but every line of it is load-bearing, so here it is in full.

1. Project once, per vertex. Before touching a triangle, all vertices are converted to the output tile's pixel space and cached in two Float32Arrays. Doing this per triangle instead would repeat the work for every triangle a vertex belongs to — six on average in a TIN.

2. Signed area as the guard and the divisor. For each triangle:

const area = (bx - ax) * (cy - ay) - (cx - ax) * (by - ay)
if (area === 0 || !Number.isFinite(area)) continue
const inv = 1 / area

Twice the signed area. Zero means the three points are collinear — a degenerate triangle covers no pixels, and dividing by it would produce NaN weights that then fail the containment test in unpredictable ways. Its reciprocal is taken once and reused, because a division per pixel is meaningfully slower than a multiply.

Note that the sign is not checked: back-facing triangles are rasterised too. A TIN has no consistent winding after projection, and discarding by winding would punch holes.

3. Bounding box, clipped to the tile. Only pixels inside the triangle's own box are visited, clamped to the tile so a triangle straddling the edge costs nothing extra:

const x0 = Math.max(0, Math.floor(Math.min(ax, bx, cx)))
const x1 = Math.min(TILE_SIZE - 1, Math.ceil(Math.max(ax, bx, cx)))

This is what keeps the cost proportional to area covered rather than triangles × pixels. A tile of a few thousand triangles scan-converts in single -digit milliseconds.

4. Barycentric weights at the pixel centre. For each candidate pixel, the same cross-product form as the area, evaluated against (col + 0.5, row + 0.5):

const w0 = ((bx - cxPix) * (cy - cyPix) - (cx - cxPix) * (by - cyPix)) * inv
const w1 = ((cx - cxPix) * (ay - cyPix) - (ax - cxPix) * (cy - cyPix)) * inv
const w2 = 1 - w0 - w1
if (w0 < -1e-6 || w1 < -1e-6 || w2 < -1e-6) continue
heights[i] = w0 * az + w1 * bz + w2 * cz

The pixel centre, not its corner, is what makes a pixel belong to exactly one triangle in the common case. w2 is derived rather than computed, since the three weights sum to 1 by definition — one fewer cross product per pixel.

All three non-negative means the point is inside. The weights are then the interpolation coefficients for free: a plane through three heights, evaluated at that pixel. That is what makes this linear interpolation across the triangle rather than nearest-vertex.

5. The -1e-6 tolerance. Two triangles sharing an edge produce weights that are exactly zero on that edge in exact arithmetic. In floating point one of them lands at -1e-9 and neither claims the pixel, leaving a one-pixel unwritten seam along every shared edge — which a hillshade turns into a visible wireframe. The tolerance is large enough to absorb that error and small enough that a pixel genuinely outside a triangle is still rejected. The cost is that edge pixels may be written twice, by both neighbours; they agree to within rounding, so the second write is harmless.

6. filled is a separate array. Height 0 is a perfectly valid elevation, so "was this pixel written" cannot be inferred from the height buffer. The parallel Uint8Array is what the hole pass reads afterwards to decide between a real 0 m and an untouched pixel.

Pixels no triangle covered are holes: 0 m with alpha 254, the app-wide convention. A tile with no coverage at all throws a 404, which MapLibre treats as an ordinary absent tile.

Cesium World Terrain: measured

The datum is the finding worth carrying. Six points across one tile read +48.1 m mean against an orthometric reference (spread 22.8 m, which is sampling noise between two models in steep ground). That is the alpine geoid separation almost exactly — so Cesium World Terrain behaves like ArcticDEM and REMA, and unlike the Esri LERC sources, which measured −0.6 m on the same points and are orthometric.

In practice: summits read roughly 50 m high in the Alps, and the offset varies by tens of metres around the world. Fine for shape and shading, wrong for absolute heights unless you correct it.

Compared with MapLibre's own example

MapLibre merged an official quantized-mesh terrain example in #7958 (August 2026, against the still-open #4493). It uses the maplibre-gl-3dtiles-terrain plugin, and it is worth knowing how it differs — it is the closest thing to a reference implementation, and it is better than this one in two specific places.

Where theirs is better, and worth taking:

1. The availability walk. Reading available means never firing a request for a level that does not exist, instead of learning by 404. (It is not a universal fix: swisstopo's layer.json carries no available array at all, so the blind fallback is still needed as a backstop.)
2. A decoded-mesh cache. At detail offset +1 four adjacent output tiles share source tiles, and this app re-decodes each one. That is a real, measurable win left on the table.

Where this one is better:

1. Forward rasterisation needs no index. Theirs builds a spatial index per tile and then queries it once per pixel; ours visits only the pixels a triangle actually covers, with no structure to build.
2. Projection per vertex, not per pixel. A tile has a few thousand vertices and 65 536 pixels.
3. Integration. ion endpoint resolution and token handling, the alpha-254 hole convention shared with every other elevation path, the shared tile cache, and returning an ImageBitmap rather than encoded bytes.

Their worker matters less than it sounds now that the PNG encode is gone — the scan-conversion itself is single-digit milliseconds, and the protocol is network-bound.

Is it worth using?

Honest answer: usually not, over lerc://. Cesium World Terrain tops out at level 15 — comparable to the z16 the keyless Esri services already reach — from a broadly similar blend of public sources, and it needs a token and a geoid correction.

Where it earns its place is ion-hosted private assets: terrain you or your organisation uploaded to ion, which nothing else in this app can reach. The protocol takes any asset id, so that is a matter of pointing a custom source at ion/<your asset>.