The goal of rsi is to address several repeated spatial infelicities, by providing utility functions that save you typing and help avoid repetitive stress injuries. Specifically, rsi provides:
The functions in rsi are designed around letting you use the tools you’re familiar with to process raster data using compute that you control – whether that means grabbing imagery with your laptop to add some context to a map, or grabbing tranches of data to a virtual server hosted near your data provider for lightning fast downloads. The outputs from rsi functions are standard objects – usually the file paths of raster files saved to your hard drive – meaning it’s easy to incorporate rsi into broader spatial data processing workflows.
You can install rsi via:
install.packages("rsi")
You can install the development version of rsi from GitHub using pak:
# install.packages("pak")
::pak("Permian-Global-Research/rsi") pak
The spectral_indices()
function provides a tibble with
data from the Awesome
Spectral Indices project:
library(rsi)
spectral_indices()
#> # A tibble: 231 × 9
#> application_domain bands contributor date_of_addition formula long_name
#> <chr> <list> <chr> <chr> <chr> <chr>
#> 1 vegetation <chr [2]> https://gith… 2021-11-17 (N - 0… Aerosol …
#> 2 vegetation <chr [2]> https://gith… 2021-11-17 (N - 0… Aerosol …
#> 3 water <chr [6]> https://gith… 2022-09-22 (B + G… Augmente…
#> 4 vegetation <chr [2]> https://gith… 2021-09-20 (1 / G… Anthocya…
#> 5 vegetation <chr [3]> https://gith… 2022-04-08 N * ((… Anthocya…
#> 6 vegetation <chr [4]> https://gith… 2021-05-11 (N - (… Atmosphe…
#> 7 vegetation <chr [4]> https://gith… 2021-05-14 sla * … Adjusted…
#> 8 vegetation <chr [2]> https://gith… 2022-04-08 (N * (… Advanced…
#> 9 water <chr [4]> https://gith… 2021-09-18 4.0 * … Automate…
#> 10 water <chr [5]> https://gith… 2021-09-18 B + 2.… Automate…
#> # ℹ 221 more rows
#> # ℹ 3 more variables: platforms <list>, reference <chr>, short_name <chr>
The first time spectral_indices()
is called it will
download the most up-to-date version of the spectral indices JSON file,
and then write the resulting table to a cache file in
tools::R_user_dir("rsi")
. After that,
spectral_indices()
will only download a new file if the
cache is older than 1 day, or if the update_cache
argument
is TRUE
, in order to provide the most up-to-date data as
quickly as possible. If offline, spectral_indices()
will
always fall back to the cache or, if no cache file exists, a (possibly
out-of-date) data file included in rsi itself.
Separately, the get_stac_data()
function provides a
generic interface for downloading composite images from any accessible
STAC catalog. For instance, we could download a cloud-masked composite
of Landsat’s visible layers using get_stac_data()
and a few
helper functions from rsi:
<- sf::st_point(c(-74.912131, 44.080410))
aoi <- sf::st_set_crs(sf::st_sfc(aoi), 4326)
aoi <- sf::st_buffer(sf::st_transform(aoi, 5070), 1000)
aoi
<- get_stac_data(
landsat_image
aoi,start_date = "2022-06-01",
end_date = "2022-06-30",
pixel_x_size = 30,
pixel_y_size = 30,
asset_names = c("red", "blue", "green"),
stac_source = "https://planetarycomputer.microsoft.com/api/stac/v1/",
collection = "landsat-c2-l2",
mask_band = "qa_pixel",
mask_function = landsat_mask_function,
output_filename = tempfile(fileext = ".tif"),
item_filter_function = landsat_platform_filter,
platforms = c("landsat-9", "landsat-8")
)
::plot(terra::rast(landsat_image)) terra
For common data sets, rsi also provides helper functions which
provide most of these arguments for you. For instance, that
get_stac_data()
call could be as simple as:
<- get_landsat_imagery(
landsat_image
aoi,start_date = "2022-06-01",
end_date = "2022-08-30",
output_filename = tempfile(fileext = ".tif")
)#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 0: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 1: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 2: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 3: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 4: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 5: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 6: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
#> ext 7: 1.66848e+06 1.67049e+06 2.52827e+06 2.53028e+06
::plot(terra::rast(landsat_image)) terra
Note that we’ve been plotting each band individually so far by
calling terra::plot()
. We could also use
terra::plotRGB()
(after terra::stretch()
ing
the band values) to see what this mosaic of images would look like to
the human eye:
|>
landsat_image ::rast(lyrs = c("R", "G", "B")) |>
terra::stretch() |>
terra::plotRGB() terra
By default, these functions download data from Microsoft’s Planetary
Computer API, using a number of configuration options set in
rsi_band_mapping
objects provided by the package. You can
see these default configuration options by printing the band mapping
objects, and can adjust them through arguments to any get_*
function in the package.
$planetary_computer_v1
landsat_band_mapping#> An rsi band mapping object with attributes:
#> names mask_band mask_function stac_source collection_name query_function download_function sign_function class
#>
#> coastal blue green red nir08 swir16 swir22 lwir lwir11
#> "A" "B" "G" "R" "N" "S1" "S2" "T" "T1"
We can put these pieces together and calculate as many spectral
indices as we can based on our downloaded Landsat imagery. The
calculate_indices()
function, well, calculates indices,
using subsets of our spectral_indices()
data frame:
<- filter_bands(
available_indices bands = names(terra::rast(landsat_image))
)
<- calculate_indices(
indices
landsat_image,
available_indices,output_filename = tempfile(fileext = ".tif")
)
# Plot the first handful of spatial indices
::plot(terra::rast(indices)) terra
And last but not least, rsi includes a utility for efficiently combining rasters containing different data about the same location into a VRT, which allows programs like GDAL to treat these separate data sources as a single file. For instance, we can combine our Landsat imagery with the derived indices:
<- stack_rasters(
raster_stack c(landsat_image, indices),
tempfile(fileext = ".vrt")
)
# The first few panels are now Landsat measurements, not indices:
::plot(terra::rast(raster_stack)) terra
This can be extremely useful as a way to create predictor bricks and other multi-band rasters from various data sources.
We love contributions! See our contribution guide for pointers on how to make your contribution as easy to accept as possible – in particular, consider opening an issue with a minimal reprex to make sure that we understand what your changes are meant to do.
Copyright 2023 Permian Global Research, Limited.
Licensed under the Apache License, Version 2.0 (the “License”); you may not use this file except in compliance with the License.
You may obtain a copy of the License at:
https://www.apache.org/licenses/LICENSE-2.0
Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an “AS IS” BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License.