Skip to content

Commit

Permalink
Merge pull request #1282 from OldLipe/feat/dev-sits
Browse files Browse the repository at this point in the history
Refactoring `sits_merge` function and tests
  • Loading branch information
gilbertocamara authored Feb 11, 2025
2 parents 9cdb1ac + 7126442 commit 03de83d
Show file tree
Hide file tree
Showing 17 changed files with 594 additions and 688 deletions.
5 changes: 4 additions & 1 deletion NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -73,6 +73,8 @@ S3method(.cube_nrows,default)
S3method(.cube_nrows,raster_cube)
S3method(.cube_paths,default)
S3method(.cube_paths,raster_cube)
S3method(.cube_period,default)
S3method(.cube_period,raster_cube)
S3method(.cube_s3class,default)
S3method(.cube_s3class,raster_cube)
S3method(.cube_source,default)
Expand Down Expand Up @@ -189,7 +191,6 @@ S3method(.source_item_get_date,cdse_cube)
S3method(.source_item_get_date,deafrica_cube)
S3method(.source_item_get_date,stac_cube)
S3method(.source_item_get_hrefs,bdc_cube)
S3method(.source_item_get_hrefs,sdc_cube)
S3method(.source_item_get_hrefs,stac_cube)
S3method(.source_item_get_hrefs,usgs_cube)
S3method(.source_items_bands_select,cdse_cube)
Expand Down Expand Up @@ -302,6 +303,8 @@ S3method(.tile_path,derived_cube)
S3method(.tile_path,raster_cube)
S3method(.tile_paths,default)
S3method(.tile_paths,raster_cube)
S3method(.tile_period,default)
S3method(.tile_period,raster_cube)
S3method(.tile_read_block,default)
S3method(.tile_read_block,derived_cube)
S3method(.tile_read_block,eo_cube)
Expand Down
7 changes: 7 additions & 0 deletions R/api_check.R
Original file line number Diff line number Diff line change
Expand Up @@ -2611,3 +2611,10 @@
)
return(invisible(NULL))
}

.check_unique_period <- function(cube) {
.check_that(
x = length(.cube_period(cube)) == 1,
msg = .conf("messages", ".check_unique_period")
)
}
9 changes: 9 additions & 0 deletions R/api_conf.R
Original file line number Diff line number Diff line change
Expand Up @@ -1287,3 +1287,12 @@ NULL
rm(leaf_map)
return(invisible(NULL))
}
#' @title Get Grid System
#' @name .conf_grid_system
#' @keywords internal
#' @noRd
#' @return Grid system name.
#'
.conf_grid_system <- function(source, collection) {
.conf("sources", source, "collections", collection, "grid_system")
}
46 changes: 44 additions & 2 deletions R/api_cube.R
Original file line number Diff line number Diff line change
Expand Up @@ -523,6 +523,26 @@ NULL
crs <- .cube_crs(cube)
return(crs)
}
#' @title Return period of a data cube
#' @keywords internal
#' @noRd
#' @name .cube_period
#' @param cube data cube
#' @return period in days associated to the cube
.cube_period <- function(cube) {
UseMethod(".cube_period", cube)
}
#' @export
.cube_period.raster_cube <- function(cube) {
.dissolve(slider::slide(cube, .tile_period))
}
#' @export
.cube_period.default <- function(cube) {
cube <- tibble::as_tibble(cube)
cube <- .cube_find_class(cube)
period <- .cube_period(cube)
return(period)
}
#' @title Adjust crs of a data cube
#' @keywords internal
#' @noRd
Expand Down Expand Up @@ -806,6 +826,16 @@ NULL
return(is_regular)
}

#' @title Check that cube has unique period
#' @name .cube_has_unique_period
#' @keywords internal
#' @noRd
#' @param cube datacube
#' @return Called for side effects.
.cube_has_unique_period <- function(cube) {
length(.cube_period(cube)) == 1
}

#' @title Check that cube is a base cube
#' @name .cube_is_base
#' @keywords internal
Expand Down Expand Up @@ -1451,8 +1481,6 @@ NULL
path <- stringr::str_replace(path, path_prefix, "")

url_parsed <- .url_parse(path)
url_parsed[["path"]] <- paste0(path_prefix, url_parsed[["path"]])

url_parsed[["query"]] <- utils::modifyList(
url_parsed[["query"]], token_parsed
)
Expand Down Expand Up @@ -1597,3 +1625,17 @@ NULL
.cube_satellite <- function(cube) {
.dissolve(slider::slide(cube, .tile_satellite))
}

#' @title Return cube grid system
#' @name .cube_grid_system
#' @keywords internal
#' @noRd
#'
#' @param cube Raster cube
#' @return Cube grid system
.cube_grid_system <- function(cube) {
.conf_grid_system(
source = .cube_source(cube),
collection = .cube_collection(cube)
)
}
166 changes: 158 additions & 8 deletions R/api_merge.R
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
}

# ---- Adjust timeline strategies strategies ----
.merge_adjust_timeline_strategy_zipper <- function(t1, t2) {
.merge_zipper_strategy <- function(t1, t2) {
# define vector to store overlapping dates
t_overlap <- c()
# define the size of the `for` - size of the reference time-series
Expand Down Expand Up @@ -107,7 +107,7 @@
if (.has(common_tiles)) {
merge_strategy <- .merge_strategy_file
} else {
# case 2: different tiles, merge cube rows
# case 2: different tiles, merge cube rows
merge_strategy <- .merge_strategy_bind
}
# merge
Expand All @@ -125,6 +125,12 @@
# extract tiles
tiles <- .merge_get_common_tiles(data1, data2)
if (!.has(tiles)) {
# It is not possible to merge non-common tiles with multiple bands using
# the same sensor
.check_that(
.cube_sensor(data1) != .cube_sensor(data2),
msg = .conf("messages", ".merge_irregular_bands")
)
# if no common tiles are available, use a global reference timeline.
# in this case, this timeline is generated by the merge of all timelines
# in the reference cube (cube 1)
Expand All @@ -134,7 +140,7 @@
# get row timeline
row_timeline <- .tile_timeline(row)
# search overlaps between the reference timeline and row timeline
t_overlap <- .merge_adjust_timeline_strategy_zipper(
t_overlap <- .merge_zipper_strategy(
t1 = reference_timeline,
t2 = row_timeline
)
Expand All @@ -156,7 +162,7 @@
ts1 <- .tile_timeline(tile1)
ts2 <- .tile_timeline(tile2)
# adjust timeline using zipper strategy
ts_overlap <- .merge_adjust_timeline_strategy_zipper(ts1, ts2)
ts_overlap <- .merge_zipper_strategy(ts1, ts2)
# filter cubes in the overlapping dates
tile1 <- .cube_filter_dates(tile1, ts_overlap)
tile2 <- .cube_filter_dates(tile2, ts_overlap)
Expand All @@ -168,8 +174,72 @@
merged_cube
}

# ---- Merge operation: Special case - DEM Cube ----
.merge_dem_cube <- function(data1, data2) {
.merge_strategy_intersects <- function(data1, data2) {
# Get data cubes timeline
t1 <- .cube_timeline(data1)[[1]]
t2 <- .cube_timeline(data2)[[1]]

# Get cubes period
t2_period <- t2[2] - t2[1]
t1_period <- t1[2] - t1[1]

# Lists to store dates
t1_date <- list()
t2_date <- list()

# Get overlapped dates
for (i in seq_len(length(t2))) {
t2_int <- lubridate::interval(
lubridate::ymd(t2[i]), lubridate::ymd(t2[i]) + t2_period - 1
)
overlapped_dates <- lapply(seq_len(length(t1)), function(j) {
t1_int <- lubridate::interval(
lubridate::ymd(t1[j]), lubridate::ymd(t1[j]) + t1_period - 1
)
lubridate::int_overlaps(t2_int, t1_int)
})

dates <- t1[unlist(overlapped_dates)]
dates <- setdiff(dates, t1_date)
if (.has(dates)) {
t1_date[[i]] <- as.Date(min(dates))
t2_date[[i]] <- as.Date(t2[i])
}
}

# Transform list to vector date
t1_date <- as.Date(unlist(t1_date))
t2_date <- as.Date(unlist(t2_date))

# Filter overlapped dates
data1 <- .cube_filter_dates(data1, t1_date)
data2 <- .cube_filter_dates(data2, t2_date)

# Change file date to match reference timeline
data2 <- slider::slide_dfr(data2, function(y) {
fi_list <- purrr::map(.tile_bands(y), function(band) {
fi_band <- .fi_filter_bands(.fi(y), bands = band)
fi_band[["date"]] <- t1_date
return(fi_band)
})
tile_fi <- dplyr::bind_rows(fi_list)
tile_fi <- dplyr::arrange(
tile_fi,
.data[["date"]],
.data[["band"]],
.data[["fid"]]
)
y[["file_info"]] <- list(tile_fi)
y
})

# Merge the cubes
data1 <- .merge_strategy_file(data1, data2)
return(data1)
}

# ---- Merge operation: DEM case ----
.merge_dem <- function(data1, data2) {
# define cubes
dem_cube <- data1
other_cube <- data2
Expand Down Expand Up @@ -200,8 +270,8 @@
.merge_strategy_file(other_cube, dem_cube)
}

# ---- Merge operation: Special case - HLS Cube ----
.merge_hls_cube <- function(data1, data2) {
# ---- Merge operation: HLS case ----
.merge_hls <- function(data1, data2) {
if ((.cube_collection(data1) == "HLSS30" ||
.cube_collection(data2) == "HLSS30")) {
data1[["collection"]] <- "HLSS30"
Expand All @@ -210,3 +280,83 @@
# merge cubes and return
.merge_strategy_file(data1, data2)
}


# ---- Merge operation: Regular case ----
.merge_regular <- function(data1, data2) {
# Rule 1: Do the cubes have same tiles?
.check_cube_tiles(data1, .cube_tiles(data2))
.check_cube_tiles(data2, .cube_tiles(data1))

# Rule 2: Do the cubes have same bands?
bands_to_merge <- setdiff(.cube_bands(data2), .cube_bands(data1))
.check_that(
length(bands_to_merge) > 0,
msg = .conf("messages", ".merge_regular_bands")
)

# Filter bands to merge
data2 <- .cube_filter_bands(data2, bands_to_merge)

# Rule 3: Do the cubes have same timeline?
if (all(.cube_timeline(data1) %in% .cube_timeline(data2)) &&
all(.cube_timeline(data2) %in% .cube_timeline(data1))) {
merged_cube <- .merge_strategy_file(data1, data2)
} else {
merged_cube <- .merge_strategy_intersects(data1, data2)
}
# Return merged cube
return(merged_cube)
}

.merge_irregular <- function(data1, data2) {
# verify if cube has the same bands
has_same_bands <- .merge_has_equal_bands(data1, data2)
# rule 1: if the bands are the same, combine cubes (`densify`)
if (has_same_bands) {
# merge!
merged_cube <- .merge_cube_densify(data1, data2)
} else {
# rule 2: if the bands are different and their timelines are
# compatible, the bands are joined. The resulting timeline is the one
# from the first cube.
merged_cube <- .merge_cube_compactify(data1, data2)
}
}

.merge_switch <- function(data1, data2, ...) {
switch(.merge_type(data1, data2),
...
)
}

.merge_type <- function(data1, data2) {
# Special cases
if (any(inherits(data1, "dem_cube"), inherits(data2, "dem_cube"))) {
"dem_case"
} else if (all(inherits(data1, "hls_cube"), inherits(data2, "hls_cube"))) {
"hls_case"
} else if (
all(
inherits(data1, "deaustralia_cube_ga_s2am_ard_3"),
inherits(data2, "deaustralia_cube_ga_s2am_ard_3")
) &&
all(
inherits(data1, "deaustralia_cube_ga_s2bm_ard_3"),
inherits(data2, "deaustralia_cube_ga_s2bm_ard_3")
)
) {
"irregular_case"
# General cases
} else if (.cube_is_regular(data1) &&
.cube_is_regular(data2) &&
.cube_has_unique_period(data1) &&
.cube_has_unique_period(data2)) {
"regular_case"
} else if (!.cube_is_regular(data1) || !.cube_is_regular(data2) ||
!.cube_has_unique_period(data1) || !.cube_has_unique_period(data2)) {
"irregular_case"
} else {
stop(.conf("messages", ".merge_type"), class(data1))
}
}
4 changes: 4 additions & 0 deletions R/api_regularize.R
Original file line number Diff line number Diff line change
Expand Up @@ -181,6 +181,10 @@
#' @noRd
#' @export
.reg_tile_convert.raster_cube <- function(cube, grid_system, roi = NULL, tiles = NULL) {
# for consistency, check if the grid is already in place
if (grid_system == .cube_grid_system(cube)) {
return(cube)
}
# if roi and tiles are not provided, use the whole cube as extent
if (!.has(roi) && !.has(tiles)) {
roi <- .cube_as_sf(cube)
Expand Down
39 changes: 0 additions & 39 deletions R/api_source_sdc.R
Original file line number Diff line number Diff line change
Expand Up @@ -70,45 +70,6 @@
)
)
}
#' @title Retrieves the paths or URLs of each file bands of an item for SDC
#' @param source Name of the STAC provider.
#' @param item \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 Returns paths to STAC item.
#' @keywords internal
#' @noRd
#' @export
.source_item_get_hrefs.sdc_cube <- function(source,
item, ...,
collection = NULL) {
hrefs <- unname(purrr::map_chr(item[["assets"]], `[[`, "href"))
asset_names <- unlist(
purrr::map(item[["assets"]], `[[`, "eo:bands"),
use.names = FALSE
)

# post-conditions
.check_chr(hrefs, allow_empty = FALSE)

# fix local images - temporary solution
is_local_images <- grepl(pattern = "^file://", x = hrefs)
if (any(is_local_images)) {
server_path <- "https://explorer.swissdatacube.org"

hrefs[is_local_images] <- gsub(
pattern = "^file://",
replacement = server_path,
x = hrefs[is_local_images]
)
}

# add gdal VSI in href urls
vsi_hrefs <- .stac_add_gdal_fs(hrefs)
vsi_hrefs <- sprintf('%s:"%s":%s', "NETCDF", vsi_hrefs, asset_names)

return(vsi_hrefs)
}
#' @title Check if roi or tiles are provided
#' @param source Data source
#' @param roi Region of interest
Expand Down
Loading

0 comments on commit 03de83d

Please sign in to comment.