Merge pull request #1482 from rolfsimoes/dev

Bug fixes and `sits_uncertainty_sampling()` improvements

Author: gilbertocamara

NAMESPACE

@@ -157,6 +157,8 @@ S3method(.source_collection_access_test,ogh_cube)
 S3method(.source_collection_access_test,stac_cube)
 S3method(.source_collection_access_test,usgs_cube)
 S3method(.source_cube,stac_cube)
+S3method(.source_filter_tiles,"bdc_cube_landsat-1y")
+S3method(.source_filter_tiles,"bdc_cube_landsat-2m")
 S3method(.source_filter_tiles,"cdse_os_cube_sentinel-1-rtc")
 S3method(.source_filter_tiles,"deafrica_cube_sentinel-1-rtc")
 S3method(.source_filter_tiles,"mpc_cube_cop-dem-glo-30")
@@ -179,6 +181,8 @@ S3method(.source_items_bands_select,stac_cube)
 S3method(.source_items_cube,stac_cube)
 S3method(.source_items_fid,stac_cube)
 S3method(.source_items_new,"aws_cube_landsat-c2-l2")
+S3method(.source_items_new,"bdc_cube_landsat-1y")
+S3method(.source_items_new,"bdc_cube_landsat-2m")
 S3method(.source_items_new,"deafrica_cube_sentinel-1-rtc")
 S3method(.source_items_new,"deafrica_cube_sentinel-2-l2a")
 S3method(.source_items_new,"mpc_cube_cop-dem-glo-30")

R/api_check.R

@@ -573,16 +573,20 @@
         msg = paste0("value should be <= ", max)
     )
     # exclusive_min and exclusive_max checks
-    .check_that(
-        all(x > exclusive_min),
-        local_msg = local_msg,
-        msg = paste0("value should be > ", exclusive_min)
-    )
-    .check_that(
-        all(x < exclusive_max),
-        local_msg = local_msg,
-        msg = paste0("value should be < ", exclusive_max)
-    )
+    if (is.finite(exclusive_min)) {
+        .check_that(
+            all(x > exclusive_min),
+            local_msg = local_msg,
+            msg = paste0("value should be > ", exclusive_min)
+        )
+    }
+    if (is.finite(exclusive_max)) {
+        .check_that(
+            all(x < exclusive_max),
+            local_msg = local_msg,
+            msg = paste0("value should be < ", exclusive_max)
+        )
+    }
 }
 #' @rdname check_functions
 #' @keywords internal

R/api_data.R

@@ -320,9 +320,16 @@
         classes <- labels[class_numbers]
         # insert classes into samples
         samples[["label"]] <- unname(classes)
+        # Preserve start_date and end_date if they exist in input samples
+        cols_to_select <- c("longitude", "latitude", "label")
+        if ("start_date" %in% names(samples)) {
+            cols_to_select <- c(cols_to_select, "start_date")
+        }
+        if ("end_date" %in% names(samples)) {
+            cols_to_select <- c(cols_to_select, "end_date")
+        }
         samples <- dplyr::select(
-            samples, dplyr::all_of("longitude"),
-            dplyr::all_of("latitude"), dplyr::all_of("label")
+            samples, dplyr::all_of(cols_to_select)
         )
         samples
     })

R/api_validate.R

@@ -27,6 +27,6 @@
     # Call caret package to the classification statistics
     acc_obj <- caret::confusionMatrix(predicted, reference)
     # Set result class and return it
-    .set_class(x = acc_obj, "sits_accuracy", class(acc_obj))
+    class(acc_obj) <- c("sits_accuracy", class(acc_obj))
     acc_obj
 }

R/sits_uncertainty.R

@@ -205,12 +205,15 @@ sits_uncertainty.default <- function(cube, ...) {
 #'                        See \code{\link[sits]{sits_uncertainty}}.
 #' @param n               Number of suggested points to be sampled per tile.
 #' @param min_uncert      Minimum uncertainty value to select a sample.
+#' @param max_uncert      Maximum uncertainty value to select a sample.
+#'                        Default is Inf (no upper limit).
 #' @param sampling_window Window size for collecting points (in pixels).
 #'                        The minimum window size is 10.
 #' @param multicores      Number of workers for parallel processing
 #'                        (integer, min = 1, max = 2048).
 #' @param memsize         Maximum overall memory (in GB) to run the
 #'                        function.
+#' @param progress        Whether to show progress bars (TRUE/FALSE).
 #'
 #' @return
 #' A tibble with longitude and latitude in WGS84 with locations
@@ -254,83 +257,182 @@ sits_uncertainty.default <- function(cube, ...) {
 sits_uncertainty_sampling <- function(uncert_cube,
                                       n = 100L,
                                       min_uncert = 0.4,
+                                      max_uncert = Inf,
                                       sampling_window = 10L,
                                       multicores = 2L,
-                                      memsize = 4L) {
+                                      memsize = 4L,
+                                      progress = FALSE) {
     .check_set_caller("sits_uncertainty_sampling")
     # Pre-conditions
     .check_is_uncert_cube(uncert_cube)
     .check_int_parameter(n, min = 1L)
     .check_num_parameter(min_uncert, min = 0.0, max = 1.0)
+    .check_num_parameter(max_uncert, min = 0.0)
     .check_int_parameter(sampling_window, min = 1L)
-    .check_int_parameter(multicores, min = 1L)
-    .check_int_parameter(memsize, min = 1L)
+    .check_int_parameter(multicores, min = 1L, max = 2048L)
+    .check_int_parameter(memsize, min = 1L, max = 16384L)
+    progress <- .message_progress(progress)
+
+    # The following functions define optimal parameters for parallel processing
+    # Get block size
+    block <- .raster_file_blocksize(.raster_open_rast(.tile_path(uncert_cube)))
+    # Overlapping pixels (no overlap for uncertainty sampling)
+    overlap <- 0L
+    # Check minimum memory needed to process one block
+    job_block_memsize <- .jobs_block_memsize(
+        block_size = .block_size(block = block, overlap = overlap),
+        npaths = 1L,
+        nbytes = 8L,
+        proc_bloat = .conf("processing_bloat_cpu")
+    )
+    # Update multicores parameter
+    multicores <- .jobs_max_multicores(
+        job_block_memsize = job_block_memsize,
+        memsize = memsize,
+        multicores = multicores
+    )
+    # Update block parameter
+    block <- .jobs_optimal_block(
+        job_block_memsize = job_block_memsize,
+        block = block,
+        image_size = .tile_size(.tile(uncert_cube)),
+        memsize = memsize,
+        multicores = multicores
+    )
+    # Prepare parallel processing
+    if (.parallel_start(workers = multicores)) {
+        on.exit(.parallel_stop(), add = TRUE)
+    }
     # Slide on cube tiles
     samples_tb <- slider::slide_dfr(uncert_cube, function(tile) {
-        # open spatial raster object
-        rast <- .raster_open_rast(.tile_path(tile))
-        # get the values
-        values <- .raster_get_values(rast)
-        # sample the maximum values
-        samples_tile <- C_max_sampling(
-            x = values,
-            nrows = nrow(rast),
-            ncols = ncol(rast),
-            window_size = sampling_window
+        # Create chunks as jobs
+        chunks <- .tile_chunks_create(
+            tile = tile,
+            overlap = overlap,
+            block = block
         )
-        # get the top most values
-        samples_tile <- samples_tile |>
-            # randomly shuffle the rows of the dataset
-            dplyr::sample_frac() |>
-            dplyr::slice_max(
-                .data[["value"]],
-                n = n,
-                with_ties = FALSE
+        # Tile path
+        tile_path <- .tile_path(tile)
+
+        # Process jobs in parallel
+        chunk_results <- .jobs_map_parallel_dfr(chunks, function(chunk) {
+            # Open tile images
+            r_obj <- .raster_open_rast(tile_path)
+
+            # Get values for this chunk only
+            values <- .raster_get_values(
+                rast = r_obj,
+                row = .block(chunk)[["row"]],
+                col = .block(chunk)[["col"]],
+                nrows = .block(chunk)[["nrows"]],
+                ncols = .block(chunk)[["ncols"]]
+            )
+
+            # Sample the maximum values in this chunk
+            samples_chunk <- C_max_sampling(
+                x = values,
+                nrows = .block(chunk)[["nrows"]],
+                ncols = .block(chunk)[["ncols"]],
+                window_size = sampling_window
             )
-        # transform to tibble
-        tb <- rast |>
-            .raster_xy_from_cell(
-                cell = samples_tile[["cell"]]
-            ) |>
-            tibble::as_tibble()
-        # find NA
-        na_rows <- which(is.na(tb))
-        # remove NA
-        if (.has(na_rows)) {
-            tb <- tb[-na_rows, ]
-            samples_tile <- samples_tile[-na_rows, ]
-        }
-        # Get the values' positions.
-        result_tile <- tb |>
-            sf::st_as_sf(
-                coords = c("x", "y"),
-                crs = .raster_crs(rast),
-                dim = "XY",
-                remove = TRUE
-            ) |>
-            sf::st_transform(crs = "EPSG:4326") |>
-            sf::st_coordinates()
 
-        colnames(result_tile) <- c("longitude", "latitude")
-        result_tile <- result_tile |>
-            dplyr::bind_cols(samples_tile) |>
-            dplyr::mutate(
-                value = .data[["value"]] *
-                    .conf("probs_cube_scale_factor")
-            ) |>
-            dplyr::filter(
-                .data[["value"]] >= min_uncert
-            ) |>
-            dplyr::select(dplyr::matches(
-                c("longitude", "latitude", "value")
-            )) |>
-            tibble::as_tibble()
+            # Skip empty chunks
+            if (nrow(samples_chunk) == 0) {
+                return(tibble::tibble(
+                    longitude = numeric(0),
+                    latitude = numeric(0),
+                    value = numeric(0)
+                ))
+            }
 
-        # All the cube's uncertainty images have the same start & end dates.
-        result_tile[["start_date"]] <- .tile_start_date(uncert_cube)
-        result_tile[["end_date"]] <- .tile_end_date(uncert_cube)
-        result_tile[["label"]] <- "NoClass"
-        result_tile
+            # Create a virtual raster object for the chunk
+            chunk_obj <- .chunks_as_raster(
+                chunk = chunk,
+                nlayers = 1L
+            )
+            
+            # transform to tibble using chunk coordinates
+            tb <- chunk_obj |>
+                .raster_xy_from_cell(
+                    cell = samples_chunk[["cell"]]
+                ) |>
+                tibble::as_tibble()
+
+            # find NA
+            na_rows <- which(is.na(tb))
+            # remove NA
+            if (.has(na_rows)) {
+                tb <- tb[-na_rows, ]
+                samples_chunk <- samples_chunk[-na_rows, ]
+            }
+
+            # Skip if all NA
+            if (nrow(tb) == 0) {
+                return(tibble::tibble(
+                    longitude = numeric(0),
+                    latitude = numeric(0),
+                    value = numeric(0)
+                ))
+            }
+
+            # Get the values' positions.
+            result_chunk <- tb |>
+                sf::st_as_sf(
+                    coords = c("x", "y"),
+                    crs = .raster_crs(chunk_obj),
+                    dim = "XY",
+                    remove = TRUE
+                ) |>
+                sf::st_transform(crs = "EPSG:4326") |>
+                sf::st_coordinates()
+
+            colnames(result_chunk) <- c("longitude", "latitude")
+            result_chunk <- result_chunk |>
+                dplyr::bind_cols(samples_chunk) |>
+                dplyr::mutate(
+                    value = .data[["value"]] *
+                        .conf("probs_cube_scale_factor")
+                ) |>
+                dplyr::filter(
+                    .data[["value"]] >= min_uncert,
+                    .data[["value"]] <= max_uncert
+                ) |>
+                dplyr::select(dplyr::matches(
+                    c("longitude", "latitude", "value")
+                )) |>
+                tibble::as_tibble()
+
+            result_chunk
+        }, progress = progress)
+
+        # Aggregate: select the top n values from all chunks in this tile
+        if (nrow(chunk_results) > 0) {
+            chunk_results |>
+                # randomly shuffle the rows of the dataset
+                dplyr::slice_sample(
+                    prop = 1
+                ) |>
+                dplyr::slice_max(
+                    .data[["value"]],
+                    n = n,
+                    with_ties = FALSE
+                ) |>
+                dplyr::mutate(
+                    start_date = .tile_start_date(uncert_cube),
+                    end_date = .tile_end_date(uncert_cube),
+                    label = "NoClass"
+                )
+        } else {
+            # Return empty result if no samples found
+            tibble::tibble(
+                longitude = numeric(0),
+                latitude = numeric(0),
+                value = numeric(0),
+                start_date = character(0),
+                end_date = character(0),
+                label = character(0)
+            )
+        }
     })
     renamed_cols <- c(uncertainty = "value")
     samples_tb <- dplyr::rename(samples_tb, dplyr::all_of(renamed_cols))