"t2future" <- function(da,ini){
# The data matrix is "da".
# Program written by Ruey S. Tsay, April 10, 2008.
#
 nr = dim(da)[1]
 nc = dim(da)[2]
 if (ini < nc) ini= nc+2
 npts= nr-ini
 T2=rep(0,npts)
 for (i in ini:nr){
 dat=da[1:(i-1),]
 cm=matrix(colMeans(dat),1,nc)
 ii=i-1
 o=matrix(rep(1,ii),ii,1)
 dev=dat-kronecker(o,cm)
 dev=as.matrix(dev)      
 s=t(dev)%*%dev/(ii-1)
 si=solve(s)
 xi=as.matrix(da[i,],1,nc)
 t2= (ii/i)*(xi-cm)%*%si%*%t(xi-cm)
 ul1=((ii-1)*nc)/(ii-nc)*qf(0.95,nc,(ii-nc))
 ul2=((ii-1)*nc)/(ii-nc)*qf(0.99,nc,(ii-nc))
 T2[(i-ini+1)]=t2
}
 yl=max(T2,ul1,ul2)
 tdx=c(ini:nr)
 plot(tdx,T2,type='l',ylim=c(0,yl))
 points(tdx,T2)
 abline(h=c(ul1),lty=2)
 abline(h=c(ul2),lty=3)
 title(main='Future Chart: the limits are 95% and 99% quantiles')
 list(Tsq=T2)

}