#
###################################################
# SCRIPT FOR TEMPLATE TO PLOT BACI DATA   
# Created by Gottfried Pestal (Solv Consulting Ltd.)
# Version 1 - August 23, 2012
###################################################

# This script illustrates how to create a pdf-based presentation of data from
# a Before/After - Control/Impact experiment.


# TEMP NOTES
# add labels on first slide, then disappear and add second plot, then add after
# add sample means across each panel


# create an array and fill it with random data for plotting (2 sets of 50 obs, with a before and after for each)
sample.data <- array(NA,dim=c(50,2,2), dimnames=list(paste("Case",1:50,sep=" "),c("Control","Impact"),c("Before","After")))
sample.data[,"Control","Before"]<- rnorm(50,mean=seq(20,25,length.out=50),sd=2)
sample.data[,"Control","After"]<- rnorm(50,mean=seq(22,27,length.out=50),sd=2)
sample.data[,"Impact","Before"]<- rnorm(50,mean=seq(21,27,length.out=50),sd=2)
sample.data[,"Impact","After"]<- rnorm(50,mean=seq(23,37,length.out=50),sd=4)



# open the pdf file as a plotting device
# note: onefile = TRUE creates a new page in the same pdf file whenever the current device is full
pdf("FigRs_Templates_BACI.pdf",width=11, height=8.5,onefile=TRUE)


# specify margins 
par(omi=c(0.5,0.5,0.5,0.5) , mai=c(0.9,0.7,0.2,0.2))

# divide the plotting device 3 into panels
layout(matrix(c(0,0,0,1,0,2,0,0,0), byrow=TRUE, nrow=3,ncol=3), heights=c(1,4,1),widths=c(2,1,2))



# loop through before/after comparison
for(set in c("Before","After")){


control.median<-median(sample.data[,"Control",set])
test.median<-median(sample.data[,"Impact",set])

# create plot in panel 1
control.bars<-barplot(sample.data[,"Control",set], ylim=c(0,40),col="lightblue",border="lightblue",
	main="Control Lake", cex=1.2,cex.lab=1.6,cex.axis=1.6,cex.main=2,names.arg="", ylab="Some Measure")
axis(1,at=c(min(control.bars),max(control.bars)), labels=c("West Shore", "East Shore"), cex.axis=2)

# To add a linear regression line throug the observations, remove the "#" in the next line
#lines(control.bars,lm(sample.data[,"Control",set]~control.bars)$fitted.values, col="darkgrey",lwd=3)

# add a median line and expand it into the middle of the page
abline(h=control.median,col="red",lwd=3)
lines(c(0,100),rep(control.median,2),col="red",lwd=1, xpd=NA)

text(max(control.bars)+1,control.median-1,labels="Median", adj=c(0,1),xpd=NA,cex=1.6)




# create plot in panel 2
test.bars<-barplot(sample.data[,"Impact",set], ylim=c(0,40),col="lightblue",border="lightblue",
	main="Test Lake",cex=1.2,cex.lab=1.6,cex.axis=1.6,cex.main=2,names.arg="", ylab="")
axis(1,at=c(min(test.bars),max(test.bars)), labels=c("West Shore", "East Shore"), cex.axis=2)

# To add a linear regression line throug the observations, remove the "#" in the next line
#lines(test.bars,lm(sample.data[,"Impact",set]~test.bars)$fitted.values, col="darkgrey",lwd=3)

# add a median line and expand it into the middle of the page
abline(h=test.median,col="red",lwd=3)
lines(c(0,-50),rep(test.median,2),col="red",lwd=1, xpd=NA)


# label the page
mtext(set,side=3,line=0.5, outer=TRUE, cex=2.5)
mtext("www.solv.ca/FigRs",side=1,line=0.5, outer=TRUE, cex=1)

} # end loop through before/after comparison


# turn off plotting device and save pdf file
dev.off()


#