Geometry & 3D Vision

Voxel Grid

A regular 3D grid of small cubes that stores occupancy, signed distance or other values per cell, and the dense, octree, hashed and VDB structures that make such grids affordable at high resolution.

intermediate

A voxel grid divides a box of 3D space into equal cubes, called voxels (volume elements), and stores one value per cube: whether it is occupied, the distance to the nearest surface, a colour, a semantic label or a feature vector. It is the 3D counterpart of the pixel grid of a digital image, and it is the form in which robots keep maps of free and occupied space, depth cameras fuse surfaces, and 3D networks consume and predict shape. Its weakness is memory, which grows with the cube of the resolution, so most practical voxel grids are stored sparsely: in an octree such as OctoMap, in blocks found through a hash table, as in voxel hashing, or in a shallow tree of blocks, as in VDB. The 3D representations article compares voxel grids with the other ways of storing shape; this article covers the grid itself, what its cells hold and how it is stored.

Definition

A voxel grid is a set of cubic cells of side hh, the voxel size or resolution, aligned with the axes of a coordinate frame and indexed by integer triples (i,j,k)(i, j, k), together with a value for each cell. A dense grid stores a value for every cell of a bounding box. A sparse grid stores values only for some cells and answers queries about the rest with a default, such as “unknown”, “free” or a constant distance.

The word occupancy grid names a voxel (or, in 2D, pixel) grid whose cells hold the belief that the cell is occupied. A grid whose cells hold truncated signed distances is a TSDF volume, the map of TSDF fusion. An octree is a tree in which each node covers a cube and has up to eight children covering its octants. Voxel hashing stores small dense blocks of voxels and finds them through a hash table keyed by integer block coordinates. VDB is a sparse grid organised as a shallow, wide tree of dense blocks (Museth, 2013).

Intuition

Think of a room filled with sugar cubes, each labelled with a number. To ask what is at a point, find the cube that contains it and read its label: no search is needed, because the cube’s position follows from the point’s coordinates by division and rounding. That constant-time lookup, and the regular neighbourhood of every cell (six face neighbours, 26 in all), make grids easy to update, query and convolve.

The cost is that the room is mostly air and the inside of solid things. A surface is two-dimensional, so it passes through a vanishing fraction of the cubes as they get smaller, yet a dense grid stores every cube. Sparse structures keep the convenience of the grid for the cubes that matter and summarise the rest: an octree merges large uniform regions into single nodes, and block-based structures allocate small dense blocks only near surfaces.

Formal Definition

Let the grid have origin o∈R3\mathbf{o} \in \mathbb{R}^3 and voxel size hh. A point x\mathbf{x} falls in the voxel with index

(i,j,k)=⌊x−oh⌋,(i, j, k) = \left\lfloor \frac{\mathbf{x} - \mathbf{o}}{h} \right\rfloor ,

taken componentwise, whose centre is o+h (i+12,j+12,k+12)\mathbf{o} + h\,(i + \tfrac12, j + \tfrac12, k + \tfrac12). A dense grid of nx×ny×nzn_x \times n_y \times n_z voxels stores its values in an array, usually in row-major order, so the voxel (i,j,k)(i, j, k) is at linear position (i ny+j) nz+k(i\,n_y + j)\,n_z + k. With bb bytes per voxel its memory is

Mdense=nxnynz b=Vh3 bM_{\text{dense}} = n_x n_y n_z \, b = \frac{V}{h^3}\, b

for a box of volume VV: halving hh multiplies it by eight. A surface of area AA meets on the order of A/h2A / h^2 voxels, so a structure whose size follows the surface grows by about four per halving (reasoning from the surface’s dimension, as in the 3D representations article).

Octree. A cube of side 2Dh2^D h is the root; a node at depth dd covers a cube of side 2D−dh2^{D-d} h and is either a leaf or has children in its octants, down to depth DD, where cells have side hh. The index bits of (i,j,k)(i, j, k) read from the most significant end choose the child at each level, so a lookup costs O(D)O(D). OctoMap’s authors bound DD at 16, which at h=1h = 1 cm covers a cube of 216×0.01=655.362^{16} \times 0.01 = 655.36 m on a side (Hornung et al., 2013, §3.1, §4.2).

Hashed blocks. Space is divided into blocks of B3B^3 voxels; block (x,y,z)=⌊(i,j,k)/B⌋(x, y, z) = \lfloor (i, j, k) / B \rfloor is stored only if allocated, and is found through a hash table. Voxel hashing uses B=8B = 8 and the hash

H(x,y,z)=(x p1⊕y p2⊕z p3) mod n,H(x, y, z) = (x\,p_1 \oplus y\,p_2 \oplus z\,p_3) \bmod n ,

with ⊕\oplus the bitwise exclusive or, nn the table size and the primes p1=73,856,093p_1 = 73{,}856{,}093, p2=19,349,669p_2 = 19{,}349{,}669, p3=83,492,791p_3 = 83{,}492{,}791, which the authors take from Teschner et al. (2003) (Nießner et al., 2013, §4). Museth (2013, §3.1) prints the second constant as 19,349,663 for a hash he also attributes to Teschner et al.; by our check that number is not prime (it is divisible by 41), while 19,349,669 is, although Museth too calls his constants large primes.

VDB. A tree of fixed depth whose node sizes are powers of two: in the configuration Museth recommends, a sparse root, two levels of internal nodes with 32332^3 and 16316^3 children, and leaves of 838^3 voxels, so each child of the root spans 409634096^3 voxels (§2.4). A voxel’s index is split into bit fields, one per level, and each level below the root is a direct array lookup.

Properties

What a voxel stores

Occupancy. A binary grid needs one bit per voxel, but mapping from noisy sensors needs a belief. OctoMap stores for each leaf the log-odds L(n)=log⁡P(n)1−P(n)L(n) = \log\frac{P(n)}{1 - P(n)} of being occupied and adds the log-odds of each new measurement (Eq. 2), the update the 3D occupancy prediction article derives from Bayes’ rule. The authors attribute this occupancy grid mapping to Moravec and Elfes (1985); Elfes (1989) is another early paper on such grids. Neither could be read for this article, and nothing here is attributed to them beyond the attribution OctoMap makes.

Two details of OctoMap matter for the grid itself. Its inverse sensor model adds locc=0.85l_{\text{occ}} = 0.85 to the voxel containing a beam’s endpoint and lfree=−0.4l_{\text{free}} = -0.4 to the voxels the beam traverses, found by a 3D variant of Bresenham’s line algorithm after Amanatides and Woo (§5.1). And it clamps the log-odds between lmin⁡=−2l_{\min} = -2 and lmax⁡=3.5l_{\max} = 3.5, probabilities 0.12 and 0.97 (Eq. 4, after Yguel et al.; values from §5.1): without a bound, a voxel seen free kk times needs kk occupied readings before it changes state, so the map adapts slowly to change (§3.2). With these values five agreeing measurements make an unknown voxel stable (§3.4). The same defaults (hit probability 0.7, miss 0.4, clamping 0.1192 and 0.971, occupancy threshold 0.5) are set in the release v1.10.1 source. OctoMap also represents three states: occupied, free, and unknown for space never initialised, which a planner must treat differently from free space (§1, §3.1).

Distances and weights. A TSDF volume stores in each voxel a truncated signed distance and an accumulated weight, fused from depth maps by a running average (TSDF fusion); the surface is the zero level set of the signed distance function, extracted by marching cubes or ray casting. Voxel hashing packs a voxel into 8 bytes: a 32-bit float distance, three bytes of colour and one byte of weight (Nießner et al., 2013, §4).

Other values. OctoMap’s nodes can be extended with other data, such as an average colour per voxel (§3.5.1); semantic occupancy grids store a class per voxel; networks store feature vectors. In VDB every voxel or tile also carries an active bit, whose meaning the application chooses: for a narrow-band level set, voxels within the band are active, and the rest hold a constant distance whose sign says inside or outside (Museth, 2013, §2.2).

Resolution and memory

The voxel size sets the smallest structure the grid can hold, and memory sets the voxel size. OctoMap’s Table 2 shows the scale on a 32-bit system: a minimal dense grid of 4-byte cells for the Freiburg campus scans, 292 × 167 × 28 m at 10 cm, would take 5,162.90 MB, against 990.66 MB for OctoMap’s pruned octree, 504.76 MB after maximum-likelihood compression, and a 13.82 MB lossy file; the authors note that such a grid at 10 cm would not fit in the addressable memory of a 32-bit machine (§5.4). For TSDF volumes, Nießner et al. (2013, §9) report on average less than 300 MB of surface data (less than 600 MB with colour) on their test scenes at 8 mm voxels, against well over 5 GB, including colour, for a regular grid of the same voxel size and an extent of 8 m in depth. By our arithmetic, their average of about 140,000 blocks of 838^3 voxels takes some 290 MB at 4 bytes per voxel and 570 MB at the paper’s 8-byte voxel, which already holds colour, so the lower figure fits a voxel without colour and weight, which the paper does not describe.

Resolution beyond the data buys nothing: Nießner et al. saw no visible improvement with voxels below 2 mm with their depth cameras (§9). Coarse resolution loses thin structures and, for occupancy, places surfaces only to within a voxel; the 3D representations article measures this against a stored distance, which places the surface inside the voxel.

Dense grids

A dense array is the simplest structure and has the fastest access: one multiply-add per lookup, regular neighbours for finite differences and convolution, and trivially parallel updates. KinectFusion used a dense TSDF volume on the GPU, which limits the mapped space to a room. OctoMap’s authors list two further drawbacks for mapping: the grid must be allocated for the whole bounding box before mapping, so the extent must be known in advance or the grid copied when it grows (§2). Dense grids remain the norm for learned 3D models at modest resolution, such as the 32332^3 grids of 3D-R2N2 on ShapeNet and the semantic grids of occupancy prediction benchmarks.

Octrees and OctoMap

An octree subdivides a cube into eight only where the contents are not uniform, so large free or solid regions become single leaves. OctoMap’s authors write that octrees for mapping were first proposed by Meagher (1982), whose paper was not read here, and they credit Payeur et al. with bringing probabilistic occupancy to octrees (§2).

Pruning. If all eight children of a node are leaves in the same state, they are replaced by the parent (§3.1). With probabilities, neighbouring voxels rarely hold identical values, so OctoMap prunes only children that have reached the clamping bounds, which loses information only near probabilities 0 and 1; combining clamping and pruning improved compression by up to 44% in their experiments (§3.4).

Multi-resolution queries. Inner nodes store a summary of their children, so a query can stop above the leaves and read a coarser map. OctoMap uses the maximum child log-odds (Eq. 6), which the authors call conservative and suited to navigation: a coarse cell is occupied if any part of it was measured occupied (§3.3).

Node layout. Each node stores its value and one pointer to an array of eight child pointers, allocated only for inner nodes; positions and sizes are not stored but recovered during traversal. On a 64-bit system an inner node takes 80 bytes and a leaf 16 bytes (40 and 8 on 32-bit), and in the authors’ robotics data 80–85% of nodes were leaves (§4.1). The v1.10.1 code computes memory the same way: a per-node size for every node plus an array of eight pointers for each inner node (memoryUsage() in OcTreeBaseImpl.hxx). For files, a maximum-likelihood map is written as a bit stream of two bits per child (unknown, occupied, free, or inner node), 16 bits per inner node, with no coordinates (§4.3).

Costs. A lookup walks the tree, O(D)O(D) steps. In very confined spaces a dense grid can take less memory than an uncompressed octree (Hornung et al., §5.4). On GPUs, Nießner et al. (2013, §2) argue that deep trees cause thread divergence during traversal, which motivated their flat alternative.

Hashed voxel blocks

Voxel hashing (Nießner et al., 2013) keeps the regular grid inside small blocks and makes the set of blocks sparse. Its table is divided into buckets of a few entries; each hash entry holds the block’s integer position, an offset for overflow and a pointer to the block’s voxels, 12 bytes in all (§4). On a collision the entry goes to the next free slot in the bucket; when the bucket is full, a linked list chains further entries stored elsewhere in the table (§4.1). Blocks are allocated by walking each depth sample’s ray through the truncation band with a voxel traversal and inserting every block it touches (§5). Blocks whose largest weight is zero or whose smallest distance exceeds a threshold are garbage-collected, and blocks leaving a sphere that contains the camera’s view frustum are streamed to host memory in chunks of 1 m³ and back when the camera returns (§6, §8).

The authors used a table of 2212^{21} entries with buckets of two, about one million buckets, and a 1 GB heap for 2182^{18} blocks. About 140,000 blocks were allocated at 8 mm voxels, bucket overflow was 0.1%, and no linked list was longer than three entries (§9). Hashing gives average constant-time access with no hierarchy to traverse, but no coarser levels, and neighbouring blocks are not adjacent in memory.

VDB

VDB (Museth, 2013) combines blocks and a tree. Leaves are dense 838^3 blocks; internal nodes are dense arrays of child pointers or tile values, a single value standing for a whole uniform child region; the root is a sparse hash map or ordered map, so the index space is unbounded within the integer range (§2.3). Bit masks in each node record which children exist and which values are active, which supports iteration over active voxels, topology operations such as dilation, and compression (§2.3). Because the depth is fixed, random access costs one root lookup and a constant number of array lookups, and a cache of recently visited nodes makes spatially coherent access cheaper still (§3).

Museth’s Table II compares configurations on a narrow-band level set of a sphere at an effective resolution of 409634096^3 with a band of ten voxels (51,033,829 active voxels). The recommended configuration with a hash-map root and nodes of 32332^3, 16316^3 and 838^3 took 420 MB; two octrees took 390 MB but needed 50 and 70 time units per lookup, against 1 for VDB’s cached access and 10 for its top-down traversal in the configuration with a dense root of 64364^3; a flat tiled grid with an internal node of 102431024^3 entries took 8,723 MB, because of its dense top-level table (§5.1). OpenVDB, the open-source library released at SIGGRAPH 2012 (Museth, 2013, §1.1), uses this four-level configuration for its standard tree types, with 32-bit signed indices, and its narrow-band level sets keep active voxels normally three voxels wide on either side of the surface (OpenVDB 14.0.0 overview, as of November 2026).

Voxel downsampling

A voxel grid is also a filter for point clouds: assign each point to a voxel and replace the points in each occupied voxel by one representative, as the point cloud article describes and demonstrates. Open3D’s voxel_down_sample (version 0.19.0 source) does this with a hash map from integer voxel index to an accumulator, so the grid it builds is itself a hashed sparse voxel grid. It averages the points of each voxel, and also their normals, colours and covariances when present. The grid’s origin is the cloud’s minimum bound minus half a voxel, so the voxel boundaries depend on the cloud’s own extent: by our reading of the code, adding one distant point can shift every boundary and change the output. The averaged normals are not renormalised, as a comment in the code notes.

Examples

The script stores one analytic scene, a corner of a room in a cube of side 2.56 m with a floor 8 cm thick, a wall 10 cm thick, a box, a sphere and a pole of radius 3 cm, at voxel sizes of 8, 4, 2 and 1 cm, and counts the memory of each structure from what the code builds:

  • Dense: occupancy packed one bit per voxel, and a dense array of voxel hashing’s 8-byte voxel.
  • Octree: the binary occupancy (every voxel known, no unknown space), pruned wherever eight children are leaves with the same value; bytes at OctoMap’s 64-bit node sizes (80 and 16) and in its maximum-likelihood bit stream (two bytes per inner node). Expanding the leaves must give back the dense grid.
  • Hashed blocks: 838^3 blocks of 8-byte voxels allocated wherever a voxel lies within a truncation band μ=4h\mu = 4h of the surface, inserted in a hash table with voxel hashing’s function (evaluated in Python’s unbounded integers) and buckets of two, and read back through the table. The band width and the number of buckets (the smallest power of two at least four times the number of blocks) are our choices; the sources adapt the band to sensor noise and fix the table in advance.
  • VDB: the same leaf blocks under internal nodes of 16316^3 and 32332^3, with node sizes summed from the fields in Museth’s listings (§2.3) for an 8-byte value, and 28 bytes per root entry, our estimate.

Accuracy is measured three ways: the occupied volume against the analytic volume (with the pole’s voxelized volume against its true volume), an exact read-back of the stored distances, and coarse queries on the 1 cm octree with OctoMap’s maximum rule against grids sampled directly at the coarser voxel size. The script uses NumPy 2.3.5 and runs in a few seconds.

import sys
import numpy as np

# Scene in the cube [0, 2.56]^3 m, z up: floor, wall, box, sphere, thin pole.
SIDE = 2.56
BOXES = {"floor": (np.array([1.28, 1.28, 0.04]), np.array([1.28, 1.28, 0.04])),
         "wall":  (np.array([0.05, 1.28, 1.32]), np.array([0.05, 1.28, 1.24])),
         "box":   (np.array([1.50, 1.20, 0.40]), np.array([0.30, 0.40, 0.30]))}
SPH_C, SPH_R = np.array([0.80, 1.80, 0.50]), 0.30
POLE_XY, POLE_R, POLE_Z = np.array([1.90, 0.50]), 0.03, (0.08, 2.00)

def box_sdf(p, c, h):
    q = np.abs(p - c) - h
    return np.linalg.norm(np.maximum(q, 0), axis=-1) + np.minimum(q.max(-1), 0)

def pole_sdf(p):
    zc, hz = (POLE_Z[0] + POLE_Z[1]) / 2, (POLE_Z[1] - POLE_Z[0]) / 2
    q = np.stack([np.linalg.norm(p[..., :2] - POLE_XY, axis=-1) - POLE_R,
                  np.abs(p[..., 2] - zc) - hz], -1)
    return np.minimum(q.max(-1), 0) + np.linalg.norm(np.maximum(q, 0), axis=-1)

def scene_sdf(p):
    f = np.linalg.norm(p - SPH_C, axis=-1) - SPH_R
    for c, h in BOXES.values():
        f = np.minimum(f, box_sdf(p, c, h))
    return np.minimum(f, pole_sdf(p))

POLE_VOL = np.pi * POLE_R**2 * (POLE_Z[1] - POLE_Z[0])
TRUE_VOL = (sum(8 * np.prod(h) for _, h in BOXES.values())
            + 4 / 3 * np.pi * SPH_R**3 + POLE_VOL)   # primitives do not overlap

def sample(n, fn=scene_sdf):
    """fn at the n^3 voxel centres, slab by slab (float64)."""
    h = SIDE / n
    c = (np.arange(n) + 0.5) * h
    Y, Z = np.meshgrid(c, c, indexing="ij")
    out = np.empty((n, n, n))
    for i in range(n):
        out[i] = fn(np.stack([np.full_like(Y, c[i]), Y, Z], -1))
    return out

# ---- octree over a binary occupancy grid (all space known) ----
def build_octree(occ):
    """Bottom-up: a node is a leaf if its 8 children are leaves with one value.
    Returns per-level (uniform, value, any_occupied) arrays, finest first."""
    levels = [(np.ones(occ.shape, bool), occ.copy(), occ.copy())]
    while levels[-1][0].shape[0] > 1:
        u, v, a = levels[-1]
        m = u.shape[0] // 2
        blk = lambda x: x.reshape(m, 2, m, 2, m, 2).transpose(0, 2, 4, 1, 3, 5).reshape(m, m, m, 8)
        U, V, A = blk(u), blk(v), blk(a)
        uni = U.all(-1) & (V.all(-1) | ~V.any(-1))
        levels.append((uni, V[..., 0], A.any(-1)))   # inner value: max of children (OctoMap Eq. 6)
    return levels

def count_nodes(levels):
    """A node exists if it is the root or its parent is an inner (non-uniform) node."""
    inner = leaves = 0
    exists = np.ones((1, 1, 1), bool)
    for uni, _, _ in reversed(levels):
        if exists.shape != uni.shape:
            exists = parent_inner.repeat(2, 0).repeat(2, 1).repeat(2, 2)
        inner += int((exists & ~uni).sum())
        leaves += int((exists & uni).sum())
        parent_inner = exists & ~uni
    return inner, leaves

def expand_leaves(levels):
    """Paint every leaf's value back onto the finest grid."""
    n = levels[0][0].shape[0]
    out = np.zeros((n, n, n), bool)
    painted = np.zeros((n, n, n), bool)
    for k in range(len(levels) - 1, -1, -1):
        uni, val, _ = levels[k]
        s = 2**k
        up = lambda x: x.repeat(s, 0).repeat(s, 1).repeat(s, 2)
        leaf = up(uni) & ~painted
        out[leaf] = up(val)[leaf]
        painted |= leaf
    return out

# ---- voxel hashing (Niessner et al. 2013, Sec. 4): 8^3 blocks, buckets of 2 ----
P1, P2, P3 = 73856093, 19349669, 83492791
VOXEL = np.dtype([("sdf", "<f4"), ("rgb", "u1", 3), ("w", "u1")])   # 8 bytes
ENTRY_BYTES = 12                                                    # short[3], short, int

class HashTable:
    def __init__(self, n_buckets, bucket_size=2):
        self.nb, self.bs = n_buckets, bucket_size
        self.pos = np.full((n_buckets * bucket_size, 3), -1, np.int64)
        self.ptr = np.full(n_buckets * bucket_size, -1, np.int64)
        self.nxt = np.full(n_buckets * bucket_size, -1, np.int64)    # overflow chain
        self.overflow = 0

    def h(self, x, y, z):
        return ((x * P1) ^ (y * P2) ^ (z * P3)) % self.nb

    def insert(self, key, ptr):
        b = self.h(*key) * self.bs
        for s in range(b, b + self.bs):
            if self.ptr[s] < 0:
                self.pos[s], self.ptr[s] = key, ptr
                return
        tail = b + self.bs - 1                                       # list head: last slot
        while self.nxt[tail] >= 0:
            tail = self.nxt[tail]
        s = (b + self.bs) % len(self.ptr)                            # linear search for a free
        while self.ptr[s] >= 0 or s % self.bs == self.bs - 1:        # slot, skipping last slots
            s = (s + 1) % len(self.ptr)
        self.pos[s], self.ptr[s], self.nxt[tail] = key, ptr, s
        self.overflow += 1

    def find(self, key):
        """Returns (pointer or -1, entries compared)."""
        b = self.h(*key) * self.bs
        probes = 0
        for s in range(b, b + self.bs):
            probes += 1
            if self.ptr[s] >= 0 and tuple(self.pos[s]) == key:
                return self.ptr[s], probes
        s = self.nxt[b + self.bs - 1]
        while s >= 0:
            probes += 1
            if tuple(self.pos[s]) == key:
                return self.ptr[s], probes
            s = self.nxt[s]
        return -1, probes

    def chain_lengths(self):
        out = []
        for b in range(self.nb):
            k, s = 0, self.nxt[b * self.bs + self.bs - 1]
            while s >= 0:
                k, s = k + 1, self.nxt[s]
            out.append(k)
        return np.array(out)

# ---------------------------------------------------------------- experiment
NS = [int(a) for a in sys.argv[1:]] or [32, 64, 128, 256]

print(f"true occupied volume {TRUE_VOL:.4f} m^3, pole {POLE_VOL * 1e3:.2f} litres")
print("\n  h(cm)  dense 1-bit  dense 8 B/voxel   octree inner   leaves  "
      "octree bytes  ML stream  vol err f<0  f<=0   pole")
fine = {}
for n in NS:
    h = SIDE / n
    f = sample(n)
    occ = f < 0
    levels = build_octree(occ)
    inner, leaves = count_nodes(levels)
    assert np.array_equal(expand_leaves(levels), occ)                # lossless
    dense_tsdf = np.zeros((n, n, n), VOXEL)
    pole = (sample(n, pole_sdf) < 0).sum() * h**3 / POLE_VOL
    vol_err = occ.sum() * h**3 / TRUE_VOL - 1
    vol_le = (f <= 0).sum() * h**3 / TRUE_VOL - 1
    print(f"{100 * h:7.0f}  {np.packbits(occ).nbytes:11,d}  {dense_tsdf.nbytes:15,d}  "
          f"{inner:13,d}  {leaves:7,d}  {80 * inner + 16 * leaves:12,d}  {2 * inner:9,d}  "
          f"{vol_err:+11.2%}  {vol_le:+6.2%}  {pole:5.0%}")
    fine[n] = (f, occ, levels, dense_tsdf)

print("\n  h(cm)  band voxels  blocks  buckets  overflow  max chain  probes  "
      "hashed bytes  VDB leaf/int1/int2  VDB bytes  max|diff|")
for n in NS:
    h = SIDE / n
    f, occ, levels, dense = fine[n]
    mu = 4 * h                                                      # truncation: 4 voxels
    band = np.abs(f) < mu
    dense["sdf"] = np.clip(f, -mu, mu)
    dense["w"] = band
    m = n // 8
    blk_band = band.reshape(m, 8, m, 8, m, 8).any((1, 3, 5))
    keys = [tuple(int(v) for v in k) for k in np.argwhere(blk_band)]
    nb = 1 << int(np.ceil(np.log2(4 * len(keys))))
    table = HashTable(nb)
    heap = np.zeros((len(keys), 8, 8, 8), VOXEL)
    for i, (x, y, z) in enumerate(keys):
        table.insert((x, y, z), i)
        heap[i] = dense[8 * x:8 * x + 8, 8 * y:8 * y + 8, 8 * z:8 * z + 8]
    diff, probes = 0.0, []
    for x, y, z in keys:                                              # read back through the table
        p, k = table.find((x, y, z))
        probes.append(k)
        got = heap[p]["sdf"]
        diff = max(diff, float(np.abs(got - dense["sdf"][8 * x:8 * x + 8, 8 * y:8 * y + 8, 8 * z:8 * z + 8]).max()))
    hashed = heap.nbytes + len(table.ptr) * ENTRY_BYTES
    kb = np.array(keys)
    l1 = len({tuple(k) for k in kb // 16})                          # 16^3 children per int1
    l2 = len({tuple(k) for k in kb // 512})                         # 32^3 children per int2
    vdb = (len(keys) * (512 * 8 + 64 + 8 + 8) + l1 * (4096 * 8 + 2 * 512 + 12)
           + l2 * (32768 * 8 + 2 * 4096 + 12) + l2 * 28)
    print(f"{100 * h:7.0f}  {band.sum():11,d}  {len(keys):6,d}  {nb:7,d}  {table.overflow:8,d}  "
          f"{table.chain_lengths().max():9d}  {np.mean(probes):6.3f}  {hashed:12,d}  "
          f"{len(keys):>8,d}/{l1}/{l2}  {vdb:10,d}  {diff:9.1e}")

# Multi-resolution queries on the finest octree (OctoMap's max rule) against
# grids built directly at the coarser resolution; volumes as % of the fine occupied volume.
n = NS[-1]
f, occ, levels, _ = fine[n]
vol = occ.sum()
print(f"\nqueries on the {100 * SIDE / n:.0f} cm octree: cell(cm)  max-rule missed  false occ  "
      f"| direct grid missed  false occ")
for c in (1, 2, 3):
    s = 2**c
    up = lambda x: x.repeat(s, 0).repeat(s, 1).repeat(s, 2)
    q = up(levels[c][2])                                            # inner-node max values
    g = up(sample(n // s) < 0)                                      # grid sampled at 2^c h
    print(f"{'':36s}{100 * SIDE / n * s:6.0f}  {(occ & ~q).sum() / vol:15.2%}  {(q & ~occ).sum() / vol:9.1%}  "
          f"| {(occ & ~g).sum() / vol:18.2%}  {(g & ~occ).sum() / vol:9.1%}")

Output:

true occupied volume 1.5657 m^3, pole 5.43 litres

  h(cm)  dense 1-bit  dense 8 B/voxel   octree inner   leaves  octree bytes  ML stream  vol err f<0  f<=0   pole
      8        4,096          262,144            756    5,293       145,168      1,512       -7.72%  -7.72%   226%
      4       32,768        2,097,152          1,189    8,324       228,304      2,378       -9.60%  +6.62%    57%
      2      262,144       16,777,216          7,361   51,528     1,413,328     14,722       -0.05%  -0.05%    57%
      1    2,097,152      134,217,728         10,845   75,916     2,082,256     21,690       +0.05%  +0.05%   113%

  h(cm)  band voxels  blocks  buckets  overflow  max chain  probes  hashed bytes  VDB leaf/int1/int2  VDB bytes  max|diff|
      8       13,397      47      256         0          0   1.149       198,656        47/1/1     500,452    0.0e+00
      4       64,220     186    1,024        13          2   1.376       786,432       186/1/1   1,080,916    0.0e+00
      2      340,792   1,027    8,192         1          1   1.042     4,403,200     1,027/1/1   4,592,932    0.0e+00
      1    1,837,528   5,162   32,768        13          2   1.067    21,929,984     5,162/7/1  22,063,516    0.0e+00

queries on the 1 cm octree: cell(cm)  max-rule missed  false occ  | direct grid missed  false occ
                                         2            0.00%       0.7%  |              0.30%       0.2%
                                         4            0.00%      11.7%  |              9.92%       0.3%
                                         8            0.00%      33.1%  |             10.09%       2.3%

Dense against sparse. Both dense columns grow eightfold per halving of the voxel size. The hashed blocks fall from about three quarters of the dense 8-byte size at 8 cm to a sixth at 1 cm, where the VDB tree, holding the same band values, is also about a sixth; at 8 cm VDB is almost twice the dense grid. The read-back difference is zero at every resolution, so the hashed table stores exactly what the dense array does. Blocks multiply by four to five and a half per halving, close to the factor of four expected for a surface; in our reading the excess comes from the coarse grids, where the band is clipped by the edge of the domain and merges between nearby objects. VDB’s leaves are the hashed blocks, and its extra cost is the internal nodes: the one 32332^3 node alone takes about 270 kB, which is why VDB is larger than the hash table at 8 cm and within 1% of it at 1 cm.

Octree. In every row, leaves equal seven times the inner nodes plus one, 87.5% of the nodes, because every inner node of a fully known tree has eight children. OctoMap’s measured 80–85% would follow if some of their children were unknown and never allocated (our reading). The octree did not grow by four per halving: its bytes rose by 1.6, 6.2 and 1.5 times. In our reading, the scene’s flat faces cost little when they lie on the boundary of a large octree cell and cost a row of finest-level leaves when they do not; at 2 cm, the wall’s inner face at 0.10 m is five voxels from the origin, an odd number. In OctoMap’s 64-bit layout the octree is about the size of the 1-bit dense grid at 1 cm (2.08 against 2.10 MB) and larger at coarser sizes, but it holds a 4-byte value per node; against a dense grid of 4-byte values (half the 8-byte column) it is 32 times smaller at 1 cm. The bit stream, 21.7 kB, is about one hundredth of the 1-bit grid.

Volume. At 2 and 1 cm the occupied volume is within 0.05% of the true volume. At 4 cm the result depends on a rule the sources leave open: the voxel centres there sit at 0.02 + 0.04k m, so the wall’s face at 0.10 m and the box’s faces at 0.10 and 0.70 m pass exactly through centres. At the wall’s face the computed distance is exactly zero, so the strict test f<0f < 0 drops those voxels while f≤0f \le 0 keeps them, which turns a loss of about a tenth of the volume into a gain of a fifteenth. At the box’s faces, and at its face at 1.80 m in the 8 cm grid, the computed distance comes out a rounding error above zero in floating point (our check of the arithmetic), so both rules drop those voxels and agree at 8 cm. Faces aligned with the grid are common in built environments, so a volume or occupancy measured on a coarse grid can depend on where the grid’s origin happens to fall, and on rounding. The thin pole shows the same sampling effect without ties: its voxelized volume swings from more than twice its true volume at 8 cm to about half at 4 and 2 cm, depending on which voxel centres fall inside a disk 6 cm across.

Coarse queries. The maximum rule never misses occupied space at any of the three cell sizes, as OctoMap’s authors intend for navigation, but the free space it marks occupied grows quickly with the cell size, to a third of the occupied volume at 8 cm. A grid built directly at the coarse voxel size inflates much less but misses occupied space, about a tenth of it at 4 and 8 cm; the 4 cm direct grid is the 4 cm grid of the first table, so its misses include the voxels lost to ties. A conservative coarse map is therefore a different object from a coarse map: it answers “could anything be here?”, which is the question a collision check asks.

The example began with a single volume-error column using f<0f < 0. After the first run, a check of the voxel centres showed that some lie exactly on primitive faces at 8 and 4 cm, so the column under f≤0f \le 0 was added and the header relabelled; nothing else changed.

Common Misconceptions

  • “An octree is always smaller than a dense grid.” It is smaller than a dense grid of the same value type when the contents are coherent. Per node it stores a value and pointers, so against a 1-bit dense grid, or in a confined space, it can be as large or larger, as OctoMap’s authors note and the example shows.
  • “Sparse grids are lossy.” Hashed blocks, VDB and pruned octrees of identical children store exactly what the dense grid stores. Loss comes from choices made on top: clamping before pruning, maximum-likelihood compression, garbage collection and truncation of distances far from the surface.
  • “A coarser query of a map is the same as a map built at the coarser resolution.” OctoMap’s inner nodes hold the maximum of their children, which over-approximates occupied space; resampling at a coarser voxel size can miss it.
  • “Finer voxels are always more accurate.” Below the sensor’s noise and spacing they add memory without detail, as Nießner et al. found below 2 mm, and for thin or grid-aligned structures the error of a given resolution depends on alignment as much as on size.

Where It Is Used

  • Robot mapping and planning: 2D and 3D occupancy grids for navigation, exploration and collision checking (robot perception, visual navigation, obstacle detection and avoidance), with octrees such as OctoMap where free and unknown space must be kept apart; maps built from LiDAR and depth cameras. OctoMap’s authors also integrated colour data from the TUM RGB-D freiburg1_360 sequence into a 2 cm map (§3.5.1, §5.2).
  • Dense reconstruction: TSDF volumes in TSDF fusion, KinectFusion and dense SLAM, sparse through voxel hashing or blocks for large scenes, meshed by marching cubes; surface reconstruction methods such as Poisson surface reconstruction solve on an octree.
  • Learning: voxel grids as inputs and outputs of 3D convolutional networks and as targets of 3D occupancy prediction and semantic scene completion; sparse convolutional networks and octree networks for higher resolutions; voxel features in hybrid neural implicit representations.
  • Graphics and simulation: narrow-band level sets, fluids and clouds in VDB and OpenVDB, the applications for which Museth designed it.
  • Point cloud processing: voxel downsampling before registration with iterative closest point and other steps that need uniform density.

Related

  • TSDF Fusion

    Fusing many depth maps with known camera poses into one voxel grid of truncated signed distances by a per-voxel weighted average, then extracting the surface at the zero level set.

  • 3D Occupancy Prediction

    Predicting, for every voxel of the space around a vehicle or robot, whether it is free or occupied and by which class, from cameras or LiDAR; covers occupancy grids and semantic scene completion, camera-based occupancy networks, the Occ3D, SurroundOcc and OpenOccupancy labels, visibility masks, and the IoU, mIoU and RayIoU scores.

  • Signed Distance Function

    A scalar field giving the distance to the nearest surface point, negative inside and positive outside, whose zero level set is the surface and whose gradient is the surface normal.

  • Marching Cubes

    An algorithm that turns a scalar field sampled on a 3D grid, such as a signed distance volume, into a triangle mesh of one level set, cube by cube, with a lookup table and linear interpolation along edges.

  • KinectFusion

    A real-time dense RGB-D SLAM system that fuses every depth frame into a TSDF volume on the GPU and tracks the camera by point-to-plane ICP against a surface ray cast from that volume.

  • Robot Perception

    How robots use cameras, depth sensors, LiDAR, inertial and joint sensors and touch to estimate their own state and the world around them fast enough to act on it, covering the perception–action loop, state estimation, latency and robustness, and the main robot settings from mobile robots and arms to drones, legged robots and cars.

References

  1. Hornung, A., Wurm, K. M., Bennewitz, M., Stachniss, C. & Burgard, W. (2013). OctoMap: An Efficient Probabilistic 3D Mapping Framework Based on Octrees. Autonomous Robots, 34(3), 189–206 (read in the authors' preprint).
  2. Nießner, M., Zollhöfer, M., Izadi, S. & Stamminger, M. (2013). Real-time 3D Reconstruction at Scale using Voxel Hashing. ACM Transactions on Graphics, 32(6), Article 169 (read in the authors' copy).
  3. Museth, K. (2013). VDB: High-Resolution Sparse Volumes with Dynamic Topology. ACM Transactions on Graphics, 32(3), Article 27 (read in the author's copy).
  4. OpenVDB. OpenVDB Overview. OpenVDB 14.0.0 documentation.
  5. OctoMap source code (release v1.10.1): octomap/include/octomap/OcTreeDataNode.h, OcTreeBaseImpl.hxx and octomap/src/AbstractOccupancyOcTree.cpp.
  6. Open3D source code (version 0.19.0): cpp/open3d/geometry/PointCloud.cpp (PointCloud::VoxelDownSample) and cpp/pybind/geometry/pointcloud.cpp.
  7. Moravec, H. & Elfes, A. (1985). High Resolution Maps from Wide Angle Sonar. IEEE International Conference on Robotics and Automation (ICRA), 2, 116–121 (not read).
  8. Elfes, A. (1989). Using Occupancy Grids for Mobile Robot Perception and Navigation. Computer, 22(6), 46–57 (not read).
  9. Meagher, D. (1982). Geometric Modeling Using Octree Encoding. Computer Graphics and Image Processing, 19(2), 129–147 (not read).