##########Figure 3 biomass

#biomass 2020 Figure 3A
df<- read.csv("2_field.csv")

data<-df%>% filter(year== "2020")

p1 <- ggplot(data=data, aes(x=factor(Date, level=c('2020/8/6', '2020/8/27', '2020/10/3')),y=biomass,
                            fill=factor(treatment,level=c('N0', 'N20', 'N60', 'N120', 'N180')))) +
  geom_bar(position = "dodge", stat="identity",colour = "black",size=0.3)+ 
  geom_errorbar(aes(x = Date, y = biomass, ymin = biomass - biomassError, 
                    ymax = biomass + biomassError),
                width = 0.3,size=0.3,position=position_dodge(.9)) +
  xlab(expression(paste(Date)))+
  labs(y=expression(Alfalfa~shoot~biomass~(kg~ha^{-1})))+
  scale_fill_manual( values = c("#00AFBB", "#E7B800", "#FC4E07","cornflowerblue","violet"),name = NULL)+
  mytheme2+ 
  theme(legend.position = c(0.15,0.80),
        legend.text = element_text(size = 8, family = "serif"),
        legend.title = element_text( size = 7),
        legend.key.size = unit(0.3, 'cm'))+   
  scale_y_continuous(expand = c(0,0),limits = c(0,5000)) +
  geom_text(label = c("a", "a", "a", "a", "a", "a", "a", "a", "a", "a", "b", "b", "ab", "ab", "a"), family="serif",
            aes(y = biomass + biomassError, x = Date), hjust = 0.45,vjust= -0.4, 
            size = 3, position = position_dodge(width = 0.9), stat="identity")+#加显著性字母并调整位置
  theme(
    axis.title.y = element_text(vjust = 3,margin = margin(t =0, r = 0, b = 0, l = 0)),
    axis.title.x = element_text(vjust = -1,margin = margin(t =0, r = 0, b = 0, l = 0)))+
  theme(plot.margin = unit(c(1,1,0.5,1), "cm")) # ("left", "top", "bottom", "right")
p1

#biomass 2021 Figure 3B

data<-df%>% filter(year== "2021")
data$cut1 = factor(data$cut, levels=c('First regrowth','Second regrowth','Third regrowth'))
p2 <-    ggplot(data=data, aes(x=factor(Date, level=c('2021/5/26', '2021/6/3', '2021/6/15', 
                                                      '2021/7/16','2021/7/30','2021/8/7','2021/9/24')),y=biomass,
                               fill=factor(treatment,level=c('N0', 'N20', 'N60', 'N120', 'N180')))) +
  facet_grid(. ~cut1,
             scales="free_x",
             space = "free") +
  geom_bar(position = "dodge", stat="identity",colour = "black",size=0.3)+   
  geom_errorbar(aes(x = Date, y = biomass, ymin = biomass - biomassError, 
                    ymax = biomass + biomassError),
                width = 0.3,size=0.3,position=position_dodge(.9)) + 
  xlab(expression(paste(Date)))+
  labs(y=expression(Alfalfa~shoot~biomass~(kg~ha^{-1})))+
  scale_fill_manual( values = c("#00AFBB", "#E7B800", "#FC4E07","cornflowerblue","violet"),name = NULL)+
  mytheme2+ 
  theme(legend.position = c(0.08,0.83),
        legend.text = element_text(size = 8, family = "serif"),
        legend.title = element_text( size = 7),
        legend.key.size = unit(0.3, 'cm'))+ 
  scale_y_continuous(expand = c(0,0),limits = c(0,12000)) +
  labs(fill=' ')+  
  geom_text(label = c("b", "b", "b","ab", "a", "b", "b", "ab", "a", "a", "b", "ab", "ab", "ab", "a",
                      "b", "ab", "ab","b", "a", "b", "ab", "b", "b", "a", "a", "a", "a", "a", "a",
                      "a", "a", "a", "a", "a"), family="serif",
            aes(y = biomass + biomassError, x = Date), hjust = 0.45,vjust= -0.4, 
            size = 3, position = position_dodge(width = 0.9), stat="identity")+#加显著性字母并调整位置
  theme(
    axis.title.y = element_text(vjust = 3,margin = margin(t =0, r = 0, b = 0, l = 0)),
    axis.title.x = element_text(vjust = -1,margin = margin(t =0, r = 0, b = 0, l = 0)))+
  theme(plot.margin = unit(c(1,1,0.5,1), "cm")) # ("left", "top", "bottom", "right")

p2
df<- read.csv("2_field.csv")
df <- read.csv("two_yield.csv")
#two year yield
p3 <- ggplot(data=df, aes(x=factor(treatment,level=c('N0', 'N20', 'N60', 'N120', 'N180')),
                          y=yield, fill=factor(treatment,level=c('N0', 'N20', 'N60', 'N120', 'N180')))) +
  facet_wrap(~year,nrow = 1)+#分面
  geom_col()+
  geom_bar(position = "dodge", stat="identity",colour = "black",size=0.3)+ 
  geom_errorbar(aes(x = treatment, y = yield, ymin = yield - yield_error, 
                    ymax = yield + yield_error),
                width = 0.3,size=0.3)+ 
  xlab(expression(paste(Treatments)))+
  labs(y=expression(Total~yield~(t~ha^{-1})))+
  mytheme1+ 
  theme(legend.position = c(0.12,0.83),
        legend.text = element_text(size = 8, family = "serif"),
        legend.title = element_text( size = 7),
        legend.key.size = unit(0.3, 'cm'),
        strip.text = element_text(size = 14, family = "serif"))+ 
  scale_y_continuous(expand = c(0,0),limits = c(0, 25)) +
  labs(fill=' ')+  
  geom_text(label = c("b", "b", "ab","ab", "a",
                      "b", "ab", "ab","ab","a"), family="serif",
            aes(y = yield + yield_error, x=treatment), hjust = 0.45,vjust= -0.4, 
            size = 3, position = position_dodge(width = 0.9), stat="identity")+
  scale_fill_manual( values = c("#00AFBB", "#E7B800", "#FC4E07","cornflowerblue","violet"),name = NULL)+
  theme(
    axis.title.y = element_text(vjust = 3,margin = margin(t =0, r = 0, b = 0, l = 0)),
    axis.title.x = element_text(vjust = -1,margin = margin(t =0, r = 0, b = 0, l = 0)))+
  theme(plot.margin = unit(c(1,1,0.5,1), "cm")) # ("left", "top", "bottom", "right")


p3

########## alfalfa yield and biomass 2021 Figure 3C

biomass <- cowplot::plot_grid(p1, p2,labels = "auto",label_fontfamily="serif",
                              nrow = 1,rel_widths = c (1,2))

biomass
biomass <- cowplot::plot_grid(p1, p2,labels = "AUTO",label_fontfamily="serif",
                              nrow = 1,rel_widths = c (1,2))

p4<- cowplot::plot_grid(p3,NULL,labels = c("C",""),label_fontfamily="serif",
                        nrow = 1)
p4
plot_all <- cowplot::plot_grid(biomass, p4,labels = NULL,
                               nrow = 2)

plot_all
ggsave("papaer_field_biomass_c.png", plot=plot_all,
       height=22, width=25, units="cm", dpi=600)