R/Heatwave.R

In alasdairhunter/EuropeanRiskIndex: Package for computing climate extremes from the CMIP data.

#Atomic heatwave function
library(evd)
Heatwave <- function(proj, ref, threshold = NULL, quantile = 0.99, min_length = 3) {
  #If the temperature threshold is not defined by the user, calculate it from the
  #reference period
  if (is.null(threshold)) {
    u <- quantile(ref, quantile)
  } else {
    u <- threshold
  }
  # Calculate the length of each cluster of heatwave events
  ref.heatwave <- c(unlist(lapply(clusters(c(ref), u, keep.names = FALSE),
                                  function(x) length(x))))

  proj.heatwave <- c(unlist(lapply(clusters(c(proj), u, keep.names = FALSE),
                              function(x) length(x))))

  # The total number of heatwave events whose individual duration
  #exceeds the minimum number of consecutive days
  ref.hw.freq <- length(ref.heatwave[ref.heatwave > min_length])
  proj.hw.freq <- length(proj.heatwave[proj.heatwave > min_length])

  #The total duration of all heatwave events
  ref.hw.duration <- sum(ref.heatwave[ref.heatwave > min_length])
  proj.hw.duration <- sum(proj.heatwave[proj.heatwave > min_length])

  invisible(result <- list(ref.freq = ref.hw.freq, proj.freq = proj.hw.freq,
                           ref.duration = ref.hw.duration, proj.duration = proj.hw.duration, threshold = u))
}