Commit 5fe91178 authored by Jana Ulrich's avatar Jana Ulrich
Browse files

Merge branch 'covariates' of gitlab.met.fu-berlin.de:Rpackages/IDF into covariates

parents 8ecdd320 14e34e19
......@@ -18,7 +18,8 @@ Imports: stats4,
evd,
ismev,
RcppRoll,
pbapply
pbapply,
fastmatch
License: GPL (>=2)
Encoding: UTF-8
LazyData: true
......
......@@ -15,6 +15,7 @@ importFrom(evd,dgev)
importFrom(evd,pgev)
importFrom(evd,qgev)
importFrom(evd,rgev)
importFrom(fastmatch,ctapply)
importFrom(grDevices,rainbow)
importFrom(grDevices,rgb)
importFrom(graphics,abline)
......
# This file contains the functions:
# -IDF.agg for the preparing the data
# -IDF.agg for preparing the data
# -IDF.plot for plotting of IDF curves at a chosen station
#### IDF.agg ####
#' Aggregation and annual maxima for choosen durations
#' Aggregation and annual maxima for chosen durations
#' @description Aggregates several time series for chosen durations and finds annual maxima
#' (either for the whole year or chosen months). Returns data.frame that can be used for
#' the function \code{\link{gev.d.fit}}.
......@@ -16,8 +16,8 @@
#' standard date format.
#' @param ds numeric vector of aggregation durations.
#' (Must be multiples of time resolution at all stations.)
#' @param na.accept numeric in [0,1] giving maximum percentage of missing values
#' for which block max. should still be calculated
#' @param na.accept numeric in [0,1) giving maximum percentage of missing values
#' for which block max. should still be calculated.
#' @param which.stations optional, subset of stations. Either numeric vector or character vector
#' containing names of elements in data. If not given, all elements in `data` will be used.
#' @param which.mon optional, subset of months of which to calculate the annual maxima from.
......@@ -36,6 +36,7 @@
#' @export
#' @importFrom pbapply pbsapply
#' @importFrom RcppRoll roll_sum
#' @importFrom fastmatch ctapply
#'
#' @examples
#' dates <- as.Date("2019-01-01")+0:729
......@@ -69,12 +70,12 @@
# function 1: aggregate over single durations and find annual maxima:
agg.ts <- function(ds){
runsum = RcppRoll::roll_sum(data.s[,names[2]],ds/dtime,fill=NA)
runsum = RcppRoll::roll_sum(data.s[,names[2]],ds/dtime,fill=NA,align='right')
#runmean <- rollapplyr(as.zoo(data.s[,names[2]]),ds/dtime,FUN=sum,fill =NA,align='right')
runsum <- runsum/ds #intensity per hour
max.subset <- lapply(1:length(which.mon),function(m.i){
subset <- is.element(as.POSIXlt(data.s[,names[1]])$mon,which.mon[[m.i]])
max <- tapply(runsum[subset],(as.POSIXlt(data.s[,names[1]])$year+1900)[subset],
max <- ctapply(runsum[subset],(as.POSIXlt(data.s[,names[1]])$year+1900)[subset],
function(vec){
n.na <- sum(is.na(vec))
max <- ifelse(n.na <= na.accept*length(vec),max(vec,na.rm = TRUE),NA)
......@@ -86,7 +87,7 @@
return(df) # maxima for single durations
}
# call function 1 in lapply to aggregate over all durations at single station
data.agg <- pbapply::pblapply(ds,agg.ts,cl=cl) #
data.agg <- pbapply::pblapply(ds,agg.ts,cl=cl)
df <- do.call(rbind,data.agg)
return(df) # maxima for all durations at one station
}
......@@ -144,12 +145,15 @@ IDF.plot <- function(durations,fitparams,probs=c(0.5,0.9,0.99),
if(is.element('ylim',names(list(...)))){
ylim <- list(...)[['ylim']]
}else{ylim <- range(idf.array,na.rm=TRUE)}
if(is.element('ylim',names(list(...)))){
if(is.element('xlim',names(list(...)))){
xlim <- list(...)[['xlim']]
}else{xlim <- range(durations,na.rm=TRUE)}
if(is.element('main',names(list(...)))){
main <- list(...)[['main']]
}else{main <- ''}
# empty plot
plot(NA,xlim=xlim,ylim=ylim,xlab="Duration [h]",ylab="Intensity [mm/h]",log="xy")
plot(NA,xlim=xlim,ylim=ylim,xlab="Duration [h]",ylab="Intensity [mm/h]",log="xy",main=main)
}
## plot IDF curves
......
......@@ -95,6 +95,9 @@ gev.d.fit<-
# test for NA values:
if(any(is.na(xdat))) stop('xdat contains NA values. NA values need to be removed first.')
# test for finite values:
if(any(is.infinite(xdat))) stop('xdat contains non finite values. Inf and -Inf need to be removed first.')
# test if covariates matrix is given correctly
npar <- max(sapply(z$model,function(x){return(ifelse(is.null(x),0,max(x)))}))
if(any(npar>ncol(ydat),npar>0 & is.null(ydat)))stop("Not enough columns in covariates matrix 'ydat'.")
......@@ -404,7 +407,7 @@ gev.d.diag <- function(fit,subset=NULL,cols=NULL,pch=NULL,which='both',mfrow=c(1
abline(0, 1, col = 1,lwd=1)
title(title[1])
if(legend){legend('bottomright',legend = round(durs,digits = 2),pch=pch,
col = cols[1:length(durs)],title = 'Durations',ncol = 2)}
col = cols[1:length(durs)],title = 'Durations[h]',ncol = 2)}
}
if(which=='both'|which=='qq'){
# qq
......@@ -413,7 +416,7 @@ gev.d.diag <- function(fit,subset=NULL,cols=NULL,pch=NULL,which='both',mfrow=c(1
abline(0, 1, col = 1,lwd=1)
title(title[2])
if(legend){legend('bottomright',legend = round(durs,digits = 2),pch=pch,
col = cols[1:length(durs)],title = 'Durations',ncol = 2)}
col = cols[1:length(durs)],title = 'Durations [h]',ncol = 2)}
}
if(which=='both') par(mfrow=c(1,1)) # reset par
}
......
......@@ -2,7 +2,7 @@
% Please edit documentation in R/IDF.R
\name{IDF.agg}
\alias{IDF.agg}
\title{Aggregation and annual maxima for choosen durations}
\title{Aggregation and annual maxima for chosen durations}
\usage{
IDF.agg(
data,
......@@ -23,8 +23,8 @@ standard date format.}
\item{ds}{numeric vector of aggregation durations.
(Must be multiples of time resolution at all stations.)}
\item{na.accept}{numeric in [0,1] giving maximum percentage of missing values
for which block max. should still be calculated}
\item{na.accept}{numeric in [0,1) giving maximum percentage of missing values
for which block max. should still be calculated.}
\item{which.stations}{optional, subset of stations. Either numeric vector or character vector
containing names of elements in data. If not given, all elements in `data` will be used.}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment