#
#
# Written by:
# -- 
# John L. Weatherwax                2009-04-21
# 
# email: wax@alum.mit.edu
# 
# Please send comments and especially bug reports to the
# above email address.
# 
#-----

# fix the random.seed, so I'll always get the same answer:
set.seed(10131985)

nsim = 10^3
X1=runif(nsim)
X2=runif(nsim)
X3=runif(nsim)

f1 = 0.5*(X1+X2)
f2 = (X1+X2+X3)/3

cBreaks = seq(0,1,length.out=20)

# call hist to get the maximum ylim value (for common plotting axes)
hf1 = hist(f1,breaks=cBreaks,plot=FALSE)
hf2 = hist(f2,breaks=cBreaks,plot=FALSE)

cMax = max( max(hf1$density), max(hf2$density) )

# Plot everything:
# 
###postscript("../../WriteUp/Graphics/Chapter2/sect_4_prob_1_hists.eps", onefile=FALSE, horizontal=FALSE)

par(mfrow=c(1,2))
hf1 = hist(f1,freq=F,main="(x1+x2)/2",breaks=cBreaks,ylim=c(0,cMax))
hf2 = hist(f2,freq=F,main="(x1+x2+x3)/3",breaks=cBreaks,ylim=c(0,cMax))

###dev.off()