From 0529ae35e17adaca90eadd8e0904d5e9c5b23520 Mon Sep 17 00:00:00 2001 From: Nicolau Manubens Date: Thu, 11 Apr 2019 16:36:41 +0200 Subject: [PATCH] Adapting vignettes to Rmd. --- DESCRIPTION | 5 +- ...gnette.md => MultiModelSkill_vignette.Rmd} | 10 +-- ..._vignette.md => MultivarRMSE_vignette.Rmd} | 11 +-- ...FARM_vignette.md => RainFARM_vignette.Rmd} | 72 +++++++++---------- 4 files changed, 52 insertions(+), 46 deletions(-) rename vignettes/{MultiModelSkill_vignette.md => MultiModelSkill_vignette.Rmd} (99%) rename vignettes/{MultivarRMSE_vignette.md => MultivarRMSE_vignette.Rmd} (98%) rename vignettes/{RainFARM_vignette.md => RainFARM_vignette.Rmd} (93%) diff --git a/DESCRIPTION b/DESCRIPTION index f4ad9b9f..dd958d8d 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -41,8 +41,9 @@ Imports: Suggests: zeallot, testthat, - R.rsp -VignetteBuilder: R.rsp + knitr, + rmarkdown +VignetteBuilder: knitr License: Apache License 2.0 | file LICENSE Encoding: UTF-8 LazyData: true diff --git a/vignettes/MultiModelSkill_vignette.md b/vignettes/MultiModelSkill_vignette.Rmd similarity index 99% rename from vignettes/MultiModelSkill_vignette.md rename to vignettes/MultiModelSkill_vignette.Rmd index 6fc67610..414955a2 100644 --- a/vignettes/MultiModelSkill_vignette.md +++ b/vignettes/MultiModelSkill_vignette.Rmd @@ -1,9 +1,9 @@ --- author: "Nuria Perez" date: "`r Sys.Date()`" -output: html_document +output: rmarkdown::html_vignette vignette: > - %\VignetteEngine{R.rsp::md} + %\VignetteEngine{knitr::knitr} %\VignetteIndexEntry{Multi-model Skill Assessment} %\usepackage[utf8]{inputenc} --- @@ -21,7 +21,6 @@ A version of s2dverification under development should be loaded by running: devtools::install_git('https://earth.bsc.es/gitlab/es/s2dverification', branch = 'mostlikely') library(s2dverification) -library(zeallot) ``` @@ -70,8 +69,11 @@ obs <- "erainterim" Using the `CST_Load` function, the data available in our data store can be loaded. The following lines, shows how this function can be used. However, the data is loaded from a previous saved `.RData` file: Ask nuria.perez at bsc.es to achieve the data to run the recipe. -```r +```r +require(zeallot) + glosea5 <- '/esnas/exp/glosea5/specs-seasonal_i1p1/$STORE_FREQ$_mean/$VAR_NAME$-allmemb/$VAR_NAME$_$START_DATE$.nc' + c(exp, obs) %<-% CST_Load(var = clim_var, exp = list(list(name = 'glosea5', path = glosea5), list(name = 'ecmwf/system4_m1'), diff --git a/vignettes/MultivarRMSE_vignette.md b/vignettes/MultivarRMSE_vignette.Rmd similarity index 98% rename from vignettes/MultivarRMSE_vignette.md rename to vignettes/MultivarRMSE_vignette.Rmd index 78af6080..3712075d 100644 --- a/vignettes/MultivarRMSE_vignette.md +++ b/vignettes/MultivarRMSE_vignette.Rmd @@ -1,9 +1,9 @@ --- author: "Deborah Verfaillie" date: "`r Sys.Date()`" -output: html_document +output: rmarkdown::html_vignette vignette: > - %\VignetteEngine{R.rsp::md} + %\VignetteEngine{knitr::knitr} %\VignetteIndexEntry{Multivariate RMSE} %\usepackage[utf8]{inputenc} --- @@ -18,7 +18,6 @@ To run this vignette, the next R packages should be installed and loaded: ```r library(s2dverification) library(RColorBrewer) -library(zeallot) ``` @@ -68,14 +67,18 @@ obs <- "erainterim" Using the `CST_Load` function from **CSTool package**, the data available in our data store can be loaded. The following lines show how this function can be used. Here, the data is loaded from a previous saved `.RData` file: Ask nuria.perez at bsc.es for the data to run the recipe. -```r +```r +require(zeallot) + glosea5 <- list(path = '/esnas/exp/glosea5/specs-seasonal_i1p1/$STORE_FREQ$_mean/$VAR_NAME$-allmemb/$VAR_NAME$_$START_DATE$.nc') + c(exp_T, obs_T) %<-% CST_Load(var = temp, exp = list(glosea5), obs = obs, sdates = dateseq, leadtimemin = 2, leadtimemax = 4, latmin = 25, latmax = 75, lonmin = -20, lonmax = 70, output = 'lonlat', nprocs = 1, storefreq = "monthly", sampleperiod = 1, nmember = 9, method = "bilinear", grid = paste("r", grid, sep = "")) + c(exp_P, obs_P) %<-% CST_Load(var = precip, exp = list(glosea5), obs = obs, sdates = dateseq, leadtimemin = 2, leadtimemax = 4, diff --git a/vignettes/RainFARM_vignette.md b/vignettes/RainFARM_vignette.Rmd similarity index 93% rename from vignettes/RainFARM_vignette.md rename to vignettes/RainFARM_vignette.Rmd index 7b9817d1..90a882c8 100644 --- a/vignettes/RainFARM_vignette.md +++ b/vignettes/RainFARM_vignette.Rmd @@ -2,18 +2,16 @@ title: "Rainfall Filtered Autoregressive Model (RainFARM) precipitation downscaling" author: "Jost von Hardenberg (ISAC-CNR)" date: "26/03/2019" -output: - html_document: default - pdf_document: default +output: rmarkdown::html_vignette vignette: > + %\VignetteEngine{knitr::knitr} %\VignetteIndexEntry{RainFARM} %\usepackage[utf8]{inputenc} - %\VignetteEngine{R.rsp::md} --- @@ -32,32 +30,33 @@ An efficient implementation of RainFARM is provided for CSTools by the `CST_Rain In order to run the examples in this vignette, the *CSTools* package and some other support R packages need to be loaded by running: -```r +```{r} install.packages('CSTools') library(CSTools) ``` We use test data provided by CSTools to load a seasonal precipitation forecast: -```r -exp %<-% lonlat_prec +```{r} +library(CSTools) +exp <- lonlat_prec ``` This gives us a CSTools object `exp`, containing an element `exp$data` with dimensions: -```r +```{r} dim(exp$data) -dataset member sdate ftime lat lon - 1 6 3 31 4 4 +#dataset member sdate ftime lat lon +# 1 6 3 31 4 4 ``` There are 6 ensemble members available in the data set, 3 starting dates and 31 forecast times, which refer to daily values in the month of March following starting dates on November 1st in the years 2010, 2011, 2012. Please notice that RainFARM (in this version) only accepts square domains, possibly with an even number of pixels on each side, so we always need to select an appropriate cutout. Also, there are time and memory limitations when a large ensemble of downscaled realizations is generated with RainFARM, so that selecting a smaller target area is advised. On the other hand, if spectral slopes are to be determined from the large scales we will still need enough resolution to allow this estimation. In this example we have preselected a 4x4 pixel cutout at resolution 1 degree in a smaller area lon=[6,9], lat=[44,47] covering Northern Italy. -#```r -#ilon <- which ( exp$lon %in% 5:12 ) -#ilat <- which ( exp$lat %in% 40:47 ) -#exp$data <- exp$data[ , , , , ilon, ilat, drop=FALSE] -#exp$lon <- exp$lon[ilon] -#exp$lat <- exp$lat[ilat] -#``` +```{r} +ilon <- which ( exp$lon %in% 5:12 ) +ilat <- which ( exp$lat %in% 40:47 ) +exp$data <- exp$data[ , , , , ilon, ilat, drop=FALSE] +exp$lon <- exp$lon[ilon] +exp$lat <- exp$lat[ilat] +``` ### Standard downscaling without climatological weights @@ -65,24 +64,25 @@ Our goal is to downscale with RainFARM these data from the resolution of 1 degre RainFARM can compute automatically its only free parameter, i.e. the spatial spectral slope, from the large-scale field (here only with size 4x4 pixel, but in general we reccomend selecting at least 8x8 pixels). In this example we would like to compute this slope as an average over the _member_ and _ftime_ dimensions, while we will use different slopes for the remaining _dataset_ and _sdate_ dimensions (a different choice may be more appropriate in a real application). To obtain this we specify the parameter `time_dim = c("member", "ftime")`. The slope is computed starting from the wavenumber corresponding to the box, `kmin=1`. We create 3 stochastic realizations for each dataset, member, starting date and forecast time with `nens=5`. The command to donwscale and the resulting fields are: -```r +```{r} exp_down <- CST_RainFARM(exp, nf=20, kmin = 1, nens = 3, time_dim = c("member", "ftime")) dim(exp_down$data) - dataset member realization sdate ftime lat lon - 1 6 3 3 31 80 80 +# dataset member realization sdate ftime lat lon +# 1 6 3 3 31 80 80 str(exp_down$lon) - num [1:80] 5.53 5.58 5.62 5.67 5.72 ... +# num [1:80] 5.53 5.58 5.62 5.67 5.72 ... str(exp_down$lat) - num [1:80] 47.5 47.4 47.4 47.3 47.3 ... +# num [1:80] 47.5 47.4 47.4 47.3 47.3 ... ``` The function returns an array `exp_down$data` with the additional "realization" dimension for the stochastic ensemble with 3 members. The longitudes and latitudes have been correspondingly interpolated to the finer resolution. Alternatively we could have used the "reduced" function `RainFARM` which accepts directly a data array (with arbitrary dimensions, provided a longitude, a latitude and a "time" dimension exist) and two arrays to describe longitudes and latitudes: -```r -downscaled <- RainFARM(exp$data, exp$lon, exp$lat, nf = 20, kmin = 1, nens = 3, +```{r} +downscaled <- RainFARM(exp$data, exp$lon, exp$lat, + nf = 20, kmin = 1, nens = 3, time_dim = c("member", "ftime")) ``` The resulting array has the same dimensions as the `$data` element in the result of `CST_RainFARM`. @@ -93,7 +93,7 @@ Each instant and each realization will of course be different, but let's plot an png("Figures/RainFARM_fig1.png", width = 640, height = 365) par(mfrow=c(1,2)) --> -```r +```{r} a <- exp$data[1, 1, 1, 17, , ] * 86400 * 1000 a[a > 60] <- 60 image(exp$lon, rev(exp$lat), t(apply(a, 2, rev)), xlab = "lon", ylab = "lat", @@ -119,13 +119,13 @@ The area of interest in our example presents a complex orography, but the basic Suitable climatology files could be for example a fine-scale precipitation climatology from a high-resolution regional climate model (see e.g. Terzago et al. 2018), a local high-resolution gridded climatology from observations, or a reconstruction such as those which can be downloaded from the WORLDCLIM (http://www.worldclim.org) or CHELSA (http://chelsa-climate.org) websites. The latter data will need to be converted to NetCDF format before being used (see for example the GDAL tools (https://www.gdal.org). We will assume that a copy of the WORLDCLIM precipitation climatology at 30 arcseconds (about 1km resolution) is available in the local file `medscope.nc`. From this file we can derive suitable weights to be used with RainFARM using the `CST_RFWeights` functions as follows: -```r +```{r} ww <- CST_RFWeights("./worldclim.nc", nf = 20, lon = exp$lon, lat = exp$lat) ``` The result is a two-dimensional weights matrix with the same `lon`and `lat` dimensions as requested. The weights (varying around an average value of 1) encode how to distribute differently precipitation in each stochastic realization of RainFARM. We call again `CST_RainFARM()`, this time passing climatological weights: -```r +```{r} exp_down_weights <- CST_RainFARM(exp, nf = 20, kmin = 1, nens = 3, weights = ww, time_dim = c("member", "ftime")) ``` @@ -134,8 +134,8 @@ We plot again the same realization as before in Fig. 2 (left panel). The WORLDCL From a single realization and time it is not possible to see that a more realistic precipitation has been achieved, but this becomes clear if we compare the climatological average over all 3 stochastic ensemble realizations, over all forecast times and over all forecast times. The resulting plot in panels Fig2b,c show that the RainFARM climatology downscaled without weights presents on average a very coarse structure, comparable with that of the original fields, while when using the weights a much more realistic distribution is achieved. -