# install.packages("rstac")
# install.packages("tidyverse")
# install.packages("stars")
# install.packages("terra")More examples analysing EOPF STAC Zarr data with R
By: @sharlagelfand
Introduction
This tutorial expands on the previous tutorials (Access the EOPF Zarr STAC API with R and Access and analyse EOPF STAC Zarr data with R), going into further details on analysing and visualising Zarr data from the EOPF Sample Service STAC catalog programmatically using R. We recommend reviewing the previous tutorials if you have not done so already.
What we will learn
- đź—‚ How to extract measurements, such as Ocean Wind field and GIFAPAR, along with latitude and longitude
- 🔎 How to open and visualise quicklook composite images
- đź“Š How to format, scale, and visualise satellite data on a curvillinear grid
Prerequisites
An R environment is required to follow this tutorial, with R version >= 4.5.0. We recommend using either RStudio or Positron (or a cloud computing environment) and making use of RStudio projects for a self-contained coding environment.
Dependencies
We will use the following packages in this tutorial:
rstac(for accessing the STAC catalog)tidyverse(for data manipulation)stars) (for working with spatiotemporal data)terra(for working with spatial - data in raster format)
You can install them directly from CRAN:
We will also use the Rarr package (version >= 1.10.1) to read Zarr data. It must be installed from Bioconductor, so first install the BiocManager package:
# install.packages("BiocManager")Then, use this package to install Rarr:
# BiocManager::install("Rarr")Finally, load the packages into your environment:
library(rstac)
library(tidyverse)
library(Rarr)
library(stars)
library(terra)Sentinel-1
The first example looks at Sentinel-1 Level 2 Ocean (OCN) data, which consists of data for oceanographic study, such as monitoring sea surface conditions, detecting oil spills, and studying ocean currents. This example will show how to access and plot Wind Direction data.
First, select the relevant collection and item from STAC:
l2_ocn <- stac("https://stac.core.eopf.eodc.eu/") |>
collections(collection_id = "sentinel-1-l2-ocn") |>
items(feature_id = "S1C_IW_OCN__2SDH_20250921T194014_20250921T194043_004227_008634_1302") |>
get_request()
l2_ocn###Item
- id: S1C_IW_OCN__2SDH_20250921T194014_20250921T194043_004227_008634_1302
- collection: sentinel-1-l2-ocn
- bbox: xmin: -36.82763, ymin: 66.72476, xmax: -29.99337, ymax: 68.93933
- datetime: 2025-09-21T19:40:29.291309Z
- assets: osw, owi, rvl, product, zipped_product, product_metadata
- item's fields:
assets, bbox, collection, geometry, id, links, properties, stac_extensions, stac_version, typeWe can look at each of the assets’ titles to understand what the item contains:
l2_ocn |>
pluck("assets") |>
map("title")$osw
[1] "Ocean Swell spectra"
$owi
[1] "Ocean Wind field"
$rvl
[1] "Surface Radial Velocity"
$product
[1] "EOPF Product"
$zipped_product
[1] "Zipped EOPF Product"
$product_metadata
[1] "Consolidated Metadata"We are interested in the “Ocean Wind field” data, and will hold onto the owi key for now.
To access all of the owi data, we get the “product” asset and then the full Zarr store, again using our helper function from the previous tutorial to extract array information from the full array path:
derive_store_array <- function(store, product_url) {
store |>
mutate(array = str_remove(path, product_url)) |>
relocate(array, .before = path)
}
l2_ocn_url <- l2_ocn |>
assets_select(asset_names = "product") |>
assets_url()
l2_ocn_store <- l2_ocn_url |>
zarr_overview(as_data_frame = TRUE) |>
derive_store_array(l2_ocn_url)
l2_ocn_store# A tibble: 114 Ă— 8
array path data_type endianness compressor dim chunk_dim nchunks
<chr> <chr> <chr> <chr> <chr> <lis> <list> <list>
1 /osw/S01SIWOCN… http… float32 little blosc <int> <int [3]> <dbl>
2 /osw/S01SIWOCN… http… float32 little blosc <int> <int [2]> <dbl>
3 /osw/S01SIWOCN… http… float32 little blosc <int> <int [2]> <dbl>
4 /osw/S01SIWOCN… http… float32 little blosc <int> <int [2]> <dbl>
5 /osw/S01SIWOCN… http… float32 little blosc <int> <int [2]> <dbl>
6 /osw/S01SIWOCN… http… float32 little blosc <int> <int [2]> <dbl>
7 /osw/S01SIWOCN… http… float32 little blosc <int> <int [5]> <dbl>
8 /osw/S01SIWOCN… http… float32 little blosc <int> <int [5]> <dbl>
9 /osw/S01SIWOCN… http… float32 little blosc <int> <int [3]> <dbl>
10 /osw/S01SIWOCN… http… float32 little blosc <int> <int [3]> <dbl>
# â„ą 104 more rowsNext, we filter to access owi measurement data only:
l2_ocn_store |>
filter(str_starts(array, "/owi"), str_detect(array, "measurements"))# A tibble: 4 Ă— 8
array path data_type endianness compressor dim chunk_dim nchunks
<chr> <chr> <chr> <chr> <chr> <lis> <list> <list>
1 /owi/S01SIWOCN_… http… float32 little blosc <int> <int [2]> <dbl>
2 /owi/S01SIWOCN_… http… float32 little blosc <int> <int [2]> <dbl>
3 /owi/S01SIWOCN_… http… float32 little blosc <int> <int [2]> <dbl>
4 /owi/S01SIWOCN_… http… float32 little blosc <int> <int [2]> <dbl> Since all of these arrays start with a long ID, we can remove that to get a clearer idea of what each array is:
owi <- l2_ocn_store |>
filter(str_starts(array, "/owi"), str_detect(array, "measurements"))
array_id_prefix <- str_split(owi[["array"]], "measurements", simplify = TRUE)[, 1] |>
unique()
array_id_prefix[1] "/owi/S01SIWOCN_20250921T194014_0029_C024_1302_008634_HH/"array_id_prefix <- paste0(array_id_prefix, "measurements/")
owi <- owi |>
mutate(array = str_remove(array, array_id_prefix))
owi# A tibble: 4 Ă— 8
array path data_type endianness compressor dim chunk_dim nchunks
<chr> <chr> <chr> <chr> <chr> <lis> <list> <list>
1 latitude https:… float32 little blosc <int> <int [2]> <dbl>
2 longitude https:… float32 little blosc <int> <int [2]> <dbl>
3 wind_direction https:… float32 little blosc <int> <int [2]> <dbl>
4 wind_speed https:… float32 little blosc <int> <int [2]> <dbl> We are interested in wind_direction, as well as the coordinate arrays (latitude and longitude). We can get an overview of the arrays’ dimensions and structures:
owi |>
filter(array == "wind_direction") |>
pull(path) |>
zarr_overview()Type: Array
Path: https://objects.eodc.eu:443/e05ab01a9d56408d82ac32d69a5aae2a:202509-s01siwocn-eu/21/products/cpm_v256/S1C_IW_OCN__2SDH_20250921T194014_20250921T194043_004227_008634_1302.zarr/owi/S01SIWOCN_20250921T194014_0029_C024_1302_008634_HH/measurements/wind_direction
Shape: 194 x 254
Chunk Shape: 194 x 254
No. of Chunks: 1 (1 x 1)
Data Type: float32
Endianness: little
Compressor: bloscowi |>
filter(array == "latitude") |>
pull(path) |>
zarr_overview()Type: Array
Path: https://objects.eodc.eu:443/e05ab01a9d56408d82ac32d69a5aae2a:202509-s01siwocn-eu/21/products/cpm_v256/S1C_IW_OCN__2SDH_20250921T194014_20250921T194043_004227_008634_1302.zarr/owi/S01SIWOCN_20250921T194014_0029_C024_1302_008634_HH/measurements/latitude
Shape: 194 x 254
Chunk Shape: 194 x 254
No. of Chunks: 1 (1 x 1)
Data Type: float32
Endianness: little
Compressor: bloscowi |>
filter(array == "longitude") |>
pull(path) |>
zarr_overview()Type: Array
Path: https://objects.eodc.eu:443/e05ab01a9d56408d82ac32d69a5aae2a:202509-s01siwocn-eu/21/products/cpm_v256/S1C_IW_OCN__2SDH_20250921T194014_20250921T194043_004227_008634_1302.zarr/owi/S01SIWOCN_20250921T194014_0029_C024_1302_008634_HH/measurements/longitude
Shape: 194 x 254
Chunk Shape: 194 x 254
No. of Chunks: 1 (1 x 1)
Data Type: float32
Endianness: little
Compressor: bloscHere, we can see that all of the arrays are of the same shape: 194 x 254, with only one chunk. Since these are small, we can read all of the data in at once.
owi_wind_direction <- owi |>
filter(array == "wind_direction") |>
pull(path) |>
read_zarr_array()
owi_wind_direction[1:5, 1:5] [,1] [,2] [,3] [,4] [,5]
[1,] NaN NaN NaN NaN NaN
[2,] NaN NaN NaN NaN NaN
[3,] NaN NaN NaN NaN NaN
[4,] NaN NaN NaN NaN NaN
[5,] NaN NaN NaN NaN NaNowi_lat <- owi |>
filter(array == "latitude") |>
pull(path) |>
read_zarr_array()
owi_lat[1:5, 1:5] [,1] [,2] [,3] [,4] [,5]
[1,] 66.73324 66.73544 66.73763 66.73981 66.74200
[2,] 66.74212 66.74432 66.74651 66.74870 66.75089
[3,] 66.75100 66.75320 66.75540 66.75759 66.75978
[4,] 66.75988 66.76208 66.76428 66.76648 66.76867
[5,] 66.76875 66.77097 66.77317 66.77537 66.77757owi_long <- owi |>
filter(array == "longitude") |>
pull(path) |>
read_zarr_array()
owi_lat[1:5, 1:5] [,1] [,2] [,3] [,4] [,5]
[1,] 66.73324 66.73544 66.73763 66.73981 66.74200
[2,] 66.74212 66.74432 66.74651 66.74870 66.75089
[3,] 66.75100 66.75320 66.75540 66.75759 66.75978
[4,] 66.75988 66.76208 66.76428 66.76648 66.76867
[5,] 66.76875 66.77097 66.77317 66.77537 66.77757As described in the previous R tutorial, Zarr data arrays are often packed or compressed in order to limit space, and may need to be scaled or offset to their actual physical units or meaningful values.
This information is contained in the metadata associated with the Zarr store. We created a helper function to obtain these values, setting the offset to 0 and scale to 1 if they do not need to be offset or scaled.
get_scale_and_offset <- function(zarr_url, array) {
metadata <- Rarr:::.read_zmetadata(
zarr_url,
s3_client = Rarr:::.create_s3_client(zarr_url)
)
metadata <- metadata[["metadata"]]
array_metadata <- metadata[[paste0(array, "/.zattrs")]]
scale <- array_metadata[["scale_factor"]]
scale <- ifelse(is.null(scale), 1, scale)
offset <- array_metadata[["add_offset"]]
offset <- ifelse(is.null(offset), 0, offset)
list(
scale = scale,
offset = offset
)
}
get_scale_and_offset(l2_ocn_url, paste0(array_id_prefix, "wind_direction"))$scale
[1] 1
$offset
[1] 0get_scale_and_offset(l2_ocn_url, paste0(array_id_prefix, "latitude"))$scale
[1] 1
$offset
[1] 0get_scale_and_offset(l2_ocn_url, paste0(array_id_prefix, "longitude"))$scale
[1] 1
$offset
[1] 0None of this data need to be scaled or offset.
Note that both longitude and latitude are 2-dimensional arrays, and they are not evenly spaced. Rather, the data grid is curvilinear — it has grid lines that are not straight, and there is a longitude and latitude for every pixel of the other layers (i.e., wind_direction). This format is very common in satellite data.
We use functions from the stars package, loaded earlier, to format the data for visualisation. stars is specifically designed for reading, manipulating, and plotting spatiotemporal data, such as satellite data.
The function st_as_stars() is used to get our data into the correct format for visualisation:
owi_stars <- st_as_stars(wind_direction = owi_wind_direction) |>
st_as_stars(curvilinear = list(X1 = owi_long, X2 = owi_lat))Getting the data into this format is also beneficial because it allows for a quick summary of the data and its attributes, providing information such as the median and mean wind_direction, the number of NAs, and information on the grid:
owi_starsstars object with 2 dimensions and 1 attribute
attribute(s):
Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
wind_direction 0.6547565 49.15513 59.34981 59.91768 67.66366 359.6461 37952
dimension(s):
from to refsys point values x/y
X1 1 194 WGS 84 (CRS84) FALSE [194x254] -36.8,...,-30.03 [x]
X2 1 254 WGS 84 (CRS84) FALSE [194x254] 66.73,...,68.92 [y]
curvilinear gridFinally, we can plot this object:
plot(owi_stars, main = "Wind Direction", as_points = FALSE, axes = TRUE, breaks = "equal", col = hcl.colors)
Sentinel-2
For this example, we return to the Sentinel-2 Level-2A Collection. The Sentinel-2 mission is based on two satellites with 13 spectral bands, with four bands at 10-metre resolution, six bands at 20-metres resolution, and three bands at 60-metre resolution. The mission supports applications for land services, including the monitoring of vegetation, soil and water cover, as well as the observation of inland waterways and coastal areas.
EOPF Zarr assets include quicklook RGB composites, which are readily viewable representations of the satellite image. We will open the 10-metre resolution quicklook and visualise it. This is available as an asset, so we can access it directly from the STAC item.
s2_l2a_item <- stac("https://stac.core.eopf.eodc.eu/") |>
collections(collection_id = "sentinel-2-l2a") |>
items(feature_id = "S2B_MSIL2A_20250530T101559_N0511_R065_T32TPT_20250530T130924") |>
get_request()
tci_10m_asset <- s2_l2a_item |>
assets_select(asset_names = "TCI_10m")
tci_10m_url <- tci_10m_asset |>
assets_url()
tci_10m_url |>
zarr_overview()Type: Array
Path: https://objects.eodc.eu:443/e05ab01a9d56408d82ac32d69a5aae2a:202505-s02msil2a/30/products/cpm_v256/S2B_MSIL2A_20250530T101559_N0511_R065_T32TPT_20250530T130924.zarr/quality/l2a_quicklook/r10m/tci
Shape: 3 x 10980 x 10980
Chunk Shape: 1 x 1830 x 1830
No. of Chunks: 108 (3 x 6 x 6)
Data Type: uint8
Endianness: NA
Compressor: bloscFrom the overview, we can see that the quicklook array has three dimensions to it, each of size 10980 x 10980. The three dimensions correspond to red, green, and blue spectral bands (B04, B03, and B02, respectively), since this is an RGB composite. This information is also available by looking at the assets’ bands:
s2_l2a_item[["assets"]][["TCI_10m"]][["bands"]] |>
map_dfr(as_tibble)# A tibble: 3 Ă— 5
name description `eo:common_name` `eo:center_wavelength`
<chr> <chr> <chr> <dbl>
1 B04 Red (band 4) red 0.665
2 B03 Green (band 3) green 0.56
3 B02 Blue (band 2) blue 0.49
# ℹ 1 more variable: `eo:full_width_half_max` <dbl>We can read in a small chunk of the array to get an idea of its shape, using the same indexing process we’ve used before. Note that we want to select all of the bands (the first dimension listed). Rather than writing 1:3, we can simply use NULL as the first dimension, indicating to get all data at this dimension. To preview the data, we will just get the first 2 entries (in each dimension) along the three bands.
tci_10m_preview <- tci_10m_url |>
read_zarr_array(list(NULL, 1:2, 1:2))
tci_10m_preview, , 1
[,1] [,2]
[1,] 16 20
[2,] 31 34
[3,] 15 19
, , 2
[,1] [,2]
[1,] 16 18
[2,] 30 34
[3,] 13 20For visualisation purposes, we need the data in a different configuration — note the dimensions of the data:
dim(tci_10m_preview)[1] 3 2 2Instead, we need to get it into e.g. 2 x 2 x 3, with the third dimension reflecting the number of bands (or layers) To do this, we use the aperm() function to transpose an array, with argument c(2, 3, 1) – moving the second dimension to the first, the third to the second, and the first to the third. Then, we can see that the dimensions of the array are correct:
tci_10m_preview_perm <- tci_10m_preview |>
aperm(c(2, 3, 1))
tci_10m_preview_perm, , 1
[,1] [,2]
[1,] 16 16
[2,] 20 18
, , 2
[,1] [,2]
[1,] 31 30
[2,] 34 34
, , 3
[,1] [,2]
[1,] 15 13
[2,] 19 20dim(tci_10m_preview_perm)[1] 2 2 3Let’s read in the full TCI array to visualise it.
tci_10m <- tci_10m_url |>
read_zarr_array()
tci_10m_perm <- tci_10m |>
aperm(c(2, 3, 1))
dim(tci_10m)[1] 3 10980 10980For visualisation, we use terra’s plotRGB() function, first converting the array into a raster object with rast():
tci_10m_perm |>
rast() |>
plotRGB()
We can do the same with the quicklook at the 60-metre resolution, first accessing the full Zarr store for the product.
s2_l2a_url <- s2_l2a_item |>
assets_select(asset_names = "product") |>
assets_url()
s2_l2a_zarr_store <- s2_l2a_url |>
zarr_overview(as_data_frame = TRUE) |>
derive_store_array(s2_l2a_url)Then, the full visualisation process can be done in a single step:
s2_l2a_zarr_store |>
filter(array == "/quality/l2a_quicklook/r60m/tci") |>
pull(path) |>
read_zarr_array() |>
aperm(c(2, 3, 1)) |>
rast() |>
plotRGB()
Sentinel-3
Finally, we look at an example from the Sentinel-3 mission. The Sentinel-3 mission measures sea-surface topography and land- and sea-surface temperature and colour, in support of environmental and climate monitoring. The Sentinel-3 OLCI L2 LFR product provides this data, computed for full resolution.
Again, we will access a specific item from this collection:
l2_lfr <- stac("https://stac.core.eopf.eodc.eu/") |>
collections(collection_id = "sentinel-3-olci-l2-lfr") |>
items(feature_id = "S3B_OL_2_LFR____20260105T103813_20260105T104113_20260106T120044_0179_115_165_2160_ESA_O_NT_003") |>
get_request()
l2_lfr###Item
- id:
S3B_OL_2_LFR____20260105T103813_20260105T104113_20260106T120044_0179_115_165_2160_ESA_O_NT_003
- collection: sentinel-3-olci-l2-lfr
- bbox: xmin: -14.13710, ymin: 39.53910, xmax: 5.86104, ymax: 52.44410
- datetime: 2026-01-05T10:39:42.550682Z
- assets:
iwv, lagp, lqsf, otci, rc681, rc865, gifapar, product, zipped_product, product_metadata
- item's fields:
assets, bbox, collection, geometry, id, links, properties, stac_extensions, stac_version, typeTo access all of the data, we get the “product” asset and then the full Zarr store, again using our helper function to extract array information from the full array path:
l2_lfr_url <- l2_lfr |>
assets_select(asset_names = "product") |>
assets_url()
l2_lfr_store <- l2_lfr_url |>
zarr_overview(as_data_frame = TRUE) |>
derive_store_array(l2_lfr_url)
l2_lfr_store# A tibble: 46 Ă— 8
array path data_type endianness compressor dim chunk_dim nchunks
<chr> <chr> <chr> <chr> <chr> <lis> <list> <list>
1 /conditions/ge… http… int32 little blosc <int> <int [2]> <dbl>
2 /conditions/ge… http… int32 little blosc <int> <int [2]> <dbl>
3 /conditions/ge… http… int32 little blosc <int> <int [2]> <dbl>
4 /conditions/ge… http… uint32 little blosc <int> <int [2]> <dbl>
5 /conditions/ge… http… int32 little blosc <int> <int [2]> <dbl>
6 /conditions/ge… http… uint32 little blosc <int> <int [2]> <dbl>
7 /conditions/im… http… int16 little blosc <int> <int [2]> <dbl>
8 /conditions/im… http… int16 little blosc <int> <int [2]> <dbl>
9 /conditions/im… http… int8 <NA> blosc <int> <int [2]> <dbl>
10 /conditions/im… http… int32 little blosc <int> <int [2]> <dbl>
# â„ą 36 more rowsNext, we filter to access measurement data only:
l2_lfr_measurements <- l2_lfr_store |>
filter(str_starts(array, "/measurements")) |>
mutate(array = str_remove(array, "/measurements/"))
l2_lfr_measurements# A tibble: 9 Ă— 8
array path data_type endianness compressor dim chunk_dim nchunks
<chr> <chr> <chr> <chr> <chr> <lis> <list> <list>
1 altitude https://ob… int16 little blosc <int> <int [2]> <dbl>
2 gifapar https://ob… uint8 <NA> blosc <int> <int [2]> <dbl>
3 iwv https://ob… uint8 <NA> blosc <int> <int [2]> <dbl>
4 latitude https://ob… int32 little blosc <int> <int [2]> <dbl>
5 longitude https://ob… int32 little blosc <int> <int [2]> <dbl>
6 otci https://ob… uint8 <NA> blosc <int> <int [2]> <dbl>
7 rc681 https://ob… uint16 little blosc <int> <int [2]> <dbl>
8 rc865 https://ob… uint16 little blosc <int> <int [2]> <dbl>
9 time_stamp https://ob… int64 little blosc <int> <int [1]> <dbl> Of these, we are interested in Green Instantaneous FAPAR (GIFAPAR). FAPAR is the fraction of absorbed photosynthetically active radiation in the plant canopy. We extract gifapar as well as longitude and latitude. We can get an overview of the arrays’ dimensions and structures:
l2_lfr_measurements |>
filter(array == "gifapar") |>
pull(path) |>
zarr_overview()Type: Array
Path: https://objects.eodc.eu:443/e05ab01a9d56408d82ac32d69a5aae2a:202601-s03olclfr-eu/05/products/cpm_v262/S3B_OL_2_LFR____20260105T103813_20260105T104113_20260106T120044_0179_115_165_2160_ESA_O_NT_003.zarr/measurements/gifapar
Shape: 4090 x 4865
Chunk Shape: 1024 x 1024
No. of Chunks: 20 (4 x 5)
Data Type: uint8
Endianness: NA
Compressor: bloscl2_lfr_measurements |>
filter(array == "longitude") |>
pull(path) |>
zarr_overview()Type: Array
Path: https://objects.eodc.eu:443/e05ab01a9d56408d82ac32d69a5aae2a:202601-s03olclfr-eu/05/products/cpm_v262/S3B_OL_2_LFR____20260105T103813_20260105T104113_20260106T120044_0179_115_165_2160_ESA_O_NT_003.zarr/measurements/longitude
Shape: 4090 x 4865
Chunk Shape: 1024 x 1024
No. of Chunks: 20 (4 x 5)
Data Type: int32
Endianness: little
Compressor: bloscl2_lfr_measurements |>
filter(array == "latitude") |>
pull(path) |>
zarr_overview()Type: Array
Path: https://objects.eodc.eu:443/e05ab01a9d56408d82ac32d69a5aae2a:202601-s03olclfr-eu/05/products/cpm_v262/S3B_OL_2_LFR____20260105T103813_20260105T104113_20260106T120044_0179_115_165_2160_ESA_O_NT_003.zarr/measurements/latitude
Shape: 4090 x 4865
Chunk Shape: 1024 x 1024
No. of Chunks: 20 (4 x 5)
Data Type: int32
Endianness: little
Compressor: bloscSimilar to the previous example, we can see that all of the arrays are of the same shape: 4090 x 4865. We read in all of the arrays:
gifapar <- l2_lfr_measurements |>
filter(array == "gifapar") |>
pull(path) |>
read_zarr_array()
gifapar_long <- l2_lfr_measurements |>
filter(array == "longitude") |>
pull(path) |>
read_zarr_array()
gifapar_long[1:5, 1:5] [,1] [,2] [,3] [,4] [,5]
[1,] -12519564 -12515610 -12511656 -12507702 -12503748
[2,] -12519948 -12515994 -12512040 -12508087 -12504133
[3,] -12520333 -12516379 -12512425 -12508471 -12504518
[4,] -12520717 -12516763 -12512810 -12508856 -12504903
[5,] -12521101 -12517147 -12513194 -12509241 -12505288gifapar_lat <- l2_lfr_measurements |>
filter(array == "latitude") |>
pull(path) |>
read_zarr_array()
gifapar_lat[1:5, 1:5] [,1] [,2] [,3] [,4] [,5]
[1,] 52444051 52443830 52443608 52443386 52443165
[2,] 52441491 52441270 52441048 52440827 52440605
[3,] 52438931 52438710 52438488 52438267 52438045
[4,] 52436371 52436150 52435928 52435707 52435485
[5,] 52433811 52433590 52433368 52433147 52432925We can immediately tell that the longitude and latitude will need to be scaled, since they are not typical values. We again find the scale and offset values for this data:
gifapar_scale_offset <- get_scale_and_offset(l2_lfr_url, "measurements/gifapar")
gifapar_scale_offset$scale
[1] 0.003937008
$offset
[1] 0long_scale_offset <- get_scale_and_offset(l2_lfr_url, "measurements/longitude")
long_scale_offset$scale
[1] 0.000001
$offset
[1] 0lat_scale_offset <- get_scale_and_offset(l2_lfr_url, "measurements/latitude")
lat_scale_offset$scale
[1] 0.000001
$offset
[1] 0Again, both longitude and latitude are unevenly spaced 2-dimensional arrays. This tells us that the data grid is curvilinear, and we use st_as_stars() to get our data into the correct format for visualisation, and scale it:
gifapar_stars <- st_as_stars(gifapar = gifapar) |>
st_as_stars(curvilinear = list(
X1 = gifapar_long * long_scale_offset[["scale"]],
X2 = gifapar_lat * lat_scale_offset[["scale"]]
)) |>
mutate(gifapar = gifapar * gifapar_scale_offset[["scale"]])
gifapar_starsstars object with 2 dimensions and 1 attribute
attribute(s), summary of first 100000 cells:
Min. 1st Qu. Median Mean 3rd Qu. Max.
gifapar 1.003937 1.003937 1.003937 1.003937 1.003937 1.003937
dimension(s):
from to refsys point values x/y
X1 1 4090 WGS 84 (CRS84) FALSE [4090x4865] -14.14,...,5.861 [x]
X2 1 4865 WGS 84 (CRS84) FALSE [4090x4865] 39.54,...,52.44 [y]
curvilinear gridFinally, we plot the GIFAPAR:
plot(gifapar_stars, as_points = FALSE, axes = TRUE, breaks = "equal", col = hcl.colors)
đź’Ş Now it is your turn
The following exercises will help you understand how to analyse and visualise different measurements.
Task 1: Visualise Ocean Swell spectra
Following the steps from the Ocean Wind field analysis, extract the Ocean Swell spectra data and visualise it. Check if the data needs to be scaled or offset, and do so if necessary.
Task 2: Visualise a quicklook composite for another item or mission
Following the steps for the quicklook composite visualisation, choose data from another item, or another mission altogether, and visualize its quicklook image. Hint: you can tell if a mission has a quicklook image if it has an entry for “True colour image” under its assets on the EOPF Sentinel Zarr Sample Service STAC Catalog!
Task 2: Visualise another measurement from Sentinel-3 OLCI Level-2 LFR
Review the Sentinel-3 OLCI Level-2 LFR page from the EOPF Sentinel Zarr Sample Service STAC Catalog and choose another measurement of interest, then visualise it.
Conclusion
In this section, we accessed additional Zarr data from the EOPF Sentinel Zarr Sample Service STAC Catalog using rstac and Rarr. We learned how to visualise a quicklook image at different resolutions. We also extracted measurement variables and curvillinear coordinates, formatted and scaled them using the stars package’s st_as_stars() function, and visualised the data.