Skip to content

Commit

Permalink
Merge branch 'devel' of github.com:SantanderMetGroup/drought4R into d…
Browse files Browse the repository at this point in the history
…evel
  • Loading branch information
miturbide committed Jun 27, 2019
2 parents cb296cc + 2bbe393 commit a98b4c7
Show file tree
Hide file tree
Showing 21 changed files with 485 additions and 39 deletions.
11 changes: 6 additions & 5 deletions DESCRIPTION
Original file line number Diff line number Diff line change
Expand Up @@ -6,22 +6,23 @@ Imports:
transformeR(>= 1.4.4),
utils,
abind,
magrittr
magrittr,
geosphere
Suggests:
loadeR,
visualizeR,
convertR
Type: Package
Title: A climate4R package for drought assessment <http://www.meteo.unican.es/climate4R>
Version: 0.2.0
Date: 2019-04-24
Title: A climate4R package for PET estimation and drought assessment <http://www.meteo.unican.es/climate4R>
Version: 0.3.0
Date: 2019-06-17
Authors@R: as.person(c(
"Santander Meteorology Group <http://meteo.unican.es> [cph]",
"Joaquín Bedia <[email protected]> [aut, cre]",
"Maialen Iturbide <[email protected]> [aut]"))
BugReports: https://github.com/SantanderMetGroup/drought4R/issues
URL: https://github.com/SantanderMetGroup/transformeR/wiki
Description: Evapotranspiration routines and drought index calculation using the climate4R framework <http://www.meteo.unican.es/climate4R>.
Description: Potential evapotranspiration methods and drought index calculation using the climate4R framework <http://www.meteo.unican.es/climate4R>.
License: file LICENSE
LazyData: true
Encoding: UTF-8
Expand Down
6 changes: 6 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
@@ -1,18 +1,23 @@
# Generated by roxygen2: do not edit by hand

export(effectiveTempGrid)
export(petGrid)
export(photoperiodGrid)
export(speiGrid)
import(transformeR)
importFrom(SPEI,hargreaves)
importFrom(SPEI,spei)
importFrom(SPEI,thornthwaite)
importFrom(abind,abind)
importFrom(geosphere,daylength)
importFrom(magrittr,"%<>%")
importFrom(magrittr,"%>%")
importFrom(magrittr,extract)
importFrom(magrittr,extract2)
importFrom(stats,ts)
importFrom(transformeR,aggregateGrid)
importFrom(transformeR,array3Dto2Dmat)
importFrom(transformeR,bindGrid)
importFrom(transformeR,checkDim)
importFrom(transformeR,checkTemporalConsistency)
importFrom(transformeR,getCoordinates)
Expand All @@ -21,6 +26,7 @@ importFrom(transformeR,getRefDates)
importFrom(transformeR,getSeason)
importFrom(transformeR,getShape)
importFrom(transformeR,getTimeResolution)
importFrom(transformeR,gridArithmetics)
importFrom(transformeR,mat2Dto3Darray)
importFrom(transformeR,redim)
importFrom(transformeR,subsetGrid)
Expand Down
8 changes: 8 additions & 0 deletions NEWS
Original file line number Diff line number Diff line change
Expand Up @@ -20,3 +20,11 @@
* Implements the new PET method Hargreaves-Samani for daily data
* Other minor improvements
* Removed vignette (now linked to paper supplementary information http://www.meteo.unican.es/work/climate4r/drought4R/drought4R_notebook.html)

## v0.3.0 (17 Jun 2019)
* Implementation of the calibrated Thorthwaite's method (Pereira and Pruitt 2004, DOI:10.1016/j.agwat.2003.11.003)
* New auxiliary functions:
* `effectiveTempGrid`, for effective temperature calculation
* `photoperiodGrid`, for daylength estimation
* New built-in daily datasets for the Iberian peninsula of tmax, tmin and precip (E-OBSv17 0.25 reg, year 2000)
* Other minor changes and documentation updates
50 changes: 50 additions & 0 deletions R/effectiveTempGrid.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,50 @@
## effectiveTemp.R Calculation of the effective temperature for improving the Thornthwaite's ET0 estimation method
##
## Copyright (C) 2019 Santander Meteorology Group (http://www.meteo.unican.es)
##
## This program is free software: you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program. If not, see <http://www.gnu.org/licenses/>.

#' @title Effective temperature
#' @description Calculation of the \dQuote{effective temperature} for improving the Thornthwaite's ET0 estimation method
#' @param tasmax Grid of maximum monthly temperature (degC)
#' @param tasmin Grid of minimum monthly temperature (degC)
#' @param k Calibration coefficient. Default to 0.69, as proposed by Pereira and Pruitt (2004)
#' @references \itemize{
#' \item de Camargo, A.P., Marin, F.R., Sentelhas, P. C., Picini, A.G., 1999. Adjust of the Thornthwaite’s method to estimate the potential evapotranspiration for arid and superhumid climates, based on daily temperature amplitude. Revista Brasileira de Agrometeorologia 7.
#' \item Pereira, A.R., Pruitt, W.O., 2004. Adaptation of the Thornthwaite scheme for estimating daily reference evapotranspiration. Agricultural Water Management 66, 251–257. https://doi.org/10.1016/j.agwat.2003.11.003
#' }
#' @details The function is internally used by \code{\link{petGrid}} when the Thornthwaite's method with calibration coefficient k is chosen.
#' @seealso \code{\link{petGrid}}
#' @author J. Bedia
#' @export

effectiveTempGrid <- function(tasmin = NULL, tasmax = NULL, k = 0.69) {
if (is.null(tasmin) | is.null(tasmax)) {
stop("Both \'tasmin\' and \'tasmax\' arguments are required for effective temperature computation", call. = FALSE)
}
if (is.null(k)) stop("A calibration factor \'k\' is required", call. = FALSE)
checkDim(tasmin, tasmax, dimensions = c("time", "lat", "lon"))
checkTemporalConsistency(tasmin, tasmax)
if ((getTimeResolution(tasmax) != "DD" | getTimeResolution(tasmin) != "DD")) stop("Daily data are required for effective temperature computation", call. = FALSE)
message("[", Sys.time(), "] Calculating effective temperature...")
tasmin <- gridArithmetics(tasmax, 3, tasmin, operator = c("*", "-")) %>% gridArithmetics(., 0.5 * k, operator = "*")
tasmax <- NULL
tasmin$Variable$varName <- "Effective_temperature"
attr(tasmin$Variable, "description") <- "effective temperature for Thorthwaite's method calibration"
attr(tasmin$Variable, "longname") <- "effective_temperature"
attr(tasmin, "R_package_desc") <- paste0("drought4R-v", packageVersion("drought4R"))
attr(tasmin, "R_package_URL") <- "https://github.com/SantanderMetGroup/drought4R"
attr(tasmin, "R_package_ref") <- "http://dx.doi.org/10.1016/j.envsoft.2018.09.009"
return(tasmin)
}
83 changes: 69 additions & 14 deletions R/petGrid.R
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,9 @@
#' @param method Potential evapotranspiration method. Currently \code{"thornthwaite"} and
#' \code{"hargreaves"} methods are available (monthly), using the implementation of package \pkg{SPEI}.
#' In addition, \code{"hargreaves-samani"} is available for daily data. See details.
#' @param k Optional calibration coefficient for the Thornthwaite method. Unused by default. See Details.
#' @param photoperiod.corr Correction with the day-night ratio applied to effective temperature (default to \code{TRUE}).
#' Only applies to the calibrated Thorthwaite method.
#' @param what Optional character string, only applied for the Hargreaves-Samani method.
#' If set to \code{what = "rad"}, it returns the estimated radiation (it is a function of latitude and date).
#' Otherwise, by default, returns the estimated daily potential evapotranspiration.
Expand All @@ -33,12 +36,33 @@
#' from the \pkg{SPEI} package (Begueria and Vicente-Serrano). Monthly input data are thus required.
#' The latitude of the sites/grid points is internally used.
#' In case of multimember grids (e.g. seasonal forecast data), the PET is calculated for each member sepparately.
#'
#' \strong{Calibration coefficient for Thornthwaite}
#'
#' The use of a calibration coefficient (\code{k.th} argument) can provide better PET estimates under certain conditions. For instance,
#' Camargo \emph{et al.} (1999) found that a value of k=0.72 is the best for estimating monthly ET0, while Pereira and Pruitt (2004)
#' recommended k=0.69 for daily ET0. Trajkovic et al. (2019) propose an optimal calibration factor of 0.65 after an intercomparison of
#' several PET estimation methods and calibration coefficients in northern Serbia.
#'
#' In addition, the correction proposed by Pereira and Pruitt (equation 7) to account for the effect of differing photoperiods
#' can be omitted through the logical argument \code{photoperiod.corr = FALSE} (it is activated by default).
#'
#' \strong{Note:} the calibration factor for the Thorthwaite method requires minimum and maximum temperatures instead of mean temperature,
#' as it also includes temperature amplitude in its formulation.
#'
#' @importFrom transformeR array3Dto2Dmat mat2Dto3Darray checkDim getCoordinates getDim
#' @importFrom utils packageVersion
#' @importFrom abind abind
#' @importFrom magrittr %>%
#' @author J Bedia
#' @export
#' @references
#' \itemize{
#' \item Camargo AP, Marin FR, Sentelhas PC, Picini AG (1999) Adjust of the Thornthwaite’s method to estimate the potential evapotranspiration for
#' arid and superhumid climates, based on daily temperature amplitude. Rev Bras Agrometeorol 7(2):251–257
#' \item Pereira, A.R., Pruitt, W.O., 2004. Adaptation of the Thornthwaite scheme for estimating daily reference evapotranspiration. Agricultural Water Management 66, 251–257. https://doi.org/10.1016/j.agwat.2003.11.003
#' \item Trajkovic, S., Gocic, M., Pongracz, R., Bartholy, J., 2019. Adjustment of Thornthwaite equation for estimating evapotranspiration in Vojvodina. Theor Appl Climatol. https://doi.org/10.1007/s00704-019-02873-1
#' }
#' @examples \donttest{
#' # Thorthwaite requires monthly mean temperature data as input:
#' data("tas.cru.iberia")
Expand Down Expand Up @@ -87,11 +111,13 @@ petGrid <- function(tasmin = NULL,
pr = NULL,
method = c("thornthwaite", "hargreaves", "hargreaves-samani"),
what = c("PET", "rad"),
k = NULL,
photoperiod.corr = TRUE,
...) {
method <- match.arg(method, choices = c("thornthwaite", "hargreaves", "hargreaves-samani"))
message("[", Sys.time(), "] Computing PET-", method, " ...")
out <- switch(method,
"thornthwaite" = petGrid.th(tas, ...),
"thornthwaite" = petGrid.th(tas, tasmin, tasmax, k, phc = photoperiod.corr, ...),
"hargreaves" = petGrid.har(tasmin, tasmax, pr, ...),
"hargreaves-samani" = petGrid.hs(tasmin, tasmax, what))
message("[", Sys.time(), "] Done")
Expand All @@ -112,31 +138,60 @@ petGrid <- function(tasmin = NULL,
} else if (tres == "DD") {
attr(pet.grid$Variable, "units") <- "mm.day-1"
}
if (method == "thornthwaite" | method == "hargreaves") {
attr(pet.grid, "origin") <- paste0("Calculated with R package 'SPEI' v",
packageVersion("SPEI"), " using 'drought4R' v",
packageVersion("drought4R"))
} else {
attr(pet.grid, "origin") <- paste0("Calculated with 'drought4R' v",
packageVersion("drought4R"))
}
attr(pet.grid, "URL") <- "https://github.com/SantanderMetGroup/drought4R"
attr(pet.grid, "R_package_desc") <- paste0("drought4R-v", packageVersion("drought4R"))
attr(pet.grid, "R_package_URL") <- "https://github.com/SantanderMetGroup/drought4R"
attr(pet.grid, "R_package_ref") <- "http://dx.doi.org/10.1016/j.envsoft.2018.09.009"
pet.grid <- redim(pet.grid, drop = TRUE)
invisible(pet.grid)
}



#' @importFrom SPEI thornthwaite
#' @importFrom transformeR getCoordinates getSeason array3Dto2Dmat getTimeResolution redim getShape subsetGrid
#' @importFrom transformeR getCoordinates getSeason array3Dto2Dmat getTimeResolution redim getShape subsetGrid gridArithmetics checkTemporalConsistency checkDim aggregateGrid
#' @importFrom magrittr %>% %<>%
#' @keywords internal
#' @author J Bedia

petGrid.th <- function(tas, ...) {
petGrid.th <- function(tas, tasmin, tasmax, k, phc, ...) {
if (is.null(tas)) {
stop("Mean temperature grid is required by Thornthwaite method", call. = FALSE)
if (is.null(tasmin) | is.null(tasmax)) {
stop("\'tas\' has been set to NULL. Therefore, both \'tasmin\' and \'tasmax\' arguments are required by Thornthwaite method", call. = FALSE)
}
checkTemporalConsistency(tasmin, tasmax)
checkDim(tasmin, tasmax, dimensions = c("time", "lat", "lon"))
if (getTimeResolution(tasmin) == "MM") {
if (!is.null(k)) message("Input data is monthly. Argument \'k\' will be ignored (no calibration will be applied)")
tas <- gridArithmetics(tasmax, tasmin, 2, operator = c("+", "/"))
} else if (getTimeResolution(tasmin) == "DD") {
tasmean <- gridArithmetics(tasmax, tasmin, 2, operator = c("+", "/"))
if (!is.null(k)) {
tas <- effectiveTempGrid(tasmin = tasmin, tasmax = tasmax, k)
if (isTRUE(phc)) {
ph <- photoperiodGrid(tas)
dn.ratio <- gridArithmetics(ph, ph, 24, ph, operator = c("-","+","-"))
tas <- gridArithmetics(ph, dn.ratio, tas, operator = c("/", "*"))
ind.sup <- which(tas$Data > tasmax$Data)
ind.inf <- which(tas$Data < tasmean$Data)
tasmean <- NULL
tas$Data[ind.sup] <- tasmax$Data[ind.sup]
tas$Data[ind.inf] <- tasmin$Data[ind.inf]
tasmax <- tasmin <- NULL
}
suppressMessages({tas %<>% aggregateGrid(aggr.m = list(FUN = "mean", na.rm = TRUE))})
} else {
suppressMessages({tas <- aggregateGrid(tasmean, aggr.m = list(FUN = "mean", na.rm = TRUE))})
}
} else {
stop("Invalid temporal resolution of input data", call. = FALSE)
}
} else {
if (!is.null(tasmin) | !is.null(tasmax)) {
message("\'tas\' argument has been provided. Therefore, \'tasmin\' and \'tasmax\' arguments will be ignored.")
}
}
if (!identical(1:12, getSeason(tas))) stop("The input grid must encompass a whole-year season")
if (getTimeResolution(tas) != "MM") stop("A monthly input grid is required by the Thornthwaite method")
# if (getTimeResolution(tas) != "MM") stop("A monthly input grid is required by the Thornthwaite method", call. = FALSE)
tas <- redim(tas, member = TRUE)
ref.grid <- tas
coords <- getCoordinates(tas)
Expand Down
68 changes: 68 additions & 0 deletions R/photoperiodGrid.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,68 @@
## photoperiod.R Compute photoperiod (daylength) as a funtion of latitude and date
##
## Copyright (C) 2019 Santander Meteorology Group (http://www.meteo.unican.es)
##
## This program is free software: you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program. If not, see <http://www.gnu.org/licenses/>.

#' @title Photoperiod
#' @description Compute photoperiod (daylength) as a funtion of latitude and date
#' @param c4r.obj climate4R object. Note that this is independent of the variable stored, as it will use only the calendar
#' an geographical information of the object to compute the photoperiod
#' @importFrom geosphere daylength
#' @importFrom transformeR getRefDates getCoordinates getShape redim subsetGrid mat2Dto3Darray bindGrid
#' @importFrom magrittr %>% %<>%
#' @details This is just a wrapper from the \code{\link[geosphere]{daylength}} function of package \pkg{geosphere}.
#' @export
#' @author J. Bedia
#' @examples
#' require(transformeR)
#' data("tasmin.eobs.iberia.daily")
#' daylength <- photoperiodGrid(tasmin.eobs.iberia.daily)
#' lat <- range(getCoordinates(daylength)[["y"]])
#' lat.ranges <- seq(from = lat[1], to = tail(lat, 1), length.out = 5)
#' lat.list <- lapply(1:(length(lat.ranges) - 1), function(x) {
#' subsetGrid(daylength, latLim = c(lat.ranges[x], lat.ranges[x + 1]), outside = TRUE)
#' })
#' require(visualizeR)
#' temporalPlot(lat.list, xyplot.custom = list(main = "Photoperiod as a function of latitude",
#' ylab = "Daylight hours",
#' xlab = "Date",
#' key = list(corner = c(1,.5),
#' lines = list(col = 1:4),
#' horizontal = TRUE,
#' text = list(c("37ºN","39ºN","41ºN","43ºN")))))

photoperiodGrid <- function(c4r.obj) {
message("[", Sys.time(), "] Calculating photoperiod...")
j <- getRefDates(c4r.obj) %>% as.POSIXlt() %>% format(format = "%j") %>% as.integer()
coords <- getCoordinates(c4r.obj)
lats <- expand.grid(coords$y, coords$x)[2:1][,2]
c4r.obj %<>% redim(member = TRUE)
nmem <- getShape(c4r.obj, "member")
aux.list <- lapply(1:nmem, function(x) {
ref <- subsetGrid(c4r.obj, members = x) %>% redim(member = FALSE)
ref$Data <- vapply(lats, FUN = geosphere::daylength, j, FUN.VALUE = numeric(length(j))) %>% mat2Dto3Darray(x = coords$x, y = coords$y)
return(ref)
})
c4r.obj <- suppressWarnings(bindGrid(aux.list, dimension = "member"))
aux.list <- NULL
c4r.obj$Variable$varName <- "photoperiod"
attr(c4r.obj$Variable, "description") <- "daylength"
attr(c4r.obj$Variable, "longname") <- "photoperiod_length"
attr(c4r.obj$Variable, "units") <- "hours.day-1"
attr(c4r.obj, "R_package_desc") <- paste0("drought4R-v", packageVersion("drought4R"))
attr(c4r.obj, "R_package_URL") <- "https://github.com/SantanderMetGroup/drought4R"
attr(c4r.obj, "R_package_ref") <- "http://dx.doi.org/10.1016/j.envsoft.2018.09.009"
return(c4r.obj)
}
Loading

0 comments on commit a98b4c7

Please sign in to comment.