| Title: | DuckDB-Backed Spatial Data Structures |
|---|---|
| Description: | Provides DuckDB-backed implementations of spatial data structures for efficient out-of-memory operations on large spatial datasets. The package includes 68 sf-compatible spatial methods (st_intersects, st_contains, st_within, etc.) and GeoParquet 1.0 I/O support via writeGeoParquet(). Ideal for spatial analysis, genomic interval analysis, and integration with BiocDuckDB for MultiAssaySpatialExperiment workflows requiring memory-efficient operations on large datasets. |
| Authors: | Patrick Aboyoun [aut, cre], Genentech, Inc. [cph] |
| Maintainer: | Patrick Aboyoun <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.99.00 |
| Built: | 2026-07-23 11:09:57 UTC |
| Source: | https://github.com/BiocStaging/DuckDBSpatial |
Turns an RFC-5 transform dict (a named list with a type field, as read
from JSON) into a coordinate-transforms object, recursing
through sequence / byDimension / bijection. A
CoordinateTransform passes through unchanged.
asCoordinateTransform(x)asCoordinateTransform(x)
x |
A |
A CoordinateTransform.
Build the RFC-5 affine-family coordinate transforms. Each returns a
CoordinateTransform (a classed list mirroring the RFC-5 / SpatialData
on-disk dict), consumed by ctAffineMatrix /
ctApply and the coordinate-transform graph
(ctGraph). Transforms are positional: they act on the
coordinate vector in the coordinate system's axis order.
ctIdentity() ctScale(scale) ctTranslation(translation) ctRotation(rotation) ctAffine(affine) ctSequence(transformations) ctMapAxis(map) ctByDimension(parts) ctBijection(forward, inverse)ctIdentity() ctScale(scale) ctTranslation(translation) ctRotation(rotation) ctAffine(affine) ctSequence(transformations) ctMapAxis(map) ctByDimension(parts) ctBijection(forward, inverse)
scale, translation
|
Numeric per-axis vectors. |
rotation, affine
|
A rotation (square) or affine ( |
transformations |
A list of |
map |
Integer 0-based axis permutation: output axis |
parts |
For |
forward, inverse
|
|
A CoordinateTransform.
ctScale(c(10, 20)) ctSequence(list(ctScale(c(10, 20)), ctTranslation(c(1, 1))))ctScale(c(10, 20)) ctSequence(list(ctScale(c(10, 20)), ctTranslation(c(1, 1))))
A coordinate system is a name plus an ordered set of axes (RFC-5
{name, type, unit}). Coordinate transforms map points between
coordinate systems; see ctGraph.
coordinateSystem(name, axes)coordinateSystem(name, axes)
name |
Coordinate-system name. |
axes |
A character vector of axis names, or a list of axis specs each
with |
A CoordinateSystem.
coordinateSystem("global", c("y", "x"))coordinateSystem("global", c("y", "x"))
Returns the homogeneous matrix of a
coordinate-transforms transform, so that
. Every RFC-5 type reduces to such a matrix;
sequence composes by matrix multiplication and byDimension
scatters its per-axis sub-transforms into the block. bijection uses
its forward transform.
ctAffineMatrix(transform, input_axes, output_axes)ctAffineMatrix(transform, input_axes, output_axes)
transform |
A |
input_axes, output_axes
|
Character vectors of the source / target axis names (only their lengths, the dimensionalities, are used). |
A numeric homogeneous matrix.
ctAffineMatrix(ctScale(c(10, 20)), c("y", "x"), c("y", "x"))ctAffineMatrix(ctScale(c(10, 20)), c("y", "x"), c("y", "x"))
Transforms an matrix of coordinates to via the transform's homogeneous affine matrix
(ctAffineMatrix). This is the RFC-5 conformance oracle.
ctApply(transform, points, input_axes, output_axes)ctApply(transform, points, input_axes, output_axes)
transform |
A |
points |
A numeric matrix (rows = points) or a vector (one point). |
input_axes, output_axes
|
Character vectors of source / target axis names. |
A numeric matrix of transformed points ().
ctApply(ctTranslation(c(10, 20)), rbind(c(1, 2)), c("y", "x"), c("y", "x"))ctApply(ctTranslation(c(10, 20)), rbind(c(1, 2)), c("y", "x"), c("y", "x"))
Resolves the transform with ctPath and applies it to a matrix
of points with ctApply, using the axis orders recorded for the
two coordinate systems in the graph.
ctBetween(graph, from, to, points)ctBetween(graph, from, to, points)
graph |
A |
from, to
|
Coordinate-system names. |
points |
A numeric matrix (rows = points) or a vector (one point). |
A numeric matrix of transformed points.
g <- ctGraph(list(list(input = "a", output = "b", transform = ctScale(c(2, 3)))), systems = list(a = c("y", "x"), b = c("y", "x"))) ctBetween(g, "a", "b", rbind(c(1, 1)))g <- ctGraph(list(list(input = "a", output = "b", transform = ctScale(c(2, 3)))), systems = list(a = c("y", "x"), b = c("y", "x"))) ctBetween(g, "a", "b", rbind(c(1, 1)))
Assembles the directed graph of coordinate systems and the transforms between
them. Each invertible edge gets an automatically added reverse edge (via
ctInvert), so a target reachable by inversion is found by
ctPath. Mirrors SpatialData's transformation graph.
ctGraph(edges, systems = NULL)ctGraph(edges, systems = NULL)
edges |
A list of edges, each |
systems |
Optional named list of coordinate systems
( |
A CTgraph.
g <- ctGraph(list(list(input = "a", output = "b", transform = ctScale(c(2, 2)))), systems = list(a = c("y", "x"), b = c("y", "x"))) ctPath(g, "a", "b")g <- ctGraph(list(list(input = "a", output = "b", transform = ctScale(c(2, 2)))), systems = list(a = c("y", "x"), b = c("y", "x"))) ctPath(g, "a", "b")
Returns the inverse coordinate-transforms transform: identity
is its own inverse; scale inverts elementwise; translation negates; rotation
/ affine invert by solve (errors for a non-invertible or
dimensionality-changing affine); mapAxis inverts the permutation; sequence
reverses a list of inverses; bijection swaps its declared forward / inverse.
Used by ctGraph to add reverse edges.
ctInvert(transform)ctInvert(transform)
transform |
A |
The inverse CoordinateTransform.
ctInvert(ctScale(c(10, 20)))ctInvert(ctScale(c(10, 20)))
Finds the shortest path in a ctGraph from from to
to and returns the composed coordinate-transforms
transform (a ctSequence of the edge transforms, or the single edge
transform for a one-hop path, or ctIdentity() when from == to).
Errors if either coordinate system is unknown, if no path exists, or if two
distinct shortest paths tie (an ambiguous resolution).
ctPath(graph, from, to)ctPath(graph, from, to)
graph |
A |
from, to
|
Coordinate-system names. |
A CoordinateTransform.
g <- ctGraph(list(list(input = "a", output = "b", transform = ctScale(c(2, 2)))), systems = list(a = c("y", "x"), b = c("y", "x"))) ctPath(g, "b", "a") # traverses the auto-added inverse edgeg <- ctGraph(list(list(input = "a", output = "b", transform = ctScale(c(2, 2)))), systems = list(a = c("y", "x"), b = c("y", "x"))) ctPath(g, "b", "a") # traverses the auto-added inverse edge
Spatial operations on DuckDBColumn objects.
Method return types are documented in the sections above.
In the code snippets below, x is a DuckDBColumn object.
st_as_sfc(x, ..., crs = NA_integer_, GeoJSON = FALSE, WKB = FALSE):Parses WKT, GeoJSON, or WKB into GEOMETRY.
st_as_binary(x, hex = FALSE):Serialises to WKB or hex-encoded WKB.
st_as_text(x, geojson = FALSE):Serialises to WKT or GeoJSON.
Scalar properties of each geometry.
st_coordinates(x):Returns a DuckDBDataFrame with X, Y (and Z, M when present) coordinate columns.
st_dimension(x):Topological dimension (0, 1, or 2).
st_end_point(x):End point of a linestring.
st_geometry_type(x):Geometry type name as character.
st_is_closed(x):Logical: is the linestring closed (first point = last point)?
st_is_empty(x):Logical: is the geometry empty?
st_is_ring(x):Logical: is the linestring both closed and simple?
st_is_simple(x):Logical: is the geometry simple (e.g. not self-intersecting)?
st_is_valid(x):Logical: is the geometry valid?
st_num_geometries(x):Number of geometries in a geometry collection.
st_num_interior_rings(x):Number of interior rings in a polygon.
st_num_points(x):Number of points within a geometry.
st_start_point(x):Start point of a linestring.
st_area(x):Area of each geometry.
st_distance(x, y):Euclidean distance from each geometry to y.
st_length(x):Length of linestring or multilinestring; zero for points and polygons.
st_perimeter(x):Perimeter of polygon or multipolygon; zero for points and lines.
Each returns a DuckDBColumn with a transformed geometry per row.
st_boundary(x):Geometry boundary.
st_buffer(x, dist):Buffer by dist units.
st_build_area(x):Build area of a geometry.
st_centroid(x):Centroid point.
st_collection_extract(x, type):Extract geometries of given type from collections.
st_concave_hull(x, ratio, allow_holes = FALSE):Concave hull.
st_convex_hull(x):Convex hull polygon.
st_envelope(x):Axis-aligned bounding-box polygon.
st_exterior_ring(x):Exterior ring of a polygon.
st_flip_coordinates(x):Flip the coordinates of a geometry.
st_inscribed_circle(x, dTolerance, ..., nQuadSegs = 30):Maximum inscribed circle of polygon.
st_line_interpolate(line, dist, normalized = FALSE):Point at distance dist along line.
st_line_merge(x, directed = FALSE):Merge connected line segments.
st_line_project(line, point, normalized = FALSE):Distance or fraction along line to nearest point.
st_line_substring(line, start, end):Substring of line between start and end fractions.
st_make_valid(x):Repair invalid geometry.
st_minimum_rotated_rectangle(x):Minimum-area rotated bounding rectangle.
st_node(x):Add nodes at intersection points for line geometry.
st_normalize(x):Normalise vertex order.
st_point_on_surface(x):A point guaranteed to lie on the surface.
st_reverse(x):Reverse vertex order.
st_reduce_precision(x, precision):Reduce precision of a geometry.
st_remove_repeated_points(x):Remove repeated points from a geometry.
st_simplify(x, preserveTopology = FALSE, dTolerance = 0.0):Simplify geometry.
st_voronoi(x):Voronoi diagram.
Each takes a second geometry argument y and returns a BOOLEAN
DuckDBColumn. See DuckDBTable-spatial for accepted
y types.
st_contains(x, y):Does x contain y?
st_contains_properly(x, y):Does x properly contain y (boundary excluded)?
st_covered_by(x, y):Is x covered by y?
st_covers(x, y):Does x cover y?
st_crosses(x, y):Does x cross y?
st_disjoint(x, y):Are x and y disjoint?
st_equals(x, y):Are x and y geometrically equal?
st_intersects(x, y):Do x and y intersect?
st_is_within_distance(x, y, dist):Is x within dist of y?
st_overlaps(x, y):Does x overlap y?
st_touches(x, y):Does x touch y?
st_within(x, y):Is x within y?
st_within_properly(x, y):Is x properly within y (boundary excluded)?
st_difference(x, y):Portion of x that does not intersect y.
st_intersection(x, y):Portion of x that intersects y.
st_nearest_points(x, y):Shortest line between x and y.
st_union(x, y):Binary: returns a DuckDBColumn. Aggregate (no y):
Apply ST_Union_Agg and materialise to sfc.
These methods materialise and return concrete R objects.
st_bbox(obj):Returns a bbox named numeric vector
c(xmin, ymin, xmax, ymax) via MIN/MAX
of ST_XMin/ST_YMin/ST_XMax/ST_YMax.
st_collect(x):Returns an sfc GEOMETRYCOLLECTION via ST_Collect.
st_union(x) (no y):Returns an sfc via ST_Union_Agg.
Patrick Aboyoun
DuckDBColumn-class for the main class
Vector for the base class
spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") df <- DuckDBDataFrame(spatial_path) df <- df[which(!is.na(df$type)), ] geom <- df[["geometry"]] st_geometry_type(geom)[1:3] st_area(geom)[1:3] st_intersects(geom, sf::st_sfc(sf::st_point(c(30, 10))))[1:3]spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") df <- DuckDBDataFrame(spatial_path) df <- df[which(!is.na(df$type)), ] geom <- df[["geometry"]] st_geometry_type(geom)[1:3] st_area(geom)[1:3] st_intersects(geom, sf::st_sfc(sf::st_point(c(30, 10))))[1:3]
Spatial operations on DuckDBTable objects.
st_make_envelope_sql(xmin, ymin, xmax, ymax) st_geomfromtext_sql(wkt) st_point_sql(x_col, y_col) st_intersects_table(x, y, sparse = TRUE, ...) ## S3 method for class 'DuckDBTable' st_filter(x, y, ..., .predicate = st_intersects) ## S3 method for class 'DuckDBTable' st_join(x, y, join = st_intersects, ...) st_union_agg(x, by = NULL, ...) st_collect_agg(x, by = NULL, ...) st_envelope_agg(x, by = NULL, ...)st_make_envelope_sql(xmin, ymin, xmax, ymax) st_geomfromtext_sql(wkt) st_point_sql(x_col, y_col) st_intersects_table(x, y, sparse = TRUE, ...) ## S3 method for class 'DuckDBTable' st_filter(x, y, ..., .predicate = st_intersects) ## S3 method for class 'DuckDBTable' st_join(x, y, join = st_intersects, ...) st_union_agg(x, by = NULL, ...) st_collect_agg(x, by = NULL, ...) st_envelope_agg(x, by = NULL, ...)
xmin, ymin, xmax, ymax
|
Bounding box limits. |
wkt |
A WKT character string. |
x_col, y_col
|
Column names for x and y coordinates. |
x |
A |
y |
A |
sparse |
Logical; return a lazy table when |
... |
Ignored. |
.predicate |
Spatial predicate function (default |
join |
Spatial predicate function (default |
by |
Optional grouping column name(s). |
Method return types are documented in the sections above.
In the code snippets below, x is a DuckDBTable object.
st_as_sfc(x, ..., crs = NA_integer_, GeoJSON = FALSE, WKB = FALSE):Parses WKT (default), GeoJSON, or WKB into GEOMETRY via
ST_GeomFromText, ST_GeomFromGeoJSON, or
ST_GeomFromWKB.
st_as_binary(x, hex = FALSE):Converts to WKB (ST_AsWKB) or hex-encoded WKB
(ST_AsHEXWKB).
st_as_text(x, geojson = FALSE):Converts to WKT (ST_AsText) or GeoJSON
(ST_AsGeoJSON).
Scalar properties of each geometry.
st_dimension(x):Topological dimension: 0 (point), 1 (line), or 2 (polygon).
st_end_point(x):End point of a linestring.
st_geometry_type(x):Geometry type name as character (e.g. "POINT",
"POLYGON").
st_is_closed(x):Logical: is the linestring closed (first point = last point)?
st_is_empty(x):Logical: is the geometry empty?
st_is_ring(x):Logical: is the linestring both closed and simple?
st_is_simple(x):Logical: is the geometry simple (e.g. not self-intersecting)?
st_is_valid(x):Logical: is the geometry valid?
st_num_geometries(x):Number of geometries in a geometry collection.
st_num_interior_rings(x):Number of interior rings in a polygon.
st_num_points(x):Number of points within a geometry.
st_start_point(x):Start point of a linestring.
st_area(x):Area of each geometry.
st_distance(x, y):Euclidean distance from each geometry in x to y.
Argument y accepts WKT character, sfg, sfc
(length 1), or call.
st_length(x):Length of linestring or multilinestring; zero for points and polygons.
st_perimeter(x):Perimeter of polygon or multipolygon; zero for points and lines.
Each returns a DuckDBTable with a transformed geometry per row.
st_boundary(x):Geometry boundary.
st_buffer(x, dist):Buffer by dist units.
st_build_area(x):Build area of a geometry.
st_centroid(x):Centroid point.
st_collection_extract(x, type = c("POLYGON", "POINT", "LINESTRING")):Extract geometries of given type from collections.
st_concave_hull(x, ratio, allow_holes = FALSE):Concave hull; ratio in [0,1] (1 = convex).
st_convex_hull(x):Convex hull polygon.
st_envelope(x):Axis-aligned bounding-box polygon.
st_exterior_ring(x):Exterior ring of a polygon.
st_flip_coordinates(x):Flip the coordinates of a geometry.
st_inscribed_circle(x, dTolerance, ..., nQuadSegs = 30):Maximum inscribed circle of polygon; dTolerance required.
st_line_interpolate(line, dist, normalized = FALSE):Point at distance dist along line; when
normalized = FALSE, dist is in line units; when
TRUE, dist is fraction in [0,1].
st_line_merge(x, directed = FALSE):Merge connected line segments.
st_line_project(line, point, normalized = FALSE):Distance or fraction along line to nearest point to
point; normalized controls units vs fraction.
st_line_substring(line, start, end):Substring of line between start and end fractions.
st_make_valid(x):Repair invalid geometry.
st_minimum_rotated_rectangle(x):Minimum-area rotated bounding rectangle.
st_node(x):Add nodes at intersection points (noding) for line geometry.
st_normalize(x):Normalise vertex order.
st_point_on_surface(x):A point guaranteed to lie on the surface.
st_reduce_precision(x, precision):Reduce precision of a geometry.
st_remove_repeated_points(x):Remove repeated points from a geometry.
st_reverse(x):Reverse vertex order.
st_simplify(x, preserveTopology = FALSE, dTolerance = 0.0):Simplify geometry. Uses ST_SimplifyPreserveTopology when
preserveTopology = TRUE.
st_voronoi(x):Voronoi diagram of points (DuckDB: no envelope/tolerance).
Each takes a second geometry argument y and returns a BOOLEAN
DuckDBTable. Argument y accepts a WKT character string, an
sfc/sfg object from sf, or a call (raw SQL
expression). Results are usable for row filtering via
extractROWS or subset.
st_contains(x, y):Does x contain y?
st_contains_properly(x, y):Does x properly contain y (boundary excluded)?
st_covered_by(x, y):Is x covered by y?
st_covers(x, y):Does x cover y?
st_crosses(x, y):Does x cross y?
st_disjoint(x, y):Are x and y disjoint?
st_equals(x, y):Are x and y geometrically equal?
st_intersects(x, y):Do x and y intersect?
st_is_within_distance(x, y, dist):Is x within dist of y?
st_overlaps(x, y):Does x overlap y?
st_touches(x, y):Does x touch y?
st_within(x, y):Is x within y?
st_within_properly(x, y):Is x properly within y (boundary excluded)?
Each takes a second geometry argument y (same accepted types as
predicates) and returns a GEOMETRY DuckDBTable.
st_difference(x, y):Portion of x that does not intersect y.
st_intersection(x, y):Portion of x that intersects y.
st_nearest_points(x, y):Shortest line between x and y (DuckDB ST_ShortestLine).
st_union(x, y):Union of x and y. Binary form only at the
DuckDBTable level; aggregate union is available on DuckDBColumn.
st_collect(x):Combine all geometries into a single collection.
st_union_agg(x, by = NULL):Aggregate union over a geometry column via ST_Union_Agg.
st_collect_agg(x, by = NULL):Aggregate collect over a geometry column via ST_Collect.
st_envelope_agg(x, by = NULL):Aggregate envelope over a geometry column via ST_Envelope_Agg.
st_filter(x, y, .predicate = st_intersects):Return rows of x that spatially match y.
st_join(x, y, join = st_intersects):Spatial inner join of two lazy tables on geometry.
st_intersects_table(x, y, sparse = TRUE):Return a sparse index table of intersecting row pairs.
Low-level helpers that return dplyr::sql expressions for use in
lazy queries.
st_make_envelope_sql(xmin, ymin, xmax, ymax):Bounding-box envelope expression via ST_MakeEnvelope.
st_point_sql(x_col, y_col):Point expression from column names via ST_Point.
st_geomfromtext_sql(wkt):Geometry expression from WKT via ST_GeomFromText.
Patrick Aboyoun
DuckDBTable-class for the main class
RectangularData for the base class
spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") df <- DuckDBDataFrame(spatial_path) df <- df[which(!is.na(df$type)), ] query_pt <- sf::st_sfc(sf::st_point(c(30, 10))) filtered <- st_filter(df, query_pt) nrow(filtered) st_make_envelope_sql(0, 0, 100, 100)spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") df <- DuckDBDataFrame(spatial_path) df <- df[which(!is.na(df$type)), ] query_pt <- sf::st_sfc(sf::st_point(c(30, 10))) filtered <- st_filter(df, query_pt) nrow(filtered) st_make_envelope_sql(0, 0, 100, 100)
Sets enable_geoparquet_conversion = true so geometry columns in
Parquet files are read as native GEOMETRY types.
enableGeoParquetConversion(conn = NULL)enableGeoParquetConversion(conn = NULL)
conn |
A DuckDB DBI connection. When |
The connection, invisibly.
enableGeoParquetConversion()enableGeoParquetConversion()
Returns the rows of point layer x inside [xmin, xmax] x [ymin, ymax]
using a plain coordinate range predicate
(x_col >= xmin & x_col <= xmax & y_col >= ymin & y_col <= ymax). Unlike
layerSubsetByBbox (which builds a per-row ST_Point and tests
ST_Intersects against an envelope, and returns row indices), the range
predicate pushes into DuckDB's Parquet reader, so on a Morton-sorted (coord-indexed)
points file (see writeSpatialPointsParquet) it prunes row groups
outside the box. The result is a lazy DuckDBDataFrame (nothing is
materialized until collected), the R counterpart of scibis' query_points.
layerBboxRange( x, xmin, xmax, ymin, ymax, x_col = "x", y_col = "y", zmin = NULL, zmax = NULL, z_col = "z" )layerBboxRange( x, xmin, xmax, ymin, ymax, x_col = "x", y_col = "y", zmin = NULL, zmax = NULL, z_col = "z" )
x |
A DuckDBDataFrame point layer. |
xmin, xmax, ymin, ymax
|
Bounding-box limits (inclusive). |
x_col, y_col
|
Coordinate column names (default |
zmin, zmax
|
Optional inclusive depth range for a 3-D viewport; when both
are non- |
z_col |
Depth column name (default |
A lazy DuckDBDataFrame filtered to the box.
p <- tempfile(fileext = ".parquet") writeSpatialPointsParquet( data.frame(x = c(1, 5, 50), y = c(1, 5, 50), gene = c("A", "B", "C")), p) pts <- DuckDBDataFrame::DuckDBDataFrame(p) as.data.frame(layerBboxRange(pts, 0, 10, 0, 10)) unlink(p)p <- tempfile(fileext = ".parquet") writeSpatialPointsParquet( data.frame(x = c(1, 5, 50), y = c(1, 5, 50), gene = c("A", "B", "C")), p) pts <- DuckDBDataFrame::DuckDBDataFrame(p) as.data.frame(layerBboxRange(pts, 0, 10, 0, 10)) unlink(p)
Match points in a lazy table layer to a geometry table
layerSpatialMatch(x, table, coords, geom = "geometry", join = NULL)layerSpatialMatch(x, table, coords, geom = "geometry", join = NULL)
x |
A |
table |
A |
coords |
Length-two character vector of x/y column names in |
geom |
Geometry column name in |
join |
Optional |
Integer vector of match indices (NA where unmatched).
spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") shapes <- DuckDBDataFrame(spatial_path) shapes <- shapes[which(!is.na(shapes$type)), ] pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 30), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) layerSpatialMatch(pts, shapes, coords = c("x", "y"))spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") shapes <- DuckDBDataFrame(spatial_path) shapes <- shapes[which(!is.na(shapes$type)), ] pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 30), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) layerSpatialMatch(pts, shapes, coords = c("x", "y"))
Layer-level intersection predicate for DuckDBDataFrame objects. Accepts coordinate columns or a geometry column.
layerSpatialOverlaps(x, y, coords = NULL, geom = "geometry")layerSpatialOverlaps(x, y, coords = NULL, geom = "geometry")
x |
A |
y |
A geometry ( |
coords |
Optional length-two character vector of x/y column names. |
geom |
Geometry column name when |
Logical vector of length nrow(x).
pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 10), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) poly <- sf::st_as_sfc("POLYGON((0 0, 6 0, 6 6, 0 6, 0 0))") layerSpatialOverlaps(pts, poly, coords = c("x", "y"))pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 10), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) poly <- sf::st_as_sfc("POLYGON((0 0, 6 0, 6 6, 0 6, 0 0))") layerSpatialOverlaps(pts, poly, coords = c("x", "y"))
Row indices of a layer inside a bounding box
layerSubsetByBbox( x, xmin, xmax, ymin, ymax, x_col = "x", y_col = "y", geom = "geometry" )layerSubsetByBbox( x, xmin, xmax, ymin, ymax, x_col = "x", y_col = "y", geom = "geometry" )
x |
A |
xmin, xmax, ymin, ymax
|
Bounding box limits. |
x_col, y_col
|
Coordinate column names for point layers. |
geom |
Geometry column name for polygon layers. |
Integer row indices.
pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 10), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) layerSubsetByBbox(pts, 0, 10, 0, 10, x_col = "x", y_col = "y")pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 10), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) layerSubsetByBbox(pts, 0, 10, 0, 10, x_col = "x", y_col = "y")
Row indices of a layer intersecting a geometry
layerSubsetByGeometry(x, y, coords = NULL, geom = "geometry")layerSubsetByGeometry(x, y, coords = NULL, geom = "geometry")
x |
A |
y |
A geometry to test against. |
coords |
Optional x/y column names for point layers. |
geom |
Geometry column name for polygon layers. |
Integer row indices.
pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 10), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) poly <- sf::st_as_sfc("POLYGON((0 0, 6 0, 6 6, 0 6, 0 0))") layerSubsetByGeometry(pts, poly, coords = c("x", "y"))pts_path <- tempfile(fileext = ".csv") on.exit(unlink(pts_path), add = TRUE) write.csv(data.frame(x = c(1, 5, 10), y = c(1, 5, 10)), pts_path, row.names = FALSE) pts <- DuckDBDataFrame(pts_path, datacols = c("x", "y")) poly <- sf::st_as_sfc("POLYGON((0 0, 6 0, 6 6, 0 6, 0 0))") layerSubsetByGeometry(pts, poly, coords = c("x", "y"))
DuckDBDataFrame
Read a GeoParquet file as a lazy DuckDBDataFrame
readGeoParquet(path, ...)readGeoParquet(path, ...)
path |
Path to a GeoParquet file or directory. |
... |
Passed to |
spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") ddb <- readGeoParquet(spatial_path) nrow(ddb) colnames(ddb)spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial") ddb <- readGeoParquet(spatial_path) nrow(ddb) colnames(ddb)
Returns df with its rows reordered by the Morton code of
(x_col, y_col) (or (x_col, y_col, z_col) when z_col is
given and present), so that spatially neighboring points become contiguous.
Writing the reordered table to Parquet gives each row group a tight
coordinate zonemap, which lets DuckDB prune groups outside a viewport query
(see writeSpatialPointsParquet and
layerBboxRange). A pure row permutation: no column is added, so
it round-trips identically. A no-op (returns df unchanged) when
df is empty or lacks a required coordinate column.
spatialSortPoints(df, x_col = "x", y_col = "y", z_col = NULL, bits = 16L)spatialSortPoints(df, x_col = "x", y_col = "y", z_col = NULL, bits = 16L)
df |
A |
x_col, y_col
|
Coordinate column names (default |
z_col |
Optional third coordinate column for 3-D clustering; used only
when non- |
bits |
Grid resolution per axis (default 16, a |
This is the low-dimensional spatial convenience wrapper over
clusterSort(df, zorder(c(x_col, y_col)));
the Morton math lives in DuckDBDataFrame and is already N-D, so
(x, y, z) is a genuine 3-D Morton curve, not a 2-D curve with z
appended.
df reordered by Morton code (same class, same columns).
df <- data.frame(x = c(9, 1, 5), y = c(9, 1, 5), gene = c("C", "A", "B")) spatialSortPoints(df)df <- data.frame(x = c(9, 1, 5), y = c(9, 1, 5), gene = c("C", "A", "B")) spatialSortPoints(df)
Applies a 2-D affine map to a geometry column via DuckDB ST_Affine:
, . The affine may be a
coordinate-transforms object (e.g. ctScale,
ctRotation, ctSequence), lowered to its 2-D coefficients, or a
numeric matrix (2x2 linear, 2x3 linear+offset, or 3x3
homogeneous). Errors if the transform is not 2-D-expressible (e.g. a
dimensionality change). For point layers stored as x/y columns
rather than a geometry column, use transformLayer.
st_affine(x, affine, ...)st_affine(x, affine, ...)
x |
A |
affine |
A |
... |
Ignored. |
An object of the same class as x with transformed geometries.
coordinate-transforms, transformLayer
Transforms a spatial layer by a 2-D affine coordinate transform, in place and
lazily (nothing is materialized until collected). transform is a
coordinate-transforms object (ctScale,
ctTranslation, ctRotation, ctAffine, ctSequence,
…) or an affine matrix, lowered to a 2-D affine (,
). A point layer (columns x_col/
y_col) is transformed by arithmetic on those columns; a geometry layer
(column geom) by st_affine (DuckDB ST_Affine). To
move a layer between named coordinate systems, resolve the transform with
ctPath first and pass it here.
transformLayer(x, transform, x_col = "x", y_col = "y", geom = "geometry")transformLayer(x, transform, x_col = "x", y_col = "y", geom = "geometry")
x |
A DuckDBDataFrame layer. |
transform |
A |
x_col, y_col
|
Coordinate column names for point layers (default
|
geom |
Geometry column name for geometry layers (default
|
A lazy DuckDBDataFrame with transformed coordinates.
coordinate-transforms, ctGraph,
st_affine
p <- tempfile(fileext = ".parquet") writeSpatialPointsParquet(data.frame(x = c(1, 2, 3), y = c(0, 0, 0)), p) pts <- DuckDBDataFrame::DuckDBDataFrame(p) # rotate 90 degrees about the origin: (x, y) -> (-y, x) as.data.frame(transformLayer(pts, ctRotation(rbind(c(0, -1), c(1, 0))))) unlink(p)p <- tempfile(fileext = ".parquet") writeSpatialPointsParquet(data.frame(x = c(1, 2, 3), y = c(0, 0, 0)), p) pts <- DuckDBDataFrame::DuckDBDataFrame(p) # rotate 90 degrees about the origin: (x, y) -> (-y, x) as.data.frame(transformLayer(pts, ctRotation(rbind(c(0, -1), c(1, 0))))) unlink(p)
Writes sf objects or data frames with geometry columns as GeoParquet
files conforming to the GeoParquet 1.0.0 specification. The geometry column
is encoded as Well-Known Binary (WKB) and appropriate metadata is added to
ensure compatibility with DuckDB's spatial extension, GeoPandas, and other
GeoParquet readers.
writeGeoParquet( x, path, geom = "geometry", compression = "zstd", options = nanoparquet::parquet_options(compression_level = 3L), ... )writeGeoParquet( x, path, geom = "geometry", compression = "zstd", options = nanoparquet::parquet_options(compression_level = 3L), ... )
x |
An |
path |
Character string specifying the output file path. Should end
in |
geom |
Character string specifying the name of the geometry column.
Defaults to |
compression |
Character string specifying the compression algorithm.
Defaults to |
options |
Optional list of parquet options created by
|
... |
Additional arguments (currently unused). |
The function performs the following steps:
Converts the geometry column to Well-Known Binary (WKB) format
Generates GeoParquet 1.0.0 metadata including:
Geometry type(s) (Point, LineString, Polygon, etc.)
Bounding box (xmin, ymin, xmax, ymax)
Encoding format (WKB)
Primary geometry column name
Writes the data to a Parquet file using nanoparquet
The resulting files can be read by:
DuckDB: SELECT ST_AsText(geometry) FROM 'file.parquet'
GeoPandas: gpd.read_parquet('file.parquet')
GDAL/OGR: ogr2ogr -f GeoJSON output.json file.parquet
sf: sf::st_read('file.parquet')
Invisibly returns NULL. Called for its side effect of writing
a GeoParquet file to disk.
https://geoparquet.org/ - GeoParquet specification
write_parquet - Underlying Parquet writer
st_write - Alternative sf-based writer
## Not run: library(sf) # Write a simple features object nc <- st_read(system.file("shape/nc.shp", package="sf")) writeGeoParquet(nc, "nc.parquet") # Write with custom geometry column name pts <- st_sf(id = 1:3, geom = st_sfc(st_point(c(0,0)), st_point(c(1,1)), st_point(c(2,2)))) writeGeoParquet(pts, "points.parquet", geom = "geom") # Read back in DuckDB library(DuckDBDataFrame) ddb <- DuckDBDataFrame("nc.parquet") ddb$geometry # Native GEOMETRY column ## End(Not run)## Not run: library(sf) # Write a simple features object nc <- st_read(system.file("shape/nc.shp", package="sf")) writeGeoParquet(nc, "nc.parquet") # Write with custom geometry column name pts <- st_sf(id = 1:3, geom = st_sfc(st_point(c(0,0)), st_point(c(1,1)), st_point(c(2,2)))) writeGeoParquet(pts, "points.parquet", geom = "geom") # Read back in DuckDB library(DuckDBDataFrame) ddb <- DuckDBDataFrame("nc.parquet") ddb$geometry # Native GEOMETRY column ## End(Not run)
Writes a points table to a single Parquet file, Morton-sorting it first (by
default) so viewport queries prune row groups. The write goes through the
shared DuckDBDataFrame connection and
buildParquetCopySQL with a bounded
ROW_GROUP_SIZE (so there are groups to prune). With
spatial_sort = FALSE the acquisition-order baseline is written
instead.
writeSpatialPointsParquet( df, path, x_col = "x", y_col = "y", z_col = NULL, spatial_sort = TRUE, row_group_size = 131072L, bits = 16L )writeSpatialPointsParquet( df, path, x_col = "x", y_col = "y", z_col = NULL, spatial_sort = TRUE, row_group_size = 131072L, bits = 16L )
df |
A |
path |
Output Parquet file path. |
x_col, y_col
|
Coordinate column names (default |
z_col |
Optional third coordinate column for 3-D Morton clustering.
Default |
spatial_sort |
If |
row_group_size |
Parquet |
bits |
Morton grid resolution per axis (default 16). |
path, invisibly.
df <- data.frame(x = runif(1000, 0, 100), y = runif(1000, 0, 100)) p <- tempfile(fileext = ".parquet") writeSpatialPointsParquet(df, p) unlink(p)df <- data.frame(x = runif(1000, 0, 100), y = runif(1000, 0, 100)) p <- tempfile(fileext = ".parquet") writeSpatialPointsParquet(df, p) unlink(p)