R/forecast-extractPS.R

Defines functions getPS

Documented in getPS

##' Extraction des pluies spatiales 
##'
##' \command{getPS} extrait les précipitations spatiales prévues 
##' par un modèle climatique sur un ensemble de bassins versants.
##' 
##' @param forecast un jeu de prévision de précipiations sous forme
##' d'objet de classe \sQuote{\samp{Forecast*}}.
##' @param bassin les limites des bassins versant sous forme d'objet de
##' classe \sQuote{\samp{SpatialPolygon*}}.
##' @param logique indiquant si les précipitations spatailes doivent 
##' être regroupées dans par scénario météorologique ou par bassin 
##' versant
##' 
##' @return Retourne un liste contenant la table des précipitations 
##' journalières sur chaque bassin versant ainsi que la table des
##' précipitations au pas de temps du modèle  
##' 
##' @author Thomas Esclaffer \email{[email protected]@edf.fr}
##' @importFrom raster extract
##' @export
getPS <- 
  function(forecast, bassin, label, by.sce=TRUE) 
  {
    if (!inherits(forecast,c("Forecast","ForecastEnsemble"))) 
      stop("Un objet de classe 'Forecast*' est requis pour 'forecast' !")
    ## Préparation des arguments 
    if (missing(bassin)) {
<<<<<<< HEAD
      data(bvPrevi, package = "DtgRecup")
=======
      data(bvPrevi)
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
      bassin <- bv.previ.ll
      label <- bv.previ.ll$mordor_name
    }
    if (missing(label)){
      lcl <- sapply(bassin@data, class)
      id <- match("character", lcl, nomatch = 1)
      label <- bassin@data[,id]
    }
<<<<<<< HEAD
    is.tot <- grepl("tot", forecast$variable, ignore.case = TRUE)
=======
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
    ## Extraction des pluies spatiales pour 
    ## Un jeu de prévision déterministe : Forecast
    if (inherits(forecast,"Forecast")) {
      return(.getPS(forecast$data, dates=forecast$dates, 
                    ech=forecast$echeances, drun=forecast$run, 
<<<<<<< HEAD
                    bassin=bassin, label=label, is.tot=is.tot))
    }
    
=======
                    bassin=bassin, label=label))
    }
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
    ## Un jeu de prévision d'ensemble : Forecast.ensemble
    if (inherits(forecast,"Forecast.ensemble")) {
      tmp <- lapply(forecast$data, .getPS, dates=forecast$dates, 
                    ech=forecast$echeances, drun=forecast$run,
<<<<<<< HEAD
                    bassin=bassin, label=label, is.tot=is.tot)
      return(.arrangePS(tmp, by.sce = by.sce))
=======
                    bassin=bassin, label=label)
      return(.arangePS(tmp, by.sce=by.sce))
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
    }
  }

## Fonction interne de calcul recomposition des pluies spatiales
## au pas de temps horaire et journlaier
.getPS <- 
<<<<<<< HEAD
  function(data, dates, ech, drun, bassin, label, is.tot=FALSE) 
=======
  function(data, dates, ech, drun, bassin, label) 
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
  {
    ##
    ## Calcul des pluies spatiales
    ## Extraction
    tmp <- raster::extract(data, bassin, mean, small=TRUE, 
                           df=TRUE, factors=FALSE)[,-1]
<<<<<<< HEAD
    tmp <- round(t(tmp), 1)
=======
    tmp <- round(t(tmp),1)
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
    tmp <- data.frame(date=dates, tmp)
    colnames(tmp) <- c("date",label)
    rownames(tmp) <- NULL
    ## On rajoute la ligne de cumul nul si la première écheance ne 
    ## correspond à la date du run 
<<<<<<< HEAD
    if (tmp$date[1] > drun) {
=======
    if (tmp$date[1]>drun) {
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
      cini <- data.frame(drun, as.list(rep(0, ncol(tmp)-1)))
      names(cini) <- names(tmp)
      tmp <- rbind(cini, tmp)
    }
<<<<<<< HEAD
    ##
    ## Pour le calcul des valeurs horaires et journalières 
    ## on s'appuie sur les précipitations cumulés que l'on 
    ## interpole aux dates d'intérêt
    ##
    ## Calcul des précipitations cumulées si besoin
    if (!is.tot) {
      onames <- names(tmp)
      ptot <- data.frame(tmp$date, lapply(tmp[-1], cumsum))
      names(ptot) <- onames
    } else ptot <- tmp
    ## Fonction d'interpolation
    dPS <-  function(ps, dates, dout) {
      pc <- approx(x=dates, y=ps, xout=dout, rule = 1)$y
=======
    ## 
    ## Pour le calcul des valeurs horaires et journalières 
    ## on s'appuie sur les précipitations cumulés que l'on 
    ## interpole aux dates d'intérêt
    ## 
    dPS <-  function(ps, dates, dout) {
      pc <- approx(x=dates, y=cumsum(ps), xout=dout)$y
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
      return(c(0, round(diff(pc),2)))
    }
    ##
    ## Désagregation horaire
    ## 
    h.dates <- seq(from=drun, to=max(dates), by="hour")
<<<<<<< HEAD
    tmp2  <- lapply(ptot[-1], dPS, dates=ptot$date, dout=h.dates)
    ps.hor <- data.frame(h.dates, tmp2)
    names(ps.hor) <- names(ptot)
=======
    tmp2  <- lapply(tmp[-1], dPS, dates=tmp$date, dout=h.dates)
    ps.hor <- data.frame(h.dates, tmp2)
    names(ps.hor) <- names(tmp)
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
    ##
    ## Agregation journalière
    ##
    d.deb <- as.POSIXlt(drun, tz="UTC") ; d.deb$hour <- 6
    d.dates <- seq(from=d.deb, to=max(dates), by="day")
<<<<<<< HEAD
    tmp2  <- lapply(ptot[-1], dPS, dates=ptot$date, dout=d.dates)
    ps.day <- data.frame(as.Date(d.dates)-1, tmp2)
    names(ps.day) <- names(ptot)
    if (d.deb > drun)
      ps.day <- subset(ps.day, date >= as.Date(drun))
    else 
      ps.day <- subset(ps.day, date > as.Date(drun))
=======
    tmp2  <- lapply(tmp[-1], dPS, dates=tmp$date, dout=d.dates)
    ps.day <- data.frame(as.Date(d.dates - 1), tmp2)
    names(ps.day) <- names(tmp)
    ps.day <- subset(ps.day, date > as.Date(drun))
>>>>>>> 53ae716ce8775c9fd36538db2ce00c6b918e4778
    ## Résultats
    invisible(list(raw=tmp, day=ps.day, hour=ps.hor))
  } 


## Ré-arrrange la structure des données des pluies spatiales afin 
## d'obtenir pour les prévisions d'ensemble : 
## soit une liste contenant pour chaque pas de temps les tableaux 
## chronologiques des pluies spatiales de chaque scenario.
## soit une liste contenant pour chaque bassins les tableaux
## chronologiques des scenario de pluies pour chaque pas de temps.
.arrangePS <- 
  function(x, by.sce=TRUE)
  {
    if (!by.sce) {
      return(list(
        raw=lapply(x, with, "raw"), 
        day=lapply(x, with, "day"), 
        hour=lapply(x, with, "hour"))
      )
    } else {
      extract.bv <- function(bv, x) {
        list(raw=data.frame(date=x[[1]]$raw$date, 
                            lapply(x, with, raw[[bv]])),
             day=data.frame(date=x[[1]]$day$date, 
                            lapply(x, with, day[[bv]])),
             hour=data.frame(date=x[[1]]$hour$date, 
                             lapply(x, with, hour[[bv]]))
             )
      }
      lbvs <- colnames(x[[1]][-1])
      return(lapply(lbvs, extract.bv, x=x))
    }
  }
coolTot/DtgRecup documentation built on May 12, 2017, 9:45 a.m.