SpatioX

An in-memory spatio-temporal index written in C++17 with Python bindings (pybind11). It stores (lat, lon, t) points and answers "what happened near here, during this time window?" queries: radius, bounding box and k-nearest-neighbour, each with or without a time filter.

On a Windows laptop, radius+time queries over 50,000 points run about 7x faster than a pure-Python linear scan, and box+time queries about 10x faster (see Benchmarks).

Features

• KD-tree spatial index with per-node spatial and temporal subtree bounds, so radius+time queries prune on distance and on time in the same traversal.
• Balanced rebuild: build() rebuilds the tree with median splits (std::nth_element) over Morton (Z-order) sorted input.
• Streaming inserts after build() go into an LSM-style delta buffer (four levels, 2k / 8k / 32k / 128k records) that compacts upward and finally rebuilds the tree. Every query reads both the tree and the buffer.
• Bloom filters over time buckets (3,600 time units wide, one hour if t is in seconds) let time-filtered queries skip whole delta levels that cannot contain the requested window.
• Query instrumentation: query_radius_time_instrumented returns nodes visited, subtrees pruned by distance and by time, and Bloom skips.
• Differential tests: results are checked against independent brute-force oracles, not just against pruning counters.

Installation

Requirements: Python 3.8+, a C++17 compiler (GCC, Clang, or MSVC with the Windows SDK). pybind11 is pulled in automatically as a build dependency.

git clone https://github.com/coderconnoisseur/SpatioX.git
cd SpatioX
pip install -e .

pip install -e . compiles the extension spatiox._spatio_core through setup.py. A CMakeLists.txt is also included for building the module with CMake directly.

Quick start

The package is imported as spatiox. SpatioIndex is a small Python wrapper that also stores an arbitrary payload per record:

from spatiox import SpatioIndex

index = SpatioIndex()

# insert(lat, lon, t, payload=None) -> record id
index.insert(40.7128, -74.0060, 1000.0, payload={"pickup": "Financial District"})
index.insert(40.7589, -73.9851, 1500.0, payload={"pickup": "Times Square"})
index.insert(40.7614, -73.9776, 2000.0, payload={"pickup": "Midtown East"})

# Records within 2 km of Times Square with 1000 <= t <= 1600
ids = index.query_radius_time(40.7589, -73.9851, radius_km=2.0, t_start=1000.0, t_end=1600.0)
for rid in ids:
    rec = index.get_record(rid)
    print(rid, rec.lat, rec.lon, rec.t, index.get_payload(rid))
# 2 40.7588996887207 -73.98509979248047 1500.0 {'pickup': 'Times Square'}

# Records inside a bounding box during a time window
print(index.query_box_time(40.70, -74.01, 40.80, -73.90, t_start=500.0, t_end=2500.0))
# [1, 2, 3]

SpatioIndexCore is the C++ class itself and exposes everything, including bulk loading, k-NN and statistics:

import random
from spatiox import SpatioIndexCore

core = SpatioIndexCore()
random.seed(0)
points = [(random.uniform(40.70, 40.88), random.uniform(-74.02, -73.91),
           random.uniform(0, 86_400)) for _ in range(10_000)]
core.bulk_insert(points)   # list of (lat, lon, t) tuples -> list of ids
core.build()               # rebuild as a balanced KD-tree

# 5 nearest records to Times Square, nearest first
print(core.query_knn(40.7589, -73.9851, k=5))

# Nearest records with t in the first 12 hours.
# May return fewer than k -- see Known limitations.
print(core.query_knn_time(40.7589, -73.9851, k=5, t_start=0, t_end=12 * 3600))

# Inserts after build() land in the delta buffer and are still visible to queries
new_id = core.insert(40.7589, -73.9851, 9 * 3600)
print(new_id in core.query_radius(40.7589, -73.9851, radius_km=0.1))   # True

ids, stats = core.query_radius_time_instrumented(40.7589, -73.9851, 1.0, 8 * 3600, 10 * 3600)
print(stats)                    # QueryStats(nodes=..., bbox_prunes=..., temporal_prunes=..., ...)
print(core.get_index_stats())   # IndexStats(records=10001, ..., built=True, delta_records=1, ...)

API reference

All coordinates are degrees; latitude/longitude are stored as float32. Distances use the haversine formula. Time ranges are inclusive on both ends. Record ids start at 1.

spatiox.SpatioIndex (Python wrapper)

Method	Returns	Notes
insert(lat, lon, t, payload=None)	int	Payload is kept in a Python dict.
query_radius_time(center_lat, center_lon, radius_km, t_start, t_end)	list[int]
query_box_time(lat_min, lon_min, lat_max, lon_max, t_start, t_end)	list[int]
get_record(record_id)	Record or None
get_payload(record_id)	payload or None
size()	int
clear()	None	Clears records and payloads.

The wrapper does not expose build(), bulk loading or k-NN; use SpatioIndexCore for those.

spatiox.SpatioIndexCore (C++ core)

Method	Returns	Notes
insert(lat, lon, t)	int	Before build(): inserted into the KD-tree. After: appended to the delta buffer.
bulk_insert(records)	list[int]	records is a list of (lat, lon, t) tuples.
build()	None	Rebuilds a balanced KD-tree from all records and empties the delta buffer.
query_radius(center_lat, center_lon, radius_km)	list[int]	No time filter.
query_box(lat_min, lon_min, lat_max, lon_max)	list[int]	No time filter.
query_knn(lat, lon, k)	list[int]	Exact; nearest first, ties broken by id.
query_radius_time(center_lat, center_lon, radius_km, t_start, t_end)	list[int]	Distance and time pruning in one traversal.
query_box_time(lat_min, lon_min, lat_max, lon_max, t_start, t_end)	list[int]	Box query, then time filter.
query_knn_time(lat, lon, k, t_start, t_end)	list[int]	Approximate; see Known limitations.
query_radius_time_instrumented(center_lat, center_lon, radius_km, t_start, t_end)	(list[int], QueryStats)	Same results as query_radius_time.
get_record(id)	Record or None
size()	int
clear()	None
get_index_stats()	IndexStats

Other exported types: Record (lat, lon, t, id, morton), RecordInput, BuildPoint, QueryStats, IndexStats (record/node counts, time range, is_built, delta-buffer and compaction counters) and QueryStrategy.

Design

Why spatial first, then time

Each query walks the KD-tree on (lat, lon) and treats time as a filter. Every tree node also stores the min/max timestamp of its subtree, so query_radius_time drops a subtree when either its bounding box is farther than the radius or its time range misses the query window. In practice the time check is a cheap extra prune in the spatial walk.

The other option is a time-first plan: range-scan the sorted time index, then check distance for each hit. That is only better when the time window is very narrow relative to the data. QueryStats.strategy_used reports which plan would be chosen (time-first when the window covers less than 10% of the data's time span), but only the spatial-first plan is implemented; the field is for diagnostics only.

Complexity

N = records, k = results.

Operation	Cost
insert before build()	O(depth): O(log N) for random input, O(N) for sorted input (the tree is not rebalanced)
insert after build()	amortised O(1) append, plus periodic compaction
build()	O(N log N) (Morton sort + median splits)
radius / box query on a balanced tree	about O(√N + k) for a 2-D KD-tree, plus a linear scan of the delta buffer
query_knn	O(log N) typical on the tree, plus the delta buffer

Tradeoffs

• Delta buffer instead of rebalancing. Inserts after build() are appended to the buffer rather than threaded into the tree, so the tree stays balanced and inserts stay cheap. The cost is that queries scan the buffer linearly, which is why its levels are capped and per-level Bloom filters skip levels that cannot match the time window.
• float32 coordinates keep records small (float32 resolution is about 0.4 m of latitude near 40°N). Points within about a metre of a query boundary can be classified differently from a float64 calculation.
• In memory only. There is no persistence and no deletion API.

Benchmarks

Produced by running python benchmark.py three times on a Windows 10 laptop (AMD Ryzen 5 3550H, Python 3.11.5). The extension was built from this repository's sources with MinGW-w64 g++ 13.2 at -O2. Raw output is in benchmarks/benchmark2.txt.

The baseline is the pure-Python linear scan in benchmark.py (NaiveSpatioIndex), which tests every record with the same haversine and box predicates. Both engines get the same random points (uniform over Manhattan, t uniform over one day) and the same 50 random queries (radius 0.5-5 km or boxes 0.01-0.05° on a side, time windows 1-6 h). The script checks that both engines return the same number of results. SpatioX is called through the SpatioIndex wrapper with one-at-a-time inserts and no build().

Mean latency per query, median of the 3 runs:

Records	Query	SpatioX	Python linear scan	Speedup	Avg. results
1,000	radius + time	109 µs	437 µs	4.6x	18
1,000	box + time	9 µs	121 µs	13.0x	7
10,000	radius + time	745 µs	5.6 ms	7.6x	174
10,000	box + time	85 µs	1.7 ms	20.1x	69
50,000	radius + time	3.4 ms	27.2 ms	6.8x	868
50,000	box + time	704 µs	7.3 ms	10.3x	342

How to read these:

• Speedups are moderate because these queries return many records (hundreds at 50k). Output size and building the Python result list take a large share of the time. The baseline is pure Python; a NumPy or compiled scan would close much of the gap.
• Timings are noisy on a laptop. The 10k box speedup ranged from 11.7x to 25.6x across the three runs; the others varied by less than 2x.
• Inserting 50,000 records one at a time from Python took 0.19 s (run 3).

Testing

pip install pytest numpy
python -m pytest tests test_basic.py

• tests/test_correctness.py: differential tests against brute-force oracles (radius, box, k-NN, post-build() inserts, Bloom filter edge cases, deep trees).
• tests/test_pruning_regression.py, test_basic.py: check that subtree pruning actually happens on clustered data.
• tests/test_bloom_filter.py: Bloom pruning never drops a matching record.

Current local result: 33 passed, 4 expected failures (xfail) that track the known bugs listed below. CI runs the same command on Ubuntu for every push and pull request.

Known limitations

These are known bugs or approximations, each covered by a test:

• Negative timestamps with the delta buffer. bucket_id() casts a negative floor() to uint64_t. For time ranges that cross zero, the Bloom check can skip a delta level that holds matches. This affects records inserted after build() with t < 0. Tracked by three xfail tests in tests/test_correctness.py.
• Deep recursion on sorted input. Before build(), the KD-tree is not rebalanced and insertion is recursive. About 30,000 inserts in coordinate order (for example one GPS track in time order) can overflow the stack on Windows. Workaround: bulk_insert + build(). Tracked by an xfail test.
• query_knn_time is approximate. It takes the 3k spatially nearest records from the tree, filters them by time, and merges in the delta buffer. When the time window is narrow, it can return fewer than k records even if more exist (the quick-start example returns 4 for k=5). The records it returns are correct and in order. query_knn without a time filter is exact.

Project structure

SpatioX/
├── .github/workflows/ci.yml   # build + tests on Ubuntu
├── CMakeLists.txt             # optional CMake build
├── setup.py                   # pybind11 extension build used by pip
├── pyproject.toml
├── benchmark.py               # benchmark vs pure-Python linear scan
├── benchmarks/                # recorded benchmark output
├── include/
│   ├── bloom_filter.hpp       # Bloom filter over time buckets
│   ├── morton.hpp             # Z-order encoding
│   ├── record.hpp
│   ├── record_store.hpp
│   ├── spatial_index.hpp      # KD-tree with spatio-temporal node bounds
│   ├── spatio_index_core.hpp  # core engine, LSM delta levels
│   ├── temporal_index.hpp
│   └── utils.hpp              # haversine
├── src/
│   ├── bindings.cpp           # pybind11 module spatiox._spatio_core
│   ├── record.cpp
│   ├── record_store.cpp
│   ├── spatial_index.cpp
│   ├── spatio_index_core.cpp
│   └── temporal_index.cpp
├── python/spatiox/__init__.py # SpatioIndex wrapper
├── tests/
│   ├── test_bloom_filter.py
│   ├── test_correctness.py
│   └── test_pruning_regression.py
└── test_basic.py

License

MIT. See LICENSE.