Implementing polygon containment checks in C++
This page solves one concrete problem: deciding, deterministically and allocation-free, whether a streaming GPS fix falls inside a predefined geofence polygon on a constrained edge gateway — a header-only C++ predicate that runs on everything from an ARM Cortex-M7 microcontroller to an x86_64 industrial gateway. Within the Local Spatial Processing Patterns practice, and specifically as a concrete build of on-device geometry filtering, the containment check is the precise second stage that runs only on points a cheap bounding-box reject has already let through. Pushing this decision to the network edge eliminates cloud round-trip latency, spares metered backhaul, and keeps the geofence verdict available even when the uplink has been gone for hours.
The deployment context is unforgiving: field nodes cannot afford dynamic allocation, garbage-collection pauses, or unpredictable cache behavior while a modem and a sensor bus contend for the same core. The implementation below is therefore built to be allocation-free, branch-predictable, and tolerant of raw GNSS drift — the three properties that separate a containment check that survives a fanless gateway in the field from one that only ever ran on a workstation.
Algorithm selection rationale
The ray-casting algorithm (even-odd rule) remains the right choice for on-device polygon containment because its cost profile maps cleanly onto the constraint envelope: O(n) in the vertex count, with minimal state and a strictly sequential memory access pattern. Unlike winding-number approaches — which track directional crossings and carry more topological bookkeeping — ray-casting only counts edge intersections along a single horizontal projection. On embedded Linux or a bare-metal RTOS that translates directly into deterministic cache behavior and zero heap fragmentation, which is exactly what a per-packet hot path needs.
Three deployment constraints are baked into the design rather than bolted on afterward:
- Zero dynamic allocation. All vertex data is passed via raw pointers or fixed-size arrays. No STL containers, no heap allocation, no exceptions — the entire working set is the polygon, the cached bounding box, and a handful of stack locals.
- Early-exit bounding box. A cheap AABB (axis-aligned bounding box) check rejects the overwhelming majority of out-of-bounds points — typically 90%+ during real telemetry — before the intersection loop ever runs, preserving CPU cycles for high-frequency ingestion.
- Boundary tolerance. GNSS jitter rarely lands exactly on a vertex or edge. A configurable epsilon prevents false negatives caused by coordinate quantization and atmospheric multipath, the same drift-tolerance idea applied elsewhere to movement and threshold gating for sensor events.
Ray-casting containment with an AABB early-exit and vertex tolerance.
A note on coordinate systems before any code runs: ray-casting assumes a planar Cartesian space. For geofences under roughly 5 km radius, treating WGS84 lat/lon as planar introduces negligible error (under 0.1 m at mid-latitudes). Larger boundaries must be pre-projected to a local UTM zone or EPSG:3857 at ingestion — the projection trade-offs for low-memory nodes are covered under coordinate reference systems at the edge, and the Cortex-M specifics under CRS transformations on ARM Cortex-M devices.
The complete, self-contained implementation
The following header-only implementation is designed for direct inclusion in gateway firmware. Compile with -O2 -fno-exceptions -fno-rtti for deterministic latency on constrained targets. Avoid -ffast-math: it can reorder floating-point operations in ways that break the strict inequality guard that keeps the division safe.
// geofence_containment.hpp
#pragma once
#include <cmath>
#include <cstddef>
struct Point2D {
double lat;
double lon;
};
struct Polygon {
const Point2D* vertices;
std::size_t count;
};
struct AABB {
double min_lat, max_lat;
double min_lon, max_lon;
};
// Precompute AABB for early-exit filtering.
// Call once during geofence provisioning, not per-telemetry packet.
inline AABB compute_aabb(const Polygon& poly) {
AABB box{poly.vertices[0].lat, poly.vertices[0].lat,
poly.vertices[0].lon, poly.vertices[0].lon};
for (std::size_t i = 1; i < poly.count; ++i) {
if (poly.vertices[i].lat < box.min_lat) box.min_lat = poly.vertices[i].lat;
if (poly.vertices[i].lat > box.max_lat) box.max_lat = poly.vertices[i].lat;
if (poly.vertices[i].lon < box.min_lon) box.min_lon = poly.vertices[i].lon;
if (poly.vertices[i].lon > box.max_lon) box.max_lon = poly.vertices[i].lon;
}
return box;
}
// Deterministic point-in-polygon check with epsilon tolerance and AABB early-exit.
// Division by zero is mathematically impossible due to the strict inequality guard.
[[nodiscard]] inline bool point_in_polygon(const Point2D& pt,
const Polygon& poly,
const AABB& box,
double epsilon = 1e-9) {
// Early-exit AABB check (expanded by epsilon for boundary tolerance)
if (pt.lat < box.min_lat - epsilon || pt.lat > box.max_lat + epsilon ||
pt.lon < box.min_lon - epsilon || pt.lon > box.max_lon + epsilon) {
return false;
}
bool inside = false;
for (std::size_t i = 0, j = poly.count - 1; i < poly.count; j = i++) {
double yi = poly.vertices[i].lat, xi = poly.vertices[i].lon;
double yj = poly.vertices[j].lat, xj = poly.vertices[j].lon;
// Vertex proximity check (handles GPS quantization on corners)
if (std::abs(pt.lat - yi) < epsilon && std::abs(pt.lon - xi) < epsilon) {
return true;
}
// Ray-casting intersection logic
// (yi > pt.lat) != (yj > pt.lat) ensures yi != yj, preventing division by zero
if (((yi > pt.lat) != (yj > pt.lat)) &&
(pt.lon < (xj - xi) * (pt.lat - yi) / (yj - yi) + xi)) {
inside = !inside;
}
}
return inside;
}
The compute_aabb call belongs in the provisioning path, not the hot path: derive the box once when a geofence is loaded and cache it in read-only memory. Every subsequent point_in_polygon call then pays only the four float comparisons of the reject, and the loop runs solely for the survivors.
Constraint validation table
Each design choice in the code maps to a specific hardware limit and the mitigation it carries. Size the pipeline against the real budget of the target — the broader picture is in device constraints and resource limits; the table below is the subset that bears directly on the containment predicate.
| Constraint | Expected impact | Mitigation built into the code |
|---|---|---|
| RAM | A heap-backed geometry library can spike tens of MB per polygon on a 128–256 MB node | Raw-pointer Polygon, fixed AABB, only stack locals — zero allocation per call |
| CPU / thermal | Per-point topology on every fix saturates a throttling core during 50–100 Hz bursts | AABB early-exit discards 90%+ before the O(n) loop; loop body is branch-light |
| Latency | Heap churn and GC-style pauses push tail latency past the polling interval | Allocation-free, exception-free path keeps per-call latency sub-millisecond and bounded |
| Power | Software-emulated double math burns cycles (and battery) on single-precision FPUs |
Optional float quantization; AABB reject keeps the expensive division off the common path |
Gotchas and edge cases
The traps here are floating-point and coordinate-system traps, not algorithmic ones. The crossing test (yi > pt.lat) != (yj > pt.lat) is doing double duty: it both selects edges that straddle the ray and guarantees yi != yj, so the /(yj - yi) division can never see a zero denominator. Re-deriving that guard incorrectly — for example switching > to >= — reintroduces the division-by-zero on horizontal edges that the strict inequality is there to prevent.
- GNSS drift tolerance. A receiver parked next to a boundary oscillates by a few meters on multipath alone. Pick epsilon from your precision contract, not arbitrarily:
1.5e-5degrees is roughly 1.5 m at the equator. Quantizing to integer microdegrees first, per spatial data precision standards, turns the boundary comparison into exact integer math and sidesteps float-equality entirely. - Coordinate system assumptions. The planar approximation only holds for small geofences. Mixing a WGS84 polygon with already-projected fixes silently produces wrong verdicts; pin one reference system at ingestion and reject anything that does not match it.
- Compiler flag pitfalls.
-ffast-mathis the classic foot-gun — it lets the optimizer assume no special values and reassociate the multiply/divide, which can move a point across the boundary by an ULP at exactly the edge case epsilon exists to cover. Keep it off and prefer explicit quantization for reproducibility. - RTOS quirks. On a single-precision FPU (Cortex-M4,
-mfpu=fpv4-sp-d16) thedoublefields fall back to slow software emulation. Either target a double-precision FPU such as the Cortex-M7’s-mfpu=fpv5-d16, or quantize coordinates tofloat. Align vertex arrays to a cache line (alignas(64)) to avoid false sharing and pipeline stalls during the loop. - Concave boundaries. Highly concave geofences make
inside = !insidetoggle often, raising branch-misprediction rates. Watch them with hardware performance counters and, during provisioning, pre-flatten complex boundaries into convex pieces or a triangulated mesh to lower instruction-cache pressure.
Integrating with the gateway pipeline
In the telemetry loop, the predicate sits behind the parser and ahead of the alerting path. Feed raw NMEA or binary fixes directly into Point2D structs — no intermediate string parsing — and reuse the cached AABB:
// High-frequency telemetry loop integration
void process_geofence_stream(const Point2D* telemetry_buffer,
std::size_t buffer_len,
const Polygon& geofence,
const AABB& geofence_aabb) {
for (std::size_t i = 0; i < buffer_len; ++i) {
if (point_in_polygon(telemetry_buffer[i], geofence, geofence_aabb, 1.5e-5)) {
// ~1.5e-5 epsilon ≈ 1.5 m at equator
trigger_edge_alert(telemetry_buffer[i]);
}
}
}
For Python-fronted gateways, do not serialize coordinates as JSON across the boundary. Pack lat/lon pairs into fixed-width binary frames (struct.pack('<dd', lat, lon)) and hand the buffer to the compiled predicate through ctypes or a memory view, which keeps the ingress path allocation-free and the GIL out of the hot loop. The same fixed-width framing is what makes the downstream delta sync for GPS coordinate streams cheap once a fix has cleared the geofence and earned its place on the uplink.
Validate against field conditions before trusting the epsilon: feed coordinates oscillating ±2 m around polygon edges to confirm the tolerance behaves, and cross-reference verdicts against the OGC Simple Features Specification test suites when migrating to multi-polygon or hole-containing geometries.
Related
- On-device geometry filtering — the parent pattern this predicate plugs into, including the two-stage reject model.
- Spatial data precision standards — fixing a quantization grid so boundary comparisons become exact integer math.
- Configuring spatial thresholds for sensor event triggers — the complementary movement/hysteresis gate for drift-prone fixes.
- Reducing RAM usage for GeoJSON parsing on Raspberry Pi — loading geofence geometry without an OOM kill on a constrained node.