Merge pull request #1144 from M3nin0/feature/dem-regularization

Support for DEM from MPC

NAMESPACE

@@ -126,6 +126,7 @@ S3method(.raster_xres,terra)
 S3method(.raster_ymax,terra)
 S3method(.raster_ymin,terra)
 S3method(.raster_yres,terra)
+S3method(.reg_s2tile_convert,dem_cube)
 S3method(.reg_s2tile_convert,grd_cube)
 S3method(.reg_s2tile_convert,rtc_cube)
 S3method(.slice_dfr,numeric)
@@ -137,6 +138,7 @@ S3method(.source_collection_access_test,usgs_cube)
 S3method(.source_cube,stac_cube)
 S3method(.source_filter_tiles,"cdse_cube_sentinel-1-rtc")
 S3method(.source_filter_tiles,"deafrica_cube_sentinel-1-rtc")
+S3method(.source_filter_tiles,"mpc_cube_cop-dem-glo-30")
 S3method(.source_filter_tiles,"mpc_cube_sentinel-1-grd")
 S3method(.source_filter_tiles,stac_cube)
 S3method(.source_item_get_bands,stac_cube)
@@ -157,6 +159,7 @@ S3method(.source_items_fid,stac_cube)
 S3method(.source_items_new,"aws_cube_landsat-c2-l2")
 S3method(.source_items_new,"deafrica_cube_sentinel-1-rtc")
 S3method(.source_items_new,"deafrica_cube_sentinel-2-l2a")
+S3method(.source_items_new,"mpc_cube_cop-dem-glo-30")
 S3method(.source_items_new,"mpc_cube_landsat-c2-l2")
 S3method(.source_items_new,"mpc_cube_sentinel-1-grd")
 S3method(.source_items_new,"mpc_cube_sentinel-1-rtc")
@@ -172,6 +175,7 @@ S3method(.source_items_tile,"aws_cube_landsat-c2-l2")
 S3method(.source_items_tile,"cdse_cube_sentinel-1-rtc")
 S3method(.source_items_tile,"deafrica_cube_rainfall-chirps-daily")
 S3method(.source_items_tile,"deafrica_cube_rainfall-chirps-monthly")
+S3method(.source_items_tile,"mpc_cube_cop-dem-glo-30")
 S3method(.source_items_tile,"mpc_cube_landsat-c2-l2")
 S3method(.source_items_tile,"mpc_cube_sentinel-1-grd")
 S3method(.source_items_tile,"mpc_cube_sentinel-1-rtc")
@@ -187,6 +191,7 @@ S3method(.source_roi_tiles,"mpc_cube_landsat-c2-l2")
 S3method(.source_roi_tiles,sdc_cube)
 S3method(.source_roi_tiles,stac_cube)
 S3method(.source_tile_get_bbox,"cdse_cube_sentinel-1-rtc")
+S3method(.source_tile_get_bbox,"mpc_cube_cop-dem-glo-30")
 S3method(.source_tile_get_bbox,"mpc_cube_sentinel-1-grd")
 S3method(.source_tile_get_bbox,"mpc_cube_sentinel-1-rtc")
 S3method(.source_tile_get_bbox,stac_cube)
@@ -373,6 +378,7 @@ S3method(sits_reclassify,default)
 S3method(sits_reduce,raster_cube)
 S3method(sits_reduce,sits)
 S3method(sits_regularize,default)
+S3method(sits_regularize,dem_cube)
 S3method(sits_regularize,derived_cube)
 S3method(sits_regularize,raster_cube)
 S3method(sits_regularize,sar_cube)

R/api_cube.R

@@ -350,12 +350,25 @@ NULL
  },
  .default = FALSE
  )
+
+ dem_cube <- .try({
+ .conf("sources", source, "collections", collection, "dem_cube")
+ },
+ .default = FALSE
+ )
+
  if (sar_cube) {
  if (grepl("rtc", col_class, fixed = TRUE))
  unique(c(col_class, "rtc_cube", "sar_cube", s3_class, class(cube)))
  else
  unique(c(col_class, "grd_cube", "sar_cube", s3_class, class(cube)))
- } else {
+ }
+
+ else if (dem_cube) {
+ unique(c(col_class, "dem_cube", s3_class, class(cube)))
+ }
+
+ else {
  unique(c(col_class, s3_class, class(cube)))
  }
 }

R/api_plot_raster.R

@@ -310,7 +310,7 @@
  # verifies if stars package is installed
  .check_require_packages("stars")
  # verifies if tmap package is installed
- .check_require_packages("plot")
+ .check_require_packages("tmap")
 
  # deal with color palette
  .check_palette(palette)
@@ -394,7 +394,7 @@
  # verifies if stars package is installed
  .check_require_packages("stars")
  # verifies if tmap package is installed
- .check_require_packages("plot")
+ .check_require_packages("tmap")
  # precondition - check color palette
  .check_palette(palette)
  # revert the palette

R/api_regularize.R

@@ -236,3 +236,51 @@
  cube_class <- c(cube_class[[1]], "sar_cube", cube_class[-1])
  .cube_set_class(cube, cube_class)
 }
+#' @noRd
+#' @export
+#'
+.reg_s2tile_convert.dem_cube<- function(cube, roi = NULL, tiles = NULL) {
+ # generate Sentinel-2 tiles and intersects it with doi
+ tiles_mgrs <- .s2tile_open(roi, tiles)
+
+ # create a new cube according to Sentinel-2 MGRS
+ cube_class <- .cube_s3class(cube)
+
+ cube <- tiles_mgrs |>
+ dplyr::rowwise() |>
+ dplyr::group_map(~{
+ # prepare a sf object representing the bbox of each image in
+ # file_info
+ cube_crs <- dplyr::filter(cube, .data[["crs"]] == .x[["crs"]])
+ if (nrow(cube_crs) == 0) {
+ cube_crs <- cube
+ }
+ fi_bbox <- .bbox_as_sf(.bbox(
+ x = .fi(cube_crs),
+ default_crs = .crs(cube_crs),
+ by_feature = TRUE
+ ))
+ file_info <- .fi(cube_crs)[.intersects({{fi_bbox}}, .x), ]
+ .cube_create(
+ source = .tile_source(cube_crs),
+ collection = .tile_collection(cube_crs),
+ satellite = .tile_satellite(cube_crs),
+ sensor = .tile_sensor(cube_crs),
+ tile = .x[["tile_id"]],
+ xmin = .xmin(.x),
+ xmax = .xmax(.x),
+ ymin = .ymin(.x),
+ ymax = .ymax(.x),
+ crs = paste0("EPSG:", .x[["epsg"]]),
+ file_info = file_info
+ )
+ }) |>
+ dplyr::bind_rows()
+
+ # Filter non-empty file info
+ cube <- .cube_filter_nonempty(cube)
+
+ # Finalize customizing cube class
+ cube_class <- c(cube_class[[1]], "dem_cube", cube_class[-1])
+ .cube_set_class(cube, cube_class)
+}

R/api_source.R

@@ -65,9 +65,18 @@ NULL
  },
  .default = FALSE
  )
+ # is this a collection of DEM data ?
+ dem_cube <- .try({
+ .conf("sources", source, "collections", collection, "dem_cube")
+ },
+ .default = FALSE
+ )
  # if this is a SAR collection, add "sar_cube" to the class
  if (sar_cube)
  class(source) <- c("sar_cube", class(source))
+ # if this is a DEM collection, add "dem_cube" to the class
+ if (dem_cube)
+ class(source) <- c("dem_cube", class(source))
  # add a class combining source and collection
  class_source_col <- paste(classes[[1]], tolower(collection), sep = "_")
  class(source) <- unique(c(class_source_col, class(source)))

R/api_source_mpc.R

@@ -233,6 +233,35 @@
  collection = collection
  )
 }
+#' @title Get bbox from file info for COP-DEM-GLO-30
+#' @keywords internal
+#' @noRd
+#' @param source Data source
+#' @param file_info File info
+#' @param ... Additional parameters.
+#' @param collection Image collection
+#' @return vector (xmin, ymin, xmax, ymax).
+#' @export
+`.source_tile_get_bbox.mpc_cube_cop-dem-glo-30` <- function(source,
+ file_info, ...,
+ collection = NULL) {
+ .check_set_caller(".source_tile_get_bbox_mpc_dem_30")
+
+ # pre-condition
+ .check_num(nrow(file_info), min = 1)
+
+ # get bbox based on file_info
+ xmin <- min(file_info[["xmin"]])
+ ymin <- min(file_info[["ymin"]])
+ xmax <- max(file_info[["xmax"]])
+ ymax <- max(file_info[["ymax"]])
+
+ # post-condition
+ .check_that(xmin < xmax && ymin < ymax)
+ # create a bbox
+ bbox <- c(xmin = xmin, ymin = ymin, xmax = xmax, ymax = ymax)
+ return(bbox)
+}
 #' @keywords internal
 #' @noRd
 #' @export
@@ -429,6 +458,53 @@
 #' @keywords internal
 #' @noRd
 #' @export
+`.source_items_new.mpc_cube_cop-dem-glo-30` <- function(source,
+ collection,
+ stac_query, ...,
+ tiles = NULL) {
+ .check_set_caller(".source_items_new_mpc_cube_cop-dem-glo-30")
+
+ # COP-DEM-GLO-30 does not support tiles - convert to ROI
+ if (!is.null(tiles)) {
+ roi <- .s2_mgrs_to_roi(tiles)
+ stac_query[["params"]][["intersects"]] <- NULL
+ stac_query[["params"]][["bbox"]] <- c(roi[["lon_min"]],
+ roi[["lat_min"]],
+ roi[["lon_max"]],
+ roi[["lat_max"]]
+ )
+ }
+
+ # Fix temporal interval (All data available in the same date)
+ stac_query[["params"]][["datetime"]] <- "2021-04-21/2021-04-23"
+
+ # Search content
+ items_info <- rstac::post_request(q = stac_query, ...)
+ .check_stac_items(items_info)
+ # fetching all the metadata
+ items_info <- suppressWarnings(
+ rstac::items_fetch(items = items_info, progress = FALSE)
+ )
+
+ # assign href
+ access_key <- Sys.getenv("MPC_TOKEN")
+ if (!nzchar(access_key)) {
+ access_key <- NULL
+ }
+ # Clean old tokens cached in rstac
+ .mpc_clean_token_cache()
+ items_info <- suppressWarnings(
+ rstac::items_sign(
+ items_info, sign_fn = rstac::sign_planetary_computer(
+ httr::add_headers("Ocp-Apim-Subscription-Key" = access_key)
+ )
+ )
+ )
+ return(items_info)
+}
+#' @keywords internal
+#' @noRd
+#' @export
 `.source_items_tile.mpc_cube_sentinel-1-grd` <- function(source,
  items, ...,
  collection = NULL) {
@@ -485,6 +561,20 @@
  })
  rstac::items_reap(items, field = c("properties", "tile"))
 }
+#' @title Organizes items for MPC COP-DEM-GLO-30 collections
+#' @param source Name of the STAC provider.
+#' @param items \code{STACItemcollection} object from rstac package.
+#' @param ... Other parameters to be passed for specific types.
+#' @param collection Collection to be searched in the data source.
+#' @return A list of items.
+#' @keywords internal
+#' @noRd
+#' @export
+`.source_items_tile.mpc_cube_cop-dem-glo-30` <- function(source,
+ items, ...,
+ collection = NULL) {
+ rep("NoTilingSystem", rstac::items_length(items))
+}
 #' @title Filter S1 GRD tiles
 #' @noRd
 #' @param source Data source
@@ -509,6 +599,19 @@
  tiles = tiles)
 
 }
+#' @title Filter COP-DEM-GLO-30 tiles
+#' @noRd
+#' @param source Data source
+#' @param cube Cube to be filtered
+#' @param tiles Tiles to be selected
+#' @return Filtered cube
+#' @export
+`.source_filter_tiles.mpc_cube_cop-dem-glo-30` <- function(source,
+ collection,
+ cube,
+ tiles) {
+ return(cube)
+}
 #' @title Check if roi or tiles are provided
 #' @param source Data source
 #' @param roi Region of interest

R/sits_merge.R

@@ -207,16 +207,13 @@ sits_merge.raster_cube <- function(data1, data2, ...) {
  # Get data1 timeline.
  d1_tl <- unique(as.Date(.cube_timeline(data1)[[1]]))
  # Create new `file_info` using dates from `data1` timeline.
- fi_new <- purrr::map(seq_len(nrow(data2)), function(row) {
- data_row <- data2[row,]
-
- fi <- .fi(data_row)
- fi[["date"]] <- as.Date(d1_tl[1:nrow(data_row)])
-
- return(fi)
+ fi_new <- purrr::map(sits_timeline(data1), function(date_row) {
+ fi <- .fi(data2)
+ fi[["date"]] <- as.Date(date_row)
+ fi
  })
  # Assign the new `file_into` into `data2`
- data2[["file_info"]] <- fi_new
+ data2[["file_info"]] <- list(dplyr::bind_rows(fi_new))
  # Merge cubes and return
  .cube_merge(data1, data2)
 }

R/sits_regularize.R

@@ -205,6 +205,52 @@ sits_regularize.sar_cube <- function(cube, ...,
  )
  return(cube)
 }
+
+#' @rdname sits_regularize
+#' @export
+sits_regularize.dem_cube <- function(cube, ...,
+ res,
+ output_dir,
+ roi = NULL,
+ tiles = NULL,
+ multicores = 2L,
+ progress = TRUE) {
+ # Preconditions
+ .check_raster_cube_files(cube)
+ .check_num_parameter(res, exclusive_min = 0)
+ output_dir <- .file_normalize(output_dir)
+ .check_output_dir(output_dir)
+ .check_num_parameter(multicores, min = 1, max = 2048)
+ .check_progress(progress)
+ # Check for ROI and tiles
+ .check_roi_tiles(roi, tiles)
+ # Display warning message in case STAC cube
+ # Prepare parallel processing
+ .parallel_start(workers = multicores)
+ on.exit(.parallel_stop(), add = TRUE)
+ # Convert input sentinel1 cube to sentinel2 grid
+ cube <- .reg_s2tile_convert(cube = cube, roi = roi, tiles = tiles)
+ .check_that(nrow(cube) > 0,
+ msg = .conf("messages", "sits_regularize_roi")
+ )
+ # Filter tiles
+ if (is.character(tiles)) {
+ cube <- .cube_filter_tiles(cube, tiles)
+ }
+ # DEMs don't have the temporal dimension, so the period is fixed in 1 day.
+ period <- "P1D"
+ # Call regularize in parallel
+ cube <- .reg_cube(
+ cube = cube,
+ res = res,
+ roi = roi,
+ period = period,
+ output_dir = output_dir,
+ progress = progress
+ )
+ return(cube)
+}
+
 #' @rdname sits_regularize
 #' @export
 sits_regularize.derived_cube <- function(cube, ...) {

inst/extdata/config_messages.yml

@@ -278,6 +278,7 @@
 .source_url: "invalid URL for requested source provider"
 .source_tile_get_bbox: "unable to retrieve images given a bounding box"
 .source_tile_get_bbox_cdse_s1_rtc: "unable to retrieve file information for S1 RTC cubes"
+.source_tile_get_bbox_mpc_dem_30: "unable to retrieve file information for COP DEM GLO 30m"
 .source_tile_get_bbox_mpc_s1_grd: "unable to retrieve file information for S1 GRD cubes"
 .stac_format_platform: "platform name should be unique (e.g. Landsat-8 or Sentinel-1A)"
 .stac_select_bands: "some bands for this product are not pre-configured in sits\n please include them in you user configuration file."