From 0b9ec69ea37bc33f2c5d5b719f1c1854cbd605f4 Mon Sep 17 00:00:00 2001 From: Felipe Date: Mon, 7 Oct 2024 14:05:20 +0000 Subject: [PATCH 01/19] update roi API --- R/api_mosaic.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/R/api_mosaic.R b/R/api_mosaic.R index 361acdab1..856ec4bee 100644 --- a/R/api_mosaic.R +++ b/R/api_mosaic.R @@ -225,7 +225,7 @@ quiet = TRUE ) # Delete temporary roi file - on.exit(.mosaic_del_roi(roi)) + on.exit(.roi_delete(roi)) } # Crop and reproject tile image out_file <- .gdal_crop_image( @@ -252,7 +252,7 @@ #' @noRd #' @param roi Region of interest #' @return Called for side effects -.mosaic_del_roi <- function(roi) { +.roi_delete <- function(roi) { if (is.null(roi)) { return(roi) } From 2ad6f136397161b6985540d7cfdbb0e41309fedd Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Fri, 22 Nov 2024 07:05:59 -0300 Subject: [PATCH 02/19] updates in sits_merge --- R/api_merge.R | 352 ++++++++++++++++++++++++++++++------ R/api_timeline.R | 24 +++ R/sits_merge.R | 104 +++++------ man/sits_merge.Rd | 6 +- tests/testthat/test-merge.R | 102 +++++++++++ 5 files changed, 480 insertions(+), 108 deletions(-) create mode 100644 tests/testthat/test-merge.R diff --git a/R/api_merge.R b/R/api_merge.R index 85ff66922..ff902999b 100644 --- a/R/api_merge.R +++ b/R/api_merge.R @@ -1,50 +1,169 @@ -.merge_diff_timelines <- function(t1, t2) { - abs(as.Date(t1) - as.Date(t2)) +# ---- General utilities ---- +.merge_bands_intersects <- function(data1, data2) { + # Extract bands + d1_bands <- .cube_bands(data1) + d2_bands <- .cube_bands(data2) + # Extract overlaps + intersect(d1_bands, d2_bands) } -.cube_merge <- function(data1, data2) { - data1 <- slider::slide2_dfr(data1, data2, function(x, y) { - .fi(x) <- dplyr::arrange( - dplyr::bind_rows(.fi(x), .fi(y)), +.merge_tiles_overlaps <- function(data1, data2) { + # Extract common tiles + d1_tiles <- .cube_tiles(data1) + d2_tiles <- .cube_tiles(data2) + # Extract overlaps + intersect(d1_tiles, d2_tiles) +} + +.merge_adjust_timeline <- function(data1, data2) { + # reference timeline + reference_tl <- .cube_timeline(data1) + # Adjust dates / bands + 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"]] <- reference_tl + return(fi_band) + }) + tile_fi <- dplyr::bind_rows(fi_list) + tile_fi <- dplyr::arrange( + tile_fi, .data[["date"]], .data[["band"]], .data[["fid"]] ) - # remove duplicates - .fi(x) <- dplyr::distinct( - .fi(x), - .data[["band"]], - .data[["date"]], - .keep_all = TRUE - ) - - return(x) + y[["file_info"]] <- list(tile_fi) + y }) - return(data1) } -.merge_irregular_cube <- function(data1, data2) { - merged_cube <- dplyr::bind_rows(data1, data2) - class(merged_cube) <- c("combined_cube", class(data1)) - return(merged_cube) +.merge_get_ts_within <- function(data1, data2) { + # extract timelines + d1_tl <- .cube_timeline(data1) + d2_tl <- .cube_timeline(data2) + + # Check if all dates in A are in B + if (all(datesA %in% datesB)) { + return(datesA) # A is contained in B + } + + # Check if all dates in B are in A + if (all(datesB %in% datesA)) { + return(datesB) # B is contained in A + } + + return(NULL) # Neither is contained } -.merge_equal_cube <- function(data1, data2) { - if (inherits(data1, "hls_cube") && inherits(data2, "hls_cube") && - (.cube_collection(data1) == "HLSS30" || - .cube_collection(data2) == "HLSS30")) { - data1[["collection"]] <- "HLSS30" +.merge_check_bands_intersects <- function(data1, data2) { + # Extract band intersects + bands_intersects <- .merge_bands_intersects(data1, data2) + # Check if there are intersects + .check_that(length(bands_intersects) >= 1) +} + +.merge_check_band_sensor <- function(data1, data2) { + # Extract band intersects + bands_intersects <- .merge_bands_intersects(data1, data2) + # If has overlaps, the sensor must be the same + if (length(bands_intersects) >= 1) { + .check_that(data1[["sensor"]] == data2[["sensor"]]) } +} - data1 <- .cube_merge(data1, data2) - return(data1) +# ---- Merge strategies ---- +.merge_strategy_file <- function(data1, data2, adjust_timeline) { + # adjust second cube timeline, based on the first cube + if (adjust_timeline) { + data2 <- .merge_adjust_timeline(data1, data2) + } + # extract tiles + tiles <- .merge_tiles_overlaps(data1, data2) + # merge cubes + .map_dfr(tiles, function(tile) { + # select data in the selected tile + data1_in_tile <- .select_raster_tiles(data1, tile) + data2_in_tile <- .select_raster_tiles(data2, tile) + # change file name to match reference timeline + slider::slide2_dfr(data1_in_tile, data2_in_tile, function(x, y) { + # arrange by `date`, `band` and `fid` + .fi(x) <- dplyr::arrange( + dplyr::bind_rows(.fi(x), .fi(y)), + .data[["date"]], + .data[["band"]], + .data[["fid"]] + ) + # remove duplicates + .fi(x) <- dplyr::distinct( + .fi(x), + .data[["band"]], + .data[["date"]], + .keep_all = TRUE + ) + # return + return(x) + }) + }) +} + +.merge_strategy_bind <- function(data1, data2, adjust_timeline) { + # Adjust second cube timeline, based on the first cube + if (adjust_timeline) { + data2 <- .merge_adjust_timeline(data1, data2) + } + # Extract band intersects + bands_intersects <- .merge_bands_intersects(data1, data2) + # Use only intersect bands + data1 <- .select_raster_bands(data1, bands_intersects) + data2 <- .select_raster_bands(data2, bands_intersects) + # Merge + dplyr::bind_rows(data1, data2) } -.merge_distinct_cube <- function(data1, data2) { - # Get cubes timeline +# ---- Regular cubes ---- +# .merge_regular_cut_timeline <- function(data1, data2) { +# # extract timelines +# d1_tl <- .cube_timeline(data1) +# d2_tl <- .cube_timeline(data2) +# # extract tiles +# tiles <- .merge_tiles_overlaps(data1, data2) +# # merge cubes +# .map_dfr(tiles, function(tile) { +# # select data in the selected tile +# data1_in_tile <- .select_raster_tiles(data1, tile) +# data2_in_tile <- .select_raster_tiles(data2, tile) +# # extract timelines +# d1_tl <- .tile_timeline(data1_in_tile) +# d2_tl <- .tile_timeline(data2_in_tile) +# # get min/max dates +# min_tl <- min( +# min(d1_tl), min(d2_tl) +# ) +# max_tl <- max( +# max(d1_tl), max(d2_tl) +# ) +# # cut timeline +# .tile_filter_interval(tile, min_tl, max_tl) +# }) +# } + +.merge_regular_check_timeline_overlaps <- function(data1, data2) { + # extract timelines + d1_tl <- .cube_timeline(data1) + d2_tl <- .cube_timeline(data2) + # check overlaps + slider::slide2_vec(d1_tl, d2_tl, function(x, y) { + x <- .dissolve(x) + y <- .dissolve(y) + + .check_that(length(.timeline_has_overlap(x, y)) >= 1) + }) +} + +.merge_regular_check_periods <- function(data1, data2) { + # get cubes timeline d1_tl <- unique(as.Date(.cube_timeline(data1)[[1]])) d2_tl <- unique(as.Date(.cube_timeline(data2)[[1]])) - # get intervals d1_period <- as.integer( lubridate::as.period(lubridate::int_diff(d1_tl)), "days" @@ -60,34 +179,105 @@ .check_that( unique(d1_period) == unique(d2_period) ) - # pre-condition - are the cubes start date less than period timeline? - .check_that( - abs(d1_period[[1]] - d2_period[[2]]) <= unique(d2_period) - ) +} + +.merge_regular_cube <- function(data1, data2) { + # pre-condition - timelines overlaps + # in case of regular cube it is assumed the timeline must overlap + # to avoid the creation of inconsistent / irregular cubes + .merge_regular_check_timeline_overlaps(data1, data2) + # pre-condition - timelines with same period + .merge_regular_check_periods(data1, data2) + # pre-condition - equal bands must be from the same sensor + # bands with the same name, must be from the same sensor to avoid confusion + .merge_check_band_sensor(data1, data2) + # ToDo: Cut timeline at overlapping intervals when length(ts1) != length(ts2) + # get tile overlaps + tiles_overlaps <- .merge_tiles_overlaps(data1, data2) + # define the strategy (default - merge tiles) + merge_strategy <- NULL + # case: same tiles, merge file info + tiles_overlaps <- .merge_tiles_overlaps(data1, data2) + if (.has(tiles_overlaps)) { + merge_strategy <- .merge_strategy_file + # case 2: different tiles, merge tile rows + } else { + merge_strategy <- .merge_strategy_bind + } + # merge + merge_strategy(data1, data2, TRUE) +} + +# ---- Irregular cubes ---- +.merge_irregular_cube <- function(data1, data2) { + # pre-condition - equal bands from the same sensor + # bands with the same name, must be from the same sensor to avoid confusion + .merge_check_band_sensor(data1, data2) + # get tile overlaps + tiles_overlaps <- .merge_tiles_overlaps(data1, data2) + # define the strategy (default - merge tiles) + merge_strategy <- NULL + # case: same tiles, merge file info + if (.has(tiles_overlaps)) { + merged_cube <- .merge_strategy_file(data1, data2, FALSE) + # case 2: different tiles, merge tile rows + } else { + merged_cube <- .merge_strategy_bind(data1, data2, FALSE) + class(merged_cube) <- c("combined_cube", class(data1)) + } + # return + return(merged_cube) +} + +# .merge_irregular_cube <- function(data1, data2) { +# # pre-condition - intersecting bands +# .merge_check_bands_intersects(data1, data2) +# # pre-condition - equal bands from the same sensor +# # bands with the same name, must be from the same sensor to avoid confusion +# .merge_regular_check_band_sensor(data1, data2) +# # merge +# merged_cube <- .merge_strategy_tile(data1, data2, FALSE) +# # assign a new class, meaning the cube must be regularized to be used +# class(merged_cube) <- c("combined_cube", class(data1)) +# # return +# return(merged_cube) +# } + +# Already incorporated +# .cube_merge <- function(data1, data2) { +# data1 <- slider::slide2_dfr(data1, data2, function(x, y) { +# .fi(x) <- dplyr::arrange( +# dplyr::bind_rows(.fi(x), .fi(y)), +# .data[["date"]], +# .data[["band"]], +# .data[["fid"]] +# ) +# # remove duplicates +# .fi(x) <- dplyr::distinct( +# .fi(x), +# .data[["band"]], +# .data[["date"]], +# .keep_all = TRUE +# ) +# +# return(x) +# }) +# return(data1) +# } + +# ToDo: Special case +.merge_equal_cube <- function(data1, data2) { + if (inherits(data1, "hls_cube") && inherits(data2, "hls_cube") && + (.cube_collection(data1) == "HLSS30" || + .cube_collection(data2) == "HLSS30")) { + data1[["collection"]] <- "HLSS30" + } - # Change file name 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"]] <- d1_tl - 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 <- .cube_merge(data1, data2) - # Return cubes merged return(data1) } +# ToDo: Special case .merge_single_timeline <- function(data1, data2) { tiles <- .cube_tiles(data1) # update the timeline of the cube with single time step (`data2`) @@ -109,3 +299,59 @@ # Merge cubes and return .cube_merge(data1, data2) } + +# Generalized case +# .merge_distinct_cube <- function(data1, data2) { +# # Get cubes timeline +# d1_tl <- unique(as.Date(.cube_timeline(data1)[[1]])) +# d2_tl <- unique(as.Date(.cube_timeline(data2)[[1]])) +# +# # get intervals +# d1_period <- as.integer( +# lubridate::as.period(lubridate::int_diff(d1_tl)), "days" +# ) +# d2_period <- as.integer( +# lubridate::as.period(lubridate::int_diff(d2_tl)), "days" +# ) +# # pre-condition - are periods regular? +# .check_that( +# length(unique(d1_period)) == 1 && length(unique(d2_period)) == 1 +# ) +# # pre-condition - Do cubes have the same periods? +# .check_that( +# unique(d1_period) == unique(d2_period) +# ) +# # pre-condition - are the cubes start date less than period timeline? +# .check_that( +# abs(d1_period[[1]] - d2_period[[2]]) <= unique(d2_period) +# ) +# +# # Change file name 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"]] <- d1_tl +# 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 <- .cube_merge(data1, data2) +# # Return cubes merged +# return(data1) +# } + +# not used anywhere! +# .merge_diff_timelines <- function(t1, t2) { +# abs(as.Date(t1) - as.Date(t2)) +# } + + diff --git a/R/api_timeline.R b/R/api_timeline.R index cc0798b3b..2cf78acdd 100644 --- a/R/api_timeline.R +++ b/R/api_timeline.R @@ -311,3 +311,27 @@ return(FALSE) } } + +#' @title Check if two timelines overlaps. +#' @name .timeline_has_overlap +#' @keywords internal +#' @noRd +#' +#' @description This function checks if the given two timeline overlaps. +#' +#' @param timeline1 First timeline +#' @param timeline2 Second timeline. +#' @return TRUE if first and second timeline overlaps. +#' +.timeline_has_overlap <- function(timeline1, timeline2) { + start1 <- min(timeline1) + end1 <- max(timeline1) + start2 <- min(timeline2) + end2 <- max(timeline2) + + if (start1 <= end2 && start2 <= end1) { + return(TRUE) + } else { + return(FALSE) + } +} diff --git a/R/sits_merge.R b/R/sits_merge.R index efb151bbe..a2711f68c 100644 --- a/R/sits_merge.R +++ b/R/sits_merge.R @@ -30,8 +30,6 @@ #' these suffixes will be added #' (character vector). #' -#' @param irregular Are those irregular data cubes? -#' #' @return merged data sets (tibble of class "sits" or #' tibble of class "raster_cube") #' @examples @@ -99,72 +97,76 @@ sits_merge.sits <- function(data1, data2, ..., suffix = c(".1", ".2")) { #' @rdname sits_merge #' @export -sits_merge.sar_cube <- function(data1, data2, ..., irregular = FALSE) { +sits_merge.sar_cube <- function(data1, data2, ...) { .check_set_caller("sits_merge_sar_cube") # pre-condition - check cube type .check_is_raster_cube(data1) .check_is_raster_cube(data2) - if (any(!.cube_is_regular(data1), !.cube_is_regular(data2))) { - .check_that( - irregular, msg = .conf("messages", "sits_merge_sar_cube_irregular") - ) - return(.merge_irregular_cube(data1, data2)) - } - # pre-condition for merge is having the same tiles - common_tiles <- intersect(data1[["tile"]], data2[["tile"]]) - .check_that(length(common_tiles) > 0) - # filter cubes by common tiles and arrange them - data1 <- dplyr::arrange( - dplyr::filter(data1, .data[["tile"]] %in% common_tiles), - .data[["tile"]] - ) - data2 <- dplyr::arrange( - dplyr::filter(data2, .data[["tile"]] %in% common_tiles), - .data[["tile"]] - ) - if (length(.cube_timeline(data2)[[1]]) == 1) { - return(.merge_single_timeline(data1, data2)) - } - if (inherits(data2, "sar_cube")) { - return(.merge_equal_cube(data1, data2)) + # Define merged cube + merged_cube <- NULL + # Check if cube is regular + is_regular <- all(.cube_is_regular(data1), .cube_is_regular(data2)) + if (is_regular) { + # Regular cube case + merged_cube <- .merge_regular_cube(data1, data2) } else { - return(.merge_distinct_cube(data1, data2)) + # Irregular cube case + merged_cube <- .merge_irregular_cube(data1, data2) } + merged_cube + + # # Irregular cube case + # if (all(!.cube_is_regular(data1), !.cube_is_regular(data2))) { + # return(.merge_irregular_cube(data1, data2)) + # } + # # Regular cube case + # Check if timelines has overlaps + # if (.timeline_has_overlap(d1_tl, )) + # # pre-condition for merge is having the same tiles + # common_tiles <- intersect(data1[["tile"]], data2[["tile"]]) + # .check_that(length(common_tiles) > 0) + # # filter cubes by common tiles and arrange them + # data1 <- dplyr::arrange( + # dplyr::filter(data1, .data[["tile"]] %in% common_tiles), + # .data[["tile"]] + # ) + # data2 <- dplyr::arrange( + # dplyr::filter(data2, .data[["tile"]] %in% common_tiles), + # .data[["tile"]] + # ) + # if (length(.cube_timeline(data2)[[1]]) == 1) { + # return(.merge_single_timeline(data1, data2)) + # } + # if (inherits(data2, "sar_cube")) { + # return(.merge_equal_cube(data1, data2)) + # } else { + # return(.merge_distinct_cube(data1, data2)) + # } } #' @rdname sits_merge #' @export -sits_merge.raster_cube <- function(data1, data2, ..., irregular = FALSE) { +sits_merge.raster_cube <- function(data1, data2, ...) { .check_set_caller("sits_merge_raster_cube") # pre-condition - check cube type .check_is_raster_cube(data1) .check_is_raster_cube(data2) - if (any(!.cube_is_regular(data1), !.cube_is_regular(data2))) { - .check_that( - irregular, msg = .conf("messages", "sits_merge_raster_cube_irregular") - ) - return(.merge_irregular_cube(data1, data2)) - } - # pre-condition for merge is having the same tiles - common_tiles <- intersect(data1[["tile"]], data2[["tile"]]) - .check_that(length(common_tiles) > 0) - # filter cubes by common tiles and arrange them - data1 <- dplyr::arrange( - dplyr::filter(data1, .data[["tile"]] %in% common_tiles), - .data[["tile"]] - ) - data2 <- dplyr::arrange( - dplyr::filter(data2, .data[["tile"]] %in% common_tiles), - .data[["tile"]] - ) - if (length(.cube_timeline(data2)[[1]]) == 1) { - return(.merge_single_timeline(data1, data2)) + # Define merged cube + merged_cube <- NULL + # Special case: Unique timeline cubes + if (...) { + merged_cube <- .merge_unique_timeline() } - if (inherits(data2, "sar_cube")) { - return(.merge_distinct_cube(data1, data2)) + # Check if cube is regular + is_regular <- all(.cube_is_regular(data1), .cube_is_regular(data2)) + if (is_regular) { + # Regular cube case + merged_cube <- .merge_regular_cube(data1, data2) } else { - return(.merge_equal_cube(data1, data2)) + # Irregular cube case + merged_cube <- .merge_irregular_cube(data1, data2) } + merged_cube } #' @rdname sits_merge diff --git a/man/sits_merge.Rd b/man/sits_merge.Rd index 42e039970..933e4451e 100644 --- a/man/sits_merge.Rd +++ b/man/sits_merge.Rd @@ -12,9 +12,9 @@ sits_merge(data1, data2, ...) \method{sits_merge}{sits}(data1, data2, ..., suffix = c(".1", ".2")) -\method{sits_merge}{sar_cube}(data1, data2, ..., irregular = FALSE) +\method{sits_merge}{sar_cube}(data1, data2, ...) -\method{sits_merge}{raster_cube}(data1, data2, ..., irregular = FALSE) +\method{sits_merge}{raster_cube}(data1, data2, ...) \method{sits_merge}{default}(data1, data2, ...) } @@ -30,8 +30,6 @@ or data cube (tibble of class "raster_cube") .} \item{suffix}{If there are duplicate bands in data1 and data2 these suffixes will be added (character vector).} - -\item{irregular}{Are those irregular data cubes?} } \value{ merged data sets (tibble of class "sits" or diff --git a/tests/testthat/test-merge.R b/tests/testthat/test-merge.R new file mode 100644 index 000000000..b540ad6f0 --- /dev/null +++ b/tests/testthat/test-merge.R @@ -0,0 +1,102 @@ +test_that("sits_merge - same sensor, same bands, same tiles, compatible timeline", { + s2a_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("BLUE"), + tiles = c("53HQE","53HPE"), + start_date = "2019-01-01", + end_date = "2019-07-10", + progress = FALSE + ) + }, + .default = NULL + ) + + s2b_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53HQE","53HPE"), + start_date = "2019-01-01", + end_date = "2019-07-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), + message = "DEAustralia is not accessible" + ) + + sentinel_cube <- sits_merge(s2a_cube, s2b_cube) + + expect_true(all(sits_bands(sentinel_cube) %in% c("BLUE"))) + expect_equal(nrow(sentinel_cube), 2) + +}) + +test_that("sits_merge - same sensor, different bands, same tiles, compatible timeline", { + +}) + +test_that("sits_merge - different sensor, different bands, same tiles, compatible timeline", { + +}) + +test_that("sits_merge - different sensor, same bands, same tiles, compatible timeline", { + +}) + +test_that("sits_merge - different sensor, same bands, different tiles, compatible timeline", { + +}) + +test_that("sits_merge - different sensor, different bands, different tiles, compatible timeline", { + +}) + +test_that("sits_merge - different sensor, different bands, different tiles, different timeline", { + +}) + +test_that("sits_merge - same sensor, same bands, same tiles, different timeline", { + s2a_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-02-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s2b_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-02-01", + progress = FALSE + ) + }, + .default = NULL + ) + + expect_error( + .cube_is_regular(sits_merge(s2a_cube, s2b_cube)) + ) +}) + From c17abae6b3523399e2710d6815ae911b4e30d73c Mon Sep 17 00:00:00 2001 From: Felipe Date: Fri, 22 Nov 2024 11:59:54 +0000 Subject: [PATCH 03/19] update sits_merge function --- R/api_merge.R | 60 +++++++++++++++++++++++++++++++++++-- R/sits_merge.R | 11 ++++--- tests/testthat/test-merge.R | 2 +- 3 files changed, 66 insertions(+), 7 deletions(-) diff --git a/R/api_merge.R b/R/api_merge.R index ff902999b..b18cec99e 100644 --- a/R/api_merge.R +++ b/R/api_merge.R @@ -160,6 +160,27 @@ }) } +.tile_has_unique_period <- function(tile) { + # get cubes timeline + d1_tl <- unique(as.Date(.cube_timeline(tile)[[1]])) + # get unique period + period_count <- length(unique(as.integer( + lubridate::as.period(lubridate::int_diff(d1_tl)), "days" + ))) + if (inherits(tile, "bdc_cube") && period_count > 1) { + .check_that( + length(unique(lubridate::year(.cube_timeline(tile)[[1]]))) > 1, + msg = "Cube has different lengths in the same year." + ) + period_count <- 1 + } + period_count == 1 +} + +.cube_has_unique_period <- function(cube) { + all(slider::slide_lgl(cube, .tile_has_unique_period)) +} + .merge_regular_check_periods <- function(data1, data2) { # get cubes timeline d1_tl <- unique(as.Date(.cube_timeline(data1)[[1]])) @@ -181,6 +202,10 @@ ) } +.merge_timeline_same_length <- function(data1, data2) { + length(.cube_timeline(data1)) == length(.cube_timeline(data2)) +} + .merge_regular_cube <- function(data1, data2) { # pre-condition - timelines overlaps # in case of regular cube it is assumed the timeline must overlap @@ -191,8 +216,39 @@ # pre-condition - equal bands must be from the same sensor # bands with the same name, must be from the same sensor to avoid confusion .merge_check_band_sensor(data1, data2) - # ToDo: Cut timeline at overlapping intervals when length(ts1) != length(ts2) - # get tile overlaps + if (!.merge_timeline_same_length(data1, data2)) { + # TODO: warning avisando o usuário que os cubos tem timelines + # com lengths diferentes + t1 <- .cube_timeline(data1)[[1]] + t2 <- .cube_timeline(data2)[[1]] + + if (length(t1) > length(t2)) { + ref <- t1[t1 >= min(t2) & t1 <= max(t2)] + } else { + ref <- t2[t2 >= min(t1) & t2 <= max(t1)] + } + + data1 <- .cube_filter_interval( + data1, start_date = min(ref), end_date = max(ref) + ) + + data2 <- .cube_filter_interval( + data2, start_date = min(ref), end_date = max(ref) + ) + + if (length(.cube_timeline(data1)) != length(.cube_timeline(data2))) { + min_length <- min(c(length(.cube_timeline(data1)), + length(.cube_timeline(data2)))) + + data1 <- .cube_filter_dates( + data1, .cube_timeline(data1)[[1]][seq_len(min_length)] + ) + data2 <- .cube_filter_dates( + data2, .cube_timeline(data2)[[1]][seq_len(min_length)] + ) + } + } + tiles_overlaps <- .merge_tiles_overlaps(data1, data2) # define the strategy (default - merge tiles) merge_strategy <- NULL diff --git a/R/sits_merge.R b/R/sits_merge.R index a2711f68c..9460cea5c 100644 --- a/R/sits_merge.R +++ b/R/sits_merge.R @@ -154,12 +154,15 @@ sits_merge.raster_cube <- function(data1, data2, ...) { # Define merged cube merged_cube <- NULL # Special case: Unique timeline cubes - if (...) { - merged_cube <- .merge_unique_timeline() - } + # if (...) { + # merged_cube <- .merge_unique_timeline() + # } # Check if cube is regular is_regular <- all(.cube_is_regular(data1), .cube_is_regular(data2)) - if (is_regular) { + has_unique_period <- all( + .cube_has_unique_period(data1), .cube_has_unique_period(data2) + ) + if (is_regular && has_unique_period) { # Regular cube case merged_cube <- .merge_regular_cube(data1, data2) } else { diff --git a/tests/testthat/test-merge.R b/tests/testthat/test-merge.R index b540ad6f0..f49c261df 100644 --- a/tests/testthat/test-merge.R +++ b/tests/testthat/test-merge.R @@ -88,7 +88,7 @@ test_that("sits_merge - same sensor, same bands, same tiles, different timeline" bands = c("BLUE"), tiles = c("53HQE"), start_date = "2019-01-01", - end_date = "2019-02-01", + end_date = "2019-02-10", progress = FALSE ) }, From 230c12976907a9fcbe7e1c25a1e8afb96d95cd21 Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Fri, 22 Nov 2024 14:26:05 -0300 Subject: [PATCH 04/19] enhance sits_merge implementation, including dem and hls cases --- NAMESPACE | 1 - R/api_cube.R | 5 + R/api_merge.R | 257 +++------- R/api_regularize.R | 4 + R/api_tile.R | 19 + R/sits_merge.R | 63 +-- man/sits_merge.Rd | 3 - tests/testthat/test-cube-deaustralia.R | 14 +- tests/testthat/test-merge.R | 642 ++++++++++++++++++++++++- 9 files changed, 727 insertions(+), 281 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 3b52a6476..f38f97c2e 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -455,7 +455,6 @@ S3method(sits_labels,sits_model) S3method(sits_labels_summary,sits) S3method(sits_merge,default) S3method(sits_merge,raster_cube) -S3method(sits_merge,sar_cube) S3method(sits_merge,sits) S3method(sits_mixture_model,default) S3method(sits_mixture_model,derived_cube) diff --git a/R/api_cube.R b/R/api_cube.R index bbd2e6edd..4e03c6958 100644 --- a/R/api_cube.R +++ b/R/api_cube.R @@ -1596,3 +1596,8 @@ NULL .cube_has_base_info <- function(cube) { return(.has(cube[["base_info"]])) } + +.cube_has_unique_period <- function(cube) { + all(slider::slide_lgl(cube, .tile_has_unique_period)) +} + diff --git a/R/api_merge.R b/R/api_merge.R index b18cec99e..66a492b2b 100644 --- a/R/api_merge.R +++ b/R/api_merge.R @@ -17,7 +17,7 @@ .merge_adjust_timeline <- function(data1, data2) { # reference timeline - reference_tl <- .cube_timeline(data1) + reference_tl <- .cube_timeline(data1)[[1]] # Adjust dates / bands slider::slide_dfr(data2, function(y) { fi_list <- purrr::map(.tile_bands(y), function(band) { @@ -37,24 +37,6 @@ }) } -.merge_get_ts_within <- function(data1, data2) { - # extract timelines - d1_tl <- .cube_timeline(data1) - d2_tl <- .cube_timeline(data2) - - # Check if all dates in A are in B - if (all(datesA %in% datesB)) { - return(datesA) # A is contained in B - } - - # Check if all dates in B are in A - if (all(datesB %in% datesA)) { - return(datesB) # B is contained in A - } - - return(NULL) # Neither is contained -} - .merge_check_bands_intersects <- function(data1, data2) { # Extract band intersects bands_intersects <- .merge_bands_intersects(data1, data2) @@ -71,6 +53,10 @@ } } +.merge_timeline_has_same_length <- function(data1, data2) { + length(.cube_timeline(data1)[[1]]) == length(.cube_timeline(data2)[[1]]) +} + # ---- Merge strategies ---- .merge_strategy_file <- function(data1, data2, adjust_timeline) { # adjust second cube timeline, based on the first cube @@ -111,42 +97,11 @@ if (adjust_timeline) { data2 <- .merge_adjust_timeline(data1, data2) } - # Extract band intersects - bands_intersects <- .merge_bands_intersects(data1, data2) - # Use only intersect bands - data1 <- .select_raster_bands(data1, bands_intersects) - data2 <- .select_raster_bands(data2, bands_intersects) # Merge dplyr::bind_rows(data1, data2) } # ---- Regular cubes ---- -# .merge_regular_cut_timeline <- function(data1, data2) { -# # extract timelines -# d1_tl <- .cube_timeline(data1) -# d2_tl <- .cube_timeline(data2) -# # extract tiles -# tiles <- .merge_tiles_overlaps(data1, data2) -# # merge cubes -# .map_dfr(tiles, function(tile) { -# # select data in the selected tile -# data1_in_tile <- .select_raster_tiles(data1, tile) -# data2_in_tile <- .select_raster_tiles(data2, tile) -# # extract timelines -# d1_tl <- .tile_timeline(data1_in_tile) -# d2_tl <- .tile_timeline(data2_in_tile) -# # get min/max dates -# min_tl <- min( -# min(d1_tl), min(d2_tl) -# ) -# max_tl <- max( -# max(d1_tl), max(d2_tl) -# ) -# # cut timeline -# .tile_filter_interval(tile, min_tl, max_tl) -# }) -# } - .merge_regular_check_timeline_overlaps <- function(data1, data2) { # extract timelines d1_tl <- .cube_timeline(data1) @@ -160,27 +115,6 @@ }) } -.tile_has_unique_period <- function(tile) { - # get cubes timeline - d1_tl <- unique(as.Date(.cube_timeline(tile)[[1]])) - # get unique period - period_count <- length(unique(as.integer( - lubridate::as.period(lubridate::int_diff(d1_tl)), "days" - ))) - if (inherits(tile, "bdc_cube") && period_count > 1) { - .check_that( - length(unique(lubridate::year(.cube_timeline(tile)[[1]]))) > 1, - msg = "Cube has different lengths in the same year." - ) - period_count <- 1 - } - period_count == 1 -} - -.cube_has_unique_period <- function(cube) { - all(slider::slide_lgl(cube, .tile_has_unique_period)) -} - .merge_regular_check_periods <- function(data1, data2) { # get cubes timeline d1_tl <- unique(as.Date(.cube_timeline(data1)[[1]])) @@ -202,11 +136,9 @@ ) } -.merge_timeline_same_length <- function(data1, data2) { - length(.cube_timeline(data1)) == length(.cube_timeline(data2)) -} - .merge_regular_cube <- function(data1, data2) { + # get tile overlaps + tiles_overlaps <- .merge_tiles_overlaps(data1, data2) # pre-condition - timelines overlaps # in case of regular cube it is assumed the timeline must overlap # to avoid the creation of inconsistent / irregular cubes @@ -216,7 +148,17 @@ # pre-condition - equal bands must be from the same sensor # bands with the same name, must be from the same sensor to avoid confusion .merge_check_band_sensor(data1, data2) - if (!.merge_timeline_same_length(data1, data2)) { + # pre-condition - bands must intersect if more then 1 tile is provided + if (length(tiles_overlaps) > 1) { + .merge_check_bands_intersects(data1, data2) + } + # Extract band intersects + bands_intersects <- .merge_bands_intersects(data1, data2) + # Use only intersect bands + data1 <- .select_raster_bands(data1, bands_intersects) + data2 <- .select_raster_bands(data2, bands_intersects) + # Check timeline consistency + if (!.merge_timeline_has_same_length(data1, data2)) { # TODO: warning avisando o usuário que os cubos tem timelines # com lengths diferentes t1 <- .cube_timeline(data1)[[1]] @@ -228,6 +170,8 @@ ref <- t2[t2 >= min(t1) & t2 <= max(t1)] } + .check_that(.has(ref)) + data1 <- .cube_filter_interval( data1, start_date = min(ref), end_date = max(ref) ) @@ -248,12 +192,9 @@ ) } } - - tiles_overlaps <- .merge_tiles_overlaps(data1, data2) # define the strategy (default - merge tiles) merge_strategy <- NULL # case: same tiles, merge file info - tiles_overlaps <- .merge_tiles_overlaps(data1, data2) if (.has(tiles_overlaps)) { merge_strategy <- .merge_strategy_file # case 2: different tiles, merge tile rows @@ -266,11 +207,20 @@ # ---- Irregular cubes ---- .merge_irregular_cube <- function(data1, data2) { + # get tile overlaps + tiles_overlaps <- .merge_tiles_overlaps(data1, data2) # pre-condition - equal bands from the same sensor # bands with the same name, must be from the same sensor to avoid confusion .merge_check_band_sensor(data1, data2) - # get tile overlaps - tiles_overlaps <- .merge_tiles_overlaps(data1, data2) + # pre-condition - bands must intersect if more then 1 tile is provided + if (length(tiles_overlaps) > 1) { + .merge_check_bands_intersects(data1, data2) + } + # Extract band intersects + bands_intersects <- .merge_bands_intersects(data1, data2) + # Use only intersect bands + data1 <- .select_raster_bands(data1, bands_intersects) + data2 <- .select_raster_bands(data2, bands_intersects) # define the strategy (default - merge tiles) merge_strategy <- NULL # case: same tiles, merge file info @@ -285,129 +235,44 @@ return(merged_cube) } -# .merge_irregular_cube <- function(data1, data2) { -# # pre-condition - intersecting bands -# .merge_check_bands_intersects(data1, data2) -# # pre-condition - equal bands from the same sensor -# # bands with the same name, must be from the same sensor to avoid confusion -# .merge_regular_check_band_sensor(data1, data2) -# # merge -# merged_cube <- .merge_strategy_tile(data1, data2, FALSE) -# # assign a new class, meaning the cube must be regularized to be used -# class(merged_cube) <- c("combined_cube", class(data1)) -# # return -# return(merged_cube) -# } - -# Already incorporated -# .cube_merge <- function(data1, data2) { -# data1 <- slider::slide2_dfr(data1, data2, function(x, y) { -# .fi(x) <- dplyr::arrange( -# dplyr::bind_rows(.fi(x), .fi(y)), -# .data[["date"]], -# .data[["band"]], -# .data[["fid"]] -# ) -# # remove duplicates -# .fi(x) <- dplyr::distinct( -# .fi(x), -# .data[["band"]], -# .data[["date"]], -# .keep_all = TRUE -# ) -# -# return(x) -# }) -# return(data1) -# } - -# ToDo: Special case -.merge_equal_cube <- function(data1, data2) { - if (inherits(data1, "hls_cube") && inherits(data2, "hls_cube") && - (.cube_collection(data1) == "HLSS30" || - .cube_collection(data2) == "HLSS30")) { - data1[["collection"]] <- "HLSS30" +# ---- Special case: DEM Cube ---- +.merge_dem_cube <- function(data1, data2) { + # define cubes + dem_cube <- data1 + other_cube <- data2 + # check which cube is the DEM + if (inherits(data2, "dem_cube")) { + # move DEM cube (de) + dem_cube <- data2 + other_cube <- data1 } - - data1 <- .cube_merge(data1, data2) - return(data1) -} - -# ToDo: Special case -.merge_single_timeline <- function(data1, data2) { - tiles <- .cube_tiles(data1) - # update the timeline of the cube with single time step (`data2`) - data2 <- .map_dfr(tiles, function(tile_name) { - tile_data1 <- .cube_filter_tiles(data1, tile_name) - tile_data2 <- .cube_filter_tiles(data2, tile_name) + tiles <- .cube_tiles(other_cube) + # update the timeline of the cube with single time step + dem_cube <- .map_dfr(tiles, function(tile_name) { + tile_other <- .cube_filter_tiles(other_cube, tile_name) + tile_dem <- .cube_filter_tiles(dem_cube, tile_name) # Get data1 timeline. - d1_tl <- unique(as.Date(.cube_timeline(tile_data1)[[1]])) + d1_tl <- unique(as.Date(.cube_timeline(tile_other)[[1]])) # Create new `file_info` using dates from `data1` timeline. - fi_new <- purrr::map(.tile_timeline(tile_data1), function(date_row) { - fi <- .fi(tile_data2) + fi_new <- purrr::map(.tile_timeline(tile_other), function(date_row) { + fi <- .fi(tile_dem) fi[["date"]] <- as.Date(date_row) fi }) # Assign the new `file_into` into `data2` - tile_data2[["file_info"]] <- list(dplyr::bind_rows(fi_new)) - tile_data2 + tile_dem[["file_info"]] <- list(dplyr::bind_rows(fi_new)) + tile_dem }) - # Merge cubes and return - .cube_merge(data1, data2) + # merge cubes and return + .merge_strategy_file(other_cube, dem_cube, FALSE) } -# Generalized case -# .merge_distinct_cube <- function(data1, data2) { -# # Get cubes timeline -# d1_tl <- unique(as.Date(.cube_timeline(data1)[[1]])) -# d2_tl <- unique(as.Date(.cube_timeline(data2)[[1]])) -# -# # get intervals -# d1_period <- as.integer( -# lubridate::as.period(lubridate::int_diff(d1_tl)), "days" -# ) -# d2_period <- as.integer( -# lubridate::as.period(lubridate::int_diff(d2_tl)), "days" -# ) -# # pre-condition - are periods regular? -# .check_that( -# length(unique(d1_period)) == 1 && length(unique(d2_period)) == 1 -# ) -# # pre-condition - Do cubes have the same periods? -# .check_that( -# unique(d1_period) == unique(d2_period) -# ) -# # pre-condition - are the cubes start date less than period timeline? -# .check_that( -# abs(d1_period[[1]] - d2_period[[2]]) <= unique(d2_period) -# ) -# -# # Change file name 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"]] <- d1_tl -# 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 <- .cube_merge(data1, data2) -# # Return cubes merged -# return(data1) -# } - -# not used anywhere! -# .merge_diff_timelines <- function(t1, t2) { -# abs(as.Date(t1) - as.Date(t2)) -# } - +.merge_hls_cube <- function(data1, data2) { + if ((.cube_collection(data1) == "HLSS30" || + .cube_collection(data2) == "HLSS30")) { + data1[["collection"]] <- "HLSS30" + } + # merge cubes and return + .merge_strategy_file(data1, data2, FALSE) +} diff --git a/R/api_regularize.R b/R/api_regularize.R index 114ac527e..6d491b4e0 100644 --- a/R/api_regularize.R +++ b/R/api_regularize.R @@ -396,3 +396,7 @@ cube_class <- c(cube_class[[1]], "rainfall_cube", cube_class[-1]) .cube_set_class(cube, cube_class) } + +.reg_tile_convert.default <- function(cube, grid_system, roi = NULL, tiles = NULL) { + return(cube) +} diff --git a/R/api_tile.R b/R/api_tile.R index 55372b145..df41654f3 100644 --- a/R/api_tile.R +++ b/R/api_tile.R @@ -1682,3 +1682,22 @@ NULL .tile_base_info <- function(tile) { return(tile[["base_info"]][[1]]) } + +.tile_has_unique_period <- function(tile) { + # get cubes timeline + d1_tl <- unique(as.Date(.cube_timeline(tile)[[1]])) + # get unique period + period_count <- length(unique(as.integer( + lubridate::as.period(lubridate::int_diff(d1_tl)), "days" + ))) + if (inherits(tile, "bdc_cube") && period_count > 1) { + .check_that( + length(unique(lubridate::year(.cube_timeline(tile)[[1]]))) > 1, + msg = "Cube has different lengths in the same year." + ) + period_count <- 1 + } + period_count == 1 +} + + diff --git a/R/sits_merge.R b/R/sits_merge.R index 9460cea5c..73c4f1738 100644 --- a/R/sits_merge.R +++ b/R/sits_merge.R @@ -95,55 +95,6 @@ sits_merge.sits <- function(data1, data2, ..., suffix = c(".1", ".2")) { return(result) } -#' @rdname sits_merge -#' @export -sits_merge.sar_cube <- function(data1, data2, ...) { - .check_set_caller("sits_merge_sar_cube") - # pre-condition - check cube type - .check_is_raster_cube(data1) - .check_is_raster_cube(data2) - # Define merged cube - merged_cube <- NULL - # Check if cube is regular - is_regular <- all(.cube_is_regular(data1), .cube_is_regular(data2)) - if (is_regular) { - # Regular cube case - merged_cube <- .merge_regular_cube(data1, data2) - } else { - # Irregular cube case - merged_cube <- .merge_irregular_cube(data1, data2) - } - merged_cube - - # # Irregular cube case - # if (all(!.cube_is_regular(data1), !.cube_is_regular(data2))) { - # return(.merge_irregular_cube(data1, data2)) - # } - # # Regular cube case - # Check if timelines has overlaps - # if (.timeline_has_overlap(d1_tl, )) - # # pre-condition for merge is having the same tiles - # common_tiles <- intersect(data1[["tile"]], data2[["tile"]]) - # .check_that(length(common_tiles) > 0) - # # filter cubes by common tiles and arrange them - # data1 <- dplyr::arrange( - # dplyr::filter(data1, .data[["tile"]] %in% common_tiles), - # .data[["tile"]] - # ) - # data2 <- dplyr::arrange( - # dplyr::filter(data2, .data[["tile"]] %in% common_tiles), - # .data[["tile"]] - # ) - # if (length(.cube_timeline(data2)[[1]]) == 1) { - # return(.merge_single_timeline(data1, data2)) - # } - # if (inherits(data2, "sar_cube")) { - # return(.merge_equal_cube(data1, data2)) - # } else { - # return(.merge_distinct_cube(data1, data2)) - # } -} - #' @rdname sits_merge #' @export sits_merge.raster_cube <- function(data1, data2, ...) { @@ -153,10 +104,16 @@ sits_merge.raster_cube <- function(data1, data2, ...) { .check_is_raster_cube(data2) # Define merged cube merged_cube <- NULL - # Special case: Unique timeline cubes - # if (...) { - # merged_cube <- .merge_unique_timeline() - # } + # Special case: DEM cube + is_dem_cube <- any(inherits(data1, "dem_cube"), inherits(data2, "dem_cube")) + if (is_dem_cube) { + return(.merge_dem_cube(data1, data2)) + } + # Special case: HLS cube + is_hls_cube <- all(inherits(data1, "hls_cube"), inherits(data2, "hls_cube")) + if (is_hls_cube) { + return(.merge_hls_cube(data1, data2)) + } # Check if cube is regular is_regular <- all(.cube_is_regular(data1), .cube_is_regular(data2)) has_unique_period <- all( diff --git a/man/sits_merge.Rd b/man/sits_merge.Rd index 933e4451e..ccacd7dc7 100644 --- a/man/sits_merge.Rd +++ b/man/sits_merge.Rd @@ -3,7 +3,6 @@ \name{sits_merge} \alias{sits_merge} \alias{sits_merge.sits} -\alias{sits_merge.sar_cube} \alias{sits_merge.raster_cube} \alias{sits_merge.default} \title{Merge two data sets (time series or cubes)} @@ -12,8 +11,6 @@ sits_merge(data1, data2, ...) \method{sits_merge}{sits}(data1, data2, ..., suffix = c(".1", ".2")) -\method{sits_merge}{sar_cube}(data1, data2, ...) - \method{sits_merge}{raster_cube}(data1, data2, ...) \method{sits_merge}{default}(data1, data2, ...) diff --git a/tests/testthat/test-cube-deaustralia.R b/tests/testthat/test-cube-deaustralia.R index 8898707a3..b66119a29 100644 --- a/tests/testthat/test-cube-deaustralia.R +++ b/tests/testthat/test-cube-deaustralia.R @@ -417,7 +417,7 @@ test_that( sits_cube( source = "DEAUSTRALIA", collection = "GA_S2BM_ARD_3", - bands = c("SWIR-2"), + bands = c("BLUE", "RED"), tiles = c("53HQE","53HPE"), start_date = "2019-01-01", end_date = "2019-08-28", @@ -431,16 +431,14 @@ test_that( message = "DEAustralia is not accessible" ) - sentinel_cube <- sits_merge(s2a_cube, s2b_cube, irregular = TRUE) + sentinel_cube <- sits_merge(s2a_cube, s2b_cube) - expect_true(all(sits_bands(sentinel_cube) %in% c( - "BLUE", "NIR-2", "SWIR-2" - ))) + expect_true(all(sits_bands(sentinel_cube) %in% c("BLUE"))) expect_equal(nrow(sentinel_cube), 2) r <- .raster_open_rast(.tile_path(sentinel_cube)) - expect_equal(sentinel_cube$xmax[[1]], .raster_xmax(r), tolerance = 1) - expect_equal(sentinel_cube$xmin[[1]], .raster_xmin(r), tolerance = 1) - expect_true(all(sentinel_cube$tile %in% c("53HQE","53HPE"))) + expect_equal(sentinel_cube[["xmax"]][[1]], .raster_xmax(r), tolerance = 1) + expect_equal(sentinel_cube[["xmin"]][[1]], .raster_xmin(r), tolerance = 1) + expect_true(all(sentinel_cube[["tile"]] %in% c("53HQE","53HPE"))) }) test_that("Creating GA_LS_FC_3 cubes from DEAustralia", { diff --git a/tests/testthat/test-merge.R b/tests/testthat/test-merge.R index f49c261df..d73c458bc 100644 --- a/tests/testthat/test-merge.R +++ b/tests/testthat/test-merge.R @@ -1,11 +1,113 @@ -test_that("sits_merge - same sensor, same bands, same tiles, compatible timeline", { +test_that("sits_merge - irregular cubes with same bands and tile", { + # Test case: If the bands are the same, the cube will have the combined + # timeline of both cubes. This is useful to merge data from the same sensors + # from different satellites (e.g, Sentinel-2A with Sentinel-2B). + # For irregular cubes, all dates are returned. + + # Test 1a: Single tile with different time period (irregular cube) s2a_cube <- .try( { sits_cube( source = "DEAUSTRALIA", collection = "ga_s2am_ard_3", bands = c("BLUE"), - tiles = c("53HQE","53HPE"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s2b_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), + message = "DEAustralia is not accessible" + ) + + merged_cube <- sits_merge(s2a_cube, s2b_cube) + + expect_equal(nrow(merged_cube), 1) + expect_equal(sits_bands(merged_cube), "BLUE") + expect_equal( + length(sits_timeline(merged_cube)), + length(sits_timeline(s2a_cube)) + length(sits_timeline(s2b_cube)) + ) + + r <- .raster_open_rast(.tile_path(merged_cube)) + expect_equal(merged_cube[["xmax"]][[1]], .raster_xmax(r), tolerance = 1) + expect_equal(merged_cube[["xmin"]][[1]], .raster_xmin(r), tolerance = 1) + + # Test 1b: Single tile with different time period (irregular cube) + s2_cube <- .try( + { + sits_cube( + source = "DEAFRICA", + collection = "SENTINEL-2-L2A", + bands = c("B02"), + tiles = c("36NWJ"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s1_cube <- .try( + { + sits_cube( + source = "DEAFRICA", + collection = "SENTINEL-1-RTC", + bands = c("VV"), + orbit = "ascending", + tiles = c("36NWJ"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(s2_cube, s1_cube)), + message = "DEAFRICA is not accessible" + ) + + merged_cube <- sits_merge(s2_cube, s1_cube) + + expect_true(inherits(merged_cube, "combined_cube")) + expect_equal( + length(merged_cube[["tile"]]), + length(s2_cube[["tile"]]) + length(s1_cube[["tile"]]) + ) + expect_equal(suppressWarnings(length(sits_timeline(merged_cube))), 5) + expect_equal( + unique(slider::slide_chr(merged_cube, .tile_bands)), c("B02", "VV") + ) + + # Test 2: Multiple tiles with different time period (irregular cube) + s2a_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("BLUE"), + tiles = c("53HQE", "53HPE"), start_date = "2019-01-01", end_date = "2019-07-10", progress = FALSE @@ -20,7 +122,7 @@ test_that("sits_merge - same sensor, same bands, same tiles, compatible timeline source = "DEAUSTRALIA", collection = "GA_S2BM_ARD_3", bands = c("BLUE"), - tiles = c("53HQE","53HPE"), + tiles = c("53HQE", "53HPE"), start_date = "2019-01-01", end_date = "2019-07-10", progress = FALSE @@ -33,38 +135,336 @@ test_that("sits_merge - same sensor, same bands, same tiles, compatible timeline message = "DEAustralia is not accessible" ) - sentinel_cube <- sits_merge(s2a_cube, s2b_cube) - - expect_true(all(sits_bands(sentinel_cube) %in% c("BLUE"))) - expect_equal(nrow(sentinel_cube), 2) + merged_cube <- sits_merge(s2a_cube, s2b_cube) + expect_equal(nrow(merged_cube), 2) + expect_equal(sits_bands(merged_cube), "BLUE") + expect_equal( + length(sits_timeline(merged_cube)), + length(sits_timeline(s2a_cube)) + length(sits_timeline(s2b_cube)) + ) + r <- .raster_open_rast(.tile_path(merged_cube)) + expect_equal(merged_cube[["xmax"]][[1]], .raster_xmax(r), tolerance = 1) + expect_equal(merged_cube[["xmin"]][[1]], .raster_xmin(r), tolerance = 1) }) -test_that("sits_merge - same sensor, different bands, same tiles, compatible timeline", { +test_that("sits_merge - irregular cubes with same bands and different tile", { + # Test case: If the bands are the same, the cube will have the combined + # timeline of both cubes. This is useful to merge data from the same sensors + # from different satellites (e.g, Sentinel-2A with Sentinel-2B). + # For irregular cubes, all dates are returned. -}) + # Test 1: Different tiles with different time period (irregular cube) + s2a_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s2b_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53JQF"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), + message = "DEAustralia is not accessible" + ) -test_that("sits_merge - different sensor, different bands, same tiles, compatible timeline", { + merged_cube <- sits_merge(s2a_cube, s2b_cube) + expect_true(inherits(merged_cube, "combined_cube")) + expect_equal(suppressWarnings(length(sits_timeline(merged_cube))), 2) }) -test_that("sits_merge - different sensor, same bands, same tiles, compatible timeline", { +test_that("sits_merge - irregular cubes with different bands and tile", { + # Test case: if the bands are different and their timelines should be + # compatible, the bands are joined. The resulting timeline is the one from + # the first cube. This is useful to merge data from different sensors + # (e.g, Sentinel-1 with Sentinel-2). + # For irregular cubes, all dates are returned. + + s2a_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("BLUE", "RED"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s2b_cube <- .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53HQE", "53JQF"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), + message = "DEAustralia is not accessible" + ) + merged_cube <- sits_merge(s2a_cube, s2b_cube) + expect_equal(sits_bands(merged_cube), "BLUE") + expect_equal(merged_cube[["tile"]], "53HQE") }) -test_that("sits_merge - different sensor, same bands, different tiles, compatible timeline", { +test_that("sits_merge - regular cubes with same bands and tile", { + # Test case: If the bands are the same, the cube will have the combined + # timeline of both cubes. This is useful to merge data from the same sensors + # from different satellites (e.g, Sentinel-2A with Sentinel-2B). + # For regular cubes, when timeline has the same length, use them. Otherwise, + # use as timeline the intersect between timelines. + + # Test 1: Tiles with same time period (regular cube) + modis_cube_a <- .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + modis_cube_b <- .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-03-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), + message = "BDC is not accessible" + ) + + merged_cube <- sits_merge(modis_cube_a, modis_cube_b) + + expect_equal(length(sits_timeline(merged_cube)), 2) + expect_equal(sits_bands(merged_cube), "NDVI") + expect_equal(merged_cube[["tile"]], "013011") + + # Test 2: no time-series overlaps (regular cube) + modis_cube_a <- .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + modis_cube_b <- .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), + message = "BDC is not accessible" + ) + expect_error(sits_merge(modis_cube_a, modis_cube_b)) }) -test_that("sits_merge - different sensor, different bands, different tiles, compatible timeline", { +test_that("sits_merge - regular cubes with same bands and different tile", { + # Test case: If the bands are the same, the cube will have the combined + # timeline of both cubes. This is useful to merge data from the same sensors + # from different satellites (e.g, Sentinel-2A with Sentinel-2B). + # For regular cubes, then timeline has the same length, use them. Otherwise, + # use as timeline the intersect between timelines. + + # Test 1: Different tiles + modis_cube_a <- .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22LBH"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + modis_cube_b <- .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), + message = "BDC is not accessible" + ) + + merged_cube <- sits_merge(modis_cube_a, modis_cube_b) + + expect_equal(length(sits_timeline(merged_cube)), 4) + expect_equal( + sits_timeline(modis_cube_b)[seq_len(4)], sits_timeline(merged_cube) + ) + expect_equal(sits_bands(merged_cube), "NDVI") + expect_equal(merged_cube[["tile"]], c("012010", "013011")) + # Test 2: Tile variation in one of the cubes + s2_cube_a <- .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LMR", "20LNR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s2_cube_b <- .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LMR")), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + merged_cube <- sits_merge(s2_cube_a, s2_cube_b) + + expect_equal(sits_timeline(merged_cube), sits_timeline(s2_cube_a)) + expect_equal(nrow(merged_cube), 4) }) -test_that("sits_merge - different sensor, different bands, different tiles, different timeline", { +test_that("sits_merge - regular cubes with different bands and tile", { + # Test case: if the bands are different and their timelines should be + # compatible, the bands are joined. The resulting timeline is the one from + # the first cube. This is useful to merge data from different sensors + # (e.g, Sentinel-1 with Sentinel-2). + # For regular cubes, then timeline has the same length, use them. Otherwise, + # use as timeline the intersect between timelines. + + s2_cube_a <- .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LMR", "20LNR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s2_cube_b <- .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02", "B03"), + roi = sits_tiles_to_roi(c("20LMR")), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + merged_cube <- sits_merge(s2_cube_a, s2_cube_b) + expect_equal(sits_timeline(merged_cube), sits_timeline(s2_cube_a)) + expect_equal(nrow(merged_cube), 4) + expect_equal(sits_bands(merged_cube), "B02") }) -test_that("sits_merge - same sensor, same bands, same tiles, different timeline", { +test_that("sits_merge - regularize combined cubes", { + # Test 1: Same sensor + output_dir <- paste0(tempdir(), "/merge-reg-1") + dir.create(output_dir, showWarnings = FALSE) + s2a_cube <- .try( { sits_cube( @@ -73,7 +473,7 @@ test_that("sits_merge - same sensor, same bands, same tiles, different timeline" bands = c("BLUE"), tiles = c("53HQE"), start_date = "2019-01-01", - end_date = "2019-02-01", + end_date = "2019-04-01", progress = FALSE ) }, @@ -86,17 +486,219 @@ test_that("sits_merge - same sensor, same bands, same tiles, different timeline" source = "DEAUSTRALIA", collection = "GA_S2BM_ARD_3", bands = c("BLUE"), - tiles = c("53HQE"), + tiles = c("53JQF"), + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), + message = "DEAustralia is not accessible" + ) + + # merge + merged_cube <- sits_merge(s2a_cube, s2b_cube) + + # regularize + regularized_cube <- sits_regularize( + cube = merged_cube, + period = "P8D", + res = 720, + output_dir = output_dir + ) + + # test + expect_equal(nrow(regularized_cube), 2) + expect_equal(length(sits_timeline(regularized_cube)), 7) + expect_equal(sits_bands(regularized_cube), "BLUE") + expect_equal(.cube_xres(regularized_cube), 720) + + # Test 2: Different sensor + output_dir <- paste0(tempdir(), "/merge-reg-2") + dir.create(output_dir, showWarnings = FALSE) + + s2_cube <- .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + bands = c("B02"), + tiles = c("19LEF"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s1_cube <- .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-1-RTC", + bands = c("VV"), + tiles = c("19LEF"), + orbit = "descending", + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(s2_cube, s1_cube)), + message = "MPC is not accessible" + ) + + # merge + merged_cube <- sits_merge(s2_cube, s1_cube) + + # regularize + regularized_cube <- sits_regularize( + cube = merged_cube, + period = "P8D", + res = 720, + output_dir = output_dir + ) + + # test + expect_equal(regularized_cube[["tile"]], "19LEF") + expect_equal(length(sits_timeline(regularized_cube)), 7) + expect_equal(sits_bands(regularized_cube), c("B02", "VV")) + expect_equal(.cube_xres(regularized_cube), 720) +}) + +test_that("sits_merge - cubes with different classes", { + s2_cube <- .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + bands = c("B02"), + tiles = c("19LEF"), start_date = "2019-01-01", - end_date = "2019-02-10", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + + s1_cube <- .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-1-RTC", + bands = c("VV"), + tiles = c("19LEF"), + orbit = "descending", + start_date = "2019-02-01", + end_date = "2019-06-10", progress = FALSE ) }, .default = NULL ) - expect_error( - .cube_is_regular(sits_merge(s2a_cube, s2b_cube)) + testthat::skip_if(purrr::is_null(c(s2_cube, s1_cube)), + message = "MPC is not accessible" + ) + + # merge + merged_cube_1 <- sits_merge(s2_cube, s1_cube) + merged_cube_2 <- sits_merge(s1_cube, s2_cube) + + # test + expect_equal(nrow(merged_cube_1), nrow(merged_cube_2)) + expect_equal(sort(merged_cube_1[["tile"]]), sort(merged_cube_2[["tile"]])) +}) + +test_that("sits_merge - special case - dem cube", { + # create S2 cube + s2_dir <- paste0(tempdir(), "/s2") + dir.create(s2_dir, showWarnings = FALSE) + s2_cube <- sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + tiles = "19HBA", + bands = c("B04", "B8A", "B12", "CLOUD"), + start_date = "2021-01-01", + end_date = "2021-03-31", + progress = FALSE + ) + + s2_cube_reg <- sits_regularize( + cube = s2_cube, + period = "P16D", + res = 720, + output_dir = s2_dir, + progress = FALSE ) + + # create DEM cube + dem_dir <- paste0(tempdir(), "/dem") + dir.create(dem_dir, showWarnings = FALSE) + dem_cube <- sits_cube( + source = "MPC", + collection = "COP-DEM-GLO-30", + bands = "ELEVATION", + tiles = "19HBA", + progress = FALSE + ) + + dem_cube_reg <- sits_regularize( + cube = dem_cube, + res = 720, + bands = "ELEVATION", + tiles = "19HBA", + output_dir = dem_dir, + progress = FALSE + ) + + # merge + merged_cube <- sits_merge(s2_cube_reg, dem_cube_reg) + + # test + expect_equal(nrow(merged_cube[["file_info"]][[1]]), 24) + expect_equal(sits_bands(merged_cube), c("B04", "B12", "B8A", "ELEVATION")) }) +test_that("sits_merge - special case - hls cube", { + # define roi + roi <- c( + lon_min = -45.6422, lat_min = -24.0335, + lon_max = -45.0840, lat_max = -23.6178 + ) + + hls_cube_s2 <- sits_cube( + source = "HLS", + collection = "HLSS30", + roi = roi, + bands = c("BLUE", "GREEN", "RED", "CLOUD"), + start_date = as.Date("2020-06-01"), + end_date = as.Date("2020-09-01"), + progress = FALSE + ) + + hls_cube_l8 <- sits_cube( + source = "HLS", + collection = "HLSL30", + roi = roi, + bands = c("BLUE", "GREEN", "RED", "CLOUD"), + start_date = as.Date("2020-06-01"), + end_date = as.Date("2020-09-01"), + progress = FALSE + ) + + # merge + merged_cube <- sits_merge(hls_cube_s2, hls_cube_l8) + + # test + expect_equal(length(sits_timeline(merged_cube)), 19) + expect_equal(sits_bands(merged_cube), c("BLUE", "GREEN", "RED", "CLOUD")) +}) From 437800cde6b466e1fc3f72ae35fe49849303d0b5 Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Fri, 22 Nov 2024 14:26:46 -0300 Subject: [PATCH 05/19] s1-rtc does not require token anymore --- inst/extdata/sources/config_source_mpc.yml | 2 -- 1 file changed, 2 deletions(-) diff --git a/inst/extdata/sources/config_source_mpc.yml b/inst/extdata/sources/config_source_mpc.yml index caaf294a6..9be2e7d78 100644 --- a/inst/extdata/sources/config_source_mpc.yml +++ b/inst/extdata/sources/config_source_mpc.yml @@ -355,5 +355,3 @@ sources: ext_tolerance: 0 grid_system : "MGRS" dates : "2014 to now" - token_vars : ["MPC_TOKEN"] - From 546d44154ed3243bdc86f882f669e72ba3da3334 Mon Sep 17 00:00:00 2001 From: Felipe Date: Tue, 26 Nov 2024 20:06:07 +0000 Subject: [PATCH 06/19] remove gdalUtilities package --- DESCRIPTION | 1 - R/api_raster.R | 62 +++++++++++++++++++++++++++----------------------- 2 files changed, 34 insertions(+), 29 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index f009354c0..d7ab6f1da 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -49,7 +49,6 @@ LazyData: true Imports: yaml, dplyr (>= 1.0.0), - gdalUtilities, grDevices, graphics, lubridate, diff --git a/R/api_raster.R b/R/api_raster.R index fdbae12a0..243a4ec1c 100644 --- a/R/api_raster.R +++ b/R/api_raster.R @@ -784,16 +784,18 @@ .raster_template <- function(base_file, out_file, nlayers, data_type, missing_value) { # Create an empty image template - gdalUtilities::gdal_translate( - src_dataset = .file_path_expand(base_file), - dst_dataset = .file_path_expand(out_file), - ot = .raster_gdal_datatype(data_type), - of = "GTiff", - b = rep(1, nlayers), - scale = c(0, 1, missing_value, missing_value), - a_nodata = missing_value, - co = .conf("gdal_creation_options"), - q = TRUE + .gdal_translate( + file = .file_path_expand(out_file), + base_file = .file_path_expand(base_file), + params = list( + "-ot" = .raster_gdal_datatype(data_type), + "-of" = .conf("gdal_presets", "image", "of"), + "-b" = rep(1, nlayers), + "-scale" = c(0, 1, missing_value, missing_value), + "-a_nodata" = missing_value, + "-co" = .conf("gdal_creation_options") + ), + quiet = TRUE ) # Delete auxiliary files on.exit(unlink(paste0(out_file, ".aux.xml")), add = TRUE) @@ -860,14 +862,16 @@ { # merge using gdal warp suppressWarnings( - gdalUtilities::gdalwarp( - srcfile = merge_files, - dstfile = out_file, - wo = paste0("NUM_THREADS=", multicores), - co = .conf("gdal_creation_options"), - multi = TRUE, - q = TRUE, - overwrite = TRUE + .gdal_warp( + file = out_file, + base_files = merge_files, + params = list( + "-wo" = paste0("NUM_THREADS=", multicores), + "-co" = .conf("gdal_creation_options"), + "-multi" = TRUE, + "-overwrite" = TRUE + ), + quiet = TRUE ) ) }, @@ -881,16 +885,18 @@ { # merge using gdal warp suppressWarnings( - gdalUtilities::gdalwarp( - srcfile = merge_files, - dstfile = out_file, - wo = paste0("NUM_THREADS=", multicores), - ot = .raster_gdal_datatype(data_type), - multi = TRUE, - of = "GTiff", - q = TRUE, - co = .conf("gdal_creation_options"), - overwrite = FALSE + .gdal_warp( + file = out_file, + base_files = merge_files, + params = list( + "-wo" = paste0("NUM_THREADS=", multicores), + "-ot" = .raster_gdal_datatype(data_type), + "-multi" = TRUE, + "-of" = .conf("gdal_presets", "image", "of"), + "-co" = .conf("gdal_creation_options"), + "-overwrite" = FALSE + ), + quiet = TRUE ) ) }, From efb5cd6454f36011d8539c345e0b19baafc401de Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Wed, 27 Nov 2024 22:19:51 -0300 Subject: [PATCH 07/19] fix gdal_translate parameters --- R/api_raster.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/api_raster.R b/R/api_raster.R index 243a4ec1c..69421dcbb 100644 --- a/R/api_raster.R +++ b/R/api_raster.R @@ -791,7 +791,7 @@ "-ot" = .raster_gdal_datatype(data_type), "-of" = .conf("gdal_presets", "image", "of"), "-b" = rep(1, nlayers), - "-scale" = c(0, 1, missing_value, missing_value), + "-scale" = list(0, 1, missing_value, missing_value), "-a_nodata" = missing_value, "-co" = .conf("gdal_creation_options") ), From e83ce4b36c5fe82ffdd212759e4d71b58d6fdab7 Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Wed, 27 Nov 2024 22:20:05 -0300 Subject: [PATCH 08/19] review sits_merge implementation --- R/api_check.R | 12 + R/api_cube.R | 7 +- R/api_merge.R | 292 ++++------- R/sits_merge.R | 32 +- tests/testthat/test-merge.R | 999 ++++++++++++++++++++++-------------- 5 files changed, 765 insertions(+), 577 deletions(-) diff --git a/R/api_check.R b/R/api_check.R index 0f6a3d7ce..ceb9e28cf 100644 --- a/R/api_check.R +++ b/R/api_check.R @@ -1923,6 +1923,18 @@ .check_that(all(bands %in% cube_bands)) return(invisible(cube)) } +#' @title Check if all rows in a cube has the same bands +#' @name .check_cube_row_same_bands +#' @param cube Data cube +#' @return Called for side effects. +#' @keywords internal +#' @noRd +.check_cube_row_same_bands <- function(cube) { + bands <- purrr::map(.compact(slider::slide(cube, .tile_bands)), length) + bands <- .dissolve(bands) + + .check_that(length(unique(bands)) == 1) +} #' @title Check if cubes have the same bbox #' @name .check_cubes_same_bbox #' @keywords internal diff --git a/R/api_cube.R b/R/api_cube.R index 4e03c6958..2a4553717 100644 --- a/R/api_cube.R +++ b/R/api_cube.R @@ -1597,7 +1597,10 @@ NULL return(.has(cube[["base_info"]])) } -.cube_has_unique_period <- function(cube) { - all(slider::slide_lgl(cube, .tile_has_unique_period)) +.cube_sensor <- function(cube) { + .dissolve(slider::slide(cube, .tile_sensor)) } +.cube_satellite <- function(cube) { + .dissolve(slider::slide(cube, .tile_satellite)) +} diff --git a/R/api_merge.R b/R/api_merge.R index 66a492b2b..2a61dc9f7 100644 --- a/R/api_merge.R +++ b/R/api_merge.R @@ -1,13 +1,24 @@ # ---- General utilities ---- -.merge_bands_intersects <- function(data1, data2) { - # Extract bands - d1_bands <- .cube_bands(data1) - d2_bands <- .cube_bands(data2) - # Extract overlaps - intersect(d1_bands, d2_bands) +.merge_has_equal_bands <- function(data1, data2) { + # get cube bands + data1_bands <- .cube_bands(data1) + data2_bands <- .cube_bands(data2) + # verify if both cubes have the same bands + has_same_bands <- all(data1_bands %in% data2_bands) + # if has the same bands, do check for consistency + if (has_same_bands) { + # get bands intersects + bands_intersects <- setdiff(data1_bands, data2_bands) + # no extra bands are allowed when the same bands are defined + .check_that(length(bands_intersects) == 0) + # same sensor is required when bands with the same names are defined + .check_that(all(.cube_sensor(data1) %in% .cube_sensor(data2))) + } + # return + has_same_bands } -.merge_tiles_overlaps <- function(data1, data2) { +.merge_get_common_tiles <- function(data1, data2) { # Extract common tiles d1_tiles <- .cube_tiles(data1) d2_tiles <- .cube_tiles(data2) @@ -15,56 +26,45 @@ intersect(d1_tiles, d2_tiles) } -.merge_adjust_timeline <- function(data1, data2) { - # reference timeline - reference_tl <- .cube_timeline(data1)[[1]] - # Adjust dates / bands - 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"]] <- reference_tl - 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_check_bands_intersects <- function(data1, data2) { - # Extract band intersects - bands_intersects <- .merge_bands_intersects(data1, data2) - # Check if there are intersects - .check_that(length(bands_intersects) >= 1) -} - -.merge_check_band_sensor <- function(data1, data2) { - # Extract band intersects - bands_intersects <- .merge_bands_intersects(data1, data2) - # If has overlaps, the sensor must be the same - if (length(bands_intersects) >= 1) { - .check_that(data1[["sensor"]] == data2[["sensor"]]) +# ---- Adjust timeline strategies strategies ---- +.merge_adjust_timeline_strategy_zipper <- function(t1, t2) { + # define vector to store overlapping dates + t_overlap <- c() + # define the size of the `for` - size of the reference time-series + ts_reference_len <- length(t1) - 1 + # search the overlapping dates + for (idx in seq_len(ts_reference_len)) { + # reference interval (`t1`) + reference_interval <- t1[idx: (idx + 1)] + # verify which dates are in the reference interval + t2_in_interval <- t2 >= t1[idx] & t2 <= t1[idx + 1] + # get the interval dates + t2_interval_dates <- t2[t2_in_interval] + # if have interval, process them + if (.has(t2_interval_dates)) { + # if all t2 dates are in the interval, just save them + if (all(reference_interval %in% t2_interval_dates)) { + t_overlap <- c(t_overlap, t2_interval_dates) + } else { + # if not, save the reference interval and the min value of + # the t2 interval dates. + # this ensure there are not two dates in the same interval + t_overlap <- c( + t_overlap, # dates storage + reference_interval, # current interval + min(t2_interval_dates) # min t2 interval date + ) + } + } } -} - -.merge_timeline_has_same_length <- function(data1, data2) { - length(.cube_timeline(data1)[[1]]) == length(.cube_timeline(data2)[[1]]) + # sort and remove duplicated values + sort(unique(as.Date(t_overlap))) } # ---- Merge strategies ---- -.merge_strategy_file <- function(data1, data2, adjust_timeline) { - # adjust second cube timeline, based on the first cube - if (adjust_timeline) { - data2 <- .merge_adjust_timeline(data1, data2) - } +.merge_strategy_file <- function(data1, data2) { # extract tiles - tiles <- .merge_tiles_overlaps(data1, data2) + tiles <- .merge_get_common_tiles(data1, data2) # merge cubes .map_dfr(tiles, function(tile) { # select data in the selected tile @@ -92,150 +92,83 @@ }) } -.merge_strategy_bind <- function(data1, data2, adjust_timeline) { - # Adjust second cube timeline, based on the first cube - if (adjust_timeline) { - data2 <- .merge_adjust_timeline(data1, data2) - } +.merge_strategy_bind <- function(data1, data2) { # Merge dplyr::bind_rows(data1, data2) } -# ---- Regular cubes ---- -.merge_regular_check_timeline_overlaps <- function(data1, data2) { - # extract timelines - d1_tl <- .cube_timeline(data1) - d2_tl <- .cube_timeline(data2) - # check overlaps - slider::slide2_vec(d1_tl, d2_tl, function(x, y) { - x <- .dissolve(x) - y <- .dissolve(y) - - .check_that(length(.timeline_has_overlap(x, y)) >= 1) - }) -} - -.merge_regular_check_periods <- function(data1, data2) { - # get cubes timeline - d1_tl <- unique(as.Date(.cube_timeline(data1)[[1]])) - d2_tl <- unique(as.Date(.cube_timeline(data2)[[1]])) - # get intervals - d1_period <- as.integer( - lubridate::as.period(lubridate::int_diff(d1_tl)), "days" - ) - d2_period <- as.integer( - lubridate::as.period(lubridate::int_diff(d2_tl)), "days" - ) - # pre-condition - are periods regular? - .check_that( - length(unique(d1_period)) == 1 && length(unique(d2_period)) == 1 - ) - # pre-condition - Do cubes have the same periods? - .check_that( - unique(d1_period) == unique(d2_period) - ) -} - -.merge_regular_cube <- function(data1, data2) { +# ---- Merge operations - Densify cube ---- +.merge_cube_densify <- function(data1, data2) { # get tile overlaps - tiles_overlaps <- .merge_tiles_overlaps(data1, data2) - # pre-condition - timelines overlaps - # in case of regular cube it is assumed the timeline must overlap - # to avoid the creation of inconsistent / irregular cubes - .merge_regular_check_timeline_overlaps(data1, data2) - # pre-condition - timelines with same period - .merge_regular_check_periods(data1, data2) - # pre-condition - equal bands must be from the same sensor - # bands with the same name, must be from the same sensor to avoid confusion - .merge_check_band_sensor(data1, data2) - # pre-condition - bands must intersect if more then 1 tile is provided - if (length(tiles_overlaps) > 1) { - .merge_check_bands_intersects(data1, data2) - } - # Extract band intersects - bands_intersects <- .merge_bands_intersects(data1, data2) - # Use only intersect bands - data1 <- .select_raster_bands(data1, bands_intersects) - data2 <- .select_raster_bands(data2, bands_intersects) - # Check timeline consistency - if (!.merge_timeline_has_same_length(data1, data2)) { - # TODO: warning avisando o usuário que os cubos tem timelines - # com lengths diferentes - t1 <- .cube_timeline(data1)[[1]] - t2 <- .cube_timeline(data2)[[1]] - - if (length(t1) > length(t2)) { - ref <- t1[t1 >= min(t2) & t1 <= max(t2)] - } else { - ref <- t2[t2 >= min(t1) & t2 <= max(t1)] - } - - .check_that(.has(ref)) - - data1 <- .cube_filter_interval( - data1, start_date = min(ref), end_date = max(ref) - ) - - data2 <- .cube_filter_interval( - data2, start_date = min(ref), end_date = max(ref) - ) - - if (length(.cube_timeline(data1)) != length(.cube_timeline(data2))) { - min_length <- min(c(length(.cube_timeline(data1)), - length(.cube_timeline(data2)))) - - data1 <- .cube_filter_dates( - data1, .cube_timeline(data1)[[1]][seq_len(min_length)] - ) - data2 <- .cube_filter_dates( - data2, .cube_timeline(data2)[[1]][seq_len(min_length)] - ) - } - } + common_tiles <- .merge_get_common_tiles(data1, data2) # define the strategy (default - merge tiles) merge_strategy <- NULL - # case: same tiles, merge file info - if (.has(tiles_overlaps)) { + # case 1: same tiles, merge file info + if (.has(common_tiles)) { merge_strategy <- .merge_strategy_file - # case 2: different tiles, merge tile rows } else { + # case 2: different tiles, merge cube rows merge_strategy <- .merge_strategy_bind } # merge - merge_strategy(data1, data2, TRUE) + merged_cube <- merge_strategy(data1, data2) + # include `combined` in cubes merged with bind + if (!.has(common_tiles)) { + class(merged_cube) <- c("combined_cube", class(data1)) + } + # return + merged_cube } -# ---- Irregular cubes ---- -.merge_irregular_cube <- function(data1, data2) { - # get tile overlaps - tiles_overlaps <- .merge_tiles_overlaps(data1, data2) - # pre-condition - equal bands from the same sensor - # bands with the same name, must be from the same sensor to avoid confusion - .merge_check_band_sensor(data1, data2) - # pre-condition - bands must intersect if more then 1 tile is provided - if (length(tiles_overlaps) > 1) { - .merge_check_bands_intersects(data1, data2) - } - # Extract band intersects - bands_intersects <- .merge_bands_intersects(data1, data2) - # Use only intersect bands - data1 <- .select_raster_bands(data1, bands_intersects) - data2 <- .select_raster_bands(data2, bands_intersects) - # define the strategy (default - merge tiles) - merge_strategy <- NULL - # case: same tiles, merge file info - if (.has(tiles_overlaps)) { - merged_cube <- .merge_strategy_file(data1, data2, FALSE) - # case 2: different tiles, merge tile rows - } else { - merged_cube <- .merge_strategy_bind(data1, data2, FALSE) +# ---- Merge operations - Temporal overlaps ---- +.merge_cube_compactify <- function(data1, data2) { + # extract tiles + tiles <- .merge_get_common_tiles(data1, data2) + if (!.has(tiles)) { + # 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) + reference_timeline <- as.Date(unlist(.cube_timeline(data1))) + # based on the global timeline, cut the timeline of all tiles in cube 2 + merged_cube <- .cube_foreach_tile(data2, function(row) { + # get row timeline + row_timeline <- .tile_timeline(row) + # search overlaps between the reference timeline and row timeline + t_overlap <- .merge_adjust_timeline_strategy_zipper( + t1 = reference_timeline, + t2 = row_timeline + ) + # cut the timeline + .cube_filter_dates(row, t_overlap) + }) + # as there is no tile reference, merge using `bind` strategy (cube row) + merged_cube <- .merge_strategy_bind(data1, merged_cube) + # assign `combined cube` class, meaning the cube is a combination of + # cubes that contains different timelines in different tiles class(merged_cube) <- c("combined_cube", class(data1)) + } else { + # align timeline tile by tile. + merged_cube <- .map_dfr(tiles, function(tile) { + # get tiles + tile1 <- .cube_filter_tiles(data1, tile) + tile2 <- .cube_filter_tiles(data2, tile) + # get tile timelines + ts1 <- .tile_timeline(tile1) + ts2 <- .tile_timeline(tile2) + # adjust timeline using zipper strategy + ts_overlap <- .merge_adjust_timeline_strategy_zipper(ts1, ts2) + # filter cubes in the overlapping dates + tile1 <- .cube_filter_dates(tile1, ts_overlap) + tile2 <- .cube_filter_dates(tile2, ts_overlap) + # merge by file + .merge_strategy_file(tile1, tile2) + }) } # return - return(merged_cube) + merged_cube } -# ---- Special case: DEM Cube ---- +# ---- Merge operation: Special case - DEM Cube ---- .merge_dem_cube <- function(data1, data2) { # define cubes dem_cube <- data1 @@ -264,9 +197,10 @@ tile_dem }) # merge cubes and return - .merge_strategy_file(other_cube, dem_cube, FALSE) + .merge_strategy_file(other_cube, dem_cube) } +# ---- Merge operation: Special case - HLS Cube ---- .merge_hls_cube <- function(data1, data2) { if ((.cube_collection(data1) == "HLSS30" || .cube_collection(data2) == "HLSS30")) { @@ -274,5 +208,5 @@ } # merge cubes and return - .merge_strategy_file(data1, data2, FALSE) + .merge_strategy_file(data1, data2) } diff --git a/R/sits_merge.R b/R/sits_merge.R index 73c4f1738..825f710dd 100644 --- a/R/sits_merge.R +++ b/R/sits_merge.R @@ -102,30 +102,36 @@ sits_merge.raster_cube <- function(data1, data2, ...) { # pre-condition - check cube type .check_is_raster_cube(data1) .check_is_raster_cube(data2) - # Define merged cube + # pre-condition - cube rows has same bands + .check_cube_row_same_bands(data1) + .check_cube_row_same_bands(data2) + # define merged cube merged_cube <- NULL - # Special case: DEM cube + # special case: DEM cube is_dem_cube <- any(inherits(data1, "dem_cube"), inherits(data2, "dem_cube")) if (is_dem_cube) { return(.merge_dem_cube(data1, data2)) } - # Special case: HLS cube + # special case: HLS cube is_hls_cube <- all(inherits(data1, "hls_cube"), inherits(data2, "hls_cube")) if (is_hls_cube) { return(.merge_hls_cube(data1, data2)) } - # Check if cube is regular - is_regular <- all(.cube_is_regular(data1), .cube_is_regular(data2)) - has_unique_period <- all( - .cube_has_unique_period(data1), .cube_has_unique_period(data2) - ) - if (is_regular && has_unique_period) { - # Regular cube case - merged_cube <- .merge_regular_cube(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 { - # Irregular cube case - merged_cube <- .merge_irregular_cube(data1, data2) + # rule 2: if the bands are different and their timelines should be + # compatible, the bands are joined. The resulting timeline is the one + # from the first cube. + merged_cube <- .merge_cube_compactify(data1, data2) } + # empty results are not possible, meaning the input data is wrong + .check_that(nrow(merged_cube) > 0) + # return merged_cube } diff --git a/tests/testthat/test-merge.R b/tests/testthat/test-merge.R index d73c458bc..0035f6e08 100644 --- a/tests/testthat/test-merge.R +++ b/tests/testthat/test-merge.R @@ -1,10 +1,9 @@ -test_that("sits_merge - irregular cubes with same bands and tile", { +test_that("sits_merge - same bands case - equal tiles", { # Test case: If the bands are the same, the cube will have the combined # timeline of both cubes. This is useful to merge data from the same sensors # from different satellites (e.g, Sentinel-2A with Sentinel-2B). - # For irregular cubes, all dates are returned. - # Test 1a: Single tile with different time period (irregular cube) + # Test 1: Single tile with different time period s2a_cube <- .try( { sits_cube( @@ -27,7 +26,7 @@ test_that("sits_merge - irregular cubes with same bands and tile", { collection = "GA_S2BM_ARD_3", bands = c("BLUE"), tiles = c("53HQE"), - start_date = "2019-04-01", + start_date = "2019-03-01", end_date = "2019-06-10", progress = FALSE ) @@ -52,55 +51,7 @@ test_that("sits_merge - irregular cubes with same bands and tile", { expect_equal(merged_cube[["xmax"]][[1]], .raster_xmax(r), tolerance = 1) expect_equal(merged_cube[["xmin"]][[1]], .raster_xmin(r), tolerance = 1) - # Test 1b: Single tile with different time period (irregular cube) - s2_cube <- .try( - { - sits_cube( - source = "DEAFRICA", - collection = "SENTINEL-2-L2A", - bands = c("B02"), - tiles = c("36NWJ"), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL - ) - - s1_cube <- .try( - { - sits_cube( - source = "DEAFRICA", - collection = "SENTINEL-1-RTC", - bands = c("VV"), - orbit = "ascending", - tiles = c("36NWJ"), - start_date = "2019-04-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL - ) - - testthat::skip_if(purrr::is_null(c(s2_cube, s1_cube)), - message = "DEAFRICA is not accessible" - ) - - merged_cube <- sits_merge(s2_cube, s1_cube) - - expect_true(inherits(merged_cube, "combined_cube")) - expect_equal( - length(merged_cube[["tile"]]), - length(s2_cube[["tile"]]) + length(s1_cube[["tile"]]) - ) - expect_equal(suppressWarnings(length(sits_timeline(merged_cube))), 5) - expect_equal( - unique(slider::slide_chr(merged_cube, .tile_bands)), c("B02", "VV") - ) - - # Test 2: Multiple tiles with different time period (irregular cube) + # Test 2: Multiple tiles with different time period s2a_cube <- .try( { sits_cube( @@ -146,43 +97,91 @@ test_that("sits_merge - irregular cubes with same bands and tile", { r <- .raster_open_rast(.tile_path(merged_cube)) expect_equal(merged_cube[["xmax"]][[1]], .raster_xmax(r), tolerance = 1) expect_equal(merged_cube[["xmin"]][[1]], .raster_xmin(r), tolerance = 1) + + # Test 3: Tiles with same time period + modis_cube_a <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + modis_cube_b <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-03-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), + message = "BDC is not accessible" + ) + + merged_cube <- sits_merge(modis_cube_a, modis_cube_b) + + expect_equal(length(sits_timeline(merged_cube)), 11) + expect_equal(sits_bands(merged_cube), "NDVI") + expect_equal(merged_cube[["tile"]], "013011") }) -test_that("sits_merge - irregular cubes with same bands and different tile", { +test_that("sits_merge - same bands case - different tiles", { # Test case: If the bands are the same, the cube will have the combined # timeline of both cubes. This is useful to merge data from the same sensors # from different satellites (e.g, Sentinel-2A with Sentinel-2B). - # For irregular cubes, all dates are returned. - # Test 1: Different tiles with different time period (irregular cube) - s2a_cube <- .try( - { - sits_cube( - source = "DEAUSTRALIA", - collection = "ga_s2am_ard_3", - bands = c("BLUE"), - tiles = c("53HQE"), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL - ) - - s2b_cube <- .try( - { - sits_cube( - source = "DEAUSTRALIA", - collection = "GA_S2BM_ARD_3", - bands = c("BLUE"), - tiles = c("53JQF"), - start_date = "2019-04-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL + # Test 1: Aligned timelines + s2a_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2b_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53JQF"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) ) testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), @@ -193,69 +192,111 @@ test_that("sits_merge - irregular cubes with same bands and different tile", { expect_true(inherits(merged_cube, "combined_cube")) expect_equal(suppressWarnings(length(sits_timeline(merged_cube))), 2) + + # Test 2: Overlapping timelines + modis_cube_a <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22LBH"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + modis_cube_b <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + bands = c("NDVI"), + roi = sits_tiles_to_roi("22KGA"), + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), + message = "BDC is not accessible" + ) + + merged_cube <- sits_merge(modis_cube_a, modis_cube_b) + expect_equal(length(sits_timeline(merged_cube)), 2) + expect_equal(sits_bands(merged_cube), "NDVI") + expect_equal(merged_cube[["tile"]], c("012010", "013011")) }) -test_that("sits_merge - irregular cubes with different bands and tile", { +test_that("sits_merge - different bands case - equal tiles", { # Test case: if the bands are different and their timelines should be # compatible, the bands are joined. The resulting timeline is the one from # the first cube. This is useful to merge data from different sensors # (e.g, Sentinel-1 with Sentinel-2). - # For irregular cubes, all dates are returned. - s2a_cube <- .try( - { - sits_cube( - source = "DEAUSTRALIA", - collection = "ga_s2am_ard_3", - bands = c("BLUE", "RED"), - tiles = c("53HQE"), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL - ) - - s2b_cube <- .try( - { - sits_cube( - source = "DEAUSTRALIA", - collection = "GA_S2BM_ARD_3", - bands = c("BLUE"), - tiles = c("53HQE", "53JQF"), - start_date = "2019-04-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL + # Test 1a: Aligned timelines + s2a_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("RED"), + tiles = c("53HQE"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2b_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) ) testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), message = "DEAustralia is not accessible" ) + # timeline created with the zipper algorithm merged_cube <- sits_merge(s2a_cube, s2b_cube) - expect_equal(sits_bands(merged_cube), "BLUE") + expect_equal(length(sits_timeline(merged_cube)), 21) + expect_equal(sits_bands(merged_cube), c("BLUE", "RED")) expect_equal(merged_cube[["tile"]], "53HQE") -}) -test_that("sits_merge - regular cubes with same bands and tile", { - # Test case: If the bands are the same, the cube will have the combined - # timeline of both cubes. This is useful to merge data from the same sensors - # from different satellites (e.g, Sentinel-2A with Sentinel-2B). - # For regular cubes, when timeline has the same length, use them. Otherwise, - # use as timeline the intersect between timelines. - - # Test 1: Tiles with same time period (regular cube) - modis_cube_a <- .try( + # Test 1b: Aligned timelines + s2_cube_a <- .try( { sits_cube( source = "BDC", - collection = "MOD13Q1-6.1", - bands = c("NDVI"), - roi = sits_tiles_to_roi("22KGA"), + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LMR")), start_date = "2019-01-01", end_date = "2019-04-01", progress = FALSE @@ -264,200 +305,348 @@ test_that("sits_merge - regular cubes with same bands and tile", { .default = NULL ) - modis_cube_b <- .try( - { - sits_cube( - source = "BDC", - collection = "MOD13Q1-6.1", - bands = c("NDVI"), - roi = sits_tiles_to_roi("22KGA"), - start_date = "2019-03-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL - ) - - testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), - message = "BDC is not accessible" + s2_cube_b <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B03"), + roi = sits_tiles_to_roi(c("20LMR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) ) - merged_cube <- sits_merge(modis_cube_a, modis_cube_b) - - expect_equal(length(sits_timeline(merged_cube)), 2) - expect_equal(sits_bands(merged_cube), "NDVI") - expect_equal(merged_cube[["tile"]], "013011") - - # Test 2: no time-series overlaps (regular cube) - modis_cube_a <- .try( - { - sits_cube( - source = "BDC", - collection = "MOD13Q1-6.1", - bands = c("NDVI"), - roi = sits_tiles_to_roi("22KGA"), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL - ) + merged_cube <- sits_merge(s2_cube_a, s2_cube_b) + expect_equal(sits_timeline(merged_cube), sits_timeline(s2_cube_a)) + expect_equal(nrow(merged_cube), 4) - modis_cube_b <- .try( - { - sits_cube( - source = "BDC", - collection = "MOD13Q1-6.1", - bands = c("NDVI"), - roi = sits_tiles_to_roi("22KGA"), - start_date = "2019-04-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL + # Test 2a: Overlapping timelines + s2a_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("RED"), + tiles = c("53HQE"), + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2b_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-03-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) ) - testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), - message = "BDC is not accessible" + testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), + message = "DEAustralia is not accessible" ) - expect_error(sits_merge(modis_cube_a, modis_cube_b)) -}) - -test_that("sits_merge - regular cubes with same bands and different tile", { - # Test case: If the bands are the same, the cube will have the combined - # timeline of both cubes. This is useful to merge data from the same sensors - # from different satellites (e.g, Sentinel-2A with Sentinel-2B). - # For regular cubes, then timeline has the same length, use them. Otherwise, - # use as timeline the intersect between timelines. + merged_cube <- sits_merge(s2a_cube, s2b_cube) + # timeline created with the zipper algorithm + expect_equal(length(sits_timeline(merged_cube)), 30) + expect_equal(sits_bands(merged_cube), c("BLUE", "RED")) + expect_equal(merged_cube[["tile"]], "53HQE") - # Test 1: Different tiles - modis_cube_a <- .try( - { - sits_cube( - source = "BDC", - collection = "MOD13Q1-6.1", - bands = c("NDVI"), - roi = sits_tiles_to_roi("22LBH"), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL + # Test 2b: Overlapping timelines + rainfall <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAFRICA", + collection = "RAINFALL-CHIRPS-MONTHLY", + roi = sits_tiles_to_roi("38LQK"), + start_date = "2022-01-01", + end_date = "2022-06-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2b_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAFRICA", + collection = "SENTINEL-2-L2A", + bands = c("B02"), + tiles = c("38LQK"), + start_date = "2022-01-01", + end_date = "2022-06-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + testthat::skip_if(purrr::is_null(c(rainfall, s2b_cube)), + message = "DEAustralia is not accessible" ) - modis_cube_b <- .try( - { - sits_cube( - source = "BDC", - collection = "MOD13Q1-6.1", - bands = c("NDVI"), - roi = sits_tiles_to_roi("22KGA"), - start_date = "2019-02-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL + # merge + merged_cube <- sits_merge(rainfall, s2b_cube) + # test + expect_true("combined_cube" %in% class(merged_cube)) + # test timeline compatibility + merged_tl <- suppressWarnings(unname(sits_timeline(merged_cube))) + # result timeline must be compatible (cube 1 is the reference in this case) + expect_true( + min(merged_tl[[2]]) >= min(merged_tl[[1]]) & + max(merged_tl[[2]]) <= max(merged_tl[[2]]) + ) + + # Test 3: Different timelines + s2a_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("RED"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2b_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-04-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) ) - testthat::skip_if(purrr::is_null(c(modis_cube_a, modis_cube_b)), - message = "BDC is not accessible" + testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), + message = "DEAustralia is not accessible" ) - merged_cube <- sits_merge(modis_cube_a, modis_cube_b) - - expect_equal(length(sits_timeline(merged_cube)), 4) - expect_equal( - sits_timeline(modis_cube_b)[seq_len(4)], sits_timeline(merged_cube) + merged_cube <- expect_error(sits_merge(s2a_cube, s2b_cube)) + + # Test 4: Different sensor with same timeline + s2_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + bands = c("B02"), + tiles = c("19LEF"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s1_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-1-RTC", + bands = c("VV"), + tiles = c("19LEF"), + orbit = "descending", + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) ) - expect_equal(sits_bands(merged_cube), "NDVI") - expect_equal(merged_cube[["tile"]], c("012010", "013011")) - # Test 2: Tile variation in one of the cubes - s2_cube_a <- .try( - { - sits_cube( - source = "BDC", - collection = "SENTINEL-2-16D", - bands = c("B02"), - roi = sits_tiles_to_roi(c("20LMR", "20LNR")), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL + testthat::skip_if(purrr::is_null(c(s2_cube, s1_cube)), + message = "MPC is not accessible" ) - s2_cube_b <- .try( - { - sits_cube( - source = "BDC", - collection = "SENTINEL-2-16D", - bands = c("B02"), - roi = sits_tiles_to_roi(c("20LMR")), - start_date = "2019-04-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL + # merge + merged_cube <- sits_merge(s2_cube, s1_cube) + expect_equal(sits_bands(merged_cube[1,]), "B02") + expect_equal(sits_bands(merged_cube[2,]), "VV") + expect_equal(merged_cube[["tile"]], c("19LEF", "NoTilingSystem")) + expect_true("combined_cube" %in% class(merged_cube)) + # test timeline compatibility + merged_tl <- suppressWarnings(unname(sits_timeline(merged_cube))) + # result timeline must be compatible (cube 1 is the reference in this case) + expect_true( + min(merged_tl[[2]]) >= min(merged_tl[[1]]) & + max(merged_tl[[2]]) <= max(merged_tl[[2]]) ) - - merged_cube <- sits_merge(s2_cube_a, s2_cube_b) - - expect_equal(sits_timeline(merged_cube), sits_timeline(s2_cube_a)) - expect_equal(nrow(merged_cube), 4) }) -test_that("sits_merge - regular cubes with different bands and tile", { +test_that("sits_merge - different bands case - different tiles", { # Test case: if the bands are different and their timelines should be # compatible, the bands are joined. The resulting timeline is the one from # the first cube. This is useful to merge data from different sensors # (e.g, Sentinel-1 with Sentinel-2). - # For regular cubes, then timeline has the same length, use them. Otherwise, - # use as timeline the intersect between timelines. - s2_cube_a <- .try( - { - sits_cube( - source = "BDC", - collection = "SENTINEL-2-16D", - bands = c("B02"), - roi = sits_tiles_to_roi(c("20LMR", "20LNR")), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL + # Test 1: Aligned timelines + s2_cube_a <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LNR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2_cube_b <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B03"), + roi = sits_tiles_to_roi(c("20LMR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) ) - - s2_cube_b <- .try( - { - sits_cube( - source = "BDC", - collection = "SENTINEL-2-16D", - bands = c("B02", "B03"), - roi = sits_tiles_to_roi(c("20LMR")), - start_date = "2019-04-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL - ) - + # merge merged_cube <- sits_merge(s2_cube_a, s2_cube_b) - + # test expect_equal(sits_timeline(merged_cube), sits_timeline(s2_cube_a)) - expect_equal(nrow(merged_cube), 4) - expect_equal(sits_bands(merged_cube), "B02") + expect_equal(nrow(merged_cube), 2) + expect_equal(sits_bands(merged_cube), c("B02", "B03")) + # as we have intersecting tiles with the same bands, they are merged! + expect_equal(sits_bands(merged_cube[1,]), c("B02", "B03")) + expect_equal(sits_bands(merged_cube[2,]), c("B02", "B03")) + + # Test 2: Overlapping timelines + s2_cube_a <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LNR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2_cube_b <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B03"), + roi = sits_tiles_to_roi(c("20LMR")), + start_date = "2019-02-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + # merge + merged_cube <- sits_merge(s2_cube_a, s2_cube_b) + # test + expect_equal(nrow(merged_cube), 2) + expect_equal(merged_cube[["tile"]], c("013014", "013015")) + expect_equal(sits_bands(merged_cube), c("B02", "B03")) + # as we have intersecting tiles with the same bands, they are merged! + expect_equal(sits_bands(merged_cube[1,]), c("B02", "B03")) + expect_equal(sits_bands(merged_cube[2,]), c("B02", "B03")) + + # Test 3: Different timelines + s2_cube_a <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LNR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2_cube_b <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B03"), + roi = sits_tiles_to_roi(c("20LMR")), + start_date = "2019-05-01", + end_date = "2019-06-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + # merge and test + expect_error(sits_merge(s2_cube_a, s2_cube_b)) }) test_that("sits_merge - regularize combined cubes", { @@ -465,34 +654,38 @@ test_that("sits_merge - regularize combined cubes", { output_dir <- paste0(tempdir(), "/merge-reg-1") dir.create(output_dir, showWarnings = FALSE) - s2a_cube <- .try( - { - sits_cube( - source = "DEAUSTRALIA", - collection = "ga_s2am_ard_3", - bands = c("BLUE"), - tiles = c("53HQE"), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL - ) - - s2b_cube <- .try( - { - sits_cube( - source = "DEAUSTRALIA", - collection = "GA_S2BM_ARD_3", - bands = c("BLUE"), - tiles = c("53JQF"), - start_date = "2019-02-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL + s2a_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "ga_s2am_ard_3", + bands = c("BLUE"), + tiles = c("53HQE"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s2b_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "DEAUSTRALIA", + collection = "GA_S2BM_ARD_3", + bands = c("BLUE"), + tiles = c("53JQF"), + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) ) testthat::skip_if(purrr::is_null(c(s2a_cube, s2b_cube)), @@ -507,7 +700,8 @@ test_that("sits_merge - regularize combined cubes", { cube = merged_cube, period = "P8D", res = 720, - output_dir = output_dir + output_dir = output_dir, + progress = FALSE ) # test @@ -520,35 +714,39 @@ test_that("sits_merge - regularize combined cubes", { output_dir <- paste0(tempdir(), "/merge-reg-2") dir.create(output_dir, showWarnings = FALSE) - s2_cube <- .try( - { - sits_cube( - source = "MPC", - collection = "SENTINEL-2-L2A", - bands = c("B02"), - tiles = c("19LEF"), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL - ) - - s1_cube <- .try( - { - sits_cube( - source = "MPC", - collection = "SENTINEL-1-RTC", - bands = c("VV"), - tiles = c("19LEF"), - orbit = "descending", - start_date = "2019-02-01", - end_date = "2019-06-10", - progress = FALSE - ) - }, - .default = NULL + s2_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + bands = c("B02"), + tiles = c("19LEF"), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + s1_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-1-RTC", + bands = c("VV"), + tiles = c("19LEF"), + orbit = "descending", + start_date = "2019-02-01", + end_date = "2019-06-10", + progress = FALSE + ) + }, + .default = NULL + ) ) testthat::skip_if(purrr::is_null(c(s2_cube, s1_cube)), @@ -563,7 +761,8 @@ test_that("sits_merge - regularize combined cubes", { cube = merged_cube, period = "P8D", res = 720, - output_dir = output_dir + output_dir = output_dir, + progress = FALSE ) # test @@ -622,14 +821,25 @@ test_that("sits_merge - special case - dem cube", { # create S2 cube s2_dir <- paste0(tempdir(), "/s2") dir.create(s2_dir, showWarnings = FALSE) - s2_cube <- sits_cube( - source = "MPC", - collection = "SENTINEL-2-L2A", - tiles = "19HBA", - bands = c("B04", "B8A", "B12", "CLOUD"), - start_date = "2021-01-01", - end_date = "2021-03-31", - progress = FALSE + s2_cube <- suppressWarnings( + .try( + { + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + tiles = "19HBA", + bands = c("B04", "B8A", "B12", "CLOUD"), + start_date = "2021-01-01", + end_date = "2021-03-31", + progress = FALSE + ) + }, + .default = NULL + ) + ) + + testthat::skip_if(purrr::is_null(s2_cube), + message = "MPC is not accessible" ) s2_cube_reg <- sits_regularize( @@ -643,12 +853,21 @@ test_that("sits_merge - special case - dem cube", { # create DEM cube dem_dir <- paste0(tempdir(), "/dem") dir.create(dem_dir, showWarnings = FALSE) - dem_cube <- sits_cube( - source = "MPC", - collection = "COP-DEM-GLO-30", - bands = "ELEVATION", - tiles = "19HBA", - progress = FALSE + dem_cube <- .try( + { + sits_cube( + source = "MPC", + collection = "COP-DEM-GLO-30", + bands = "ELEVATION", + tiles = "19HBA", + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(dem_cube), + message = "MPC is not accessible" ) dem_cube_reg <- sits_regularize( @@ -675,24 +894,38 @@ test_that("sits_merge - special case - hls cube", { lon_max = -45.0840, lat_max = -23.6178 ) - hls_cube_s2 <- sits_cube( - source = "HLS", - collection = "HLSS30", - roi = roi, - bands = c("BLUE", "GREEN", "RED", "CLOUD"), - start_date = as.Date("2020-06-01"), - end_date = as.Date("2020-09-01"), - progress = FALSE + hls_cube_s2 <- .try( + { + sits_cube( + source = "HLS", + collection = "HLSS30", + roi = roi, + bands = c("BLUE", "GREEN", "RED", "CLOUD"), + start_date = as.Date("2020-06-01"), + end_date = as.Date("2020-09-01"), + progress = FALSE + ) + }, + .default = NULL ) - hls_cube_l8 <- sits_cube( - source = "HLS", - collection = "HLSL30", - roi = roi, - bands = c("BLUE", "GREEN", "RED", "CLOUD"), - start_date = as.Date("2020-06-01"), - end_date = as.Date("2020-09-01"), - progress = FALSE + hls_cube_l8 <- .try( + { + sits_cube( + source = "HLS", + collection = "HLSL30", + roi = roi, + bands = c("BLUE", "GREEN", "RED", "CLOUD"), + start_date = as.Date("2020-06-01"), + end_date = as.Date("2020-09-01"), + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(c(hls_cube_s2, hls_cube_l8)), + message = "HLS is not accessible" ) # merge @@ -700,5 +933,5 @@ test_that("sits_merge - special case - hls cube", { # test expect_equal(length(sits_timeline(merged_cube)), 19) - expect_equal(sits_bands(merged_cube), c("BLUE", "GREEN", "RED", "CLOUD")) + expect_equal(sits_bands(merged_cube), c("BLUE", "CLOUD", "GREEN", "RED")) }) From afac3130f36844a3d2e2910be1cbfbfe973ed53c Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Thu, 28 Nov 2024 06:10:23 -0300 Subject: [PATCH 09/19] update merge validations --- tests/testthat/test-cube-deaustralia.R | 2 +- tests/testthat/test-cube.R | 9 ++- tests/testthat/test-merge.R | 79 +++++++++++++++----------- 3 files changed, 52 insertions(+), 38 deletions(-) diff --git a/tests/testthat/test-cube-deaustralia.R b/tests/testthat/test-cube-deaustralia.R index b66119a29..ba23f7c4b 100644 --- a/tests/testthat/test-cube-deaustralia.R +++ b/tests/testthat/test-cube-deaustralia.R @@ -433,7 +433,7 @@ test_that( sentinel_cube <- sits_merge(s2a_cube, s2b_cube) - expect_true(all(sits_bands(sentinel_cube) %in% c("BLUE"))) + expect_true(all(sits_bands(sentinel_cube) %in% c("BLUE", "NIR-2", "RED"))) expect_equal(nrow(sentinel_cube), 2) r <- .raster_open_rast(.tile_path(sentinel_cube)) expect_equal(sentinel_cube[["xmax"]][[1]], .raster_xmax(r), tolerance = 1) diff --git a/tests/testthat/test-cube.R b/tests/testthat/test-cube.R index 3ce8a5c00..6acc1b2c6 100644 --- a/tests/testthat/test-cube.R +++ b/tests/testthat/test-cube.R @@ -176,12 +176,11 @@ test_that("Combining Sentinel-1 with Sentinel-2 cubes", { sits_bands(cube_merged) %in% c(sits_bands(s2_reg), sits_bands(s1_reg))) ) - testthat::expect_error( - sits_merge( - s2_cube, - s1_cube - ) + merged_cube <- sits_merge( + s2_cube, + s1_cube ) + expect_equal(nrow(merged_cube), 2) unlink(list.files(dir_images, pattern = ".tif", full.names = TRUE)) }) diff --git a/tests/testthat/test-merge.R b/tests/testthat/test-merge.R index 0035f6e08..3f44ad9e6 100644 --- a/tests/testthat/test-merge.R +++ b/tests/testthat/test-merge.R @@ -233,7 +233,7 @@ test_that("sits_merge - same bands case - different tiles", { ) merged_cube <- sits_merge(modis_cube_a, modis_cube_b) - expect_equal(length(sits_timeline(merged_cube)), 2) + expect_equal(suppressWarnings(length(sits_timeline(merged_cube))), 2) expect_equal(sits_bands(merged_cube), "NDVI") expect_equal(merged_cube[["tile"]], c("012010", "013011")) }) @@ -290,19 +290,21 @@ test_that("sits_merge - different bands case - equal tiles", { expect_equal(merged_cube[["tile"]], "53HQE") # Test 1b: Aligned timelines - s2_cube_a <- .try( - { - sits_cube( - source = "BDC", - collection = "SENTINEL-2-16D", - bands = c("B02"), - roi = sits_tiles_to_roi(c("20LMR")), - start_date = "2019-01-01", - end_date = "2019-04-01", - progress = FALSE - ) - }, - .default = NULL + s2_cube_a <- suppressWarnings( + .try( + { + sits_cube( + source = "BDC", + collection = "SENTINEL-2-16D", + bands = c("B02"), + roi = sits_tiles_to_roi(c("20LMR")), + start_date = "2019-01-01", + end_date = "2019-04-01", + progress = FALSE + ) + }, + .default = NULL + ) ) s2_cube_b <- suppressWarnings( @@ -696,12 +698,14 @@ test_that("sits_merge - regularize combined cubes", { merged_cube <- sits_merge(s2a_cube, s2b_cube) # regularize - regularized_cube <- sits_regularize( - cube = merged_cube, - period = "P8D", - res = 720, - output_dir = output_dir, - progress = FALSE + regularized_cube <- suppressWarnings( + sits_regularize( + cube = merged_cube, + period = "P8D", + res = 720, + output_dir = output_dir, + progress = FALSE + ) ) # test @@ -710,6 +714,8 @@ test_that("sits_merge - regularize combined cubes", { expect_equal(sits_bands(regularized_cube), "BLUE") expect_equal(.cube_xres(regularized_cube), 720) + unlink(output_dir, recursive = TRUE) + # Test 2: Different sensor output_dir <- paste0(tempdir(), "/merge-reg-2") dir.create(output_dir, showWarnings = FALSE) @@ -757,12 +763,14 @@ test_that("sits_merge - regularize combined cubes", { merged_cube <- sits_merge(s2_cube, s1_cube) # regularize - regularized_cube <- sits_regularize( - cube = merged_cube, - period = "P8D", - res = 720, - output_dir = output_dir, - progress = FALSE + regularized_cube <- suppressWarnings( + sits_regularize( + cube = merged_cube, + period = "P8D", + res = 720, + output_dir = output_dir, + progress = FALSE + ) ) # test @@ -770,6 +778,8 @@ test_that("sits_merge - regularize combined cubes", { expect_equal(length(sits_timeline(regularized_cube)), 7) expect_equal(sits_bands(regularized_cube), c("B02", "VV")) expect_equal(.cube_xres(regularized_cube), 720) + + unlink(output_dir, recursive = TRUE) }) test_that("sits_merge - cubes with different classes", { @@ -842,12 +852,14 @@ test_that("sits_merge - special case - dem cube", { message = "MPC is not accessible" ) - s2_cube_reg <- sits_regularize( - cube = s2_cube, - period = "P16D", - res = 720, - output_dir = s2_dir, - progress = FALSE + s2_cube_reg <- suppressWarnings( + sits_regularize( + cube = s2_cube, + period = "P16D", + res = 720, + output_dir = s2_dir, + progress = FALSE + ) ) # create DEM cube @@ -885,6 +897,9 @@ test_that("sits_merge - special case - dem cube", { # test expect_equal(nrow(merged_cube[["file_info"]][[1]]), 24) expect_equal(sits_bands(merged_cube), c("B04", "B12", "B8A", "ELEVATION")) + + unlink(s2_dir, recursive = TRUE) + unlink(dem_dir, recursive = TRUE) }) test_that("sits_merge - special case - hls cube", { From de66d0049965bff0dd48cb2efc0e89a4a489828b Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Thu, 28 Nov 2024 08:09:10 -0300 Subject: [PATCH 10/19] fix crop result files --- R/api_smooth.R | 34 ++++++++++++++++++++++++---------- 1 file changed, 24 insertions(+), 10 deletions(-) diff --git a/R/api_smooth.R b/R/api_smooth.R index 7ae976ea0..65b612281 100644 --- a/R/api_smooth.R +++ b/R/api_smooth.R @@ -99,9 +99,21 @@ # Return block file block_file }) + # Check if there is a exclusion_mask + # If exclusion_mask exists, blocks are merged to a different directory + # than output_dir, which is used to save the final cropped version + merge_out_file <- out_file + if (.has(exclusion_mask)) { + merge_out_file <- .file_derived_name( + tile = tile, + band = out_band, + version = version, + output_dir = file.path(output_dir, ".sits") + ) + } # Merge blocks into a new probs_cube tile probs_tile <- .tile_derived_merge_blocks( - file = out_file, + file = merge_out_file, band = band, labels = .tile_labels(tile), base_tile = tile, @@ -111,19 +123,21 @@ update_bbox = FALSE ) # Exclude masked areas - probs_tile <- .crop( - cube = probs_tile, - roi = exclusion_mask, - output_dir = output_dir, - multicores = 1, - overwrite = TRUE, - progress = FALSE - ) + if (.has(exclusion_mask)) { + probs_tile <- .crop( + cube = probs_tile, + roi = exclusion_mask, + output_dir = output_dir, + multicores = 1, + overwrite = TRUE, + progress = FALSE + ) + unlink(.fi_paths(.fi(probs_tile))) + } # Return probs tile probs_tile } - #---- Bayesian smoothing ---- #' @title Smooth probability cubes with spatial predictors #' @noRd From fb3ccc1398a650dbd35ac71ab129b615149abcb5 Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Thu, 28 Nov 2024 18:19:32 -0300 Subject: [PATCH 11/19] improve exclusion_mask tests --- tests/testthat/test-classification.R | 39 ++++++++++++---- tests/testthat/test-smooth.R | 69 ++++++++++++++++++++++++++++ 2 files changed, 100 insertions(+), 8 deletions(-) create mode 100644 tests/testthat/test-smooth.R diff --git a/tests/testthat/test-classification.R b/tests/testthat/test-classification.R index 3657818f4..5d9ac1960 100644 --- a/tests/testthat/test-classification.R +++ b/tests/testthat/test-classification.R @@ -110,23 +110,46 @@ test_that("Classify with exclusion mask", { multicores = 2, progress = FALSE ) + + # preparation - create exclusion mask + exclusion_mask <- sf::st_as_sfc( + x = sf::st_bbox(c( + xmin = -55.63478, + ymin = -11.63328, + xmax = -55.54080, + ymax = -11.56978 + ), + crs = "EPSG:4326" + ) + ) + + exclusion_mask <- sf::st_transform(exclusion_mask, .cube_crs(raster_cube)) + + # preparation - calculate centroid of the exclusion mask + exclusion_mask_centroid <- sf::st_centroid(exclusion_mask) + # preparation - create a random forest model rfor_model <- sits_train(samples_modis_ndvi, sits_rfor(num_trees = 40)) + # test classification with NA - class_map <- suppressWarnings( + probs_map <- suppressWarnings( sits_classify( data = raster_cube, ml_model = rfor_model, output_dir = tempdir(), - exclusion_mask = c( - xmin = -55.63478, - ymin = -11.63328, - xmax = -55.54080, - ymax = -11.56978 - ), + exclusion_mask = exclusion_mask, progress = FALSE ) ) - class_map_rst <- terra::rast(class_map[["file_info"]][[1]][["path"]]) + + # testing original data + probs_map_rst <- terra::rast(probs_map[["file_info"]][[1]][["path"]]) expect_true(anyNA(class_map_rst[])) + + probs_map_value <- terra::extract( + x = probs_map_rst, + y = terra::vect(exclusion_mask_centroid) + ) + + expect_true(any(is.na(probs_map_value))) }) diff --git a/tests/testthat/test-smooth.R b/tests/testthat/test-smooth.R new file mode 100644 index 000000000..e67a91b0e --- /dev/null +++ b/tests/testthat/test-smooth.R @@ -0,0 +1,69 @@ +test_that("Smoothing with exclusion mask", { + # preparation - create cube + data_dir <- system.file("extdata/raster/mod13q1", package = "sits") + raster_cube <- sits_cube( + source = "BDC", + collection = "MOD13Q1-6.1", + data_dir = data_dir, + tiles = "012010", + bands = "NDVI", + start_date = "2013-09-14", + end_date = "2014-08-29", + multicores = 2, + progress = FALSE + ) + + # preparation - create exclusion mask + exclusion_mask <- sf::st_as_sfc( + x = sf::st_bbox(c( + xmin = -6057482, + ymin = -1290723, + xmax = -6055209, + ymax = -1288406 + ), + crs = .cube_crs(raster_cube) + ) + ) + + # preparation - calculate centroid of the exclusion mask + exclusion_mask_centroid <- sf::st_centroid(exclusion_mask) + + # preparation - create a random forest model + rfor_model <- sits_train(samples_modis_ndvi, sits_rfor(num_trees = 40)) + + # test classification with NA + probs_map <- suppressWarnings( + sits_classify( + data = raster_cube, + ml_model = rfor_model, + output_dir = tempdir(), + progress = FALSE + ) + ) + + # smoth with exclusion mask + smooth_map <- sits_smooth( + cube = probs_map, + exclusion_mask = exclusion_mask, + output_dir = tempdir(), + multicores = 2 + ) + + # testing original data (no na) + probs_map_rst <- terra::rast(probs_map[["file_info"]][[1]][["path"]]) + probs_map_value <- terra::extract( + x = probs_map_rst, + y = terra::vect(exclusion_mask_centroid) + ) + + expect_false(any(is.na(probs_map_value))) + + # testing smooth data (with na) + smooth_map_rst <- terra::rast(smooth_map[["file_info"]][[1]][["path"]]) + smooth_map_value <- terra::extract( + x = smooth_map_rst, + y = terra::vect(exclusion_mask_centroid) + ) + + expect_true(any(is.na(smooth_map_value))) +}) From 8dae5945d6bbdb25de037f6b400f6300415aa08a Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Thu, 28 Nov 2024 18:20:16 -0300 Subject: [PATCH 12/19] fix crop in smooth --- R/api_smooth.R | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/R/api_smooth.R b/R/api_smooth.R index 65b612281..cd5390afe 100644 --- a/R/api_smooth.R +++ b/R/api_smooth.R @@ -106,7 +106,7 @@ if (.has(exclusion_mask)) { merge_out_file <- .file_derived_name( tile = tile, - band = out_band, + band = band, version = version, output_dir = file.path(output_dir, ".sits") ) @@ -124,15 +124,21 @@ ) # Exclude masked areas if (.has(exclusion_mask)) { - probs_tile <- .crop( + # crop + probs_tile_crop <- .crop( cube = probs_tile, roi = exclusion_mask, - output_dir = output_dir, + output_dir = output_dir, multicores = 1, overwrite = TRUE, progress = FALSE ) + + # delete old files unlink(.fi_paths(.fi(probs_tile))) + + # assign new cropped value in the old probs variable + probs_tile <- probs_tile_crop } # Return probs tile probs_tile From 9d3acc10e4a410cdfead23787e6e32ff7b809679 Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Thu, 28 Nov 2024 19:03:36 -0300 Subject: [PATCH 13/19] fix test-apply --- tests/testthat/test-apply.R | 21 ++++++++++++++++++++- 1 file changed, 20 insertions(+), 1 deletion(-) diff --git a/tests/testthat/test-apply.R b/tests/testthat/test-apply.R index 9f73fc823..4c9436f30 100644 --- a/tests/testthat/test-apply.R +++ b/tests/testthat/test-apply.R @@ -25,17 +25,21 @@ test_that("Testing normalized index generation", { if (!dir.exists(dir_images)) { suppressWarnings(dir.create(dir_images)) } + unlink(list.files(dir_images, pattern = "\\.tif$", full.names = TRUE )) - gc_cube <- sits_regularize( + + gc_cube <- suppressWarnings( + sits_regularize( cube = s2_cube, output_dir = dir_images, res = 160, period = "P1M", multicores = 2, progress = FALSE + ) ) gc_cube_new <- sits_apply(gc_cube, @@ -118,6 +122,8 @@ test_that("Testing normalized index generation", { values_evi2 <- .tibble_time_series(evi_tibble_2)$EVI values_evi2_new <- .tibble_time_series(evi_tibble_2)$EVI_NEW expect_equal(values_evi2, values_evi2_new, tolerance = 0.001) + + unlink(dir_images, recursive = TRUE) }) test_that("Testing non-normalized index generation", { @@ -134,6 +140,11 @@ test_that("Testing non-normalized index generation", { if (!dir.exists(dir_images)) { suppressWarnings(dir.create(dir_images)) } + unlink(list.files(dir_images, + pattern = "\\.tif$", + full.names = TRUE + )) + gc_cube_new <- sits_apply(cube, XYZ = 1 / NDVI * 0.25, normalized = FALSE, @@ -191,6 +202,8 @@ test_that("Testing non-normalized index generation", { values_xyz2 <- .tibble_time_series(xyz_tibble)$XYZ values_xyz_new <- .tibble_time_series(xyz_tibble_2)$XYZ_NEW expect_equal(values_xyz2, values_xyz_new, tolerance = 0.001) + + unlink(dir_images, recursive = TRUE) }) test_that("Kernel functions", { @@ -326,6 +339,11 @@ test_that("Error", { if (!dir.exists(output_dir)) { dir.create(output_dir) } + unlink(list.files(output_dir, + pattern = "\\.tif$", + full.names = TRUE + )) + Sys.setenv("SITS_DOCUMENTATION_MODE" = "FALSE") expect_warning({ cube_median <- sits_apply( @@ -347,4 +365,5 @@ test_that("Error", { ) expect_error(sits_apply(sinop_probs)) + unlink(output_dir, recursive = TRUE) }) From 93a117ab876b8c5c05e41703e8502979f0f80dbe Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Thu, 28 Nov 2024 19:34:19 -0300 Subject: [PATCH 14/19] fix invalid variable in classification test --- tests/testthat/test-classification.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/testthat/test-classification.R b/tests/testthat/test-classification.R index 5d9ac1960..7f8ac20aa 100644 --- a/tests/testthat/test-classification.R +++ b/tests/testthat/test-classification.R @@ -144,7 +144,7 @@ test_that("Classify with exclusion mask", { # testing original data probs_map_rst <- terra::rast(probs_map[["file_info"]][[1]][["path"]]) - expect_true(anyNA(class_map_rst[])) + expect_true(anyNA(probs_map_rst[])) probs_map_value <- terra::extract( x = probs_map_rst, From 64448d634dc17a9987c01d94aff4561919f39680 Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Fri, 29 Nov 2024 05:07:28 -0300 Subject: [PATCH 15/19] fix file management in classification and smooth tests --- tests/testthat/test-classification.R | 20 +++++++++++--------- tests/testthat/test-smooth.R | 18 +++++++----------- 2 files changed, 18 insertions(+), 20 deletions(-) diff --git a/tests/testthat/test-classification.R b/tests/testthat/test-classification.R index 7f8ac20aa..fb18f727d 100644 --- a/tests/testthat/test-classification.R +++ b/tests/testthat/test-classification.R @@ -89,11 +89,13 @@ test_that("Classify with NA values", { class_map <- sits_classify( data = raster_cube, ml_model = rfor_model, - output_dir = tempdir(), + output_dir = data_dir, progress = FALSE ) class_map_rst <- terra::rast(class_map[["file_info"]][[1]][["path"]]) expect_true(anyNA(class_map_rst[])) + # remove test files + unlink(data_dir) }) test_that("Classify with exclusion mask", { @@ -110,7 +112,9 @@ test_that("Classify with exclusion mask", { multicores = 2, progress = FALSE ) - + # preparation - create directory to save NA + data_dir <- paste0(tempdir(), "/exclusion-mask-na") + dir.create(data_dir, recursive = TRUE, showWarnings = FALSE) # preparation - create exclusion mask exclusion_mask <- sf::st_as_sfc( x = sf::st_bbox(c( @@ -122,34 +126,32 @@ test_that("Classify with exclusion mask", { crs = "EPSG:4326" ) ) - + # transform object to cube crs exclusion_mask <- sf::st_transform(exclusion_mask, .cube_crs(raster_cube)) - # preparation - calculate centroid of the exclusion mask exclusion_mask_centroid <- sf::st_centroid(exclusion_mask) - # preparation - create a random forest model rfor_model <- sits_train(samples_modis_ndvi, sits_rfor(num_trees = 40)) - # test classification with NA probs_map <- suppressWarnings( sits_classify( data = raster_cube, ml_model = rfor_model, - output_dir = tempdir(), + output_dir = data_dir, exclusion_mask = exclusion_mask, progress = FALSE ) ) - # testing original data probs_map_rst <- terra::rast(probs_map[["file_info"]][[1]][["path"]]) expect_true(anyNA(probs_map_rst[])) - + # extract values probs_map_value <- terra::extract( x = probs_map_rst, y = terra::vect(exclusion_mask_centroid) ) expect_true(any(is.na(probs_map_value))) + # remove test files + unlink(data_dir) }) diff --git a/tests/testthat/test-smooth.R b/tests/testthat/test-smooth.R index e67a91b0e..17ad9cecb 100644 --- a/tests/testthat/test-smooth.R +++ b/tests/testthat/test-smooth.R @@ -12,7 +12,9 @@ test_that("Smoothing with exclusion mask", { multicores = 2, progress = FALSE ) - + # preparation - create directory to save NA + data_dir <- paste0(tempdir(), "/smooth-exclusion-mask-na") + dir.create(data_dir, recursive = TRUE, showWarnings = FALSE) # preparation - create exclusion mask exclusion_mask <- sf::st_as_sfc( x = sf::st_bbox(c( @@ -24,46 +26,40 @@ test_that("Smoothing with exclusion mask", { crs = .cube_crs(raster_cube) ) ) - # preparation - calculate centroid of the exclusion mask exclusion_mask_centroid <- sf::st_centroid(exclusion_mask) - # preparation - create a random forest model rfor_model <- sits_train(samples_modis_ndvi, sits_rfor(num_trees = 40)) - # test classification with NA probs_map <- suppressWarnings( sits_classify( data = raster_cube, ml_model = rfor_model, - output_dir = tempdir(), + output_dir = data_dir, progress = FALSE ) ) - # smoth with exclusion mask smooth_map <- sits_smooth( cube = probs_map, exclusion_mask = exclusion_mask, - output_dir = tempdir(), + output_dir = data_dir, multicores = 2 ) - # testing original data (no na) probs_map_rst <- terra::rast(probs_map[["file_info"]][[1]][["path"]]) probs_map_value <- terra::extract( x = probs_map_rst, y = terra::vect(exclusion_mask_centroid) ) - expect_false(any(is.na(probs_map_value))) - # testing smooth data (with na) smooth_map_rst <- terra::rast(smooth_map[["file_info"]][[1]][["path"]]) smooth_map_value <- terra::extract( x = smooth_map_rst, y = terra::vect(exclusion_mask_centroid) ) - expect_true(any(is.na(smooth_map_value))) + # remove test files + unlink(data_dir) }) From a7388b17b647e4ef187e7f61363a4403e340bf1a Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Sun, 1 Dec 2024 05:05:42 -0300 Subject: [PATCH 16/19] include cube creation test with various roi types --- tests/testthat/test-cube.R | 110 +++++++++++++++++++++++++++++++++++++ 1 file changed, 110 insertions(+) diff --git a/tests/testthat/test-cube.R b/tests/testthat/test-cube.R index 6acc1b2c6..a6135a00e 100644 --- a/tests/testthat/test-cube.R +++ b/tests/testthat/test-cube.R @@ -99,6 +99,116 @@ test_that("Reading a raster cube", { expect_true(params_2$xres >= 231.5) }) +test_that("Reading raster cube with various type of ROI", { + roi <- c( + xmin = -44.58699, + ymin = -23.12016, + xmax = -44.45059, + ymax = -22.97294 + ) + + crs <- "EPSG:4326" + expected_tile <- "23KNQ" + + # Test 1a: ROI as vector + cube <- .try({ + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + roi = roi, + crs = crs, + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(cube), message = "MPC is not accessible") + expect_equal(cube[["tile"]], expected_tile) + + # Test 1b: ROI as vector - Expect a message when no CRS is specified + expect_warning( + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + roi = roi, + progress = FALSE + ) + ) + + # Test 2: ROI as SF + roi_sf <- sf::st_as_sfc( + x = sf::st_bbox( + roi, crs = crs + ) + ) + + cube <- .try({ + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + roi = roi_sf, + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(cube), message = "MPC is not accessible") + expect_equal(cube[["tile"]], expected_tile) + + # Test 3: ROI as lon/lat + roi_lonlat <- roi + names(roi_lonlat) <- c("lon_min", "lat_min", "lon_max", "lat_max") + + cube <- .try({ + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + roi = roi_lonlat, + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(cube), message = "MPC is not accessible") + expect_equal(cube[["tile"]], expected_tile) + + # Test 4a: ROI as SpatExtent + roi_raster <- terra::rast( + extent = terra::ext(roi["xmin"], roi["xmax"], roi["ymin"], roi["ymax"]), + crs = crs + ) + + roi_raster <- terra::ext(roi_raster) + + cube <- .try({ + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + roi = roi_raster, + crs = crs, + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(cube), message = "MPC is not accessible") + expect_equal(cube[["tile"]], expected_tile) + + # Test 4b: ROI as SpatExtent - Error when no CRS is specified + expect_error( + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + roi = roi_raster, + progress = FALSE + ) + ) +}) + test_that("Combining Sentinel-1 with Sentinel-2 cubes", { s2_cube <- .try( { From ee57ea6286430231f7abf05c190e74ae67d07c0a Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Sun, 1 Dec 2024 11:40:51 -0300 Subject: [PATCH 17/19] include shapefile as roi test case --- tests/testthat/test-cube.R | 24 ++++++++++++++++++++++++ 1 file changed, 24 insertions(+) diff --git a/tests/testthat/test-cube.R b/tests/testthat/test-cube.R index a6135a00e..f280610ef 100644 --- a/tests/testthat/test-cube.R +++ b/tests/testthat/test-cube.R @@ -207,6 +207,30 @@ test_that("Reading raster cube with various type of ROI", { progress = FALSE ) ) + + # Test 5: ROI as shapefile + shp_file <- tempfile(fileext = ".shp") + + sf::st_as_sfc( + x = sf::st_bbox( + roi, crs = crs + ) + ) |> + sf::st_write(shp_file, quiet = TRUE) + + cube <- .try({ + sits_cube( + source = "MPC", + collection = "SENTINEL-2-L2A", + roi = shp_file, + progress = FALSE + ) + }, + .default = NULL + ) + + testthat::skip_if(purrr::is_null(cube), message = "MPC is not accessible") + expect_equal(cube[["tile"]], expected_tile) }) test_that("Combining Sentinel-1 with Sentinel-2 cubes", { From 1cd5a6887a2beea0556691cb8cba0e1bba37b29b Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Sun, 1 Dec 2024 11:41:12 -0300 Subject: [PATCH 18/19] enhance sits_cube documentation --- R/sits_cube.R | 294 +++++++++++++++++++++++++++++------------------ man/sits_cube.Rd | 293 ++++++++++++++++++++++++++++------------------ 2 files changed, 360 insertions(+), 227 deletions(-) diff --git a/R/sits_cube.R b/R/sits_cube.R index 120cf42be..4a8372d90 100755 --- a/R/sits_cube.R +++ b/R/sits_cube.R @@ -5,10 +5,10 @@ #' in collections available in cloud services or local repositories. #' The following cloud providers are supported, based on the STAC protocol: #' Amazon Web Services (AWS), Brazil Data Cube (BDC), -#' Digital Earth Africa (DEAFRICA), Microsoft Planetary Computer (MPC), -#' Nasa Harmonized Landsat/Sentinel (HLS), USGS Landsat (USGS), and -#' Swiss Data Cube (SDC). Data cubes can also be created using local files. -#' +#' Copernicus Data Space Ecosystem (CDSE), Digital Earth Africa (DEAFRICA), +#' Digital Earth Australia (DEAUSTRALIA), Microsoft Planetary Computer (MPC), +#' Nasa Harmonized Landsat/Sentinel (HLS), Swiss Data Cube (SDC), TERRASCOPE or +#' USGS Landsat (USGS). Data cubes can also be created using local files. #' #' @param source Data source (one of \code{"AWS"}, \code{"BDC"}, #' \code{"DEAFRICA"}, \code{"MPC"}, \code{"SDC"}, @@ -25,13 +25,13 @@ #' the cube (see details below) #' (character vector of length 1). #' @param roi Region of interest (either an sf object, shapefile, -#' SpatExtent, or a numeric vector with named XY values -#' ("xmin", "xmax", "ymin", "ymax") or +#' \code{SpatExtent}, or a numeric vector with named XY +#' values ("xmin", "xmax", "ymin", "ymax") or #' named lat/long values #' ("lon_min", "lat_min", "lon_max", "lat_max"). #' @param crs The Coordinate Reference System (CRS) of the roi. It #' must be specified when roi is named XY values -#' ("xmin", "xmax", "ymin", "ymax") and SpatExtent +#' ("xmin", "xmax", "ymin", "ymax") or \code{SpatExtent} #' @param bands Spectral bands and indices to be included #' in the cube (optional - character vector). #' Use \code{\link{sits_list_collections}()} to find out @@ -62,135 +62,201 @@ #' @note{ #' To create cubes from cloud providers, users need to inform: #' \enumerate{ -#' \item \code{source}: One of "AWS", "BDC", "DEAFRICA", "HLS", "MPC", -#' "SDC" or "USGS"; +#' \item \code{source}: One of "AWS", "BDC", "CDSE", "DEAFRICA", "DEAUSTRALIA", +#' "HLS", "MPC", "SDC", "TERRASCOPE", or "USGS"; #' \item \code{collection}: Collection available in the cloud provider. -#' Use \code{sits_list_collections()} to see which +#' Use \code{\link{sits_list_collections}()} to see which #' collections are supported; #' \item \code{tiles}: A set of tiles defined according to the collection #' tiling grid; #' \item \code{roi}: Region of interest. Either -#' a named \code{vector} (\code{"lon_min"}, \code{"lat_min"}, -#' \code{"lon_max"}, \code{"lat_max"}) in WGS84, a \code{sfc} -#' or \code{sf} object from sf package in WGS84 projection. A named -#' \code{vector} (\code{"xmin"}, \code{"xmax"}, -#' \code{"ymin"}, \code{"ymax"}) and a \code{SpatExtent} can also +#' a shapefile, a named \code{vector} (\code{"lon_min"}, +#' \code{"lat_min"}, \code{"lon_max"}, \code{"lat_max"}) in WGS84, a +#' \code{sfc} or \code{sf} object from sf package in WGS84 projection. +#' A named \code{vector} (\code{"xmin"}, \code{"xmax"}, +#' \code{"ymin"}, \code{"ymax"}) or a \code{SpatExtent} can also #' be used, requiring only the specification of the \code{crs} parameter. #' } -#' Either \code{tiles} or \code{roi} must be informed. -#' The parameters \code{bands}, \code{start_date}, and -#' \code{end_date} are optional for cubes created from cloud providers. #' -#' GeoJSON geometries (RFC 7946) and shapefiles should be converted to -#' \code{sf} objects before being used to define a region of interest. -#' This parameter does not crop a region; it only selects images that -#' intersect the \code{roi}. +#' The parameter \code{bands}, \code{start_date}, and \code{end_date} are +#' optional for cubes created from cloud providers. #' -#' To create a cube from local files, users need to inform: -#' \enumerate{ -#' \item \code{source}: Provider from where the data has been downloaded -#' (e.g, "BDC"); -#' \item \code{collection}: Collection where the data has been extracted from. -#' (e.g., "SENTINEL-2-L2A" for the Sentinel-2 MPC collection level 2A); -#' \item \code{data_dir}: Local directory where images are stored. -#' \item \code{parse_info}: Parsing information for files. -#' Default is \code{c("X1", "X2", "tile", "band", "date")}. -#' \item \code{delim}: Delimiter character for parsing files. -#' Default is \code{"_"}. -#' } +#' Either \code{tiles} or \code{roi} must be informed. The \code{roi} parameter +#' is used to select images. This parameter does not crop a region; it only +#' selects images that intersect it. #' -#' To create a cube from local files, all images should have -#' the same spatial resolution and projection and each file should contain -#' a single image band for a single date. -#' Files can belong to different tiles of a spatial reference system and -#' file names need to include tile, date, and band information. -#' For example: \code{"CBERS-4_WFI_022024_B13_2018-02-02.tif"} -#' and \code{"SENTINEL-2_MSI_20LKP_B02_2018-07-18.jp2"} are accepted names. -#' The user has to provide parsing information to allow \code{sits} -#' to extract values of tile, band, and date. In the examples above, -#' the parsing info is c("X1", "X2", "tile", "band", "date") -#' and the delimiter is "_", which are the default values. -#' -#' It is also possible to create result cubes for these are local files -#' produced by classification or post-classification algorithms. In -#' this case, more parameters that are required (see below). The -#' parameter \code{parse_info} is specified differently, as follows: +#' If you want to use GeoJSON geometries (RFC 7946) as value \code{roi}, you +#' can convert it to sf object and then use it. #' -#' \enumerate{ -#' \item \code{band}: Band name associated to the type of result. Use -#' \code{"probs"}, for probability cubes produced by \code{sits_classify()}; -#' \code{"bayes"}, for smoothed cubes produced by \code{sits_smooth()}; -#' \code{"segments"}, for vector cubes produced by \code{sits_segment()}; -#' \code{"entropy"} when using \code{sits_uncertainty()}, and \code{"class"} -#' for cubes produced by \code{sits_label_classification()}; -#' \item \code{labels}: Labels associated to the classification results; -#' \item \code{parse_info}: File name parsing information -#' to deduce the values of "tile", "start_date", "end_date" from -#' the file name. Default is c("X1", "X2", "tile", "start_date", -#' "end_date", "band"). Unlike non-classified image files, -#' cubes with results have both -#' "start_date" and "end_date". -#' } +#' \code{sits} can access data from multiple providers, including +#' \code{Amazon Web Services} (AWS), \code{Microsoft Planetary Computer} (MPC), +#' \code{Brazil Data Cube} (BDC), \code{Copernicus Data Space Ecosystem} (CDSE), +#' \code{Digital Earth Africa}, \code{Digital Earth Australia}, +#' \code{NASA EarthData}, \code{Terrascope} and more. #' -#' In MPC, sits can access are two open data collections: +#' In each provider, \code{sits} can access multiple collections. For example, +#' in MPC \code{sits} can access multiple open data collections, including #' \code{"SENTINEL-2-L2A"} for Sentinel-2/2A images, and #' \code{"LANDSAT-C2-L2"} for the Landsat-4/5/7/8/9 collection. -#' (requester-pays) and \code{"SENTINEL-S2-L2A-COGS"} (open data). -#' -#' Sentinel-2/2A level 2A files in MPC are organized by sensor -#' resolution. The bands in 10m resolution are \code{"B02"}, \code{"B03"}, -#' \code{"B04"}, and \code{"B08"}. The 20m bands are \code{"B05"}, -#' \code{"B06"}, \code{"B07"}, \code{"B8A"}, \code{"B11"}, and \code{"B12"}. -#' Bands \code{"B01"} and \code{"B09"} are available at 60m resolution. -#' The \code{"CLOUD"} band is also available. -#' -#' All Landsat-4/5/7/8/9 images in MPC have bands with 30 meter -#' resolution. To account for differences between the different sensors, -#' Landsat bands in this collection have been renamed \code{"BLUE"}, -#' \code{"GREEN"}, \code{"RED"}, \code{"NIR08"}, \code{"SWIR16"} -#' and \code{"SWIR22"}. The \code{"CLOUD"} band is also available. #' #' In AWS, there are two types of collections: open data and -#' requester-pays. Currently, \code{sits} supports collection -#' \code{"SENTINEL-2-L2A"} (open data) and LANDSAT-C2-L2 (requester-pays). -#' There is no need to provide AWS credentials to access open data -#' collections. For requester-pays data, users need to provide their -#' access codes as environment variables, as follows: +#' requester-pays. Currently, \code{sits} supports collections +#' \code{"SENTINEL-2-L2A"}, \code{"SENTINEL-S2-L2A-COGS"} (open data) and +#' \code{"LANDSAT-C2-L2"} (requester-pays). There is no need to provide AWS +#' credentials to access open data collections. For requester-pays data, you +#' need to provide your AWS access codes as environment variables, as follows: #' \code{ #' Sys.setenv( #' AWS_ACCESS_KEY_ID = , #' AWS_SECRET_ACCESS_KEY = #' )} #' -#' Sentinel-2/2A level 2A files in AWS are organized by sensor -#' resolution. The AWS bands in 10m resolution are \code{"B02"}, \code{"B03"}, -#' \code{"B04"}, and \code{"B08"}. The 20m bands are \code{"B05"}, -#' \code{"B06"}, \code{"B07"}, \code{"B8A"}, \code{"B11"}, and \code{"B12"}. -#' Bands \code{"B01"} and \code{"B09"} are available at 60m resolution. -#' -#' For DEAFRICA, sits currently works with collections \code{"S2_L2A"} -#' for Sentinel-2 level 2A and \code{"LS8_SR"} for Landsat-8 ARD collection. -#' (open data). These collections are located in Africa -#' (Capetown) for faster access to African users. No payment for access -#' is required. -#' -#' For USGS, sits currently works with collection -#' \code{"LANDSAT-C2L2-SR"}, which corresponds to Landsat -#' Collection 2 Level-2 surface reflectance data, covering -#' Landsat-8 dataset. This collection is requester-pays and -#' requires payment for accessing. -#' -#' All BDC collections are regularized. -#' BDC users need to provide their credentials using environment -#' variables. To create your credentials, please see -#' . -#' Accessing data in the BDC is free. -#' After obtaining the BDC access key, please include it as -#' an environment variable, as follows: -#' \code{ -#' Sys.setenv( -#' BDC_ACCESS_KEY = -#' )} +#' In BDC, there are many collections, including \code{"LANDSAT-OLI-16D"} +#' (Landsat-8 OLI, 30 m resolution, 16-day intervals), \code{"SENTINEL-2-16D"} +#' (Sentinel-2A and 2B MSI images at 10 m resolution, 16-day intervals), +#' \code{"CBERS-WFI-16D"} (CBERS 4 WFI, 64 m resolution, 16-day intervals), and +#' others. All BDC collections are regularized. +#' +#' To explore providers and collections \code{sits} supports, use the +#' \code{\link{sits_list_collections}()} function. +#' +#' If you want to learn more details about each provider and collection +#' available in \code{sits}, please read the online sits book +#' (e-sensing.github.io/sitsbook). The chapter +#' \code{Earth Observation data cubes} provides a detailed description of all +#' collections you can use with \code{sits} +#' (e-sensing.github.io/sitsbook/earth-observation-data-cubes.html). +#' +#' To create a cube from local files, you need to inform: +#' \enumerate{ +#' \item \code{source}: The data provider from which the data was +#' downloaded (e.g, "BDC", "MPC"); +#' +#' \item \code{collection}: The collection from which the data comes from. +#' (e.g., \code{"SENTINEL-2-L2A"} for the Sentinel-2 MPC collection level 2A); +#' +#' \item \code{data_dir}: The local directory where the image files are stored. +#' +#' \item \code{parse_info}: Defines how to extract metadata from file names +#' by specifying the order and meaning of each part, separated by the +#' \code{"delim"} character. Default value is +#' \code{c("X1", "X2", "tile", "band", "date")}. +#' +#' \item \code{delim}: The delimiter character used to separate components in +#' the file names. Default is \code{"_"}. +#' } +#' +#' Note that if you are working with local data cubes created by \code{sits}, +#' you do not need to specify \code{parse_info} and \code{delim}. These elements +#' are automatically identified. This is particularly useful when you have +#' downloaded or created data cubes using \code{sits}. +#' +#' For example, if you downloaded a data cube from the Microsoft Planetary +#' Computer (MPC) using the function \code{\link{sits_cube_copy}()}, you do +#' not need to provide \code{parse_info} and \code{delim}. +#' +#' If you are using a data cube from a source supported by \code{sits} +#' (e.g., AWS, MPC) but downloaded / managed with an external tool, you will +#' need to specify the \code{parse_info} and \code{delim} parameters manually. +#' For this case, you first need to ensure that the local files meet some +#' critical requirements: +#' +#' \itemize{ +#' \item All image files must have the same spatial resolution and projection; +#' +#' \item Each file should represent a single image band for a single date; +#' +#' \item File names must include information about the \code{"tile"}, +#' \code{"date"}, and \code{"band"} in the file. +#' } +#' +#' For example, if you are creating a Sentinel-2 data cube on your local +#' machine, and the files have the same spatial resolution and projection, with +#' each file containing a single band and date, an acceptable file name could be: +#' \itemize{ +#' \item \code{"SENTINEL-2_MSI_20LKP_B02_2018-07-18.jp2"} +#' } +#' +#' This file name works because it encodes the three key pieces of information +#' used by \code{sits}: +#' \itemize{ +#' \item Tile: "20LKP"; +#' +#' \item Band: "B02"; +#' +#' \item Date: "2018-07-18" +#' } +#' +#' Other example of supported file names are: +#' \itemize{ +#' \item \code{"CBERS-4_WFI_022024_B13_2021-05-15.tif"}; +#' +#' \item \code{"SENTINEL-1_GRD_30TXL_VV_2023-03-10.tif"}; +#' +#' \item \code{"LANDSAT-8_OLI_198030_B04_2020-09-12.tif"}. +#' } +#' +#' The \code{parse_info} parameter tells \code{sits} how to extract essential +#' metadata from file names. It defines the sequence of components in the +#' file name, assigning each part a label such as \code{"tile"}, \code{"band"}, +#' and \code{"date"}. For parts of the file name that are irrelevant to +#' \code{sits}, you can use dummy labels like \code{"X1"}, \code{"X2"}, and so +#' on. +#' +#' For example, consider the file name: +#' \itemize{ +#' \item \code{"SENTINEL-2_MSI_20LKP_B02_2018-07-18.jp2"} +#' } +#' +#' With \code{parse_info = c("X1", "X2", "tile", "band", "date")} and +#' \code{delim = "_"}, the extracted metadata would be: +#' +#' \itemize{ +#' \item X1: "SENTINEL-2" (ignored) +#' \item X2: "MSI" (ignored) +#' \item tile: "20LKP" (used) +#' \item band: "B02" (used) +#' \item date: "2018-07-18" (used) +#' } +#' +#' The \code{delim} parameter specifies the character that separates components +#' in the file name. The default delimiter is \code{"_"}. +#' +#' Note that when you load a local data cube specifying the \code{source} +#' (e.g., AWS, MPC) and \code{collection}, \code{sits} assumes that the data +#' properties (e.g., scale factor, minimum, and maximum values) match those +#' defined for the selected provider. However, if you are working with +#' custom data from an unsupported source or data that does not follow the +#' standard definitions of providers in sits, refer to the Technical Annex of +#' the \code{sits} online book for guidance on handling such cases +#' (e-sensing.github.io/sitsbook/technical-annex.html). +#' +#' It is also possible to create result cubes from local files produced by +#' classification or post-classification algorithms. In this case, the +#' \code{parse_info} is specified differently, and other additional parameters +#' are required: +#' +#' \itemize{ +#' +#' \item \code{band}: Band name associated to the type of result. Use +#' \code{"probs"}, for probability cubes produced by +#' \code{\link{sits_classify}()}; +#' \code{"bayes"}, for smoothed cubes produced by \code{\link{sits_smooth}()}; +#' \code{"segments"}, for vector cubes produced by +#' \code{\link{sits_segment}()}; +#' \code{"entropy"} when using \code{\link{sits_uncertainty}()}, and +#' \code{"class"} for cubes produced by +#' \code{\link{sits_label_classification}()}; +#' +#' \item \code{labels}: Labels associated to the classification results; +#' +#' \item \code{parse_info}: File name parsing information +#' to deduce the values of \code{"tile"}, \code{"start_date"}, +#' \code{"end_date"} from the file name. Unlike non-classified image files, +#' cubes with results have both \code{"start_date"} and \code{"end_date"}. +#' Default is c("X1", "X2", "tile", "start_date", "end_date", "band"). +#' } +#' #' } #' @examples #' if (sits_run_examples()) { diff --git a/man/sits_cube.Rd b/man/sits_cube.Rd index 32b8a8305..feccf13f5 100644 --- a/man/sits_cube.Rd +++ b/man/sits_cube.Rd @@ -83,14 +83,14 @@ the cube (see details below) (character vector of length 1).} \item{roi}{Region of interest (either an sf object, shapefile, -SpatExtent, or a numeric vector with named XY values -("xmin", "xmax", "ymin", "ymax") or +\code{SpatExtent}, or a numeric vector with named XY +values ("xmin", "xmax", "ymin", "ymax") or named lat/long values ("lon_min", "lat_min", "lon_max", "lat_max").} \item{crs}{The Coordinate Reference System (CRS) of the roi. It must be specified when roi is named XY values -("xmin", "xmax", "ymin", "ymax") and SpatExtent} +("xmin", "xmax", "ymin", "ymax") or \code{SpatExtent}} \item{start_date, end_date}{Initial and final dates to include images from the collection in the cube (optional). @@ -134,143 +134,210 @@ Creates a data cube based on spatial and temporal restrictions in collections available in cloud services or local repositories. The following cloud providers are supported, based on the STAC protocol: Amazon Web Services (AWS), Brazil Data Cube (BDC), -Digital Earth Africa (DEAFRICA), Microsoft Planetary Computer (MPC), -Nasa Harmonized Landsat/Sentinel (HLS), USGS Landsat (USGS), and -Swiss Data Cube (SDC). Data cubes can also be created using local files. +Copernicus Data Space Ecosystem (CDSE), Digital Earth Africa (DEAFRICA), +Digital Earth Australia (DEAUSTRALIA), Microsoft Planetary Computer (MPC), +Nasa Harmonized Landsat/Sentinel (HLS), Swiss Data Cube (SDC), TERRASCOPE or +USGS Landsat (USGS). Data cubes can also be created using local files. } \note{ { To create cubes from cloud providers, users need to inform: \enumerate{ - \item \code{source}: One of "AWS", "BDC", "DEAFRICA", "HLS", "MPC", -"SDC" or "USGS"; + \item \code{source}: One of "AWS", "BDC", "CDSE", "DEAFRICA", "DEAUSTRALIA", + "HLS", "MPC", "SDC", "TERRASCOPE", or "USGS"; \item \code{collection}: Collection available in the cloud provider. - Use \code{sits_list_collections()} to see which + Use \code{\link{sits_list_collections}()} to see which collections are supported; \item \code{tiles}: A set of tiles defined according to the collection tiling grid; \item \code{roi}: Region of interest. Either - a named \code{vector} (\code{"lon_min"}, \code{"lat_min"}, - \code{"lon_max"}, \code{"lat_max"}) in WGS84, a \code{sfc} - or \code{sf} object from sf package in WGS84 projection. A named - \code{vector} (\code{"xmin"}, \code{"xmax"}, - \code{"ymin"}, \code{"ymax"}) and a \code{SpatExtent} can also + a shapefile, a named \code{vector} (\code{"lon_min"}, + \code{"lat_min"}, \code{"lon_max"}, \code{"lat_max"}) in WGS84, a + \code{sfc} or \code{sf} object from sf package in WGS84 projection. + A named \code{vector} (\code{"xmin"}, \code{"xmax"}, + \code{"ymin"}, \code{"ymax"}) or a \code{SpatExtent} can also be used, requiring only the specification of the \code{crs} parameter. } -Either \code{tiles} or \code{roi} must be informed. -The parameters \code{bands}, \code{start_date}, and -\code{end_date} are optional for cubes created from cloud providers. -GeoJSON geometries (RFC 7946) and shapefiles should be converted to -\code{sf} objects before being used to define a region of interest. -This parameter does not crop a region; it only selects images that -intersect the \code{roi}. +The parameter \code{bands}, \code{start_date}, and \code{end_date} are +optional for cubes created from cloud providers. -To create a cube from local files, users need to inform: -\enumerate{ - \item \code{source}: Provider from where the data has been downloaded - (e.g, "BDC"); - \item \code{collection}: Collection where the data has been extracted from. - (e.g., "SENTINEL-2-L2A" for the Sentinel-2 MPC collection level 2A); - \item \code{data_dir}: Local directory where images are stored. - \item \code{parse_info}: Parsing information for files. - Default is \code{c("X1", "X2", "tile", "band", "date")}. - \item \code{delim}: Delimiter character for parsing files. - Default is \code{"_"}. -} +Either \code{tiles} or \code{roi} must be informed. The \code{roi} parameter +is used to select images. This parameter does not crop a region; it only +selects images that intersect it. -To create a cube from local files, all images should have -the same spatial resolution and projection and each file should contain -a single image band for a single date. -Files can belong to different tiles of a spatial reference system and -file names need to include tile, date, and band information. -For example: \code{"CBERS-4_WFI_022024_B13_2018-02-02.tif"} -and \code{"SENTINEL-2_MSI_20LKP_B02_2018-07-18.jp2"} are accepted names. -The user has to provide parsing information to allow \code{sits} -to extract values of tile, band, and date. In the examples above, -the parsing info is c("X1", "X2", "tile", "band", "date") -and the delimiter is "_", which are the default values. - -It is also possible to create result cubes for these are local files -produced by classification or post-classification algorithms. In -this case, more parameters that are required (see below). The -parameter \code{parse_info} is specified differently, as follows: +If you want to use GeoJSON geometries (RFC 7946) as value \code{roi}, you +can convert it to sf object and then use it. -\enumerate{ -\item \code{band}: Band name associated to the type of result. Use - \code{"probs"}, for probability cubes produced by \code{sits_classify()}; - \code{"bayes"}, for smoothed cubes produced by \code{sits_smooth()}; - \code{"segments"}, for vector cubes produced by \code{sits_segment()}; - \code{"entropy"} when using \code{sits_uncertainty()}, and \code{"class"} - for cubes produced by \code{sits_label_classification()}; -\item \code{labels}: Labels associated to the classification results; -\item \code{parse_info}: File name parsing information - to deduce the values of "tile", "start_date", "end_date" from - the file name. Default is c("X1", "X2", "tile", "start_date", - "end_date", "band"). Unlike non-classified image files, - cubes with results have both - "start_date" and "end_date". -} +\code{sits} can access data from multiple providers, including +\code{Amazon Web Services} (AWS), \code{Microsoft Planetary Computer} (MPC), +\code{Brazil Data Cube} (BDC), \code{Copernicus Data Space Ecosystem} (CDSE), +\code{Digital Earth Africa}, \code{Digital Earth Australia}, +\code{NASA EarthData}, \code{Terrascope} and more. -In MPC, sits can access are two open data collections: +In each provider, \code{sits} can access multiple collections. For example, +in MPC \code{sits} can access multiple open data collections, including \code{"SENTINEL-2-L2A"} for Sentinel-2/2A images, and \code{"LANDSAT-C2-L2"} for the Landsat-4/5/7/8/9 collection. -(requester-pays) and \code{"SENTINEL-S2-L2A-COGS"} (open data). - -Sentinel-2/2A level 2A files in MPC are organized by sensor -resolution. The bands in 10m resolution are \code{"B02"}, \code{"B03"}, -\code{"B04"}, and \code{"B08"}. The 20m bands are \code{"B05"}, -\code{"B06"}, \code{"B07"}, \code{"B8A"}, \code{"B11"}, and \code{"B12"}. -Bands \code{"B01"} and \code{"B09"} are available at 60m resolution. -The \code{"CLOUD"} band is also available. - -All Landsat-4/5/7/8/9 images in MPC have bands with 30 meter -resolution. To account for differences between the different sensors, -Landsat bands in this collection have been renamed \code{"BLUE"}, -\code{"GREEN"}, \code{"RED"}, \code{"NIR08"}, \code{"SWIR16"} -and \code{"SWIR22"}. The \code{"CLOUD"} band is also available. In AWS, there are two types of collections: open data and -requester-pays. Currently, \code{sits} supports collection -\code{"SENTINEL-2-L2A"} (open data) and LANDSAT-C2-L2 (requester-pays). -There is no need to provide AWS credentials to access open data -collections. For requester-pays data, users need to provide their -access codes as environment variables, as follows: +requester-pays. Currently, \code{sits} supports collections +\code{"SENTINEL-2-L2A"}, \code{"SENTINEL-S2-L2A-COGS"} (open data) and +\code{"LANDSAT-C2-L2"} (requester-pays). There is no need to provide AWS +credentials to access open data collections. For requester-pays data, you +need to provide your AWS access codes as environment variables, as follows: \code{ Sys.setenv( AWS_ACCESS_KEY_ID = , AWS_SECRET_ACCESS_KEY = )} -Sentinel-2/2A level 2A files in AWS are organized by sensor -resolution. The AWS bands in 10m resolution are \code{"B02"}, \code{"B03"}, -\code{"B04"}, and \code{"B08"}. The 20m bands are \code{"B05"}, -\code{"B06"}, \code{"B07"}, \code{"B8A"}, \code{"B11"}, and \code{"B12"}. -Bands \code{"B01"} and \code{"B09"} are available at 60m resolution. - -For DEAFRICA, sits currently works with collections \code{"S2_L2A"} -for Sentinel-2 level 2A and \code{"LS8_SR"} for Landsat-8 ARD collection. -(open data). These collections are located in Africa -(Capetown) for faster access to African users. No payment for access -is required. - -For USGS, sits currently works with collection -\code{"LANDSAT-C2L2-SR"}, which corresponds to Landsat -Collection 2 Level-2 surface reflectance data, covering -Landsat-8 dataset. This collection is requester-pays and -requires payment for accessing. - -All BDC collections are regularized. -BDC users need to provide their credentials using environment -variables. To create your credentials, please see -. -Accessing data in the BDC is free. -After obtaining the BDC access key, please include it as -an environment variable, as follows: -\code{ -Sys.setenv( - BDC_ACCESS_KEY = -)} +In BDC, there are many collections, including \code{"LANDSAT-OLI-16D"} +(Landsat-8 OLI, 30 m resolution, 16-day intervals), \code{"SENTINEL-2-16D"} +(Sentinel-2A and 2B MSI images at 10 m resolution, 16-day intervals), +\code{"CBERS-WFI-16D"} (CBERS 4 WFI, 64 m resolution, 16-day intervals), and +others. All BDC collections are regularized. + +To explore providers and collections \code{sits} supports, use the +\code{\link{sits_list_collections}()} function. + +If you want to learn more details about each provider and collection +available in \code{sits}, please read the online sits book +(e-sensing.github.io/sitsbook). The chapter +\code{Earth Observation data cubes} provides a detailed description of all +collections you can use with \code{sits} +(e-sensing.github.io/sitsbook/earth-observation-data-cubes.html). + +To create a cube from local files, you need to inform: +\enumerate{ + \item \code{source}: The data provider from which the data was + downloaded (e.g, "BDC", "MPC"); + + \item \code{collection}: The collection from which the data comes from. + (e.g., \code{"SENTINEL-2-L2A"} for the Sentinel-2 MPC collection level 2A); + + \item \code{data_dir}: The local directory where the image files are stored. + + \item \code{parse_info}: Defines how to extract metadata from file names + by specifying the order and meaning of each part, separated by the + \code{"delim"} character. Default value is + \code{c("X1", "X2", "tile", "band", "date")}. + + \item \code{delim}: The delimiter character used to separate components in + the file names. Default is \code{"_"}. +} + +Note that if you are working with local data cubes created by \code{sits}, +you do not need to specify \code{parse_info} and \code{delim}. These elements +are automatically identified. This is particularly useful when you have +downloaded or created data cubes using \code{sits}. + +For example, if you downloaded a data cube from the Microsoft Planetary +Computer (MPC) using the function \code{\link{sits_cube_copy}()}, you do +not need to provide \code{parse_info} and \code{delim}. + +If you are using a data cube from a source supported by \code{sits} +(e.g., AWS, MPC) but downloaded / managed with an external tool, you will +need to specify the \code{parse_info} and \code{delim} parameters manually. +For this case, you first need to ensure that the local files meet some +critical requirements: + +\itemize{ + \item All image files must have the same spatial resolution and projection; + + \item Each file should represent a single image band for a single date; + + \item File names must include information about the \code{"tile"}, + \code{"date"}, and \code{"band"} in the file. +} + +For example, if you are creating a Sentinel-2 data cube on your local +machine, and the files have the same spatial resolution and projection, with +each file containing a single band and date, an acceptable file name could be: +\itemize{ + \item \code{"SENTINEL-2_MSI_20LKP_B02_2018-07-18.jp2"} +} + +This file name works because it encodes the three key pieces of information +used by \code{sits}: +\itemize{ + \item Tile: "20LKP"; + + \item Band: "B02"; + + \item Date: "2018-07-18" +} + +Other example of supported file names are: +\itemize{ + \item \code{"CBERS-4_WFI_022024_B13_2021-05-15.tif"}; + + \item \code{"SENTINEL-1_GRD_30TXL_VV_2023-03-10.tif"}; + + \item \code{"LANDSAT-8_OLI_198030_B04_2020-09-12.tif"}. +} + +The \code{parse_info} parameter tells \code{sits} how to extract essential +metadata from file names. It defines the sequence of components in the +file name, assigning each part a label such as \code{"tile"}, \code{"band"}, +and \code{"date"}. For parts of the file name that are irrelevant to +\code{sits}, you can use dummy labels like \code{"X1"}, \code{"X2"}, and so +on. + +For example, consider the file name: +\itemize{ + \item \code{"SENTINEL-2_MSI_20LKP_B02_2018-07-18.jp2"} +} + +With \code{parse_info = c("X1", "X2", "tile", "band", "date")} and +\code{delim = "_"}, the extracted metadata would be: + +\itemize{ + \item X1: "SENTINEL-2" (ignored) + \item X2: "MSI" (ignored) + \item tile: "20LKP" (used) + \item band: "B02" (used) + \item date: "2018-07-18" (used) +} + +The \code{delim} parameter specifies the character that separates components +in the file name. The default delimiter is \code{"_"}. + +Note that when you load a local data cube specifying the \code{source} +(e.g., AWS, MPC) and \code{collection}, \code{sits} assumes that the data +properties (e.g., scale factor, minimum, and maximum values) match those +defined for the selected provider. However, if you are working with +custom data from an unsupported source or data that does not follow the +standard definitions of providers in sits, refer to the Technical Annex of +the \code{sits} online book for guidance on handling such cases +(e-sensing.github.io/sitsbook/technical-annex.html). + +It is also possible to create result cubes from local files produced by +classification or post-classification algorithms. In this case, the +\code{parse_info} is specified differently, and other additional parameters +are required: + +\itemize{ + +\item \code{band}: Band name associated to the type of result. Use + \code{"probs"}, for probability cubes produced by + \code{\link{sits_classify}()}; + \code{"bayes"}, for smoothed cubes produced by \code{\link{sits_smooth}()}; + \code{"segments"}, for vector cubes produced by + \code{\link{sits_segment}()}; + \code{"entropy"} when using \code{\link{sits_uncertainty}()}, and + \code{"class"} for cubes produced by + \code{\link{sits_label_classification}()}; + +\item \code{labels}: Labels associated to the classification results; + +\item \code{parse_info}: File name parsing information + to deduce the values of \code{"tile"}, \code{"start_date"}, + \code{"end_date"} from the file name. Unlike non-classified image files, + cubes with results have both \code{"start_date"} and \code{"end_date"}. + Default is c("X1", "X2", "tile", "start_date", "end_date", "band"). +} + } } \examples{ From 294814bcc356634eca1a22774e2066231ff71f8e Mon Sep 17 00:00:00 2001 From: Felipe Carlos Date: Sun, 1 Dec 2024 12:07:31 -0300 Subject: [PATCH 19/19] update sits_merge documentation --- R/api_merge.R | 2 +- R/sits_merge.R | 36 ++++++++++++++++++++---------------- man/sits_merge.Rd | 34 +++++++++++++++++++--------------- 3 files changed, 40 insertions(+), 32 deletions(-) diff --git a/R/api_merge.R b/R/api_merge.R index 2a61dc9f7..51b83c5cd 100644 --- a/R/api_merge.R +++ b/R/api_merge.R @@ -48,7 +48,7 @@ } else { # if not, save the reference interval and the min value of # the t2 interval dates. - # this ensure there are not two dates in the same interval + # this ensures there are not two dates in the same interval t_overlap <- c( t_overlap, # dates storage reference_interval, # current interval diff --git a/R/sits_merge.R b/R/sits_merge.R index 825f710dd..99e51f5ca 100644 --- a/R/sits_merge.R +++ b/R/sits_merge.R @@ -3,20 +3,25 @@ #' @author Gilberto Camara, \email{gilberto.camara@@inpe.br} #' #' @description To merge two series, we consider that they contain different -#' attributes but refer to the same data cube, and spatiotemporal location. -#' This function is useful to merge different bands of the same locations. -#' For example, one may want to put the raw and smoothed bands -#' for the same set of locations in the same tibble. +#' attributes but refer to the same data cube and spatiotemporal location. +#' This function is useful for merging different bands of the same location. +#' For example, one may want to put the raw and smoothed bands for the same set +#' of locations in the same tibble. #' -#' In case of data cubes, the function merges the images based on the following -#' conditions: +#' In the case of data cubes, the function merges the images based on the +#' following conditions: #' \enumerate{ -#' \item if the bands are different and their timelines should be compatible, -#' the bands are joined. The resulting timeline is the one from the first cube. -#' This is useful to merge data from different sensors (e.g, Sentinel-1 with Sentinel-2). -#' \item if the bands are the same, the cube will have the combined -#' timeline of both cubes. This is useful to merge data from the same sensors -#' from different satellites (e.g, Sentinel-2A with Sentinel-2B). +#' \item If the two cubes have different bands but compatible timelines, the +#' bands are combined, and the timeline is adjusted to overlap. To create the +#' overlap, we align the timelines like a "zipper": for each interval defined +#' by a pair of consecutive dates in the first timeline, we include matching +#' dates from the second timeline. If the second timeline has multiple dates +#' in the same interval, only the minimum date is kept. This ensures the final +#' timeline avoids duplicates and is consistent. This is useful when merging +#' data from different sensors (e.g., Sentinel-1 with Sentinel-2). +#' \item If the bands are the same, the cube will have the combined timeline of +#' both cubes. This is useful for merging data from the same sensors from +#' different satellites (e.g., Sentinel-2A with Sentinel-2B). #' \item otherwise, the function will produce an error. #' } #' @@ -26,9 +31,8 @@ #' or data cube (tibble of class "raster_cube") . #' #' @param ... Additional parameters -#' @param suffix If there are duplicate bands in data1 and data2 -#' these suffixes will be added -#' (character vector). +#' @param suffix If data1 and data2 are tibble with duplicate bands, this +#' suffix will be added (character vector). #' #' @return merged data sets (tibble of class "sits" or #' tibble of class "raster_cube") @@ -124,7 +128,7 @@ sits_merge.raster_cube <- function(data1, data2, ...) { # merge! merged_cube <- .merge_cube_densify(data1, data2) } else { - # rule 2: if the bands are different and their timelines should be + # 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) diff --git a/man/sits_merge.Rd b/man/sits_merge.Rd index ccacd7dc7..82a823fa7 100644 --- a/man/sits_merge.Rd +++ b/man/sits_merge.Rd @@ -24,9 +24,8 @@ or data cube (tibble of class "raster_cube") .} \item{...}{Additional parameters} -\item{suffix}{If there are duplicate bands in data1 and data2 -these suffixes will be added -(character vector).} +\item{suffix}{If data1 and data2 are tibble with duplicate bands, this +suffix will be added (character vector).} } \value{ merged data sets (tibble of class "sits" or @@ -34,20 +33,25 @@ merged data sets (tibble of class "sits" or } \description{ To merge two series, we consider that they contain different -attributes but refer to the same data cube, and spatiotemporal location. -This function is useful to merge different bands of the same locations. -For example, one may want to put the raw and smoothed bands -for the same set of locations in the same tibble. +attributes but refer to the same data cube and spatiotemporal location. +This function is useful for merging different bands of the same location. +For example, one may want to put the raw and smoothed bands for the same set +of locations in the same tibble. -In case of data cubes, the function merges the images based on the following -conditions: +In the case of data cubes, the function merges the images based on the +following conditions: \enumerate{ -\item if the bands are different and their timelines should be compatible, -the bands are joined. The resulting timeline is the one from the first cube. -This is useful to merge data from different sensors (e.g, Sentinel-1 with Sentinel-2). -\item if the bands are the same, the cube will have the combined -timeline of both cubes. This is useful to merge data from the same sensors -from different satellites (e.g, Sentinel-2A with Sentinel-2B). +\item If the two cubes have different bands but compatible timelines, the +bands are combined, and the timeline is adjusted to overlap. To create the +overlap, we align the timelines like a "zipper": for each interval defined +by a pair of consecutive dates in the first timeline, we include matching +dates from the second timeline. If the second timeline has multiple dates +in the same interval, only the minimum date is kept. This ensures the final +timeline avoids duplicates and is consistent. This is useful when merging +data from different sensors (e.g., Sentinel-1 with Sentinel-2). +\item If the bands are the same, the cube will have the combined timeline of +both cubes. This is useful for merging data from the same sensors from +different satellites (e.g., Sentinel-2A with Sentinel-2B). \item otherwise, the function will produce an error. } }