statS<-function(fname,from=0,to=2, len=101){
x<-seq(from,to,length=len);
n<-length(x);
meanS<-1:n;
stderrS<-1:n;
qup<-1:n;
qbot<-1:n;
sm<-matrix(scan(fname),ncol=n , byrow=T);
nobs<-nrow(sm);
for (i in 1:n){
  meanS[i]<-mean(sm[,i]);
  stderrS[i]<-sqrt(var(sm[,i]));
  y<-quantile(sm[,i],c(0.05,0.95));
  qup[i]<-y[2];
  qbot[i]<-y[1];
  }

lines(x, meanS,lty=2);
#lines(x, meanS+stderrS, lty=3);
#lines(x, meanS-stderrS, lty=3);
lines(x, qup,lty=3);
lines(x, qbot,lty=3);
cat("M.C.stderr=", max(stderrS)/sqrt(nobs));
}