new BphEffect function with fill - not built

R/calcBphEffect.R

@@ -13,19 +13,15 @@
 calcBphEffect <-function(){
 
   x <- readSource("BphEffect", convert="onlycorrect")
-  k <- calcOutput("ClimateClass", aggregate = FALSE)
+  k <-readSource("Koeppen", subtype="cellular", convert="onlycorrect")[,1995,]
   weight <- dimSums(calcOutput("LanduseInitialisation", aggregate=FALSE, cellular=TRUE, land="fao", input_magpie=TRUE, years="y1995", round=6), dim=3)
 
 
   #mapping to connect cell names with latitudes
   map <- toolGetMapping(type="cell", name="CountryToCellMapping.csv")
-  map <- data.table(map)
 
   #assuming 0 was NA before.
-  x[x==0] <- NA
-
-  #check: How many NA cells
-  length(x[is.na(x)])
+  x[,,"ann_bph"][x[,,"ann_bph"]==0] <- NA
 
   #extract climate class names
   cclass <- getNames(k)
@@ -38,27 +34,23 @@ calcBphEffect <-function(){
   map <- cbind(map,ccl)
 
   #global mean, used in case of NA for all cells for a give lat
-  m_glo <- mean(x,na.rm=TRUE)
+  m_glo <- mean(x[,,"ann_bph"],na.rm=TRUE)
   #m_glo <- mean(b,na.rm=FALSE)
 
   #Loop over climate classes
   for (sel in cclass) {
-    print(sel)
     #get the magpie cells corresponding to cl
-    cells <- map[ccl %in% sel,cell]
+    cells <- map[ccl %in% sel,"celliso"]
     #find cells with NA
-    cells_NA <- cells[is.na(x[cells,,])]
+    cells_NA <- cells[is.na(x[cells,,"ann_bph"])]
     #If all cells are NA, use m_glo, otherwise calc mean based on the non NA cells.
     if(identical(cells,cells_NA)) {
-      x[cells_NA,,] <- m_glo
-      print("Attention: Global mean used!")
+      x[cells_NA,,"ann_bph"] <- m_glo      
     } else {
-      m <- mean(x[cells,,],na.rm=TRUE)
-      x[cells_NA,,] <- m
+      m <- mean(x[cells,,"ann_bph"],na.rm=TRUE)
+      x[cells_NA,,"ann_bph"] <- m
     }
   }
-  #check: no more NA cells
-  length(x[is.na(x)])
 
   return(list(
     x=x,

R/fullCELLULARMAGPIE.R

@@ -111,6 +111,10 @@ fullCELLULARMAGPIE <- function(rev=0.1, dev="", ctype="c200", climatetype="HadGE
   calcOutput("NpiNdcAdAolcPol", aggregate="cluster", round=6, file=paste0("npi_ndc_ad_aolc_pol_", ctype, ".mz"))
   calcOutput("NpiNdcAffPol",    aggregate="cluster", round=6, file=paste0("npi_ndc_aff_pol_", ctype, ".mz"))
 
+  calcOutput("BphEffect", aggregate="cluster", file="f32_bph_effect_noTCRE_", ctype, ".mz")
+  calcOutput("BphTCRE", aggregate="cluster", file="f32_localTCRE_", ctype, ".mz")
+  calcOutput("BphMask", aggregate="cluster", file="f32_bph_mask_", ctype, ".mz")
+
 
 
   #34