diff --git a/.buildlibrary b/.buildlibrary index 5c4ef32..0eeed70 100644 --- a/.buildlibrary +++ b/.buildlibrary @@ -1,4 +1,4 @@ -ValidationKey: '3437592' +ValidationKey: '3458270' AutocreateReadme: yes AcceptedWarnings: - 'Warning: package ''.*'' was built under R version' diff --git a/CITATION.cff b/CITATION.cff index c9a0539..35f5713 100644 --- a/CITATION.cff +++ b/CITATION.cff @@ -2,8 +2,8 @@ cff-version: 1.2.0 message: If you use this software, please cite it using the metadata from this file. type: software title: 'lpjmlkit: Toolkit for Basic LPJmL Handling' -version: 1.7.2 -date-released: '2024-09-20' +version: 1.7.3 +date-released: '2024-09-24' abstract: A collection of basic functions to facilitate the work with the Dynamic Global Vegetation Model (DGVM) Lund-Potsdam-Jena managed Land (LPJmL) hosted at the Potsdam Institute for Climate Impact Research (PIK). It provides functions for diff --git a/DESCRIPTION b/DESCRIPTION index d4e0a54..68e98be 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Package: lpjmlkit Type: Package Title: Toolkit for Basic LPJmL Handling -Version: 1.7.2 +Version: 1.7.3 Authors@R: c( person("Jannes", "Breier", , "jannesbr@pik-potsdam.de", role = c("aut", "cre"), comment = c(ORCID = "0000-0002-9055-6904")), person("Sebastian","Ostberg", , "ostberg@pik-potsdam.de", role = "aut", comment = c(ORCID = "0000-0002-2368-7015")), @@ -55,4 +55,4 @@ Suggests: sf Config/testthat/edition: 3 VignetteBuilder: knitr -Date: 2024-09-20 +Date: 2024-09-24 diff --git a/R/get_cellindex.R b/R/get_cellindex.R index ff04990..3721508 100644 --- a/R/get_cellindex.R +++ b/R/get_cellindex.R @@ -12,10 +12,13 @@ #' @param extent A numeric vector (lonmin, lonmax, latmin, latmax) containing the #' longitude and latitude boundaries between which values included in the subset. #' @param coordinates A list of two named (lon, lat) numeric vectors representing the coordinates. +#' @param shape A terra SpatVector object in the WGS 84 coordinate reference system +#' representing the shape to subset the grid. +#' @param simplify A logical indicating whether to simplify the output to a vector. #' #' -#' @return The cell index from the grid file based on the provided extent or -#' coordinates. +#' @return Either an LPJmLData object containing the grid or a vector subsetted +#' to the provided extent, coordinates or shape. #' @export #' #' @examples @@ -43,28 +46,36 @@ #' should be a list of two character vectors representing the longitude and #' latitude values as for [`subset()`]. #' -#' If both `extent` and `coordinates` are provided, the function will stop and -#' ask for only one of them. If neither `extent` nor `coordinates` are provided, -#' the function will return the cell numbers for all cells in the grid. +#' If a shape is provided as a SpatVector object, the function will return the +#' cell index for the cells that intersect with the shape. +#' +#' If more than on of `extent`, `coordinates` `shape` are provided, the function +#' will stop and ask for only one of them. If neither `extent` nor `coordinates` +#' nor `shape` are provided, the function will return the cell numbers for all +#' cells in the grid. #' #' The function also includes checks for input types and values, and gives #' specific error messages for different error conditions. For example, it #' checks if the `grid_filename` exists, if the `extent` vector has the correct #' length, and if the `coordinates` list contains two vectors of equal length. -get_cellindex <- function(grid_filename, extent = NULL, coordinates = NULL) { +get_cellindex <- function(grid_filename, extent = NULL, coordinates = NULL, shape = NULL, simplify = TRUE) { # check if filepath is valid check_filepath(grid_filename) # check if (only) one of extent or coordinates is provided - check_extent_and_coordinates(extent, coordinates) + check_multiple(extent, coordinates, shape) - grid_lonlat <- read_io(filename = grid_filename) %>% + grid_lonlat <- read_io(filename = grid_filename) |> LPJmLGridData$new() if (!is.null(extent)) { - extent <- check_extent(extent) %>% + extent <- check_extent(extent) |> correct_extent() } else if (!is.null(coordinates)) { coordinates <- check_coordinates(coordinates) + } else if (!is.null(shape)) { + if (!("SpatVector" %in% class(shape))) { + stop("shape must be a SpatVector object.") + } } # Read the grid file and create a data frame @@ -76,7 +87,7 @@ get_cellindex <- function(grid_filename, extent = NULL, coordinates = NULL) { # Check if extent values are within the longitude and latitude range in the cells if (!is.null(extent)) { - cells$cellindex <- as.numeric(row.names(cells)) + 1 + cells$cellindex <- as.numeric(row.names(cells)) out_of_bounds_lon <- extent[c(1, 2)][extent[c(1, 2)] < lon_range[1] | extent[c(1, 2)] > lon_range[2]] @@ -94,7 +105,14 @@ get_cellindex <- function(grid_filename, extent = NULL, coordinates = NULL) { cells <- cells[cells$lon >= extent[1] & cells$lon <= extent[2] & - cells$lat >= extent[3] & cells$lat <= extent[4], ]$cellindex + cells$lat >= extent[3] & cells$lat <= extent[4], ] + + grid_cell <- transform(grid_lonlat, "lon_lat") + + cells <- grid_cell$subset(coordinates = lapply(X = list(lon = cells$lon, + lat = cells$lat), + FUN = as.character)) + } # Check if coordinates are within the longitude and latitude range in the cells @@ -120,9 +138,8 @@ get_cellindex <- function(grid_filename, extent = NULL, coordinates = NULL) { grid_cell <- transform(grid_lonlat, "lon_lat") - grid_cell$subset(coordinates = lapply(X = coordinates, FUN = as.character)) - - cells <- c(stats::na.omit(c(grid_cell$data + 1))) + cells <- grid_cell$subset(coordinates = lapply(X = coordinates, + FUN = as.character)) message( col_note( @@ -131,6 +148,22 @@ get_cellindex <- function(grid_filename, extent = NULL, coordinates = NULL) { ) } + if (!is.null(shape)) { + cell_coords <- grid_lonlat |> + as_terra() |> + terra::mask(shape) |> + terra::as.data.frame(xy = TRUE) |> + dplyr::select("x", "y") + + cells <- grid_lonlat$transform("lon_lat") |> + subset(coordinates = lapply(list(lon = cell_coords$x, lat = cell_coords$y), + FUN = as.character)) + } + + if (simplify) { + cells <- c(stats::na.omit(c(cells$data + 1))) + } + cells } @@ -167,14 +200,16 @@ check_extent <- function(extent) { # Check if both extent and coordinates are provided -check_extent_and_coordinates <- function(extent, coordinates) { - if (is.null(extent) && is.null(coordinates)) { - warning("Neither extent or coordinates provided. Full grid will be returned.") +check_multiple <- function(extent, coordinates, shape) { + if (is.null(extent) && is.null(coordinates) && is.null(shape)) { + warning("Neither extent, coordinates or shape provided. Full grid will be returned.") } - if (!is.null(extent) && !is.null(coordinates)) { + if ((!is.null(extent) && !is.null(coordinates)) || + (!is.null(extent) && !is.null(shape)) || + (!is.null(coordinates) && !is.null(shape))) { stop( - "Both extent and coordinates are provided.", - " Please provide only one of them." + "Multiple subset options provided.", + " Please provide only one of coordinates, extent and shape." ) } } diff --git a/README.md b/README.md index efcde9f..1d59b75 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ # logo Toolkit for Basic LPJmL Handling -R package **lpjmlkit**, version **1.7.2** +R package **lpjmlkit**, version **1.7.3** [![CRAN status](https://www.r-pkg.org/badges/version/lpjmlkit)](https://cran.r-project.org/package=lpjmlkit) [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.7773134.svg)](https://doi.org/10.5281/zenodo.7773134) [![R build status](https://github.com/PIK-LPJmL/lpjmlkit/workflows/check/badge.svg)](https://github.com/PIK-LPJmL/lpjmlkit/actions) [![codecov](https://codecov.io/gh/PIK-LPJmL/lpjmlkit/branch/master/graph/badge.svg)](https://app.codecov.io/gh/PIK-LPJmL/lpjmlkit) [![r-universe](https://pik-piam.r-universe.dev/badges/lpjmlkit)](https://pik-piam.r-universe.dev/builds) @@ -76,7 +76,7 @@ In case of questions / problems please contact Jannes Breier , R package version 1.7.2, . +Breier J, Ostberg S, Wirth S, Minoli S, Stenzel F, Hötten D, Müller C (2024). _lpjmlkit: Toolkit for Basic LPJmL Handling_. doi:10.5281/zenodo.7773134 , R package version 1.7.3, . A BibTeX entry for LaTeX users is @@ -85,7 +85,7 @@ A BibTeX entry for LaTeX users is title = {lpjmlkit: Toolkit for Basic LPJmL Handling}, author = {Jannes Breier and Sebastian Ostberg and Stephen Björn Wirth and Sara Minoli and Fabian Stenzel and David Hötten and Christoph Müller}, year = {2024}, - note = {R package version 1.7.2}, + note = {R package version 1.7.3}, url = {https://github.com/PIK-LPJmL/lpjmlkit}, doi = {10.5281/zenodo.7773134}, } diff --git a/man/get_cellindex.Rd b/man/get_cellindex.Rd index 65198dd..1973d85 100644 --- a/man/get_cellindex.Rd +++ b/man/get_cellindex.Rd @@ -4,7 +4,13 @@ \alias{get_cellindex} \title{Get Cell Index} \usage{ -get_cellindex(grid_filename, extent = NULL, coordinates = NULL) +get_cellindex( + grid_filename, + extent = NULL, + coordinates = NULL, + shape = NULL, + simplify = TRUE +) } \arguments{ \item{grid_filename}{A string representing the grid file name.} @@ -13,10 +19,15 @@ get_cellindex(grid_filename, extent = NULL, coordinates = NULL) longitude and latitude boundaries between which values included in the subset.} \item{coordinates}{A list of two named (lon, lat) numeric vectors representing the coordinates.} + +\item{shape}{A terra SpatVector object in the WGS 84 coordinate reference system +representing the shape to subset the grid.} + +\item{simplify}{A logical indicating whether to simplify the output to a vector.} } \value{ -The cell index from the grid file based on the provided extent or -coordinates. +Either an LPJmLData object containing the grid or a vector subsetted +to the provided extent, coordinates or shape. } \description{ This function returns the cell index from a grid file based on @@ -41,9 +52,13 @@ include only those that match the specified coordinates. The \code{coordinates} should be a list of two character vectors representing the longitude and latitude values as for \code{\link[=subset]{subset()}}. -If both \code{extent} and \code{coordinates} are provided, the function will stop and -ask for only one of them. If neither \code{extent} nor \code{coordinates} are provided, -the function will return the cell numbers for all cells in the grid. +If a shape is provided as a SpatVector object, the function will return the +cell index for the cells that intersect with the shape. + +If more than on of \code{extent}, \code{coordinates} \code{shape} are provided, the function +will stop and ask for only one of them. If neither \code{extent} nor \code{coordinates} +nor \code{shape} are provided, the function will return the cell numbers for all +cells in the grid. The function also includes checks for input types and values, and gives specific error messages for different error conditions. For example, it diff --git a/tests/testthat/test-get_cellindex.R b/tests/testthat/test-get_cellindex.R index ec3c9e5..3e581f2 100644 --- a/tests/testthat/test-get_cellindex.R +++ b/tests/testthat/test-get_cellindex.R @@ -26,7 +26,7 @@ test_that("get_cellindex handles both extent and coordinates provided", { extent = c(1.25, 2.75, 3.25, 4.75), coordinates = list(c(1.25, 2.75), c(1.25, 2.75)) ), - "Both extent and coordinates are provided. Please provide only one of them." + "Multiple subset options provided. Please provide only" ) }) @@ -62,3 +62,11 @@ test_that("get_cellindex returns correct cell index for given coordinates", { ) expect_true(length(result) == 2 && result[1] == 10001 && result[2] == 10002) }) + +test_that("get_cellindex returns correct cell index for a given shape", { + result <- get_cellindex("../testdata/output/grid.bin.json", + shape = terra::vect(terra::ext(c(-87.25, -87.25, + 55.25, 55.75)), + crs = "+proj=longlat +datum=WGS84")) + expect_true(length(result) == 2 && result[1] == 10002 && result[2] == 10001) +})