ace example

_freeze/contents/ace/execute-results/html.json

@@ -0,0 +1,16 @@
+{
+  "hash": "77107a53791af3118af7d6dc5a464fab",
+  "result": {
+    "markdown": "---\n\n---\n\n\n\n::: {.cell}\n\n```{.r .cell-code}\nlibrary(stars)\nlibrary(rstac)\nlibrary(gdalcubes)\nlibrary(tmap)\nlibrary(dplyr)\ntmap_mode(\"plot\")\n\ngdalcubes::gdalcubes_options(parallel=24*2)\n```\n:::\n\n\nCDFW makes the [ACE GIS Data](https://wildlife.ca.gov/Data/Analysis/Ace) freely available, but in one large zip archive.\n\nHere I'll use a public mirror on a Berkeley-based server for a more cloud-native access pattern so we don't have to download the whole thing.  We also index this in a simple [stac catalog](https://radiantearth.github.io/stac-browser/#/external/raw.githubusercontent.com/schmidtDSE/biodiversity-catalog/main/stac/v1/collections/cdfw_ace/summary.json?.asset=asset-ds2715)\n\n\n\n::: {.cell}\n\n```{.r .cell-code}\nurl <- \"/vsicurl/https://minio.carlboettiger.info/public-biodiversity/ACE_summary/ds2715.gdb\"\n```\n:::\n\n\nCDFW uses hexes!  Note that these pre-date the modern [H3 hex tiling system](https://h3geo.org/), so do not provide the magic of [hierarchical (zoom level) tiles](https://h3geo.org/docs/highlights/indexing)\n\n\n::: {.cell}\n\n```{.r .cell-code}\nhexes <- st_read(url) \n```\n\n::: {.cell-output .cell-output-stdout}\n```\nReading layer `ds2715' from data source \n  `/vsicurl/https://minio.carlboettiger.info/public-biodiversity/ACE_summary/ds2715.gdb' \n  using driver `OpenFileGDB'\nSimple feature collection with 63890 features and 17 fields\nGeometry type: MULTIPOLYGON\nDimension:     XY\nBounding box:  xmin: -373987.9 ymin: -604495.8 xmax: 540034.3 ymax: 450023.2\nProjected CRS: NAD83 / California Albers\n```\n:::\n\n```{.r .cell-code}\ntm_shape(hexes) + tm_polygons(fill = \"RwiRankEco\", col=NA)\n```\n\n::: {.cell-output-display}\n![](ace_files/figure-html/unnamed-chunk-4-1.png){width=672}\n:::\n:::\n\n\nFind Biodiversity Intactness Index COG tiles from Planetary Computer using STAC search:\n\n\n::: {.cell}\n\n```{.r .cell-code}\nbox <- st_bbox(st_transform(hexes, 4326))\nitems <-\n  stac(\"https://planetarycomputer.microsoft.com/api/stac/v1\") |>\n  stac_search(collections = \"io-biodiversity\",\n              bbox = c(box),\n              limit = 1000) |>\n  post_request() |>\n  items_sign(sign_fn = sign_planetary_computer())\n```\n:::\n\n\nCreated a desired cube view (Note, for simplicitly I show this at a much coarser resolution than the native data, allowing the gdalwarper to aggregate so this will run much faster)\n\n\n::: {.cell}\n\n```{.r .cell-code}\ncube <- gdalcubes::stac_image_collection(items$features, asset_names = \"data\")\nv <- cube_view(srs = \"EPSG:4326\",\n               extent = list(t0 = \"2020-01-01\", t1 = \"2020-12-31\",\n                             left = box[1], right = box[3],\n                             top = box[4], bottom = box[2]),\n               nx = 1000, ny=1000, dt = \"P1Y\")\n\nQ <- raster_cube(cube,v)\n```\n:::\n\n\nQuick plot of the data at requested cube resolution.  \n\n\n::: {.cell}\n\n```{.r .cell-code}\nQ |> plot(col=viridisLite::viridis(10))\n```\n\n::: {.cell-output-display}\n![](ace_files/figure-html/unnamed-chunk-7-1.png){width=672}\n:::\n:::\n\n\nExtract mean value of the BII for each hex.\n\n\n::: {.cell}\n\n```{.r .cell-code}\nbii_hexes <- Q |> gdalcubes::extract_geom(hexes, FUN=mean)\n```\n:::\n\n\n\nPlot\n\n\n::: {.cell}\n\n```{.r .cell-code}\nbii_hexes <- hexes |>\n  tibble::rowid_to_column(\"FID\") |> \n  left_join(bii_hexes)\n\ntm_shape(bii_hexes) + tm_polygons(fill = \"data\", col=NA)\n```\n\n::: {.cell-output-display}\n![](ace_files/figure-html/unnamed-chunk-9-1.png){width=672}\n:::\n:::\n\n::: {.cell}\n\n```{.r .cell-code}\ntmap_mode(\"plot\")\nsf = st_bbox(c(xmin=-123, ymin=39, xmax=-122, ymax=37.5), crs=4326)\ntm_shape(bii_hexes, bbox = sf) + tm_polygons(fill = \"data\", col=NA)\n```\n\n::: {.cell-output-display}\n![](ace_files/figure-html/unnamed-chunk-10-1.png){width=672}\n:::\n:::\n\n\n\n## Leaflet map\n\n\n::: {.cell}\n\n```{.r .cell-code}\ntmap_mode(\"view\")\ntm_shape(bii_hexes) + tm_polygons(fill = \"data\", col=NA)\n```\n:::\n",
+    "supporting": [
+      "ace_files"
+    ],
+    "filters": [
+      "rmarkdown/pagebreak.lua"
+    ],
+    "includes": {},
+    "engineDependencies": {},
+    "preserve": {},
+    "postProcess": true
+  }
+}

contents/ace.qmd

@@ -0,0 +1,100 @@
+---
+
+---
+
+```{r include=FALSE}
+knitr::opts_chunk$set(message = FALSE, warning = FALSE)
+```
+
+```{r message=FALSE}
+library(stars)
+library(rstac)
+library(gdalcubes)
+library(tmap)
+library(dplyr)
+tmap_mode("plot")
+
+gdalcubes::gdalcubes_options(parallel=24*2)
+```
+
+CDFW makes the [ACE GIS Data](https://wildlife.ca.gov/Data/Analysis/Ace) freely available, but in one large zip archive.
+
+Here I'll use a public mirror on a Berkeley-based server for a more cloud-native access pattern so we don't have to download the whole thing.  We also index this in a simple [stac catalog](https://radiantearth.github.io/stac-browser/#/external/raw.githubusercontent.com/schmidtDSE/biodiversity-catalog/main/stac/v1/collections/cdfw_ace/summary.json?.asset=asset-ds2715)
+
+
+```{r}
+url <- "/vsicurl/https://minio.carlboettiger.info/public-biodiversity/ACE_summary/ds2715.gdb"
+```
+
+CDFW uses hexes!  Note that these pre-date the modern [H3 hex tiling system](https://h3geo.org/), so do not provide the magic of [hierarchical (zoom level) tiles](https://h3geo.org/docs/highlights/indexing)
+
+```{r}
+hexes <- st_read(url) 
+tm_shape(hexes) + tm_polygons(fill = "RwiRankEco", col=NA)
+```
+
+Find Biodiversity Intactness Index COG tiles from Planetary Computer using STAC search:
+
+```{r}
+box <- st_bbox(st_transform(hexes, 4326))
+items <-
+  stac("https://planetarycomputer.microsoft.com/api/stac/v1") |>
+  stac_search(collections = "io-biodiversity",
+              bbox = c(box),
+              limit = 1000) |>
+  post_request() |>
+  items_sign(sign_fn = sign_planetary_computer())
+
+```
+
+Created a desired cube view (Note, for simplicitly I show this at a much coarser resolution than the native data, allowing the gdalwarper to aggregate so this will run much faster)
+
+```{r}
+cube <- gdalcubes::stac_image_collection(items$features, asset_names = "data")
+v <- cube_view(srs = "EPSG:4326",
+               extent = list(t0 = "2020-01-01", t1 = "2020-12-31",
+                             left = box[1], right = box[3],
+                             top = box[4], bottom = box[2]),
+               nx = 1000, ny=1000, dt = "P1Y")
+
+Q <- raster_cube(cube,v)
+```
+
+Quick plot of the data at requested cube resolution.  
+
+```{r}
+Q |> plot(col=viridisLite::viridis(10))
+```
+
+Extract mean value of the BII for each hex.
+
+```{r}
+bii_hexes <- Q |> gdalcubes::extract_geom(hexes, FUN=mean)
+```
+
+
+Plot
+
+```{r}
+bii_hexes <- hexes |>
+  tibble::rowid_to_column("FID") |> 
+  left_join(bii_hexes)
+
+tm_shape(bii_hexes) + tm_polygons(fill = "data", col=NA)
+```
+
+```{r}
+tmap_mode("plot")
+sf = st_bbox(c(xmin=-123, ymin=39, xmax=-122, ymax=37.5), crs=4326)
+tm_shape(bii_hexes, bbox = sf) + tm_polygons(fill = "data", col=NA)
+
+```
+
+
+## Leaflet map
+
+```{r}
+tmap_mode("view")
+tm_shape(bii_hexes) + tm_polygons(fill = "data", col=NA)
+```
+