Merge pull request #1736 from gabriel-abrahao/sccfix

Multiple fixes and improvements on MAGICC for damages

modules/51_internalizeDamages/KWlikeItr/postsolve.gms

@@ -17,10 +17,10 @@ p51_marginalDamageCumul(tall,tall2,regi2)$((tall2.val ge tall.val) and (tall.val
 ;
 
 p51_sccLastItr(tall) = p51_scc(tall);
-
+* Add an epsilon to pm_consPC to avoid division by zero in case of INFES in the reference data
 p51_sccParts(tall,tall2,regi2)$((tall.val ge 2010) and (tall.val le 2150) and (tall2.val ge tall.val) and (tall2.val le 2250)) = 
 	 (1 + pm_prtp(regi2) )**(-(tall2.val - tall.val))
-	* pm_consPC(tall,regi2)/pm_consPC(tall2,regi2) 
+	* pm_consPC(tall,regi2)/(pm_consPC(tall2,regi2) + 0.0000001) 
         * pm_damage(tall2,regi2) * pm_GDPGross(tall2,regi2) 
 	* p51_marginalDamageCumul(tall,tall2,regi2) 
 	* pm_sccIneq(tall2,regi2)

scripts/input/climate_assessment_openscm_run.py

@@ -17,9 +17,21 @@
 import json
 
 # TONN REMOVE start
-# If these are already set, it shouldn't override it. We actually may not want to have a default but just throw an error if not set
-os.environ["MAGICC_EXECUTABLE_7"]       =   "/p/projects/piam/abrahao/scratch/module_climate_tests/climate-assessment-files/magicc-v7.5.3/bin/magicc"
-os.environ["MAGICC_WORKER_ROOT_DIR"]    =   "/p/projects/piam/abrahao/scratch/methane/methane_scm/workers"
+class EnvironmentError(Exception):  
+    pass  
+
+for env_var in ["MAGICC_EXECUTABLE_7", "MAGICC_WORKER_ROOT_DIR"]:  
+    if os.environ.get(env_var, '') == '':  
+        # If clause covers both cases in which the env var is not set at all 
+        # as well as the case in which it is set to an empty string
+        # If these are already set, it shouldn't override it. We actually may not want to have a default but just throw an error if not set
+        os.environ["MAGICC_EXECUTABLE_7"]       =   "/p/projects/rd3mod/climate-assessment-files/magicc-v7.5.3/bin/magicc"
+        os.environ["MAGICC_WORKER_ROOT_DIR"]    =   os.environ["PTMP"] + "/"  
+        raise EnvironmentError(f"{env_var} does not exist")  
+
+
+    # Optional debug prints  
+    print(f"Found '{env_var}' = '{os.environ.get(env_var)}' ") 
 
 LOGGER = logging.getLogger(__name__) # We don't need this
 # TONN REMOVE end
@@ -171,6 +183,7 @@ def fix_hfc_unit(variable):
     },
     output_variables=(
         "Surface Air Temperature Change",
+        "Effective Radiative Forcing|Anthropogenic",
         "Net Atmosphere to Land Flux|CO2"
     ),
 scenarios = scmdata.ScmRun(basescen)
@@ -189,10 +202,3 @@ def fix_hfc_unit(variable):
                   ).to_excel(
                       outfilename
                       )
-
-# # # %%
-# # # Show basic results
-# # runresults.filter(
-# #     variable = "Surface Air Temperature Change",
-# #     region = "World"
-# #     ).lineplot()

scripts/input/climate_assessment_run.R

@@ -204,6 +204,12 @@ if (any(!alreadySet)) do.call(Sys.setenv, as.list(environmentVariables[!alreadyS
 #
 # BUILD climate-assessment RUN COMMANDS
 #
+
+expectedHarmInfFile <- file.path(
+  climateTempDir,
+  paste0(baseFn, "_harmonized_infilled.xlsx")
+)
+
 runHarmoniseAndInfillCmd <- paste(
   "python", file.path(scriptsDir, "run_harm_inf.py"),
   climateAssessmentEmi,
@@ -213,13 +219,32 @@ runHarmoniseAndInfillCmd <- paste(
 )
 
 runClimateEmulatorCmd <- paste(
-  "python", file.path(scriptsDir, "run_clim.py"),
+  "python climate_assessment_openscm_run.py ",
   normalizePath(file.path(climateTempDir, paste0(baseFn, "_harmonized_infilled.csv"))),
-  climateTempDir,
+  "--climatetempdir", climateTempDir,
+  # Note: Option --year-filter-last requires https://github.com/gabriel-abrahao/climate-assessment/tree/yearfilter
+  "--endyear", 2250,
   "--num-cfgs", nparsets,
   "--scenario-batch-size", 1,
   "--probabilistic-file", probabilisticFile
 )
+# runClimateEmulatorCmd <- paste(
+#   "python", file.path(scriptsDir, "run_clim.py"),
+#   normalizePath(file.path(climateTempDir, paste0(baseFn, "_harmonized_infilled.csv"))),
+#   climateTempDir,
+#   "--num-cfgs", nparsets,
+#   "--scenario-batch-size", 1,
+#   "--probabilistic-file", probabilisticFile
+# )
+
+# Get conda environment folder
+condaDir <- "/p/projects/rd3mod/python/environments/scm_magicc7"
+# Command to activate the conda environment, changes depending on the cluster
+if (file.exists("/p/system/modulefiles/defaults/piam/1.25")) {
+  condaCmd <- paste0("module load conda/2023.09; source activate ", condaDir, ";")
+} else {
+  condaCmd <- paste0("module load anaconda/2023.09; source activate ", condaDir, ";")
+}
 
 logMsg <- paste0(
   date(), "  CLIMATE-ASSESSMENT ENVIRONMENT:\n",
@@ -230,6 +255,7 @@ logMsg <- paste0(
   "  scriptsDir            = '", scriptsDir, "' exists? ", dir.exists(scriptsDir), "\n",
   "  magiccBinFile         = '", magiccBinFile, "' exists? ", file.exists(magiccBinFile), "\n",
   "  magiccWorkersDir      = '", magiccWorkersDir, "' exists? ", dir.exists(magiccWorkersDir), "\n\n",
+  "  condaCmd              = '", condaCmd, "'\n",
   "  ENVIRONMENT VARIABLES:\n",
   "  MAGICC_EXECUTABLE_7    = ", Sys.getenv("MAGICC_EXECUTABLE_7"), "\n",
   "  MAGICC_WORKER_ROOT_DIR = ", Sys.getenv("MAGICC_WORKER_ROOT_DIR"), "\n",
@@ -244,7 +270,7 @@ timeStopSetUpAssessment <- Sys.time()
 ############################# HARMONIZATION/INFILLING #############################
 
 timeStartHarmInf <- Sys.time()
-system(paste(runHarmoniseAndInfillCmd, "&>>", logFile))
+system(paste(condaCmd, runHarmoniseAndInfillCmd, "&>>", logFile))
 timeStopHarmInf <- Sys.time()
 
 ############################# RUNNING MODEL #############################
@@ -258,7 +284,7 @@ logMsg <- paste0(
 capture.output(cat(logMsg), file = logFile, append = TRUE)
 
 timeStartEmulation <- Sys.time()
-system(paste(runClimateEmulatorCmd, "&>>", logFile))
+system(paste(condaCmd, runClimateEmulatorCmd, "&>>", logFile))
 timeStopEmulation <- Sys.time()
 
 ############################# POSTPROCESS CLIMATE OUTPUT #############################
@@ -283,12 +309,13 @@ timeStartPostProcessing <- Sys.time()
 # mutliple files, add another entry to the list. TODO: This could config file...
 associateVariablesAndFiles <- as.data.frame(rbind(
   c(
-    magicc7Variable = "AR6 climate diagnostics|Surface Temperature (GSAT)|MAGICCv7.5.3|50.0th Percentile",
+    # magicc7Variable = "AR6 climate diagnostics|Surface Temperature (GSAT)|MAGICCv7.5.3|50.0th Percentile",
+    magicc7Variable = "Surface Air Temperature Change",
     gamsVariable = "pm_globalMeanTemperature",
     fileName = "p15_magicc_temp"
   ),
   c(
-    magicc7Variable = "AR6 climate diagnostics|Effective Radiative Forcing|Basket|Anthropogenic|MAGICCv7.5.3|50.0th Percentile",
+    magicc7Variable = "Effective Radiative Forcing|Anthropogenic",
     gamsVariable = "p15_forc_magicc",
     fileName = "p15_forc_magicc"
   )