require(graphics)
setwd('C://Users//dweise//Documents//SERDP//RC2640//FJ Leaf Temperatures')
rm(list=ls(all=TRUE))
library(compositions);library(MASS);
library(zCompositions);library(openair)
library(ggplot2);library(plyr);library(dplyr)
library(lattice);library(readxl);
library(magrittr);library(chron);library(lubridate)
par(family="sans",font=2,cex.lab=2)
p16d1n1<-read.csv("16Dp1 N1.csv");p16d1n1$area<-"16D";p16d1n1$plot<-"16D1";p16d1n1$leafplot<-"1N";p16d1n1$fire<-207
p16d1n2<-read.csv("16Dp1 N2.csv");p16d1n2$area<-"16D";p16d1n2$plot<-"16D1";p16d1n2$leafplot<-"2N";p16d1n2$fire<-207
p16d1s1<-read.csv("16Dp1 S1.csv");p16d1s1$area<-"16D";p16d1s1$plot<-"16D1";p16d1s1$leafplot<-"1S";p16d1s1$fire<-207
p16d1s2<-read.csv("16Dp1 S2.csv");p16d1s2$area<-"16D";p16d1s2$plot<-"16D1";p16d1s2$leafplot<-"2S";p16d1s2$fire<-207
p16d5n1<-read.csv("16Dp5 1N.csv");p16d5n1$area<-"16D";p16d5n1$plot<-"16D5";p16d5n1$leafplot<-"1N";p16d5n1$fire<-206
p16d5n2<-read.csv("16Dp5 2N.csv");p16d5n2$area<-"16D";p16d5n2$plot<-"16D5";p16d5n2$leafplot<-"2N";p16d5n2$fire<-206
p16d5s1<-read.csv("16Dp5 1S.csv");p16d5s1$area<-"16D";p16d5s1$plot<-"16D5";p16d5s1$leafplot<-"1S";p16d5s1$fire<-206
p16d5s2<-read.csv("16Dp5 2S.csv");p16d5s2$area<-"16D";p16d5s2$plot<-"16D5";p16d5s2$leafplot<-"2S";p16d5s2$fire<-206
p16d6n1<-read.csv("16Dp1 N1.csv");p16d6n1$area<-"16D";p16d6n1$plot<-"16D6";p16d6n1$leafplot<-"1N";p16d6n1$fire<-205
p16d6n2<-read.csv("16Dp1 N2.csv");p16d6n2$area<-"16D";p16d6n2$plot<-"16D6";p16d6n2$leafplot<-"2N";p16d6n2$fire<-205
p16d6s1<-read.csv("16Dp1 S1.csv");p16d6s1$area<-"16D";p16d6s1$plot<-"16D6";p16d6s1$leafplot<-"1S";p16d6s1$fire<-205
p16d6s2<-read.csv("16Dp1 S2.csv");p16d6s2$area<-"16D";p16d6s2$plot<-"16D6";p16d6s2$leafplot<-"2S";p16d6s2$fire<-205
p24Ae1<-read.csv("24A E1.csv");p24Ae1$area<-"24A";p24Ae1$plot<-"24A";p24Ae1$leafplot<-"1E";p24Ae1$fire<-202
p24Ae2<-read.csv("24A E2.csv");p24Ae2$area<-"24A";p24Ae2$plot<-"24A";p24Ae2$leafplot<-"2E";p24Ae2$fire<-202
p24Aw1<-read.csv("24A W1.csv");p24Aw1$area<-"24A";p24Aw1$plot<-"24A";p24Aw1$leafplot<-"1W";p24Aw1$fire<-202
p24Aw2<-read.csv("24A W2.csv");p24Aw2$area<-"24A";p24Aw2$plot<-"24A";p24Aw2$leafplot<-"2W";p24Aw2$fire<-202
p24Be1<-read.csv("24B E1.csv");p24Be1$area<-"24B";p24Be1$plot<-"24B";p24Be1$leafplot<-"1E";p24Be1$fire<-204
p24Be2<-read.csv("24B E2.csv");p24Be2$area<-"24B";p24Be2$plot<-"24B";p24Be2$leafplot<-"2E";p24Be2$fire<-204
p24Bw1<-read.csv("24B W1.csv");p24Bw1$area<-"24B";p24Bw1$plot<-"24B";p24Bw1$leafplot<-"1W";p24Bw1$fire<-204
p24Bw2<-read.csv("24B W2.csv");p24Bw2$area<-"24B";p24Bw2$plot<-"24B";p24Bw2$leafplot<-"2W";p24Bw2$fire<-204
leaftemp<-rbind(p16d1n1,p16d1n2,p16d1s1,p16d1s2,p16d5n1,p16d5n2,p16d5s1,p16d5s2,p16d6n1,p16d6n2,
                p16d6s1,p16d6s2,p24Ae1,p24Ae2,p24Aw1,p24Aw2,p24Be1,p24Be2,p24Bw1,p24Bw2)
leaftemp$frame<-substr(leaftemp$leafplot,2,2)
leaftemp$leafplot<-factor(leaftemp$leafplot,levels=c("1E","1W","2E","2W","1N","1S","2N","2S"))
FJtemp_rate<-leaftemp %>% group_by(fire) %>% summarise(trate=as.double((TempC-lag(TempC,1))/(reltime-lag(reltime,1))),
                                                       maxT=max(TempC,na.rm=TRUE))
FJmaxtr<-FJtemp_rate %>% group_by(fire) %>% summarise(mrate = max(trate,na.rm=TRUE),
                                                      maxT = max(maxT,na.rm=TRUE))
write.csv(FJmaxtr,file="FJ_temprate.csv")

windows(9,9)
ggplot(leaftemp[leaftemp$reltime<=2000,],aes(reltime,TempC,col=leafplot))+
#  geom_point()+
  geom_line(size=1)+
#   geom_smooth(method="loess",se=TRUE)+
    theme_classic()+
  ylim(0,600)+
  xlab("Relative time (s)")+
  ylab(expression(paste("Temperature"," (",degree,"C)")))+
  facet_wrap(~factor(plot,levels=c("24B","24A","16D5","16D6","16D1")))+
  scale_color_brewer(
    name="Leaf Plot",
    aesthetics="col",
    type="qual",
    guide="legend",
    palette="Set1") +
  theme(title=element_text(size=20,face="bold",hjust=0.5),
        axis.text.x=element_text(size=20,vjust=0,hjust=1,angle=90,face="bold",color="black"),
        axis.text.y=element_text(size=20,vjust=0.5,face="bold",color="black"),
        axis.title.y=element_text(size=25,vjust=0.5,face="bold",color="black"), 
        strip.text.x=element_text(size=20,face="bold",color="blue"),
        legend.title=element_text(size=20,face="bold"),
        legend.text=element_text(size=20,face="bold"),
        legend.position=c(0.8,0.2))

savePlot(filename = "FJ_leaf_temp",
         type = "jpg",
         device = dev.cur(),
         restoreConsole = TRUE)

#
#  Wind tunnel IR temps
#
wt_fuel_temps<-read.csv(file="C://Users//dweise//Documents//SERDP//RC2640//Wind tunnel//WindTunnelFLIRTemperatureTime.csv")
wt_temps<-subset(subset(subset(wt_fuel_temps,Fire>7),Fire!=17),Fire!=84)
DF<-wt_temps
DF$ctime<-chron(dates.=DF$Date,times.=DF$Time)
DF2<-subset(DF,select=c(Fire,ctime))
elaps <- by(DF2,DF2$Fire,function(x) apply(x,2,function(x) x-min(x)))
elaps<- do.call(rbind,elaps)
DF$reltime<-elaps[,2]*24*3600  # convert elapsed time from days to seconds
temp_rate<-DF %>% group_by(Fire) %>% summarise(trate=as.double((TempC-lag(TempC,1))/((ctime-lag(ctime,1))*24*3600)),
                                               maxT=max(TempC,na.rm=TRUE))
maxtr<-temp_rate %>% group_by(Fire) %>% summarise(mrate = max(trate,na.rm=TRUE),
                                                  maxT = max(maxT, na.rm=TRUE))
write.csv(maxtr,file="C://Users//dweise//Documents//SERDP//RC2640//Wind tunnel//wt_temprate.csv")
windows(9,9)
FireLab<-function(string){
  string<-paste("WT",string)
  string
}
DF$fireplot<-"WT"
DF2<-data.frame("Fire"=rep(DF$Fire,2),"reltime"=rep(DF$reltime,2),"TempC"=c(DF$TempC,DF$TempC9),"Atten"=c(rep("N",length(DF$reltime)),rep("Y",length(DF$reltime))))
ggplot(DF[DF$Fire<34,],aes(reltime,TempC,col=fireplot))+
#  geom_point(color="red")+
  geom_line(size=2)+
  facet_wrap(~Fire,labeller=labeller(Fire=FireLab))+
  theme_classic()+
  xlab("Relative time (s)")+
  ylab("Temperature (C)")+
  scale_color_brewer(
    name="Wind Tunnel",
    aesthetics="col",
    type="qual",
    guide="legend",
    palette="Set1") +
  theme(title=element_text(size=20,face="bold",hjust=0.5),
        axis.text.x=element_text(size=15,vjust=0.5,hjust=1,angle=90,face="bold",color="black"),
        axis.text.y=element_text(size=20,vjust=0.5,face="bold",color="black"),
        strip.text.x=element_text(size=15,face="bold",color="blue"),
        legend.title=element_text(size=20,face="bold"),
        legend.text=element_text(size=20,face="bold"),
        legend.position="none" )

savePlot(filename = "Time_IRTemp_8_33",
         type = "jpg",
         device = dev.cur(),
         restoreConsole = TRUE)

ggplot(DF[DF$Fire>33 & DF$Fire<59,],aes(x=reltime,y=TempC,col=fireplot))+
  geom_line(size=2)+
  facet_wrap(~Fire,labeller=labeller(Fire=FireLab))+
  theme_classic()+
  xlab("Relative time (s)")+
  ylab("Temperature (C)")+
  scale_color_brewer(
    name="Wind Tunnel",
    aesthetics="col",
    type="qual",
    guide="legend",
    palette="Set1") +
  theme(title=element_text(size=20,face="bold",hjust=0.5),
        axis.text.x=element_text(size=15,vjust=0.5,hjust=1,angle=90,face="bold",color="black"),
        axis.text.y=element_text(size=20,vjust=0.5,face="bold",color="black"),
        strip.text.x=element_text(size=15,face="bold",color="blue"),
        legend.title=element_text(size=20,face="bold"),
        legend.text=element_text(size=20,face="bold"),
        legend.position="none" )

savePlot(filename = "Time_IRTemp_34_58",
         type = "jpg",
         device = dev.cur(),
         restoreConsole = TRUE)

ggplot(DF[DF$Fire>58 & DF$Fire<85,],aes(x=reltime,y=TempC,col=fireplot))+
  geom_line(size=2)+
  facet_wrap(~Fire,labeller=labeller(Fire=FireLab))+
  theme_classic()+
  xlab("Relative time (s)")+
  ylab("Temperature (C)")+
  scale_color_brewer(
    name="Wind Tunnel",
    aesthetics="col",
    type="qual",
    guide="legend",
    palette="Set1") +
  theme(title=element_text(size=20,face="bold",hjust=0.5),
        axis.text.x=element_text(size=15,vjust=0.5,hjust=1,angle=90,face="bold",color="black"),
        axis.text.y=element_text(size=20,vjust=0.5,face="bold",color="black"),
        strip.text.x=element_text(size=15,face="bold",color="blue"),
        legend.title=element_text(size=20,face="bold"),
        legend.text=element_text(size=20,face="bold"),
        legend.position="none" )

savePlot(filename = "Time_IRTemp_59_83",
         type = "jpg",
         device = dev.cur(),
         restoreConsole = TRUE)

ggplot(DF,aes(x=reltime,y=TempC,col=fireplot))+
  geom_line(size=1)+
  facet_wrap(~Fire,labeller=labeller(Fire=FireLab),ncol=10)+
  theme_classic()+
  xlab("Relative time (s)")+
  ylab("Temperature (C)")+
  scale_color_brewer(
    name="Wind Tunnel",
    aesthetics="col",
    type="qual",
    guide="legend",
    palette="Set1") +
  theme(title=element_text(size=20,face="bold",hjust=0.5),
        axis.text.x=element_text(size=10,vjust=0.5,hjust=1,angle=90,face="bold",color="black"),
        axis.text.y=element_text(size=10,vjust=0.5,face="bold",color="black"),
        strip.text.x=element_text(size=10,face="bold",color="blue"),
        legend.title=element_text(size=20,face="bold"),
        legend.text=element_text(size=20,face="bold"),
        legend.position="none" )

savePlot(filename = "Time_IRTemp_8_97",
         type = "jpg",
         device = dev.cur(),
         restoreConsole = TRUE)#
