Skip to content

Commit

Permalink
Modify extent formats and names
Browse files Browse the repository at this point in the history
  • Loading branch information
ghislainv committed Aug 29, 2023
1 parent 68ad5a3 commit d11db50
Show file tree
Hide file tree
Showing 6 changed files with 196 additions and 196 deletions.
26 changes: 14 additions & 12 deletions R/create_xml_legend.R
Original file line number Diff line number Diff line change
@@ -1,14 +1,14 @@
create_xml_legend <- function(unique_values, destination, name_file){
#' Create .tif.aux.xml legend for soilgrids layer
#'
#' @description
#' internal function.
#'
#' @param unique_values int vector. with unique values in soilgrids layer.
#' @param destination character. absolute path where create .tif.aux.xml file.
#' @param name_file character. Name of the file to create, `gecevar` or `environ`.
#' @import rjson
#' @export
#' Create .tif.aux.xml legend for soilgrids layer
#'
#' @description
#' internal function.
#'
#' @param unique_values int vector. with unique values in soilgrids layer.
#' @param destination character. absolute path where create .tif.aux.xml file.
#' @param name_file character. Name of the file to create, `gecevar` or `environ`.
#' @import rjson

create_xml_legend <- function(unique_values, destination, name_file) {

json <- fromJSON(file = "https://files.isric.org/soilgrids/latest/data/wrb/MostProbable.rat.json")
soil <- unlist(json)[-1]
Expand All @@ -32,7 +32,7 @@ create_xml_legend <- function(unique_values, destination, name_file){
</PAMDataset>'


for (i in 1:length(unique_values)){
for (i in 1:length(unique_values)) {
row <- glue(' <Row index="{i - 1}">
<F>{unique_values[i]}</F>
<F>{soil[i + 1]}</F>
Expand All @@ -42,3 +42,5 @@ create_xml_legend <- function(unique_values, destination, name_file){
header <- paste(header, footer, sep = '\n')
writeLines(header, paste(destination, paste0(name_file, ".tif.aux.xml"), sep = "/"))
}

# End
55 changes: 21 additions & 34 deletions R/get_aoi_extent.R
Original file line number Diff line number Diff line change
@@ -1,21 +1,11 @@
#!/usr/bin/R

## ==============================================================================
## authors :Pierre Guillaumont, Ghislain Vieilledent
## emails :pierreguillaumont@gmail.com,
## ghislain.vieilledent@cirad.fr, ghislainv@gmail.com
## web :https://ecology.ghislainv.fr
## license :GPLv3
## ==============================================================================

#' Get the extent of an area of interest (AOI) in GDAL format (xmin, ymin, xmax, ymax).
#'
#' This function gets the extent of an area of interest (AOI) in GDAL
#' format (xmin, ymin, xmax, ymax). The extent is provided both in
#' latitude and longitude and in the projected coordinate reference
#' system (CRS) specified by the user. The area of interest can be
#' specified either by a country iso3 code, a vector file, or an
#' extent (in latlong or projected CRS). When the area of interest
#' extent (in latlon or projected CRS). When the area of interest
#' is specified with a country iso3 code, the country borders are
#' downloaded from the GADM (Global ADMinistrative areas database)
#' at \url{https://gadm.org}.
Expand Down Expand Up @@ -49,9 +39,10 @@
#' @param rm_output_dir boolean. If FALSE, does not remove the output
#' directory. Default is TRUE.
#'
#' @return list. With `e_proj`: projected extent, `e_latlong`: extent
#' in latlong, and `gpkg_file`: path to .gpkg country vector file if
#' downloaded.
#' @return list. With `extent_latlon`: extent in latitude and
#' longitude (lonmin, latmin, lonmax, latmax),
#' `extent_proj`: projected extent (xmin, ymin, xmax, ymax),
#' and `gpkg_file`: path to .gpkg country vector file if downloaded.
#'
#' @import sf
#' @importFrom utils download.file
Expand All @@ -71,7 +62,7 @@ get_aoi_extent <- function(country_iso=NULL,
as.numeric(!is.null(vector_file)) +
as.numeric(!is.null(extent_proj)) +
as.numeric(!is.null(extent_latlon))) != 1) {
stop("One attribute among 'country_iso', 'vector_file', 'extent_proj', or 'exent_latlong' must be specified")
stop("One attribute among 'country_iso', 'vector_file', 'extent_proj', or 'exent_latlon' must be specified")
}

# Set gpkg_file to NULL
Expand All @@ -89,13 +80,13 @@ get_aoi_extent <- function(country_iso=NULL,
borders <- sf::st_read(gpkg_file, layer="ADM_ADM_0", quiet=TRUE)
# Get bounding box
bbox <- sf::st_bbox(borders)
# e_latlong extending to nearest degree
# e_latlon extending to nearest degree
lonmin <- floor(as.numeric(bbox$xmin))
lonmax <- ceiling(as.numeric(bbox$xmax))
latmin <- floor(as.numeric(bbox$ymin))
latmax <- ceiling(as.numeric(bbox$ymax))
e_latlong <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlong
e_latlon <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlon
e <- terra::ext(lonmin, lonmax, latmin, latmax)
e <- terra::as.polygons(e, crs="epsg:4326")
e_proj <- sf::st_bbox(terra::project(e, paste0("epsg:", EPSG_proj)))
Expand All @@ -108,18 +99,18 @@ get_aoi_extent <- function(country_iso=NULL,

# With AOI vector file
if (!is.null(vector_file)) {
# Reproject in latlong
# Reproject in latlon
borders <- sf::st_read(vector_file, quiet=TRUE)
borders <- sf::st_transform(borders, "EPSG:4326")
# Get bounding box
bbox <- sf::st_bbox(borders)
# e_latlong extending to nearest degree
# e_latlon extending to nearest degree
lonmin <- floor(as.numeric(bbox$xmin))
lonmax <- ceiling(as.numeric(bbox$xmax))
latmin <- floor(as.numeric(bbox$ymin))
latmax <- ceiling(as.numeric(bbox$ymax))
e_latlong <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlong
e_latlon <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlon
e <- terra::ext(lonmin, lonmax, latmin, latmax)
e <- terra::as.polygons(e, crs="epsg:4326")
e_proj <- sf::st_bbox(terra::project(e, paste0("epsg:", EPSG_proj)))
Expand All @@ -128,20 +119,20 @@ get_aoi_extent <- function(country_iso=NULL,

# With projected extent
if (!is.null(extent_proj)) {
# Get the bounding box in latlong
# Get the bounding box in latlon
e <- terra::ext(as.numeric(extent_proj[1]),
as.numeric(extent_proj[3]),
as.numeric(extent_proj[2]),
as.numeric(extent_proj[4]))
e <- terra::as.polygons(e, crs=paste0("epsg:", EPSG_proj))
bbox <- sf::st_bbox(terra::project(e, "epsg:4326"))
# e_latlong extending to nearest degree
# e_latlon extending to nearest degree
lonmin <- floor(as.numeric(bbox$xmin))
lonmax <- ceiling(as.numeric(bbox$xmax))
latmin <- floor(as.numeric(bbox$ymin))
latmax <- ceiling(as.numeric(bbox$ymax))
e_latlong <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlong
e_latlon <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlon
e <- terra::ext(lonmin, lonmax, latmin, latmax)
e <- terra::as.polygons(e, crs="epsg:4326")
e_proj <- sf::st_bbox(terra::project(e, paste0("epsg:", EPSG_proj)))
Expand All @@ -150,25 +141,21 @@ get_aoi_extent <- function(country_iso=NULL,

# With extent in lat long
if (!is.null(extent_latlon)) {
# e_latlong extending to nearest degree
# e_latlon extending to nearest degree
lonmin <- floor(as.numeric(extent_latlon[1]))
lonmax <- ceiling(as.numeric(extent_latlon[3]))
latmin <- floor(as.numeric(extent_latlon[2]))
latmax <- ceiling(as.numeric(extent_latlon[4]))
e_latlong <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlong
e_latlon <- c(lonmin=lonmin, latmin=latmin, lonmax=lonmax, latmax=latmax)
# e_proj from e_latlon
e <- terra::ext(lonmin, lonmax, latmin, latmax)
e <- terra::as.polygons(e, crs="epsg:4326")
e_proj <- sf::st_bbox(terra::project(e, paste0("epsg:", EPSG_proj)))
e_proj <- c(floor(e_proj$xmin), floor(e_proj$ymin), ceiling(e_proj$xmax), ceiling(e_proj$ymax))
}

# Collapse vectors to characters
e_latlong <- paste(e_latlong, collapse=" ")
e_proj <- paste(e_proj, collapse=" ")

# Return results
return(list(e_proj=e_proj, e_latlong=e_latlong, gpkg_file=gpkg_file))
return(list(extent_latlon=e_latlon, extent_proj=e_proj, gpkg_file=gpkg_file))

}

Expand Down
Loading

0 comments on commit d11db50

Please sign in to comment.