#####   R-code for making Figure 5.9. It reads data from
#####   "example55.dat" 

dat = matrix(scan("example55.dat"),ncol=2,byrow=T)
x1 = dat[,1]
x2 = dat[,2]
N = length(x1)

la1=0.1
rho1=2.595

la2=0.1
rho2=1.580

En1 = rep(0,N)
En2 = rep(0,N)
En1[1] = la1*x1[1]^2+(1-la1)*1
En2[1] = la2*x2[1]^2+(1-la2)*1

for(i in 2:N){
   En1[i] = la1*x1[i]^2+(1-la1)*En1[i-1]
   En2[i] = la2*x2[i]^2+(1-la2)*En2[i-1]
}

U = 1+rho1*sqrt(2*la1/(2-la1))
L = 1-rho2*sqrt(2*la2/(2-la2))

ii <- seq(1,N)

postscript("fig59.ps",width=6.5,height=3.5,horizontal=F)

par(mfrow=c(1,2), mar=c(4,4,2,2))

plot(ii,En1,type="o",lty=1,pch=16,xlab="n",
     ylab=expression(E[n]),mgp=c(2,1,0),xlim=c(0,N), 
     ylim=c(0,5),cex=0.8)
lines(ii,rep(U,N),lty=2,cex=0.8)
title(xlab="(a)",cex=0.9)

plot(ii,En2,type="o",lty=1,pch=16,xlab="n",
     ylab=expression(E[n]),mgp=c(2,1,0),xlim=c(0,N), 
     ylim=c(0,2),cex=0.8)
lines(ii,rep(L,N),lty=2,cex=0.8)
title(xlab="(b)",cex=0.9)

graphics.off()