# MakeFigure.R 
# histograms of the pressure
# use from unix shell: Rscript MakeFigure.R 
sSCAN<-read.table("sigmaSCAN.dat",header=TRUE)
sPBE<-read.table("sigmaPBE.dat",header=TRUE)
# compute the pressure p = (sigma_xx+sigma_yy+sigma_zz)/3
pSCAN<-rowMeans(sSCAN[,1:3])
pPBE<-rowMeans(sPBE[,1:3])
summary(pSCAN)
summary(pPBE)

nbreaks=12
min=-2.0
max=2.0

setEPS()

postscript("pSCAN.eps", width=7, height=5)
hist(pSCAN,
main="SCAN330",
xlim=c(min,max),xlab="Stress (GPa)",
ylab="N",
col="#FF0000",
breaks=nbreaks,
cex.lab=1.5, cex.axis=1.5, cex.main=1.5, cex.sub=1.5)

postscript("pPBE.eps", width=7, height=5)
hist(pPBE,main="PBE400",
xlim=c(min,max),xlab="Stress (GPa)",
ylab="N",
col="#1DC9FF",
breaks=nbreaks,
cex.lab=1.5, cex.axis=1.5, cex.main=1.5, cex.sub=1.5)