Skip to contents

Aggregates points into H3 hexagons at zoom-appropriate resolutions and writes a PMTiles archive in which low zooms show coarse hexagons, intermediate zooms show progressively finer hexagons, and zooms at or above base_zoom show individual points. Aggregations (count, sum, mean, etc.) are computed in DuckDB via the H3 community extension; the function then assembles the per-resolution hex layers and the raw-point layer via freestile().

Usage

freestile_h3(
  input,
  output,
  agg = "count",
  hex_layer_prefix = "h3",
  point_layer_name = "points",
  min_zoom = 0L,
  max_zoom = 14L,
  base_zoom = NULL,
  h3_resolutions = NULL,
  source_crs = NULL,
  db_path = NULL,
  tile_format = "mvt",
  fade = FALSE,
  fade_overlap = 1L,
  overwrite = TRUE,
  quiet = FALSE
)

Arguments

input

An sf data frame of POINT geometry, or a character SQL query that returns a geometry column when executed against DuckDB. (DuckDB's spatial functions such as ST_Read() and read_parquet() are available.)

output

Character. Path for the output .pmtiles file.

agg

Aggregation specification:

  • "count" (default): single count = COUNT(*) property.

  • Named character vector of SQL aggregation expressions, e.g. c(n = "COUNT(*)", avg_pop = "AVG(pop)", total = "SUM(pop)").

  • Named list of prop_name = c(fn, column) for callers who don't want to write SQL, e.g. list(n = c("count", "*"), avg_pop = c("mean", "pop")). Supported fn values: "count", "sum", "mean"/"avg", "min", "max", "median".

hex_layer_prefix

Character. Prefix for the per-resolution hex MVT layer names. Default "h3" produces "h3_r01", "h3_r02", etc.

point_layer_name

Character. MVT layer name for raw points (default "points").

min_zoom, max_zoom

Integer. Global zoom range (default 0–14).

base_zoom

Integer or NULL. Zoom level at and above which raw points take over from hex aggregations. Default NULL resolves to max_zoom - 2L, clamped so min_zoom <= base_zoom <= max_zoom. If base_zoom == min_zoom, no hex layers are created.

h3_resolutions

Optional override of the zoom -> H3 resolution mapping. Accepts NULL (use built-in defaults), an unnamed integer vector with length(min_zoom:(base_zoom - 1L)) entries mapped positionally, or a named integer vector with names that parse to integer zoom levels (sparse overrides; defaults fill the rest). All resolutions must be integers in 0:15. The same resolution appearing in non-contiguous zoom runs (e.g. zooms 4–5 and 8) is rejected.

source_crs

Character or NULL. CRS of geometry returned by SQL input, e.g. "EPSG:4326" (default for sf input, since sf is auto-transformed to WGS84) or "EPSG:3857". Ignored for sf input. For SQL input, if NULL, a warning is emitted and the geometry is assumed to be "EPSG:4326".

db_path

Character or NULL. Path to an existing DuckDB database file, or NULL/"" (default) for an in-memory database.

tile_format

Character. "mvt" (default) or "mlt".

fade

Logical. If TRUE, adjacent hex layer zoom windows overlap by fade_overlap zooms so view_h3_tiles() can cross-fade between resolutions. Default FALSE (disjoint windows, clean zoom breaks).

fade_overlap

Integer. Zooms of overlap on each side (only used when fade = TRUE, default 1L).

overwrite

Logical. Whether to overwrite an existing output file (default TRUE).

quiet

Logical. Suppress progress messages (default FALSE).

Value

The output file path (invisibly).

Details

Each distinct H3 resolution becomes its own MVT source-layer (named "<hex_layer_prefix>_r<resolution>", e.g. "h3_r05"); raw points are emitted as a separate source-layer (point_layer_name). With the default fade = FALSE, per-layer zoom windows are disjoint so the rendered map swaps between resolutions cleanly. With fade = TRUE, adjacent windows overlap by fade_overlap zooms so the companion view_h3_tiles() helper can cross-fade between resolutions.

DuckDB and the H3 community extension are required. With sf input the data is written to DuckDB via a temporary Parquet (or dbWriteTable) roundtrip; with character SQL input, the user query is wrapped in a temporary view inside DuckDB.

Point volume at base_zoom

The raw-point layer is encoded without thinning: every input point is written to every tile that contains it from base_zoom up (and from base_zoom - fade_overlap when fade = TRUE). For very large inputs the tiles at the first point zoom can get heavy. Raising base_zoom defers points to higher zooms, where each tile covers less area and so holds fewer points. Per-layer feature dropping for the points layer is on the roadmap.

Antimeridian and polar cells

Hexagons that cross the antimeridian are split at +/-180 degrees so they render correctly on both sides of the dateline instead of as world-spanning slivers. The rare cells that contain a pole receive the same split and render approximately (the polygon edge follows the cell boundary vertices, so the polar cap itself is not filled).

Examples

if (FALSE) { # \dontrun{
library(sf)
pts <- st_as_sf(data.frame(
  x = runif(50000, -100, -80),
  y = runif(50000, 30, 45),
  w = rnorm(50000, 100, 10)
), coords = c("x", "y"), crs = 4326)

freestile_h3(pts, "wind.pmtiles",
  agg = c(n = "COUNT(*)", avg_w = "AVG(w)"),
  min_zoom = 2, max_zoom = 12, base_zoom = 10)

view_h3_tiles("wind.pmtiles", agg_column = "n")

# Cross-fade between resolutions
freestile_h3(pts, "wind_fade.pmtiles",
  agg = "count",
  min_zoom = 2, max_zoom = 12, base_zoom = 10,
  fade = TRUE)
} # }