Peatland considered during disaggregation

CHANGELOG.md

@@ -6,18 +6,18 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).
 ## [Unreleased]
 
 ### changed
-- **38_factor_costs** updated use of USDA cost shares 
+- **38_factor_costs** updated use of USDA cost shares
 - **inputdata** changed GDP base year from 2005USD to 2017USD
 - **config** changed default input data to use 2017USD
 - **module_documentation** all references to USD05 changed to USD17
 - **scripts** REMIND coupling reads data in US$2017, not US$2005
 - **config** updated input data to rev4.114
 - **config** SHAPE scenarios start year of dietary shift changed to 2025
+- **extra/disaggregation** Peatland now considered in disaggregation of land pools
 - **core** number of age-classes doubled from 150 to 300 years for better match of growth curves with potential natural vegetation.
 - **35_natveg** revised age-class initialization of secondary forest
 - **modules** update of scaling factors in several modules
 
-
 ### added
 - **62_material** added switch to turn off future material demand for bioplastic
 - **config** added SSP1-POP-GDP SSP2-POP-GDP and SSP5-POP-GDP
@@ -26,9 +26,9 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).
 -
 
 ### fixed
-- **11_costs** changed equation to fix bug in total water cost calculation 
-- **29_cropland** treecover age-class growth was not working properly because ac_sub was erroneously not fixed 
-- **scripts** script/output/extra/resubmit.R 
+- **11_costs** changed equation to fix bug in total water cost calculation
+- **29_cropland** treecover age-class growth was not working properly because ac_sub was erroneously not fixed
+- **scripts** script/output/extra/resubmit.R
 - **28_ageclass** bugfix input data unit and code clean-up. Renamed `feb21` realization to `oct24`
 - **70_livestock** bugfix scaling.gms file in wrong folder
 - **10_land** bugfix land transition matrix for improved feasibility (variables and parameters have different accuracy)

main.gms

@@ -147,44 +147,44 @@ $title magpie
 *'  * Always try to access model outputs through the corresponding magpie package instead of accessing them directly with readGDX. It cannot be guaranteed that your script will work in the future if you do otherwise (as only the corresponding magpie package will be continuously adapted to changes in the GAMS code).
 
 *##################### R SECTION START (VERSION INFO) ##########################
-* 
+*
 * Used data set: rev4.114_h12_magpie.tgz
 * md5sum: NA
 * Repository: https://rse.pik-potsdam.de/data/magpie/public
-* 
+*
 * Used data set: rev4.114_h12_fd712c0b_cellularmagpie_c200_MRI-ESM2-0-ssp370_lpjml-8e6c5eb1.tgz
 * md5sum: NA
 * Repository: https://rse.pik-potsdam.de/data/magpie/public
-* 
+*
 * Used data set: rev4.114_h12_validation.tgz
 * md5sum: NA
 * Repository: https://rse.pik-potsdam.de/data/magpie/public
-* 
+*
 * Used data set: additional_data_rev4.56.tgz
 * md5sum: NA
 * Repository: https://rse.pik-potsdam.de/data/magpie/public
-* 
+*
 * Used data set: calibration_H12_27Sep24.tgz
 * md5sum: NA
 * Repository: https://rse.pik-potsdam.de/data/magpie/public
-* 
+*
 * Low resolution: c200
 * High resolution: 0.5
-* 
+*
 * Total number of cells: 200
-* 
+*
 * Number of cells per region:
 *   CAZ  CHA  EUR  IND  JPN  LAM  MEA  NEU  OAS  REF  SSA  USA
 *    14   23   10    7    4   26   21    9   16   23   32   15
-* 
+*
 * Regionscode: 62eff8f7
-* 
+*
 * Regions data revision: 4.114
-* 
+*
 * lpj2magpie settings:
 * * LPJmL data: MRI-ESM2-0:ssp370
 * * Revision: 4.114
-* 
+*
 * aggregation settings:
 * * Input resolution: 0.5
 * * Output resolution: c200
@@ -193,10 +193,10 @@ $title magpie
 *   CAZ  CHA  EUR  IND  JPN  LAM  MEA  NEU  OAS  REF  SSA  USA
 *    14   23   10    7    4   26   21    9   16   23   32   15
 * * Call: withCallingHandlers(expr, message = messageHandler, warning = warningHandler,     error = errorHandler)
-* 
-* 
+*
+*
 * Last modification (input data): Sun Oct 27 00:37:36 2024
-* 
+*
 *###################### R SECTION END (VERSION INFO) ###########################
 
 $offupper

scripts/output/extra/disaggregation.R

@@ -87,7 +87,8 @@ if (length(map_file) > 1) {
 }
 
 .dissagLandConsv <- function(gdx, cfg, map_file, wdpa_hr_file, consv_prio_hr_file) {
-  land_consv_lr <- readGDX(gdx, "p22_conservation_area", react = "silent")
+  land_consv_lr <- readGDX(gdx, "pm_land_conservation", react = "silent")
+  land_consv_lr <- dimSums(land_consv_lr, dim=3.2)
   wdpa_hr <- read.magpie(wdpa_hr_file)
   map <- readRDS(map_file)
 
@@ -275,6 +276,8 @@ if (cfg$gms$urban == "exo_nov21") {
 # Prepare land conservation data
 # ----------------------------------------
 
+message("Disaggregating conservation land")
+
 land_consv_hr <- NULL
 if (file.exists(wdpa_hr_file)) {
   land_consv_hr <- .dissagLandConsv(gdx, cfg, map_file, wdpa_hr_file, consv_prio_hr_file)
@@ -300,6 +303,30 @@ avl_cropland_hr <- file.path(outputdir, "avl_cropland_0.5.mz") # available cropl
 marginal_land <- cfg$gms$c29_marginal_land # marginal land scenario
 snv_pol_fader <- readGDX(gdx, "i29_snv_scenario_fader")
 
+
+# --------------------------------
+# Disaggregate peatland
+# --------------------------------
+
+message("Disaggregating peatland")
+
+# check for peatland version
+if (cfg$gms$peatland == "v2") {
+  peat_lr <- PeatlandArea(gdx, level = "cell", sum = FALSE)
+  peat_ini_hr <- read.magpie(peatland_v2_hr_file)
+  peat_ini_hr <- add_columns(peat_ini_hr, addnm = "rewetted", dim = "d3", fill = 0)
+  peat_ini_hr <- add_columns(peat_ini_hr, addnm = "unused", dim = "d3", fill = 0)
+  peat_hr <- suppressWarnings(luscale::interpolate2(peat_lr, peat_ini_hr, map_file))
+  peat_hr <- peat_hr[, getYears(peat_hr, as.integer = T) >= cfg$gms$s58_fix_peatland, ]
+} else if (cfg$gms$peatland == "on") {
+  peat_lr <- PeatlandArea(gdx, level = "cell", sum = TRUE)
+  peat_ini_hr <- mbind(setNames(read.magpie(peatland_on_intact_hr_file), "intact"), setNames(read.magpie(peatland_on_degrad_hr_file), "degrad"))
+  peat_ini_hr <- add_columns(peat_ini_hr, addnm = "rewet", dim = "d3", fill = 0)
+  peat_hr <- suppressWarnings(luscale::interpolate2(peat_lr, peat_ini_hr, map_file))
+  peat_hr <- peat_hr[, getYears(peat_hr, as.integer = T) >= cfg$gms$s58_fix_peatland, ]
+}
+peat_hr <- .fixCoords(peat_hr)
+
 # ============================================
 # Start disaggregation
 # ============================================
@@ -319,6 +346,7 @@ land_hr <- interpolateAvlCroplandWeighted(
   marginal_land = marginal_land,
   urban_land_hr = urban_land_hr,
   land_consv_hr = land_consv_hr,
+  peat_hr = peat_hr,
   snv_pol_shr = snv_pol_shr,
   snv_pol_fader = snv_pol_fader
 )
@@ -335,6 +363,28 @@ land_hr <- .fixCoords(land_hr)
 )
 gc()
 
+# -----------------------------------
+# Write peatland outputs
+# -----------------------------------
+
+# Write output
+.writeDisagg(peat_hr, peatland_hr_out_file,
+  comment = "unit: Mha per grid-cell",
+  message = "Write outputs peatland Mha"
+)
+gc()
+
+out <- peat_hr / dimSums(land_hr[, getYears(peat_hr), ], dim = 3)
+out[is.nan(out)] <- 0
+out[is.infinite(out)] <- 0
+
+.writeDisagg(out, peatland_hr_share_out_file,
+  comment = "unit: grid-cell land area fraction",
+  message = "Write outputs peatland share"
+)
+gc()
+
+
 # ---------------------------------
 #  Split land pools
 # ---------------------------------
@@ -474,6 +524,10 @@ if (grepl("grass", cfg$gms$past)) {
     "past", "manpast",
     gsub("range", "rangeland", getNames(land_lr))
   )
+  getNames(land_consv_hr) <- gsub(
+    "past", "manpast",
+    gsub("range", "rangeland", getNames(land_consv_hr))
+  )
 } else {
   # Disaggregate pasture
   land_ini_lr <- mbind(
@@ -485,6 +539,11 @@ if (grepl("grass", cfg$gms$past)) {
     land_lr[, , c("past"), invert = TRUE],
     collapseNames(land_lr[, , "past"]) * side_layers_lr[, , c("manpast", "rangeland")]
   )
+
+  land_consv_hr <- mbind(
+    land_consv_hr[, , c("past"), invert = TRUE],
+    collapseNames(land_consv_hr[, , "past"]) * side_layers_hr[, , c("manpast", "rangeland")]
+  )
 }
 
 # Sort and rename
@@ -504,6 +563,7 @@ land_bii_hr <- interpolateAvlCroplandWeighted(
   marginal_land = marginal_land,
   urban_land_hr = urban_land_hr,
   land_consv_hr = land_consv_hr,
+  peat_hr = peat_hr,
   snv_pol_shr = snv_pol_shr,
   snv_pol_fader = snv_pol_fader,
   unit = "share"
@@ -532,46 +592,4 @@ rm(bii_hr)
 gc()
 
 
-# --------------------------------
-# Disaggregate peatland
-# --------------------------------
-
-message("Disaggregating peatland")
-
-# check for peatland version
-if (cfg$gms$peatland == "v2") {
-  peat_lr <- PeatlandArea(gdx, level = "cell", sum = FALSE)
-  peat_ini_hr <- read.magpie(peatland_v2_hr_file)
-  peat_ini_hr <- add_columns(peat_ini_hr, addnm = "rewetted", dim = "d3", fill = 0)
-  peat_ini_hr <- add_columns(peat_ini_hr, addnm = "unused", dim = "d3", fill = 0)
-  peat_hr <- suppressWarnings(luscale::interpolate2(peat_lr, peat_ini_hr, map_file))
-  peat_hr <- peat_hr[, getYears(peat_hr, as.integer = T) >= cfg$gms$s58_fix_peatland, ]
-} else if (cfg$gms$peatland == "on") {
-  peat_lr <- PeatlandArea(gdx, level = "cell", sum = TRUE)
-  peat_ini_hr <- mbind(setNames(read.magpie(peatland_on_intact_hr_file), "intact"), setNames(read.magpie(peatland_on_degrad_hr_file), "degrad"))
-  peat_ini_hr <- add_columns(peat_ini_hr, addnm = "rewet", dim = "d3", fill = 0)
-  peat_hr <- suppressWarnings(luscale::interpolate2(peat_lr, peat_ini_hr, map_file))
-  peat_hr <- peat_hr[, getYears(peat_hr, as.integer = T) >= cfg$gms$s58_fix_peatland, ]
-}
-peat_hr <- .fixCoords(peat_hr)
-
-# Write output
-.writeDisagg(peat_hr, peatland_hr_out_file,
-  comment = "unit: Mha per grid-cell",
-  message = "Write outputs peatland Mha"
-)
-gc()
-
-out <- peat_hr / dimSums(land_hr[, getYears(peat_hr), ], dim = 3)
-out[is.nan(out)] <- 0
-out[is.infinite(out)] <- 0
-
-rm(land_hr, peat_hr)
-
-.writeDisagg(out, peatland_hr_share_out_file,
-  comment = "unit: grid-cell land area fraction",
-  message = "Write outputs peatland share"
-)
-gc()
-
 message("Finished disaggregation")