Repository files navigation

The goal of rsi is to address severalrepeatedspatialinfelicities, by providing utility functions that save you typingand help avoidrepetitivestressinjuries. Specifically, rsiprovides:

• An interface to theRsome – excuse me,Awesome Spectral Indicesproject,providing the list of indices directly in R as a friendly tibble,
• A method for efficientlycalculating those awesome spectral indicesusing local rasters, enablingrapidspectralinference,
• A method for downloading STAC data – excuse me,retrivingSTACinformation – from any STAC server, with additional helpers fordownloading Landsat, Sentinel-1, and Sentinel-2 data from free andpublic STAC servers providingrapidsatelliteimagery,
• Arasterstackintegration method for combining multiplerasters containing distinct data sets into a single raster stack.

The functions in rsi are designed around letting you use the toolsyou’re familiar with to process raster data using compute that youcontrol – whether that means grabbing imagery with your laptop to addsome context to a map, or grabbing tranches of data to a virtual serverhosted near your data provider for lightning fast downloads. The outputsfrom rsi functions are standard objects – usually the file paths ofraster files saved to your hard drive – meaning it’s easy to incorporatersi into broader spatial data processing workflows.

You can install rsi via:

install.packages("rsi")

You can install the development version of rsi fromGitHub usingpak:

# install.packages("pak")pak::pak("Permian-Global-Research/rsi")

Thespectral_indices() function provides a tibble with data from theAwesome Spectral Indicesproject:

library(rsi)spectral_indices()#> # A tibble: 245 × 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…#> # ℹ 235 more rows#> # ℹ 3 more variables: platforms <list>, reference <chr>, short_name <chr>

The first timespectral_indices() is called it will download the mostup-to-date version of the spectral indices JSON file, and then write theresulting table to a cache file intools::R_user_dir("rsi"). Afterthat,spectral_indices() will only download a new file if the cache isolder than 1 day, or if theupdate_cache argument isTRUE, in orderto 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 cachefile exists, a (possibly out-of-date) data file included in rsi itself.

Separately, theget_stac_data() function provides a generic interfacefor downloading composite images from any accessible STAC catalog. Forinstance, we could download a cloud-masked composite of Landsat’svisible layers usingget_stac_data() and a few helper functions fromrsi:

aoi<-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)landsat_image<- get_stac_data(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"))terra::plot(terra::rast(landsat_image))

For common data sets, rsi also provides helper functions which providemost of these arguments for you. For instance, thatget_stac_data()call could be as simple as:

landsat_image<- get_landsat_imagery(aoi,start_date="2022-06-01",end_date="2022-08-30",output_filename= tempfile(fileext=".tif"))terra::plot(terra::rast(landsat_image))

Note that we’ve been plotting each band individually so far by callingterra::plot(). We could also useterra::plotRGB() (afterterra::stretch()ing the band values) to see what this mosaic of imageswould look like to the human eye:

landsat_image|>terra::rast(lyrs= c("R","G","B"))|>terra::stretch()|>terra::plotRGB()

By default, these functions download data from Microsoft’s PlanetaryComputer API, using a number of configuration options set inrsi_band_mapping objects provided by the package. You can see thesedefault configuration options by printing the band mapping objects, andcan adjust them through arguments to anyget_* function in thepackage.

landsat_band_mapping$planetary_computer_v1#> 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 indicesas we can based on our downloaded Landsat imagery. Thecalculate_indices() function, well, calculates indices, using subsetsof ourspectral_indices() data frame:

available_indices<- filter_bands(bands= names(terra::rast(landsat_image)))indices<- calculate_indices(landsat_image,available_indices,output_filename= tempfile(fileext=".tif"))#> |---------|---------|---------|---------|=========================================# Plot the first handful of spatial indicesterra::plot(terra::rast(indices))

And last but not least, rsi includes a utility for efficiently combiningrasters containing different data about the same location into aVRT, which allowsprograms like GDAL to treat these separate data sources as a singlefile. For instance, we can combine our Landsat imagery with the derivedindices:

raster_stack<- stack_rasters(  c(landsat_image,indices),  tempfile(fileext=".vrt"))# The first few panels are now Landsat measurements, not indices:terra::plot(terra::rast(raster_stack))

This can be extremely useful as a way to create predictor bricks andother multi-band rasters from various data sources.

We love contributions! See ourcontributionguidefor pointers on how to make your contribution as easy to accept aspossible – in particular, consideropening anissuewith aminimal reprex to makesure 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 maynot 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, softwaredistributed 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 andlimitations under the License.