Binning Features into Tiles with Atomic Counters
Grouping features by the tile they fall in does not always need a sort. When the number of tiles is small and known — a viewport is rarely more than a few hundred — a two-pass count-then-scatter with atomic counters produces exactly the same grouping in roughly a third of the time, at the cost of an arbitrary order within each tile. For rendering, aggregation and culling that order never matters, which makes this the right default and the sort the exception. This page is the two passes, the buffer layout they produce, and the honest account of what the arbitrary order costs. It is one stage of GPU sorting and prefix sums for spatial data.
Runnable reference implementation
Pass one counts how many features fall in each tile. Pass two scans those counts into starting offsets on the host or in a third small pass, then scatters each feature into its tile’s range.
struct Grid {
origin : vec2<f32>,
inv_size : vec2<f32>, // 1 / tile size, in the same units as the features
dims : vec2<u32>, // tiles across, tiles down
};
@group(0) @binding(0) var<uniform> grid : Grid;
@group(0) @binding(1) var<storage, read> centroids : array<vec2<f32>>;
@group(0) @binding(2) var<storage, read_write> counts : array<atomic<u32>>;
fn tile_of(p : vec2<f32>) -> u32 {
let c = clamp(vec2<u32>((p - grid.origin) * grid.inv_size),
vec2<u32>(0u), grid.dims - vec2<u32>(1u));
return c.y * grid.dims.x + c.x;
}
@compute @workgroup_size(256)
fn count(@builtin(global_invocation_id) gid : vec3<u32>) {
let i = gid.x;
if (i >= arrayLength(¢roids)) { return; }
atomicAdd(&counts[tile_of(centroids[i])], 1u);
}
@group(0) @binding(0) var<uniform> grid : Grid;
@group(0) @binding(1) var<storage, read> centroids : array<vec2<f32>>;
@group(0) @binding(2) var<storage, read> starts : array<u32>; // scanned counts
@group(0) @binding(3) var<storage, read_write> cursor : array<atomic<u32>>;
@group(0) @binding(4) var<storage, read_write> binned : array<u32>;
@compute @workgroup_size(256)
fn scatter(@builtin(global_invocation_id) gid : vec3<u32>) {
let i = gid.x;
if (i >= arrayLength(¢roids)) { return; }
let tile = tile_of(centroids[i]);
// The atomic hands out a slot within this tile's range. Which slot is
// arbitrary; that it is inside the right range is guaranteed.
let slot = atomicAdd(&cursor[tile], 1u);
binned[starts[tile] + slot] = i;
}
The output is one dense array plus a start-and-count table — the same shape a sort would produce, arrived at in two passes rather than four, with the elements within each tile in whatever order the hardware happened to schedule them.
Parameter reference
| Value | Setting here | Guidance |
|---|---|---|
| Tile count | ≤ ~4096 | Counts and cursors are one u32 each; the arrays stay small. Past a few thousand tiles, contention falls and a sort’s fixed cost starts to look reasonable. |
| Passes | 2 (+1 small scan) | Against four for a 32-bit radix sort. |
| Scratch | counts + cursor + output | 2 × tiles words plus one word per feature. No ping-pong pair — the scatter reads and writes different buffers. |
| Cursor reset | every frame | A stale cursor writes past a tile’s range and corrupts its neighbour. |
Clamp in tile_of |
required | A centroid outside the grid otherwise indexes past the counts array. |
What the arbitrary order actually costs
Within a tile, the order is whatever the atomic happened to hand out, and it differs run to run. Three consequences follow, and only one of them is usually a problem.
Rendering is unaffected. A tile’s features are drawn as a set, and the rasterizer does not care in which order opaque geometry arrives; where blending makes order matter, it matters globally rather than within a tile, and a per-tile sort would not have helped anyway.
Aggregation is unaffected. Sums, counts, centroids and density grids are all commutative, so the order in which contributions arrive changes nothing but the floating-point rounding — and even that is only visible at the last bit or two.
Reproducibility is affected, and this is the real cost. A binned array cannot be compared against a stored fixture, because the order will differ. Debugging becomes harder for the same reason: two runs of the same input produce two different arrays, so a diff is uninformative. Where those properties are needed, the answer is the deterministic scan rather than a sort — the scan is two passes as well, and it is the middle option between this page and a full radix sort.
Sizing the grid
The tile grid is a parameter rather than a given, and the two failure modes sit at opposite ends of it.
Too few tiles and the atomics contend. Every feature in a tile increments the same address, so a viewport divided into sixteen cells with two hundred thousand features per cell serialises heavily — the hardware handles the contention correctly and slowly. The symptom is a binning pass whose cost is wildly out of proportion to the number of features, and it gets worse as the data gets denser, which is exactly the wrong direction.
Too many tiles and the tables grow while the benefit disappears. Ten thousand tiles is ten thousand counts and ten thousand cursors, a scan over ten thousand entries, and — more importantly — an average of a few features per tile, at which point the grouping has stopped being useful for anything that wanted a tile’s worth of work.
The band that works is roughly one to four thousand tiles for a viewport, which for a typical window is a grid somewhere between 32 by 32 and 64 by 64. Within that band the contention is spread thinly enough not to matter and each tile still holds enough features to be a useful unit of work.
Two refinements are worth knowing. Privatising the counts per workgroup, exactly as a histogram does, removes contention almost entirely and makes the low end of the range viable again. And matching the binning grid to the rendering tile grid, where one exists, means the binned ranges can be consumed directly by a tile-parallel pass with no remapping.
Failure modes
- A tile’s features leak into its neighbour. The cursor was not reset between frames, so the scatter wrote past the tile’s range. This is the single most common bug here and it produces plausible-looking geometry in the wrong tile.
- Some features are missing. The count pass and the scatter pass used different
tile_ofimplementations — usually because one clamped and the other did not, so a boundary feature was counted in one tile and scattered into another. - Throughput collapses at low zoom. A viewport with few tiles means few atomic addresses and heavy contention. Fix: raise the tile count, or privatise the counts per workgroup exactly as a histogram does.
binnedhas stale entries from the previous frame. The array is written per frame but never cleared, and a tile whose count shrank leaves old indices past its new end. Fix: always read exactlycountentries, never scan to the end of the array.- Counts disagree with the total feature count. Something is being dropped — usually a missing tail guard, or a centroid that is NaN and therefore fails every comparison in the clamp.
Backend / Python interop note
The grid definition is the one thing that has to agree between client and server, and it is worth transmitting rather than deriving.
A tile grid is an origin, a cell size and a pair of dimensions. If the server computes an aggregate per tile — a count for a legend, a density for a choropleth — and the client bins the same features against a grid it derived independently, the two will disagree at boundaries whenever a floating-point comparison lands differently. Sending the grid as four numbers in the response metadata, and having the client use exactly those, removes the disagreement.
There is a related consistency point about the tile assignment itself. A feature is binned by its centroid here, which is cheap and unambiguous, but a polygon whose centroid falls in one tile can extend well into three others. Whether that matters depends entirely on what consumes the bins: a density aggregate is happy with centroids, while a per-tile draw that clips to the tile boundary will drop geometry that should have been visible. Where the extent matters, the fix is to bin by bounding box and accept that a feature appears in several tiles — which changes the counting pass from one increment per feature to one per overlapped tile, and makes the output array longer than the input.
The other server-side note is about what to send when the grid is coarse. A binning pass over four million features to produce a 32-by-32 density grid is a large amount of GPU work for 1024 numbers, and if the same aggregate is needed every frame regardless of the camera, computing it once in pandas or duckdb and shipping the grid directly is both cheaper and simpler. The GPU pass earns its place when the binning depends on the current view — which, for anything interactive, it usually does.
Related
- GPU sorting and prefix sums for spatial data — the topic this page belongs to.
- Implementing a workgroup prefix sum for stream compaction — the deterministic middle option.
- Radix sorting point clouds on the GPU — when the ordering genuinely has to be spatial.
- Parallel histogram binning for GPU heatmaps — the privatisation trick that fixes contention here too.
- Using workgroup_id for parallel tile processing — how the binned ranges are consumed.