From 081711f35a699dcb4953f69cdd3088fdafbdd646 Mon Sep 17 00:00:00 2001 From: florianh Date: Fri, 19 Apr 2024 12:21:37 +0200 Subject: [PATCH] bugfix forestry yields --- .buildlibrary | 2 +- CITATION.cff | 4 +- DESCRIPTION | 4 +- R/ForestYield.R | 175 ++++++++------------------------------------- README.md | 6 +- man/ForestYield.Rd | 6 +- 6 files changed, 41 insertions(+), 156 deletions(-) diff --git a/.buildlibrary b/.buildlibrary index 40119091..51c3c63e 100644 --- a/.buildlibrary +++ b/.buildlibrary @@ -1,4 +1,4 @@ -ValidationKey: '2391914430' +ValidationKey: '2392175504' AutocreateReadme: yes AcceptedWarnings: - 'Warning: package ''.*'' was built under R version' diff --git a/CITATION.cff b/CITATION.cff index 6f33e3db..a8496d10 100644 --- a/CITATION.cff +++ b/CITATION.cff @@ -2,8 +2,8 @@ cff-version: 1.2.0 message: If you use this software, please cite it using the metadata from this file. type: software title: 'magpie4: MAgPIE outputs R package for MAgPIE version 4.x' -version: 1.206.21 -date-released: '2024-04-17' +version: 1.206.22 +date-released: '2024-04-19' abstract: Common output routines for extracting results from the MAgPIE framework (versions 4.x). authors: diff --git a/DESCRIPTION b/DESCRIPTION index 510c0746..dd6fecf5 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,8 +1,8 @@ Type: Package Package: magpie4 Title: MAgPIE outputs R package for MAgPIE version 4.x -Version: 1.206.21 -Date: 2024-04-17 +Version: 1.206.22 +Date: 2024-04-19 Authors@R: c( person("Benjamin Leon", "Bodirsky", , "bodirsky@pik-potsdam.de", role = c("aut", "cre")), person("Florian", "Humpenoeder", , "humpenoeder@pik-potsdam.de", role = "aut"), diff --git a/R/ForestYield.R b/R/ForestYield.R index 7e6541b5..9738dffd 100644 --- a/R/ForestYield.R +++ b/R/ForestYield.R @@ -1,165 +1,50 @@ #' @title ForestYield -#' @description reads timber yield out of a MAgPIE gdx file -#' +#' @description calculates timber yield out of a MAgPIE gdx file +#' #' @export #' #' @param gdx GDX file #' @param file a file name the output should be written to using write.magpie #' @param level Level of regional aggregation; "cell", "reg" (regional), "glo" (global), "regglo" (regional and global) or any secdforest aggregation level defined in superAggregate #' @details Forest yield for timber production -#' @return Forest yield for timber production -#' @author Abhijeet Mishra +#' @return Forest yield for timber production in tDM per ha per year +#' @author Abhijeet Mishra, Florian Humpenoeder #' @importFrom gdx readGDX out #' @importFrom magclass clean_magpie dimSums collapseNames setYears write.magpie #' @importFrom luscale superAggregate #' @examples -#' +#' #' \dontrun{ #' x <- ForestYield(gdx) #' } ForestYield <- function(gdx, file=NULL, level="cell"){ a <- NULL - + if(as.numeric(readGDX(gdx, "s32_hvarea")) > 0 & as.numeric(readGDX(gdx, "s35_hvarea")) > 0) { - if (level == "cell") { - #### Production and harvest area calculations - ov73_prod_forestry <-dimSums(readGDX(gdx,"ov_prod_forestry","ov73_prod_forestry",select = list(type="level")),dim=3) - ov73_hvarea_forestry <- dimSums(readGDX(gdx,"ov32_hvarea_forestry","ov73_hvarea_forestry","ov_hvarea_forestry",react = "silent",select = list(type="level")),dim=3) - - ov73_prod_natveg_secdf <- dimSums(readGDX(gdx,"ov_prod_natveg","ov73_prod_natveg",select = list(type="level"))[,,"secdforest"],dim=3) - ov_hvarea_secdf <- dimSums(readGDX(gdx,"ov35_hvarea_secdforest","ov_hvarea_secdforest",select = list(type="level")),dim=3) - - ov73_prod_natveg_primf <- dimSums(readGDX(gdx,"ov_prod_natveg","ov73_prod_natveg",select = list(type="level"))[,,"primforest"],dim=3) - ov_hvarea_primf <- dimSums(readGDX(gdx,"ov35_hvarea_primforest","ov_hvarea_primforest",select = list(type="level")),dim=3) - - ov73_prod_natveg_other <- dimSums(readGDX(gdx,"ov_prod_natveg","ov73_prod_natveg",select = list(type="level"))[,,"other"],dim=3) - ov_hvarea_other <- dimSums(readGDX(gdx,"ov35_hvarea_other","ov73_hvarea_other","ov_hvarea_other",react = "silent",select = list(type="level")),dim=3) - - #### Yield calculations - - ## Plantations - yield_forestry <- ov73_prod_forestry/ov73_hvarea_forestry - if(any(is.na(range(yield_forestry)))){ - yield_forestry[is.na(yield_forestry)] <- 0 - } - if(any(is.infinite(range(yield_forestry)))){ - div0 <- where(ov73_prod_forestry != 0 & ov73_hvarea_forestry == 0)$true$regions - yield_forestry[is.infinite(yield_forestry)] <- 0 - } - - ## Secondary forest - yield_secdf <- ov73_prod_natveg_secdf/ov_hvarea_secdf - if(any(is.na(range(yield_secdf)))){ - yield_secdf[is.na(yield_secdf)] <- 0 - } - - if(any(is.infinite(range(yield_secdf)))){ - div0 <- where(ov73_prod_natveg_secdf != 0 & ov_hvarea_secdf == 0)$true$regions - yield_secdf[is.infinite(yield_secdf)] <- 0 - } - - ## Primary forest - yield_primf <- ov73_prod_natveg_primf/ov_hvarea_primf - if(any(is.na(range(yield_primf)))){ - yield_primf[is.na(yield_primf)] <- 0 - } - - if(any(is.infinite(range(yield_primf)))){ - div0 <- where(ov73_prod_natveg_primf != 0 & ov_hvarea_primf == 0)$true$regions - yield_primf[is.infinite(yield_primf)] <- 0 - } - - ## Other land - yield_other <- ov73_prod_natveg_other/ov_hvarea_other - if(any(is.na(range(yield_other)))){ - yield_other[is.na(yield_other)] <- 0 - } - - if(any(is.infinite(range(yield_other)))){ - div0 <- where(ov73_prod_natveg_other != 0 & ov_hvarea_other == 0)$true$regions - yield_other[is.infinite(yield_other)] <- 0 - } - - a <- mbind(setNames(yield_forestry,"Forestry"), - setNames(yield_secdf,"Secondary forest"), - setNames(yield_primf,"Primary forest"), - setNames(yield_other,"Other land")) - } else if (level == "regglo" | level == "reg"){ - #### Production and harvest area calculations - ov73_prod_forestry <- dimSums(readGDX(gdx,"ov_prod_forestry","ov73_prod_forestry",select = list(type="level")),dim=3) - ov73_prod_forestry <- superAggregate(data = ov73_prod_forestry,aggr_type = "sum",level = level) - ov73_hvarea_forestry <- dimSums(readGDX(gdx,"ov32_hvarea_forestry","ov73_hvarea_forestry","ov_hvarea_forestry",react = "silent",select = list(type="level")),dim=3) - ov73_hvarea_forestry <- superAggregate(data = ov73_hvarea_forestry ,aggr_type = "sum",level = level) - - ov73_prod_natveg_secdf <- dimSums(readGDX(gdx,"ov_prod_natveg","ov73_prod_natveg",select = list(type="level"))[,,"secdforest"],dim=3) - ov73_prod_natveg_secdf <- superAggregate(data = ov73_prod_natveg_secdf,aggr_type = "sum",level = level) - ov_hvarea_secdf <- dimSums(readGDX(gdx,"ov35_hvarea_secdforest","ov_hvarea_secdforest",select = list(type="level")),dim=3) - ov_hvarea_secdf <- superAggregate(data = ov_hvarea_secdf ,aggr_type = "sum",level = level) - - ov73_prod_natveg_primf <- dimSums(readGDX(gdx,"ov_prod_natveg","ov73_prod_natveg",select = list(type="level"))[,,"primforest"],dim=3) - ov73_prod_natveg_primf <- superAggregate(data = ov73_prod_natveg_primf,aggr_type = "sum",level = level) - ov_hvarea_primf <- dimSums(readGDX(gdx,"ov35_hvarea_primforest","ov_hvarea_primforest",select = list(type="level")),dim=3) - ov_hvarea_primf <- superAggregate(data = ov_hvarea_primf ,aggr_type = "sum",level = level) - - ov73_prod_natveg_other <- dimSums(readGDX(gdx,"ov_prod_natveg","ov73_prod_natveg",select = list(type="level"))[,,"other"],dim=3) - ov73_prod_natveg_other <- superAggregate(data = ov73_prod_natveg_other,aggr_type = "sum",level = level) - ov_hvarea_other <- dimSums(readGDX(gdx,"ov35_hvarea_other","ov73_hvarea_other","ov_hvarea_other",react = "silent",select = list(type="level")),dim=3) - ov_hvarea_other <- superAggregate(data = ov_hvarea_other ,aggr_type = "sum",level = level) - - #### Yield calculations - - ## Plantations - yield_forestry <- ov73_prod_forestry/ov73_hvarea_forestry - if(any(is.na(range(yield_forestry)))){ - yield_forestry[is.na(yield_forestry)] <- 0 - } - - if(any(is.infinite(range(yield_forestry)))){ - div0 <- where(ov73_prod_forestry != 0 & ov73_hvarea_forestry == 0)$true$regions - yield_forestry[is.infinite(yield_forestry)] <- 0 - } - - ## Secondary forest - yield_secdf <- ov73_prod_natveg_secdf/ov_hvarea_secdf - if(any(is.na(range(yield_secdf)))){ - yield_secdf[is.na(yield_secdf)] <- 0 - } - - if(any(is.infinite(range(yield_secdf)))){ - div0 <- where(ov73_prod_natveg_secdf != 0 & ov_hvarea_secdf == 0)$true$regions - yield_secdf[is.infinite(yield_secdf)] <- 0 - } - - ## Primary forest - yield_primf <- ov73_prod_natveg_primf/ov_hvarea_primf - if(any(is.na(range(yield_primf)))){ - yield_primf[is.na(yield_primf)] <- 0 - } - - if(any(is.infinite(range(yield_primf)))){ - div0 <- where(ov73_prod_natveg_primf != 0 & ov_hvarea_primf == 0)$true$regions - yield_primf[is.infinite(yield_primf)] <- 0 - } - - ## Other land - yield_other <- ov73_prod_natveg_other/ov_hvarea_other - if(any(is.na(range(yield_other)))){ - yield_other[is.na(yield_other)] <- 0 - } - - if(any(is.infinite(range(yield_other)))){ - div0 <- where(ov73_prod_natveg_other != 0 & ov_hvarea_other == 0)$true$regions - yield_other[is.infinite(yield_other)] <- 0 - } - - a <- mbind(setNames(yield_forestry,"Forestry"), - setNames(yield_secdf,"Secondary forest"), - setNames(yield_primf,"Primary forest"), - setNames(yield_other,"Other land")) - } else {stop("Resolution not recognized. Select cell or reg or regglo as level. NULL returned.")} - + + #### get annual production and annual harvested area + ov73_prod_forestry <- setNames(dimSums(readGDX(gdx,"ov_prod_forestry","ov73_prod_forestry",select = list(type="level")),dim=3), "Forestry") + ov73_prod_natveg <- dimSums(readGDX(gdx,"ov_prod_natveg","ov73_prod_natveg",select = list(type="level")),dim="kforestry") + ov73_prod_natveg <- setNames(ov73_prod_natveg,c("Primary forest","Secondary forest","Other land")) + ov73_prod <- mbind(ov73_prod_forestry, ov73_prod_natveg) + ov73_prod <- gdxAggregate(gdx, ov73_prod, to = level) + ov73_hvarea <- harvested_area_timber(gdx, level = level) + + #### Yield calculations + yield <- ov73_prod / ov73_hvarea + if(any(is.na(range(yield)))){ + yield[is.na(yield)] <- 0 + } + + if(any(is.infinite(range(yield)))){ + div0 <- where(ov73_prod != 0 & ov73_hvarea == 0)$true$regions + yield[is.infinite(yield)] <- 0 + } + + a <- yield + } else {message("Disabled (no dynamic forestry) ", appendLF = FALSE)} - + out(a,file) -} \ No newline at end of file +} diff --git a/README.md b/README.md index 97d82c96..703ab1b4 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ # MAgPIE outputs R package for MAgPIE version 4.x -R package **magpie4**, version **1.206.21** +R package **magpie4**, version **1.206.22** [![CRAN status](https://www.r-pkg.org/badges/version/magpie4)](https://cran.r-project.org/package=magpie4) [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.1158582.svg)](https://doi.org/10.5281/zenodo.1158582) [![R build status](https://github.com/pik-piam/magpie4/workflows/check/badge.svg)](https://github.com/pik-piam/magpie4/actions) [![codecov](https://codecov.io/gh/pik-piam/magpie4/branch/master/graph/badge.svg)](https://app.codecov.io/gh/pik-piam/magpie4) [![r-universe](https://pik-piam.r-universe.dev/badges/magpie4)](https://pik-piam.r-universe.dev/builds) @@ -39,7 +39,7 @@ In case of questions / problems please contact Benjamin Leon Bodirsky , R package version 1.206.21, . +Bodirsky B, Humpenoeder F, Dietrich J, Stevanovic M, Weindl I, Karstens K, Wang X, Mishra A, Beier F, Breier J, Yalew A, Chen D, Biewald A, Wirth S, von Jeetze P, Leip D, Crawford M, Alves M (2024). _magpie4: MAgPIE outputs R package for MAgPIE version 4.x_. doi:10.5281/zenodo.1158582 , R package version 1.206.22, . A BibTeX entry for LaTeX users is @@ -48,7 +48,7 @@ A BibTeX entry for LaTeX users is title = {magpie4: MAgPIE outputs R package for MAgPIE version 4.x}, author = {Benjamin Leon Bodirsky and Florian Humpenoeder and Jan Philipp Dietrich and Miodrag Stevanovic and Isabelle Weindl and Kristine Karstens and Xiaoxi Wang and Abhijeet Mishra and Felicitas Beier and Jannes Breier and Amsalu Woldie Yalew and David Chen and Anne Biewald and Stephen Wirth and Patrick {von Jeetze} and Debbora Leip and Michael Crawford and Marcos Alves}, year = {2024}, - note = {R package version 1.206.21}, + note = {R package version 1.206.22}, doi = {10.5281/zenodo.1158582}, url = {https://github.com/pik-piam/magpie4}, } diff --git a/man/ForestYield.Rd b/man/ForestYield.Rd index 90f83968..a8564035 100644 --- a/man/ForestYield.Rd +++ b/man/ForestYield.Rd @@ -14,10 +14,10 @@ ForestYield(gdx, file = NULL, level = "cell") \item{level}{Level of regional aggregation; "cell", "reg" (regional), "glo" (global), "regglo" (regional and global) or any secdforest aggregation level defined in superAggregate} } \value{ -Forest yield for timber production +Forest yield for timber production in tDM per ha per year } \description{ -reads timber yield out of a MAgPIE gdx file +calculates timber yield out of a MAgPIE gdx file } \details{ Forest yield for timber production @@ -29,5 +29,5 @@ Forest yield for timber production } } \author{ -Abhijeet Mishra +Abhijeet Mishra, Florian Humpenoeder }