version of scripts used for the response to the reviewers

scripts/fig_diff_VPD.R

@@ -91,7 +91,8 @@ tib_diff_VPD %>%
                #limits = date(ymd(c("2020-07-15","2020-09-01"))),
                sec.axis = sec_axis(trans = ~.,
                                    breaks = tick_dates),
-               expand = c(0.002,0.003))+scale_y_discrete("time (hh:mm)", breaks = function(x){x[c(TRUE,rep(F,times=5))]})-> p3
+               expand = c(0.002,0.003))+
+  scale_y_discrete("Time (UT+1)", labels = str_pad(seq(0,21,by =3),side = "left",width = "2",pad = "0"),breaks = function(x){x[c(TRUE,rep(F,times=5))]})-> p3
 #######################################
 plot_grid(p1, p2, p3+theme(legend.position = "none"), get_legend(p3),
   ncol=1, align="hv",axis = "tblr", rel_heights = c(0.4,0.4,0.4,0.1)) %>% 

scripts/figure_SA_energy_storage.R

@@ -43,7 +43,7 @@ library(InraeThemes)
           legend.key.size = unit(c(2,2),"mm"), axis.ticks.length = unit(-1,"mm"), axis.ticks = element_line(size = .2), axis.line = element_line("black",size=0),
           legend.direction = "vertical", legend.spacing = unit(1,"mm"), legend.position = c(0.825,0.88), legend.background = element_blank(),
           panel.background = element_rect(colour = "black"))+
-    labs(fill = "radiative component"),NA,
+    labs(fill = "radiative component"), NA,
   
   Storage_ecosystem %>% 
     pivot_longer(-c(TIMESTAMP, S)) %>% 
@@ -53,10 +53,10 @@ library(InraeThemes)
     ggplot() + 
     geom_col(aes(TIMESTAMP, value, fill = name), size = .2)+
     geom_text(aes(x = day + dhours(to_add), label = label),y = 70, data = dates_Denis %>% mutate(label = LETTERS[6:10], to_add = rep(1,each = 5)),size = 3)+
-    scale_x_datetime(date_breaks = "4 hours", date_labels = "%H", expand = c(0.05,0.01))+
+    scale_x_datetime("Time (UT+1)",date_breaks = "4 hours", date_labels = "%H", expand = c(0.05,0.01))+
     scale_y_continuous(expression(W~m^-2),expand = c(0.02,0.01),breaks = c(-3:3*20) )+
     scale_fill_manual("Storage terms", values = RColorBrewer::brewer.pal(5,"Set2") , guide = guide_legend(direction = "vertical",title.position = "top", title.hjust = 0))+
-    theme(axis.title.x = element_blank(), plot.margin = unit(c(1,1,0,1),"mm"), text = element_text(size=6), axis.text = element_text(size = 6),legend.text.align = .5,
+    theme(axis.title.x = element_text(size = 6), plot.margin = unit(c(1,1,0,1),"mm"), text = element_text(size=6), axis.text = element_text(size = 6),legend.text.align = .5,
           strip.background = element_blank(), strip.text = element_blank(),
           legend.key.size = unit(c(2,2),"mm"), axis.ticks.length = unit(-1,"mm"), axis.ticks = element_line(size = .2), axis.line = element_line("black",size=0),
           legend.spacing = unit(1,"mm"), legend.position = c(0.825,0.165), legend.background = element_blank(),

scripts/figure_Ts_to_gs.R

@@ -123,7 +123,6 @@ rbind(tib_env_var$donnee[[1]] %>%
                  data = fig_segments %>% filter(variable == var)) +
       geom_text(aes(x=xtext,y=ytext, label = LETTER), inherit.aes = F, data = panel, size = textsize/1.5) + 
       geom_line() + geom_point() + 
-      scale_x_datetime(date_breaks = "2 hours", date_labels = "%H", expand = c(0.01,0.01))+
       themeP+
       theme(text = element_text(size = axis_textsize*1.8),
             legend.title = element_blank(),
@@ -150,6 +149,16 @@ rbind(tib_env_var$donnee[[1]] %>%
       }else{
         scale_y_continuous(name = ytitle, labels =function(x) format(x, scientific = TRUE) )}
       } +
+      {if(var !="gc"){
+        scale_x_datetime(date_breaks = "2 hours", date_labels = "%H", expand = c(0.01,0.01))
+      }else{
+        scale_x_datetime("Time (UT+1)",date_breaks = "2 hours", date_labels = "%H", expand = c(0.01,0.01))}
+      } +
+      {if(var !="gc"){
+        theme()
+      }else{
+        theme(axis.title.x = element_text(size = axis_textsize*1.5))}
+      } +
       # colours
       {if(var !="gc"){
         scale_colour_manual(var,values = c("crown" = palette_gs[1], "canopy" = palette_gs[2]), label = labels)

scripts/figure_daily_gs.R

@@ -10,17 +10,21 @@ if (sys.nframe() == 0) {# must be run manually if this is the file on which we a
   source(file.path("scripts","general_variables.R"))
   ################################################################################
   #### daily gs (sapflux methods only) under different env conditions ####
-  {gc_compa %>% 
-      mutate(day=date(TIMESTAMP)) %>% 
-      filter(day %in% ymd(dates_kept$day)) %>% 
-      left_join(dates_kept, by ="day") %>% 
-      filter(day !=ymd("2020-08-19")) %>% 
-      filter(between(hour(TIMESTAMP),6,23)) %>% 
-      mutate(gc = ifelse(gc <0,0,gc)) %>% 
-      filter(TIMESTAMP != ymd_hm("2020-07-07 23:00") | method != "PMT") %>% 
-      filter(method !="PMEC") %>% 
-      group_by(method,day) %>% 
-      mutate(gc_mean = mean(gc, na.rm =T)) %>% 
+gc_prep <- gc_compa %>% 
+  mutate(day=date(TIMESTAMP)) %>% 
+  filter(day %in% ymd(dates_kept$day)) %>% 
+  left_join(dates_kept %>% 
+              mutate(REW_f2=fct_recode(REW_f,"2~\u2265~REW" = "low", "10~\u2265~REW~\u2265~2" = "medium", "REW~\u2265~10" = "high")) %>% 
+              mutate(VPD_f2=fct_recode(VPD_f,"D \u2264 2" = "low", "2 \u2264 D <3" = "medium", "D \u2264 3.5" = "high")), by ="day") %>% 
+  filter(day !=ymd("2020-08-19")) %>% 
+  filter(between(hour(TIMESTAMP),6,23)) %>% 
+  mutate(gc = ifelse(gc <0, 0, gc)) %>% 
+  filter(TIMESTAMP != ymd_hm("2020-07-07 23:00") | method != "PMT") %>% 
+  group_by(method,day) %>% 
+  mutate(gc_mean = mean(gc, na.rm =T))
+  
+{ gc_prep %>% 
+    filter(method !="PMEC") %>% 
       filter(gc < 3*gc_mean) %>% 
       group_by(VPD_f) %>% 
       mutate(gc_max = max(gc, na.rm = T)) %>% 
@@ -31,13 +35,14 @@ if (sys.nframe() == 0) {# must be run manually if this is the file on which we a
              REW_f = fct_relabel(REW_f,~paste(.,"~REW",sep=""))) %>% 
       mutate(condition = sprintf("atop(%s,%s)",VPD_f,REW_f)) %>% 
       mutate(xtext = ymd_hm(ifelse(day == ymd("2020-07-05"),"2020-01-01 19:00","2020-01-01 23:30"))) %>% 
-      ggplot(aes(TIMESTAMP2, gc,linetype = method, colour = method))+
-      geom_line(na.rm = T) + geom_point(size = 0.3, na.rm = T) +
-      geom_text(aes(label = letter, y = gc_max*.95 ), x = ymd_hm("2020-01-01 06:15"),
-                size = textsize/2.5, inherit.aes = F, hjust = 0) +
-      geom_text(aes(label = condition, y = gc_max*.90, x = xtext),
-                size = textsize/3, inherit.aes = F, hjust = 1, parse = T) +
-      facet_grid(VPD_f~REW_f, scales="free_y", switch = "y")+
+      ggplot(aes(TIMESTAMP2, gc,shape = method, colour = method, size = method))+
+      # geom_line(na.rm = T) + 
+      geom_point(na.rm = T) +
+      # geom_text(aes(label = letter, y = gc_max*.95 ), x = ymd_hm("2020-01-01 06:15"),
+      #           size = textsize/2.5, inherit.aes = F, hjust = 0) +
+      # geom_text(aes(label = condition, y = gc_max*.90, x = xtext),
+      #           size = textsize/3, inherit.aes = F, hjust = 1, parse = T) +
+      facet_grid(VPD_f2~REW_f2, scales="free_y", drop = T)+
       theme(axis.line = element_line("black",size=0),
             legend.position=c(0.10,0.17),
             legend.text = element_text(size = textsize*1.5, hjust = 0),
@@ -63,92 +68,94 @@ if (sys.nframe() == 0) {# must be run manually if this is the file on which we a
             #axis.title.y = element_blank(),
             axis.title.y.right = element_blank(),
             strip.placement = "outside",
-            #strip.text = element_text(size = textsize), # back to before
-            strip.text = element_blank(),
+            strip.text = element_text(size = textsize*1.1), # back to before
+            #strip.text = element_blank(),
             strip.background = element_rect(colour =NA,fill=NA)) +
       #labs(colour = parse(text="g[s]~method"), linetype = parse(text="g[s]~method") )+
-      scale_colour_manual(name = expression(g[s]~ method),
-                          values = c("GT" = palette_gs[2],"GT'" = palette_gs[1],"PMT"= palette_gs[3]),
+    scale_colour_manual(name = expression(g[s]~ method),
+                        values = c("GT" = palette_gs[2],"GT'" = palette_gs[1],"PMT"= palette_gs[3]),
+                        labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")])))+
+    scale_linetype_manual(name = expression(g[s]~ method),
+                          limits = c("GT","GT'","PMT"),values = c("solid","dashed","dotdash"),
                           labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")])))+
-      scale_linetype_manual(name = expression(g[s]~ method),
-                            limits = c("GT","GT'","PMT"),values = c("solid","dashed","dotdash"),
-                            labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")])))+
-      scale_x_datetime("Time (UT+1)", date_breaks = "3 hours", date_labels = "%HH", date_minor_breaks = "1 hour")+
-      scale_y_continuous(parse(text="canopy~conductance~(m.s^-1)"),labels = fmt_dcimals(2),
-                         sec.axis = dup_axis()) } %>% 
-    
-    ggsave2("graphes/figure_daily_gs.jpeg",., 
-            width = 190, height = 120, units = "mm", bg = "white", dpi = 1000) 
-  ################################################################################
-  {gc_compa %>% 
-      mutate(day=date(TIMESTAMP)) %>% 
-      filter(day %in% ymd(dates_kept$day)) %>% 
-      left_join(dates_kept, by ="day") %>% 
-      filter(day !=ymd("2020-08-19")) %>% 
-      filter(between(hour(TIMESTAMP),6,23)) %>% 
-      mutate(gc = ifelse(gc <0,0,gc)) %>% 
-      filter(TIMESTAMP != ymd_hm("2020-07-07 23:00") | method != "PMT") %>% 
-      group_by(method,day) %>% 
-      mutate(gc_mean = mean(gc, na.rm =T)) %>% 
-      filter(gc < 5*gc_mean) %>% 
-      group_by(VPD_f) %>% 
-      mutate(gc_max = max(gc, na.rm = T)) %>% 
-      select(-gc_mean) %>%  
-      filter(letter != "G") %>% # remove this one cause it's after August
-      mutate(TIMESTAMP2 = ymd("2020-01-01")+dminutes(minute(TIMESTAMP))+dhours(hour(TIMESTAMP))) %>% 
-      mutate(VPD_f = fct_relabel(VPD_f,~paste(.,"~D[8]",sep="")),
-             REW_f = fct_relabel(REW_f,~paste(.,"~REW",sep=""))) %>% 
-      mutate(condition = sprintf("atop(%s,%s)",VPD_f,REW_f)) %>% 
-      mutate(stresslevel = ifelse(letter =="A","least stress:",ifelse(letter == "I","maximal stress:",NA))) %>% 
-      mutate(xtext = ymd_hm(ifelse(day %in% ymd("2020-07-16","2020-07-30"),"2020-01-01 19:30","2020-01-01 23:30"))) %>% 
-      ggplot(aes(TIMESTAMP2, gc,linetype = method, colour = method))+
-      geom_line(na.rm = T) + geom_point(size = 0.3, na.rm = T) +
-      geom_text(aes(label = letter, y = gc_max*.95 ), x = ymd_hm("2020-01-01 06:15"),
-                size = textsize/2.5, inherit.aes = F, hjust = 0) +
-      geom_text(aes(label = condition, y = gc_max*.90, x = xtext),
-                size = textsize/3, inherit.aes = F, hjust = 1, parse = T) +
-      geom_text(aes(label = stresslevel, y =gc_max*.90, x=xtext - dhours(5)), size = textsize/3,inherit.aes = F, hjust = 1)+
-      facet_grid(VPD_f~REW_f, scales="free_y", switch = "y")+
-      theme(axis.line = element_line("black",size=0),
-            legend.position=c(0.10,0.17),
-            legend.text = element_text(size = textsize*1.5, hjust = 0),
-            legend.title = element_text(size = textsize*1.5, hjust = .5),
-            legend.key.height  = unit(3,"mm"),
-            legend.key = element_rect(fill = NA),
-            legend.background = element_rect(fill = NA),
-            
-            panel.border = element_rect(colour="black", fill=NA, size=1),
-            panel.background = element_rect(fill = NA),
-            
-            axis.ticks = element_line(size = .3),
-            axis.ticks.length = unit(-.15, units = "cm"),
-            
-            axis.title = element_text(size = textsize),
-            axis.title.x.bottom = element_blank(),
-            axis.text.x = element_text(size =textsize , vjust = -1),
-            axis.text.x.top = element_blank(),
-            axis.ticks.x.top = element_blank(), 
-            
-            axis.text.y = element_text(size = textsize, hjust = - 0.5),
-            axis.text.y.right = element_blank(),
-            #axis.title.y = element_blank(),
-            axis.title.y.right = element_blank(),
-            strip.placement = "outside",
-            #strip.text = element_text(size = textsize), # back to before
-            strip.text = element_blank(),
-            strip.background = element_rect(colour =NA,fill=NA)) +
-      scale_colour_manual(name = expression(g[s]~ method),
-                          values = c("GT" = palette_gs[2],"GT'" = palette_gs[1],"PMT"= palette_gs[3], "PMEC" = "dark grey"),
+    scale_shape_manual(name = expression(g[s]~ method),
+                       limits = c("GT","GT'","PMT"),values = c(1,16,2),
+                       labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")])))+
+    scale_size_manual(name = expression(g[s]~ method),
+                      limits = c("GT","GT'","PMT"),values = c(0.3,0.5,0.3),
+                      labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")])))+
+    scale_x_datetime("Time (UT+1)", date_breaks = "3 hours", date_labels = "%H", date_minor_breaks = "1 hour")+
+    scale_y_continuous(parse(text="canopy~conductance~(m.s^-1)"),labels = fmt_dcimals(2),
+                       sec.axis = dup_axis()) } %>% 
+  
+  ggsave2("graphes/figure_daily_gs.jpeg",., 
+          width = 190, height = 120, units = "mm", bg = "white", dpi = 1000) 
+################################################################################
+{gc_prep %>% 
+    filter(gc < 3*gc_mean) %>% 
+    group_by(VPD_f) %>% 
+    mutate(gc_max = max(gc, na.rm = T)) %>% 
+    select(-gc_mean) %>%  
+    filter(letter != "G") %>% # remove this one cause it's after August
+    mutate(TIMESTAMP2 = ymd("2020-01-01")+dminutes(minute(TIMESTAMP))+dhours(hour(TIMESTAMP))) %>% 
+    mutate(VPD_f = fct_relabel(VPD_f,~paste(.,"~D[8]",sep="")),
+           REW_f = fct_relabel(REW_f,~paste(.,"~REW",sep=""))) %>% 
+    mutate(condition = sprintf("atop(%s,%s)",VPD_f,REW_f)) %>% 
+    mutate(xtext = ymd_hm(ifelse(day == ymd("2020-07-05"),"2020-01-01 19:00","2020-01-01 23:30"))) %>% 
+    ggplot(aes(TIMESTAMP2, gc,shape = method, colour = method, size = method))+
+    # geom_line(na.rm = T) + 
+    geom_point(na.rm = T) +
+    # geom_text(aes(label = letter, y = gc_max*.95 ), x = ymd_hm("2020-01-01 06:15"),
+    #           size = textsize/2.5, inherit.aes = F, hjust = 0) +
+    # geom_text(aes(label = condition, y = gc_max*.90, x = xtext),
+    #           size = textsize/3, inherit.aes = F, hjust = 1, parse = T) +
+    facet_grid(VPD_f2~REW_f2, scales="free_y", drop = T)+
+    theme(axis.line = element_line("black",size=0),
+          legend.position=c(0.10,0.17),
+          legend.text = element_text(size = textsize*1.5, hjust = 0),
+          legend.title = element_text(size = textsize*1.5, hjust = .5),
+          legend.key.height  = unit(3,"mm"),
+          legend.key = element_rect(fill = NA),
+          legend.background = element_rect(fill = NA),
+          
+          panel.border = element_rect(colour="black", fill=NA, size=1),
+          panel.background = element_rect(fill = NA),
+          
+          axis.ticks = element_line(size = .3),
+          axis.ticks.length = unit(-.15, units = "cm"),
+          
+          axis.title = element_text(size = textsize*1.1),
+          axis.title.x.bottom = element_text(size = textsize*1.1, vjust = 0),
+          axis.text.x = element_text(size =textsize , vjust = -1),
+          axis.text.x.top = element_blank(),
+          axis.ticks.x.top = element_blank(), 
+          
+          axis.text.y = element_text(size = textsize, hjust = - 0.5),
+          axis.text.y.right = element_blank(),
+          #axis.title.y = element_blank(),
+          axis.title.y.right = element_blank(),
+          strip.placement = "outside",
+          strip.text = element_text(size = textsize*1.1), # back to before
+          #strip.text = element_blank(),
+          strip.background = element_rect(colour =NA,fill=NA)) +
+    scale_colour_manual(name = expression(g[s]~ method),
+                        values = c("GT" = palette_gs[2],"GT'" = palette_gs[1],"PMT"= palette_gs[3], "PMEC" = "dark grey"),
+                        labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")]),expression(g[paste("s",",PMEC")])))+
+    scale_linetype_manual(name = expression(g[s]~ method),
+                          limits = c("GT","GT'","PMT","PMEC"),values = c("solid","dashed","dotdash","twodash"),
                           labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")]),expression(g[paste("s",",PMEC")])))+
-      scale_linetype_manual(name = expression(g[s]~ method),
-                            limits = c("GT","GT'","PMT","PMEC"),values = c("solid","dashed","dotdash","twodash"),
-                            labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")]),expression(g[paste("s",",PMEC")])))+
-      scale_x_datetime( date_breaks = "3 hours", date_labels = "%HH", date_minor_breaks = "1 hour")+
-      scale_y_continuous(parse(text="canopy~conductance~(m.s^-1)"),labels = fmt_dcimals(2),
-                         sec.axis = dup_axis())} %>% 
-    
-    ggsave2("graphes/SA_daily_gs_FULL.jpeg",., 
-            width = 190, height = 120, units = "mm", bg = "white", dpi = 1000) 
+    scale_shape_manual(name = expression(g[s]~ method),
+                       limits = c("GT","GT'","PMT","PMEC"),values = c(1,16,2,2),
+                       labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")]),expression(g[paste("s",",PMEC")])))+
+    scale_size_manual(name = expression(g[s]~ method),
+                      limits = c("GT","GT'","PMT","PMEC"),values = c(0.3,0.5,0.3,0.3),
+                      labels = c(expression(g[paste("s",",GT")]),expression(g[paste("s",",GT'")]),expression(g[paste("s",",PMT")]),expression(g[paste("s",",PMEC")])))+
+    scale_x_datetime("Time (UT+1)", date_breaks = "3 hours", date_labels = "%H", date_minor_breaks = "1 hour")+
+    scale_y_continuous(parse(text="canopy~conductance~(m.s^-1)"),labels = fmt_dcimals(2),
+                       sec.axis = dup_axis())} %>% 
+  
+  ggsave2("graphes/SA_daily_gs_FULL.jpeg",., 
+          width = 190, height = 120, units = "mm", bg = "white", dpi = 1000) 
   ################################################################################
 }