7.8 Example: Olinda, Brazil

r <- read_stars(system.file("tif/L7_ETMs.tif", package = "stars"))
r
## stars object with 3 dimensions and 1 attribute
## attribute(s):
##              Min. 1st Qu. Median     Mean 3rd Qu. Max.
## L7_ETMs.tif     1      54     69 68.91242      86  255
## dimension(s):
##      from  to  offset delta                     refsys point x/y
## x       1 349  288776  28.5 SIRGAS 2000 / UTM zone 25S FALSE [x]
## y       1 352 9120761 -28.5 SIRGAS 2000 / UTM zone 25S FALSE [y]
## band    1   6      NA    NA                         NA    NA
attributes
  • from: starting index
  • to: ending index
  • offset: dimension value at the start (edge) of the first pixel
  • delta: cell size; negative delta values indicate that pixel index increases with decreasing dimension values
  • refsys: reference system
  • point: logical, indicates whether cell values have point support or cell support
  • x/y: indicates whether a dimension is associated with a spatial raster x- or y-axis

7.8.1 Plots

plot(r)

par(mfrow = c(1, 2))
plot(r, rgb = c(3,2,1), reset = FALSE, main = "RGB")    # rgb
plot(r, rgb = c(4,3,2), main = "False colour (NIR-R-G)") # false colour

7.8.2 Subsetting

Example: selects from r

  • attributes 1-2, index 101-200 for dimension 1, and index 5-10 for dimension 3
  • omitting dimension 2 means that no subsetting takes place
r[1:2, 101:200,, 5:10]

Selecting discontinuous ranges is supported only when it is a regular sequence

r[,1:100, seq(1, 250, 5), 4] |> dim() #default behavior
##    x    y band 
##  100   50    1
r[,1:100, seq(1, 250, 5), 4, drop = TRUE] |> dim()
##   x   y 
## 100  50

For selecting particular ranges of dimension values, one can use dplyr::filter

filter(r, x > 289000, x < 290000)
## stars object with 3 dimensions and 1 attribute
## attribute(s):
##              Min. 1st Qu. Median    Mean 3rd Qu. Max.
## L7_ETMs.tif     5      51     63 64.3337      75  242
## dimension(s):
##      from  to  offset delta                     refsys point x/y
## x       1  35  289004  28.5 SIRGAS 2000 / UTM zone 25S FALSE [x]
## y       1 352 9120761 -28.5 SIRGAS 2000 / UTM zone 25S FALSE [y]
## band    1   6       1     1                         NA    NA

or slice

slice(r, band, 3)
## stars object with 2 dimensions and 1 attribute
## attribute(s):
##              Min. 1st Qu. Median     Mean 3rd Qu. Max.
## L7_ETMs.tif    21      49     63 64.35886      77  255
## dimension(s):
##   from  to  offset delta                     refsys point x/y
## x    1 349  288776  28.5 SIRGAS 2000 / UTM zone 25S FALSE [x]
## y    1 352 9120761 -28.5 SIRGAS 2000 / UTM zone 25S FALSE [y]

7.8.3 Cropping

b <- st_bbox(r) |>
    st_as_sfc() |>
    st_centroid() |>
    st_buffer(units::set_units(500, m))
r[b]
## stars object with 3 dimensions and 1 attribute
## attribute(s):
##              Min. 1st Qu. Median     Mean 3rd Qu. Max. NA's
## L7_ETMs.tif    22      54     66 67.68302   78.25  174 2184
## dimension(s):
##      from  to  offset delta                     refsys point x/y
## x     157 193  288776  28.5 SIRGAS 2000 / UTM zone 25S FALSE [x]
## y     159 194 9120761 -28.5 SIRGAS 2000 / UTM zone 25S FALSE [y]
## band    1   6      NA    NA                         NA    NA
Circular centre region of the Landsat 7 scene (band 1)
Circular centre region of the Landsat 7 scene (band 1)
image code 
plot(r[b][,,,1], reset = FALSE)
plot(b, border = 'brown', lwd = 2, col = NA, add = TRUE)

By default, the resulting raster is cropped to the extent of the selection object; otherwise

r[b, crop = FALSE]
## stars object with 3 dimensions and 1 attribute
## attribute(s):
##              Min. 1st Qu. Median     Mean 3rd Qu. Max.   NA's
## L7_ETMs.tif    22      54     66 67.68302   78.25  174 731280
## dimension(s):
##      from  to  offset delta                     refsys point x/y
## x       1 349  288776  28.5 SIRGAS 2000 / UTM zone 25S FALSE [x]
## y       1 352 9120761 -28.5 SIRGAS 2000 / UTM zone 25S FALSE [y]
## band    1   6      NA    NA                         NA    NA
# we can reset dimension offsets
r[b] |> st_normalize() |> st_dimensions()
##      from to  offset delta                     refsys point x/y
## x       1 37  293222  28.5 SIRGAS 2000 / UTM zone 25S FALSE [x]
## y       1 36 9116258 -28.5 SIRGAS 2000 / UTM zone 25S FALSE [y]
## band    1  6      NA    NA                         NA    NA

or simply using the stars function made for cropping

st_crop(r, b)