-
Notifications
You must be signed in to change notification settings - Fork 5
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
Showing
5 changed files
with
39 additions
and
3 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -131,5 +131,6 @@ cython_debug/ | |
#.idea/ | ||
.Rproj.user | ||
*.hdf | ||
|
||
*.tif | ||
/.quarto/ | ||
.Rhistory |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,14 @@ | ||
{ | ||
"hash": "7ecd036ecddd03d27feb81daa6bce529", | ||
"result": { | ||
"markdown": "---\ntitle: \"Mosaic images using STAC\"\noutput: rmarkdown::html_vignette\nvignette: >\n %\\VignetteIndexEntry{Mosaic images using STAC}\n %\\VignetteEngine{knitr::rmarkdown}\n %\\VignetteEncoding{UTF-8}\n---\n\n\n## *WORK IN PROGRESS*\n\n\n\nHigh-resolution satellites generate many snapshot images each with a limited field of view or spatial extent. In order to see a larger area in space, and/or observe changes across space and time, we need to assemble these many snapshots into a mosaic or \"data cube\" that we can analyze as a cohesive whole. \n\n\n[EARTHDATA STAC CATALOGS](https://radiantearth.github.io/stac-browser/#/external/cmr.earthdata.nasa.gov/stac/LPCLOUD?.language=en)\n\n\n\n\n\n::: {.cell}\n\n```{.r .cell-code}\nlibrary(earthdatalogin)\nlibrary(rstac)\nlibrary(gdalcubes)\ngdalcubes_options(parallel = TRUE)\n```\n:::\n\n\n\n### Earth Data Authentication\n\nFirst let's get EDL authentication out of the way.\nFor cloud data from almost any other STAC catalog (NOAA, USGS, Planetary Computer, etc), authentication is either unnecessary or already provided by the STAC API, but NASA EDL is special. \n\n\n\n::: {.cell}\n\n```{.r .cell-code}\nlibrary(earthdatalogin)\n```\n:::\n\n\nWe could just use `edl_set_token()` here as usual to set the environmental variable. \nThis works fine but can problems if we do not remember to `edl_unset_token()` before accessing other non-EDL resources over the http interface.\nWhen using the `gdalcubes` package, we have support for a somewhat nicer, more localized authentication that uses the configuration options instead.\nWe tell `edl_set_token` not to set the environmental variable globally, but to return\nin the header format which we can pass to `gdalcubes_set_gdal_config()`:\n\n\n::: {.cell}\n\n```{.r .cell-code}\nheader <- edl_set_token(set_env_var = FALSE, format = \"header\")\ngdalcubes_set_gdal_config(\"GDAL_HTTP_HEADERS\", header)\n```\n:::\n\n\n`earthdatalogin` also includes optional configuration settings for GDAL which can improve performance of cloud-based data access. Set the GDAL environmental variables using `gdal_cloud_config()`\n\n\n::: {.cell}\n\n```{.r .cell-code}\ngdal_cloud_config()\n```\n:::\n\n\n## Search via STAC\n\nWe will now use the `rstac` package to search one or more NASA collections for data that falls into our desired bounding \n\nSet a search box in space & time\n\n\n::: {.cell}\n\n```{.r .cell-code}\nbbox <- c(xmin=-122.5, ymin=37.5, xmax=-122.0, ymax=38) \nstart <- \"2022-01-01\"\nend <- \"2022-06-30\"\n\n# Find all assets from the desired catalog:\nitems <- stac(\"https://cmr.earthdata.nasa.gov/stac/LPCLOUD\") |> \n stac_search(collections = \"HLSL30.v2.0\",\n bbox = bbox,\n datetime = paste(start,end, sep = \"/\")) |>\n post_request() |>\n items_fetch() |>\n items_filter(filter_fn = \\(x) {x[[\"eo:cloud_cover\"]] < 20})\n```\n:::\n\n\n\nNote that many features have matched our search criteria! Each feature represents a 'snapshot' image taken by the satellite as it passes by (this is a harmonized product so actually there's quite a lot of post-processing.) Each feature thus shares the same bounding box, projection, and timestamp, but may consist of many different 'assets', different COG files representing the different spectral bands on the satellite camera instrument. Each feature can potentially include quite extensive metadata about the feature, including details of instrument itself or from post-processing, such as cloud cover. Unfortunately, EarthData's STAC metadata tends to be quite sparse. \n\n\n\n## Building a Data Cube\n\n\n\n::: {.cell}\n\n```{.r .cell-code}\n# Desired data cube shape & resolution\nv = cube_view(srs = \"EPSG:4326\",\n extent = list(t0 = as.character(start), \n t1 = as.character(end),\n left = bbox[1], right = bbox[3],\n top = bbox[4], bottom = bbox[2]),\n nx = 512, ny = 512, dt = \"P1M\")\n```\n:::\n\n::: {.cell}\n\n```{.r .cell-code}\n# RGB bands + cloud cover mask\ncol <- stac_image_collection(items$features, \n asset_names = c(\"B02\", \"B03\", \"B04\", \"Fmask\"))\n```\n:::\n\n::: {.cell}\n\n```{.r .cell-code}\n# use a cloud mask -- not sure I have this correct\n# https://lpdaac.usgs.gov/documents/1326/HLS_User_Guide_V2.pdf\ncloud_mask <- image_mask(\"Fmask\", values=1) # mask clouds and cloud shadows\nrgb_bands <- c(\"B04\",\"B03\", \"B02\")\n\n# Here we go! note eval is lazy\nraster_cube(col, v, mask=cloud_mask) |>\n select_bands(rgb_bands) |>\n plot(rgb=1:3)\n```\n:::\n\n\n\n## Scaling up\n\nSame code with larger search box:\n\n\n::: {.cell}\n\n```{.r .cell-code}\nbbox <- c(xmin=-123, ymin=37, xmax=-121, ymax=39) \n\nstart <- \"2023-01-01\"\nend <- \"2023-01-31\"\nitems <- stac(\"https://cmr.earthdata.nasa.gov/stac/LPCLOUD\") |> \n stac_search(collections = \"HLSL30.v2.0\",\n bbox = c(bbox),\n datetime =paste(start,end, sep = \"/\")) |>\n post_request() |>\n items_fetch() |>\n items_filter(filter_fn = \\(x) {x[[\"eo:cloud_cover\"]] < 20})\n```\n:::\n\n::: {.cell}\n\n```{.r .cell-code}\nview <- cube_view(srs = \"EPSG:4326\",\n extent = list(t0 = as.character(start), \n t1 = as.character(end),\n left = bbox[1], right = bbox[3],\n top = bbox[4], bottom = bbox[2]),\n nx = 1024, ny = 1024, dt = \"P1M\")\nassets <- stac_image_collection(items$features, \n asset_names = c(\"B02\", \"B03\", \"B04\", \"Fmask\"))\n```\n:::\n\n::: {.cell}\n\n```{.r .cell-code}\nraster_cube(assets, view) |>\n select_bands(c(\"B04\",\"B03\", \"B02\")) |>\n plot(rgb=1:3)\n```\n:::\n", | ||
"supporting": [], | ||
"filters": [ | ||
"rmarkdown/pagebreak.lua" | ||
], | ||
"includes": {}, | ||
"engineDependencies": {}, | ||
"preserve": {}, | ||
"postProcess": true | ||
} | ||
} |
Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.
Oops, something went wrong.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,13 @@ | ||
Version: 1.0 | ||
|
||
RestoreWorkspace: Default | ||
SaveWorkspace: Default | ||
AlwaysSaveHistory: Default | ||
|
||
EnableCodeIndexing: Yes | ||
UseSpacesForTab: Yes | ||
NumSpacesForTab: 2 | ||
Encoding: UTF-8 | ||
|
||
RnwWeave: Sweave | ||
LaTeX: pdfLaTeX |