## ----setup, include=FALSE-----------------------------------------------------
knitr::opts_chunk$set(collapse = TRUE, comment = "#>",
                      error = FALSE, warning = FALSE, message = FALSE)

# Probe once whether DuckDB's spatial extension is available in this build
# environment (it is fetched on first use, or loaded from a provisioned cache).
# The live spatial examples below are gated on this so the vignette still builds
# where the extension cannot be obtained (e.g. an offline developer machine);
# CI / a networked build renders them for real.
.spatial_ok <- tryCatch({
    con <- DuckDBDataFrame::acquireDuckDBConn()
    DBI::dbGetQuery(con,
        "SELECT ST_Area(ST_GeomFromText('POLYGON((0 0,0 1,1 1,1 0,0 0))')) AS a")
    TRUE
}, error = function(e) FALSE)

## ----spatial-unavailable-note, echo=FALSE, results='asis', eval=!.spatial_ok----
# cat("> **Note:** DuckDB's `spatial` extension could not be loaded in this build ",
#     "environment, so the live spatial examples below are shown but not evaluated. ",
#     "See the README for how to make the extension available (it is fetched ",
#     "automatically where the DuckDB extension repository is reachable).\n", sep = "")

## ----install, eval=FALSE------------------------------------------------------
# if (!require("BiocManager"))
#     install.packages("BiocManager")
# BiocManager::install("DuckDBSpatial")

## ----load---------------------------------------------------------------------
library(DuckDBSpatial)
library(sf)

## ----columns, eval=.spatial_ok------------------------------------------------
spatial_path <- system.file("extdata", "spatial", package = "DuckDBSpatial")
df <- DuckDBDataFrame(spatial_path)
df <- df[which(!is.na(df$type)), ]

geom <- df[["geometry"]]
head(st_geometry_type(geom))  # one of POINT, LINESTRING, POLYGON, ... per row
head(st_area(geom))           # planar area in the geometry's CRS units

## ----transforms, eval=.spatial_ok---------------------------------------------
centroids <- st_centroid(geom)  # each geometry's geometric center point
class(centroids)
head(st_as_text(centroids))     # WKT: a human-readable text form of a geometry

## ----predicates, eval=.spatial_ok---------------------------------------------
query_pt <- st_sfc(st_point(c(30, 10)))
hits <- st_intersects(geom, query_pt)
head(as.vector(hits))

## ----filter, eval=.spatial_ok-------------------------------------------------
filtered <- st_filter(df, query_pt)
nrow(filtered)

## ----layers, eval=.spatial_ok-------------------------------------------------
pts_path <- tempfile(fileext = ".csv")
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"))

poly <- st_as_sfc("POLYGON((0 0, 6 0, 6 6, 0 6, 0 0))")  # a 6x6 square, from WKT
layerSpatialOverlaps(pts, poly, coords = c("x", "y"))   # which points fall in poly
layerSubsetByGeometry(pts, poly, coords = c("x", "y"))  # row indices inside poly
unlink(pts_path)

## ----geoparquet-read, eval=.spatial_ok----------------------------------------
ddb <- readGeoParquet(spatial_path)
nrow(ddb)

## ----geoparquet-write, eval=.spatial_ok && requireNamespace("nanoparquet", quietly=TRUE)----
pts_sf <- st_sf(id = 1:2,
                geometry = st_sfc(st_point(0:1), st_point(2:3)))
path <- tempfile(fileext = ".parquet")
writeGeoParquet(pts_sf, path)
readGeoParquet(path)
unlink(path)

## ----sessioninfo--------------------------------------------------------------
sessionInfo()

